Navigator

sos programmingsos optimizationReferences

sos programming

YALMIP内置了求解sos问题的模块，在sos问题中，需要对多项式$p(x)$进行分解，使得$p(x)=h^T(x)h(x)$，或者分解为$p(x)=v^T(x)Qv(x)$，其中$Q$是半正定矩阵，$v(x)$是向量.

YALMIP also supports linearly and nonlinearly parameterized problems, decomposition of matrix valued polynomials.

手动实现sos的建模求解过程

%% 定义多项式 x = sdpvar(1, 1); y = sdpvar(1, 1); p = (1+x)^4+(1-y)^2; % 分解 v = monolist([x y],degree(p)/2); Q = sdpvar(length(v)); p_sos = v'*Q*v; % 系数匹配 F = [coefficients(p-p_sos,[x y]) == 0, Q >= 0]; % optimize(F) % 设置从对偶形式求解 optimize(F, [], sdpsettings('dualize', 1))

找到多项式下界

%% 找到多项式下界 sdpvar t F = [coefficients((p-t)-p_sos,[x y]) == 0, Q >= 0]; optimize(F,-t)

整合优化流程如下

sdpvar x y P = [1+x^2 -x+y+x^2;-x+y+x^2 2*x^2-2*x*y+y^2]; m = size(P,1); v = monolist([x y],degree(P)/2); Q = sdpvar(length(v)*m); R = kron(eye(m),v)'*Q*kron(eye(m),v)-P; s = coefficients(R(find(triu(R))),[x y]); optimize([Q >= 0, s==0]); sdisplay(clean(kron(eye(m),v)'*value(Q)*kron(eye(m),v),1e-6))

sos optimization

使用命令sos定义（SOS约束），使用命令solvesos求解，sosd得到sos-decomposition结果

x = sdpvar(1,1);y = sdpvar(1,1); p = (1+x)^4 + (1-y)^2; F = sos(p); solvesos(F); h = sosd(F); sdisplay(h)

检测结果，如果分解结果正确，那么$p(x)-h^T(x)h(x)$的值应该接近于0，可以使用命令clean忽略系数很小的向量

clean(p-h'*h, 1e-6)

得到分解式$p(x)=v^T(x)Qv(x)$

x = sdpvar(1,1);y = sdpvar(1,1); p = (1+x)^4 + (1-y)^2; F = sos(p); [sol,v,Q] = solvesos(F); clean(p-v{1}'*Q{1}*v{1},1e-6)

可以使用check指令返回系数最大值或者使用solvesos的第四个参数等效

checkset(F) e = checkset(F(is(F, 'sos'))) [sol,v,Q,res] = solvesos(F); res

References

pre and post processing sum-of-squares program in practice Sum-of-squares programming