Skip to content
Snippets Groups Projects
Select Git revision
  • master default protected
1 result

rbf-approx

  • Clone with SSH
  • Clone with HTTPS
    • Name Last commit Last update
    figures
    functions
    README.md
    demo.jl
    ecc19_paper.pdf
    README.md

    Overview

    This project contains supplemental Julia code for the article:

    M. Thelander Andrén, "Using Radial Basis Functions to Appriximate the LQG-Optimal Event-Based Sampling Policy", Submitted to European Control Conference, 2019

    It demonstrates a numerical method using radial basis functions (RBF) to compute an approximation of the optimal event-based sampling policy for the continuous-time sampled-data LQG problem. The optimal policy is given by the solution to a stationary partial differential equation (PDE) with free boundary. The PDE is given by:

    -V(\tilde{x}^\intercal Q\tilde{x} - J + \tilde{x}^\intercal A^\intercal\nabla V + \frac{1}{2}\Delta V) = 0,\quad \forall \tilde{x},
    -V\geq 0, \quad \tilde{x}^\intercal Q\tilde{x} - J + \tilde{x}^\intercal A^\intercal\nabla V + \frac{1}{2}\Delta V\geq 0, \quad \forall \tilde{x}.

    The solution to this PDE is the relative value function V:\tilde{x}\rightarrow\mathbb{R}, and the optimal sampling policy is to trigger sampling whenever V(\tilde{x})=0 holds. For more details we refer to our article.

    Requirements

    The demo code has been tested for Julia v0.6, which can be downloaded here. The following packages are required:

    For each required package, simply run Pkg.add("package_name") in the Julia REPL to install it.

    Running the Demo Code

    First make sure that you have installed Julia v0.6 and the required packages. Then start a Julia REPL in the directory where the file demo.jl is located. In the Julia REPL, type include("demo.jl") to run the demo.

    Note that it may take quite some time to run the demo the first time, since all the required packages are then being compiled. It takes roughly 1 min to run the entire script with 30^2=900 collocation points.

    Printouts

    While running the demo, the current status is printed to the Julia REPL. First, the script generates the set of collocation points:

    GENERATING COLLOCATION POINTS...DONE!

    If specified, it will then proceed by doing a binary search for the smallest feasible value of the shape parameter c:

    SEARCHING FOR SMALLEST FEASIBLE C...
cᵤ = 160.0, cₗ = 80.0
cᵤ = 120.0, cₗ = 80.0
cᵤ = 100.0, cₗ = 80.0
cᵤ = 90.0, cₗ = 80.0
cᵤ = 90.0, cₗ = 85.0
cᵤ = 87.5, cₗ = 85.0
cᵤ = 86.25, cₗ = 85.0
cᵤ = 86.25, cₗ = 85.625
cᵤ = 85.9375, cₗ = 85.625
cᵤ = 85.9375, cₗ = 85.78125
cᵤ = 85.859375, cₗ = 85.78125
cᵤ = 85.8203125, cₗ = 85.78125
cᵤ = 85.8203125, cₗ = 85.80078125
cᵤ = 85.810546875, cₗ = 85.80078125
FOUND C = 85.810546875

    It then proceeds by generating a linear complementarity problem (LCP) to compute the RBF weights:

    GENERATING LCP...DONE!

    With the LCP defined, it will call the QP-solver in the package OSQP to solve it:

    SOLVING LCP...
-----------------------------------------------------------------
           OSQP v0.4.1  -  Operator Splitting QP Solver
              (c) Bartolomeo Stellato,  Goran Banjac
        University of Oxford  -  Stanford University 2018
-----------------------------------------------------------------
problem:  variables n = 900, constraints m = 1800
          nnz(P) + nnz(A) = 1215861
settings: linear system solver = qdldl,
          eps_abs = 1.0e-09, eps_rel = 1.0e-09,
          eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
          rho = 1.00e-01 (adaptive),
          sigma = 1.00e-06, alpha = 1.60, max_iter = 100000
          check_termination: on (interval 25),
          scaling: on, scaled_termination: off
          warm start: on, polish: off

iter   objective    pri res    dua res    rho        time
   1  -2.1255e+01   1.81e+01   3.25e+00   1.00e-01   1.50e+00s
 200  -1.1622e+00   2.61e-01   3.72e-03   1.00e-01   2.49e+00s
 400  -3.0957e-01   6.63e-02   1.11e-03   1.00e-01   3.39e+00s
 600  -8.8080e-02   1.72e-02   2.97e-04   1.00e-01   4.34e+00s
 800  -2.9303e-02   4.88e-03   8.66e-05   1.00e-01   5.50e+00s
1000  -1.0692e-02   1.55e-03   3.35e-05   1.00e-01   6.49e+00s
1200  -4.3451e-03   5.77e-04   1.46e-05   1.00e-01   7.66e+00s
1400  -1.7629e-03   2.28e-04   6.05e-06   1.00e-01   8.71e+00s
1600  -7.1483e-04   9.20e-05   2.47e-06   1.00e-01   9.82e+00s
1800  -2.8981e-04   3.72e-05   1.00e-06   1.00e-01   1.08e+01s
2000  -1.1749e-04   1.51e-05   4.06e-07   1.00e-01   1.17e+01s
2200  -4.7627e-05   6.12e-06   1.65e-07   1.00e-01   1.28e+01s
2400  -1.9307e-05   2.48e-06   6.68e-08   1.00e-01   1.38e+01s
2600  -7.8269e-06   1.00e-06   2.71e-08   1.00e-01   1.47e+01s
2800  -3.1729e-06   4.07e-07   1.10e-08   1.00e-01   1.57e+01s
3000  -1.2862e-06   1.65e-07   4.45e-09   1.00e-01   1.66e+01s
3200  -5.2142e-07   6.70e-08   1.80e-09   1.00e-01   1.76e+01s
3400  -2.1138e-07   2.71e-08   7.31e-10   1.00e-01   1.85e+01s
3600  -8.5688e-08   1.10e-08   2.96e-10   1.00e-01   1.94e+01s
3800  -3.4736e-08   4.46e-09   1.20e-10   1.00e-01   2.04e+01s
3850  -2.7717e-08   3.56e-09   9.59e-11   1.00e-01   2.06e+01s

status:               solved
number of iterations: 3850
optimal objective:    -0.0000
run time:             2.06e+01s
optimal rho estimate: 8.75e-02

LCP SOLVED! OBJECTIVE VALUE = -2.7716763106866438e-8
COMPUTING RBF-WEIGHTS...DONE!

    Finally, it evaluates the RBF approximation on a grid for plotting. When specified (and A=0), it will also evaluate the true analytic solution on the same grid for comparison.

    EVALUATING RBF-APPROXIMATION ON GRID...DONE!
EVALUATING TRUE SOLUTION ON GRID...DONE!
MAE ON GRID = 0.0008456961574835812
GENERATING PLOT...DONE!

    Plots

    If specified, the demo script will output plots of the resulting approximation. Below is an example: Example plot output

    When validation is specified, the plot will also feature a comparison to the true solution: Validation plot output

    Modifying the Demo Code

    We encourage you to open the file demo.jl and modify it to try out different precision of the approximation and different systems. Whenever you have changed the parameters to your liking, simply run include("demo.jl") in the Julia REPL again to see the new result.

    Changing Script Behaviour

    You can specify if the script should run a validation, if it should search for the smallest feasible shape parameter, and if it should generate plots by changing the following boolean variables in demo.jl:

    validation = true
plotting = true
search_c = true

    Changing the System

    You can change the system matrix A, the (diagonal) quadratic cost weight matrix Q and the total cost J by changing the following variables in demo.jl:

    A = [0. 15.;15. 0.]
Q = diagm([1., 1.])
J = 1.

    Note that A and Q should have elements of type Float64, while J should be of type Float64. This is most easily achieved by adding a decimal point after each number, as done above.

    Changing the RBF Approximation

    You can change the number of collocation points in each dimension, the value of the shape parameter and the boundary values of the grid for the collocation points by changing the following variables in demo.jl:

    xlimits = [[-1.5, 1.5], [-1.5, 1.5]]
Ni = [30, 30]
c = 80.0

    Note that xlimits should be of the type Array{Array{Float64,1},1}, i.e. a vector of vectors, where each element contains the lower and upper bounds of the grid for each dimension. Ni, which specifies the number of collocation points in each dimension, should only contain elements of type Int64. The shape parameter c should be of type Float64. Note that if search_c=true has been specified, then the value of c will only be used as an initial value in the binary search for the smallest feasible shape parameter.