Robust ROA calculations for the Van der Pol dynamics
Contents
This file demonstrates the problem setup for robust ROA computation with parametric uncertainty. For more details on all the options see NACodeExplanation.pdf.
Affine uncertainty
Form the vector field
pvar x1 x2; x = [x1;x2]; x1dot = -x2; x2dot = x1+(x1^2-1)*x2;
Nominal system
f = [x1dot; x2dot];
Introduce an uncertain parameter
pvar d1
Specify its range
ini_cell = [-0.1 0.1];
Form the uncertain vector field
f = f + d1*[-x2;0];
Set some options
[roaconstr,opt,sys] = GetRoaOpts(f, x); opt.sim.NumConvTraj = 40;
Solve the problem for an estimate of the robust ROA
outputs = wrapper(sys,ini_cell,roaconstr,opt);
Extract the solution
[V,beta,gamma,p,multip,betaUpper] = extractSol(outputs);
Certified beta
beta
beta = 1.3030
Upper bound on beta
betaUpper
betaUpper = 2.3544
The above code implement the two-step procedure for robust ROA calculations:
- Compute an "optimal" V for the center of the cell - Assess V from the initial step for the vertex systems
We can try to compute an optimal V directly taking all the vertices in to account in constructing V
[fNOM,fVER] = getf(sys,ini_cell);
Re-generate the options, etc.
[roaconstr,opt,sys] = GetRoaOpts(fVER, x); sys.fWithDel = [];
Run the computations
outputs = wrapper(sys,[],roaconstr,opt);
Extract the solution
[V,beta,gamma,p,multip,betaUpper] = extractSol(outputs);
beta
beta = 1.3256
Non-affine uncertainty
pvar x1 x2; x = [x1;x2]; x1dot = -x2; x2dot = x1+(x1^2-1)*x2; f = [x1dot; x2dot];
Introduce the uncertain parameters
pvar d1 d2
Specify its value ini_cell = [d1min d1max d2min d2max ... dmmin dmmax]
ini_cell = [-.2 .2 -.1 .2];
Form the uncertain vector field
f = f + d1*[-x2;0] + d2*[0;-x1] + d1^2*[-x1+x2;0] [roaconstr,opt,sys] = GetRoaOpts(f, x); opt.sim.NumConvTraj = 40;
f = [ -d1^2*x1 + d1^2*x2 - d1*x2 - x2 ] [ x1^2*x2 - d2*x1 + x1 - x2 ]
Solve
outputs = wrapper(sys,ini_cell,roaconstr,opt);
Extract the solution
[V,beta,gamma,p,multip,betaUpper] = extractSol(outputs)
V = 0.47208*x1^2 - 0.19451*x1*x2 + 0.43122*x2^2 beta = 0.7321 gamma = 0.4034 p = x1^2 + x2^2 multip = r1: [8x1 polynomial] r2: [1x1 polynomial] betaUpper = 2.1455