Avoid use of 'nonlinearreplace'

[Misha Bekker]

Hello,

I've written a script which determines the robust stability margin of an affinely modeled uncertain LPV system ('YALMIP_LPV.m', you also need 'allcomb.m' to run it). The math is based on some research I've found and as you can tell from the code, the amount of LMI's to be solved increases exponentially depending on system dimension and number of uncertain parameters. That and the fact that I want to use this test in future for synthesizing robust controllers makes me want to minimize the time it takes to solve. Which leads me to my question, on my machine the YALMIP diagnostics return the following:

sol = 

  struct with fields:

    yalmiptime: 18.3967

    solvertime: 1.2185

          info: 'Successfully solved (bisection)'

       problem: 0

The matlab profiler shows that most of the yalmiptime ( 16.391 s ) is spent in the function 'nonlinearreplace', which starts with the following comment:

function Z = nonlinearreplace(X,Y,W)

% Very slow, but easily coded. Hopefully not used much

Is there a way to pose the problem in such a way that this function is avoided? 

Otherwise, I've seen that the 'optimizer' function can be used to pre-compile an optimization function. Do you recommend just using that and taking the one-time hit of pre-compiling the problem?

Thanks in advance!

[Johan Löfberg]

bisection is just syntactic sugar over an optimizer object, so that would be just as slow

I'll look into it

[Johan Löfberg]

BTW, X_v isn't symmetric, which sounds odd in a LUapunov context

>> X_v

Bilinear matrix variable 16x16 (full, real, 49 variables)

Coeffiecient range: 0.5 to 8

Not too surprising as you create some basis

% P = P0 + P1 + P2 ...

P   = sdpvar(n_d*(n_p+1),n_d);

and the matrices there are not symmetric. You must use something along the lines of

% P = P0 + P1 + P2 ...

P   = sdpvar(n_d,n_d*(n_p+1));

as that will create symmetric slices, and then you sum up accordingly along third dimension

[Johan Löfberg]

...and dimensions in general appear weird.

L = A_v'*X_v + X_v*A_v or (A_v'*X_v + X_v*A_v +X_v) looks like the standard matrix you'd want negative semidefinite for symmetric positive definite X_v

However, you create the 64x16 matrix L, which you then constrain etc

I think you're getting lost in the kronecker algebras, A tip would be to start everything with for-loops and get that to work, and look at the objects you create to ensure they have expected properties (dimensions, symmetry), and then when you have that working, compactify using kron tricks etc if you must

[Johan Löfberg]

and btw, do you know that the problem is quasi-convex in d_p, and thus making bisection viable.

the problem is quadratic in d_p, i.e not what you typically would have in a bisection affair

>> degree(sdpvar(C(2)),d_p)

ans =

     2

[Johan Löfberg]

Noticed now that the bottle-neck in replace really isn't related to the used optimizer framework, but something that happens when yalmip normalizes the problem to a maximization problem which bisection initially was coded for. Hence, a manual implementation using optimizer would not waste this time. I will try to fix this

[Johan Löfberg]

download the develop version and it is circumvented¨

but then again, your first goal is to fix your code so you actually solve as SDP and not a nonlinear elementwise problem as you do now

[Misha Bekker]

Thanks for the help!

I did not notice that I was losing symmetry and you're right, the dimensionality is a bit of a mess. I created a weird blockdiagonal constraint which repeats along the 1st dimension, don't know why I did it like that. I did make a for-loop version at first, but did not check consistency of the problem between the 2 versions, only the results. Will try to fix these things.