yalmipR20120309.zip_MATLAB/yalmip_yalmip

function [sol,m,Q,residuals,everything] = solvesos(F,obj,options,params,candidateMonomials) %SOLVESOS Sum of squares decomposition % % [sol,m,B,residual] = solvesos(F,h,options,params,monomials) is used % for finding SOS decompositions of polynomials. % % The coefficients of the polynomials are assumed linear w.r.t a set of % decision variables params and polynomial with respect to a variable x. % % An extension with a nonlinear parameterization in params is possible. % Note though that this gives BMIs or PMIs, solvable (locally) only if % PENBMI is installed, or by specifying 'moment' as solver to try to % solve the nonconvex semidefinite programming problem using a % semidefinite relaxation based on moments. % % The SOS problem can be formulated as % % min h(params) % % subject to F(i) >(=) 0 or F(i) is SOS w.r.t x % % INPUT % F : SET object with SOS constrained polynomials and constraints on variables params % h : scalar SDPVAR object (can be []) % options : options structure obtained from SDPSETTINGS (can be []) % params : SDPVAR object defining parametric variables (can be []) % monomials : SDPVAR object with user-specified monomials for decomposition (can be []) % % OUTPUT % sol : Solution diagnostic from SDP problem % v : Cell with monomials used in decompositions % Q : Cell with Gram matrices, p = v{i}'*Q{i}*v{i}, where p is the ith SOS polynomial in your model. % residuals : Mismatch between p and decompositions. Same values (modulo numerical issue) as checkset(find(is(F,'sos'))) % Warning, these residuals are not computed on matrix sos-of-squares % % EXAMPLE % x = sdpvar(1);solvesos(set(sos(x^4+x^3+1))); % Simple decompositions % x = sdpvar(1);t = sdpvar(1);solvesos(set(sos(x^4+x^3+1-t)),-t); % Lower bound by maximizing t % % NOTES % % Variables not part of params, but part of non-SOS constraints in F % or objective h will automatically be appended to the params list. % % To extract SOS decomposition, use the command SOSD (or compute from use v and Q) % % If the 5th input argument is used, no additional monomial reduction is % performed (Newton, inconstency, congruence). It is thus assumed that % the supplied candidate monomials constitute a sufficient basis. % % The field options.sos can be used to tune the SOS-calculations. See HTML help for details % % sos.model - Kernel (1) or image (2) representation of SOS problem [0|1|2 (0, YALMIP decides)] % sos.newton - Use Newton polytope to reduce size [0|1 (1)] % sos.congruence - Block-diagonalize using congruence classes [0|1|2 (2)] % sos.scale - Scale polynomial [0|1 (1)] % sos.numblkdg - Try to perform a-posteriori block-diagonalization [real (0)] % sos.inconsistent - Remove diagonal-inconsistent monomials [0|1|2 (0)] % sos.clean - Remove monomials with coefficients < clean [real > 0 (1e-4)] % sos.traceobj - Minimize trace of Gram matrix in problems without objective function [0|1 (0)] % sos.extlp - Extract simple translated LP cones when performing dualization [0|1 (1)] % % See also SOS, SOSD, SDPSETTINGS, SOLVEMOMENT, SDPVAR, SDISPLAY %% Time YALMIP yalmip_time = clock; % ************************************************ %% Check #inputs % ************************************************ if nargin<5 candidateMonomials = []; if nargin<4 params = []; if nargin<3 options = sdpsettings; if nargin<2 obj = []; if nargin<1 help solvesos return end end end end end if isa(obj,'double') obj = []; end if isa(F,'constraint') F = set(F); end if isequal(obj,0) obj = []; end if isempty(options) options = sdpsettings; end % Lazy syntax (not official...) if nargin==1 & isa(F,'sdpvar') F = set(sos(F)); end if ~isempty(options) if options.sos.numblkdg [sol,m,Q,residuals,everything] = solvesos_find_blocks(F,obj,options,params,candidateMonomials); return end end % ************************************************************************* %% Extract all SOS constraints and candidate monomials % ************************************************************************* if ~any(is(F,'sos')) error('At-least one constraint should be an SOS constraints!'); end p = []; ranks = []; for i = 1:length(F) if is(F(i),'sos') pi = sdpvar(F(i)); p{end+1} = pi; ranks(end+1) = getsosrank(pi); % Desired rank : Experimental code end end if isempty(candidateMonomials) for i = 1:length(F) candidateMonomials{i}=[]; end elseif isa(candidateMonomials,'sdpvar') cM=candidateMonomials; candidateMonomials={}; for i = 1:length(p) candidateMonomials{i}=cM; end elseif isa(candidateMonomials,'cell') if length(p)~=length(candidateMonomials) error('Dimension mismatch between the candidate monomials and the number of SOS constraints'); end end % ************************************************************************* %% Get the parametric constraints % ************************************************************************* F_original = F; F_parametric = F(find(~is(F,'sos'))); if isempty(F_parametric) F_parametric = set([]); end % Expand removes so called marker variables generated in ismember ParametricBefore = getvariables(F_parametric); % ************************************************************************* %% Expand the parametric constraints % ************************************************************************* if ~isempty(yalmip('extvariables')) [F_parametric,failure] = expandmodel(F_parametric,obj,options); F_parametric = expanded(F_parametric,1); obj = expanded(obj,1); if failure Q{1} = [];m{1} = [];residuals = [];everything = []; sol.yalmiptime = etime(clock,yalmip_time); sol.solvertime = 0; sol.info = yalmiperror(14,'YALMIP'); sol.problem = 14; end end if ~isempty(params) if ~isa(params,'sdpvar') error('Fourth argment should be a SDPVAR variable or empty') end end % ************************************************************************* % Collect all possible parametric variables % ************************************************************************* ParametricVariables = uniquestripped([depends(obj) depends(F_parametric) depends(params) ParametricBefore]); if options.verbose>0; disp('-------------------------------------------------------------------------'); disp('YALMIP SOS module started...'); disp('-------------------------------------------------------------------------'); end % ************************************************************************* %% INITIALIZE SOS-DECOMPOSITIONS SDP CONSTRAINTS % ************************************************************************* F_sos = set([]); % ************************************************************************* %% FIGURE OUT ALL USED PARAMETRIC VARIABLES % ************************************************************************* AllVariables = uniquestripped([depends(obj) depends(F_original) depends(F_parametric)]); ParametricVariables = intersect(ParametricVariables,AllVariables); MonomVariables = setdiff(AllVariables,ParametricVariables); params = recover(ParametricVariables); if isempty(MonomVariables) error('No independent variables? Perhaps you added a constraint set(p(x)) when you meant set(sos(p(x)))'); end if options.verbose>0;disp(['Detected ' num2str(length(Pa