多項式
が与えられたとき，
となるような（半正定とは限らない）
を線形計画問題 (LP) を解くことにより求める．ただし，
である．つぎに，
を満足する 
を半正定値計画問題 (SDP) を解くことにより求める．その結果得られる
は，次式を満足する．
>> sos_decomposition3
SeDuMi had unexplained problems, maybe due to linear dependence?
YALMIP tweaks the problem (adds 1e6 magnitude bounds on all variables) and restarts...
SeDuMi 1.3 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003. %%%% 線形計画問題
Alg = 2: xz-corrector, theta = 0.250, beta = 0.500
Put 5 free variables in a quadratic cone
eqs m = 6, order n = 15, dim = 19, blocks = 2
nnz(A) = 18 + 0, nnz(ADA) = 36, nnz(L) = 21
it : b*y gap delta rate t/tP* t/tD* feas cg cg prec
0 : 1.13E+007 0.000
1 : 0.00E+000 3.13E+001 0.000 0.0000 1.0000 1.0000 1.00 1 1 1.0E+000
2 : 0.00E+000 3.20E-004 0.000 0.0000 1.0000 1.0000 1.00 2 3 3.4E-001
3 : 0.00E+000 2.82E-004 0.000 0.8806 0.9000 0.9000 1.00 5 5 3.0E-001
Run into numerical problems.
iter seconds digits c*x b*y
3 0.6 9.7 3.0033423269e-004 0.0000000000e+000
|Ax-b| = 1.7e-016, [Ay-c]_+ = 6.9E-006, |x|= 1.4e+000, |y|= 4.8e+000
Detailed timing (sec)
Pre IPM Post
2.190E-001 5.200E-001 7.799E-002
Max-norms: ||b||=0, ||c|| = 1000000,
Cholesky |add|=0, |skip| = 1, ||L.L|| = 1.
sol0 =
yalmiptime: 1.2600
solvertime: 0.9760
info: 'No problems detected (SeDuMi-1.3)'
problem: 0
dimacs: [3.8005e-011 0 0 1.6018e-006 -1.8143e-010 -2.9810e-016]
SeDuMi 1.3 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003. %%%% 半正定値計画問題
Alg = 2: xz-corrector, theta = 0.250, beta = 0.500
Put 5 free variables in a quadratic cone
eqs m = 6, order n = 6, dim = 16, blocks = 3
nnz(A) = 12 + 0, nnz(ADA) = 36, nnz(L) = 21
it : b*y gap delta rate t/tP* t/tD* feas cg cg prec
0 : 9.01E+000 0.000
1 : 0.00E+000 1.66E+000 0.000 0.1844 0.9000 0.9000 1.27 1 1 1.6E+000
2 : 0.00E+000 8.88E-002 0.000 0.0534 0.9900 0.9900 1.26 1 1 7.7E-001
3 : 0.00E+000 4.48E-006 0.000 0.0001 1.0000 1.0000 1.00 1 1 2.2E-003
4 : 0.00E+000 1.53E-012 0.000 0.0000 1.0000 1.0000 1.00 1 3 7.1E-010
iter seconds digits c*x b*y
4 0.1 12.7 9.6169903348e-013 0.0000000000e+000
|Ax-b| = 2.4e-013, [Ay-c]_+ = 4.4E-013, |x|= 1.1e+000, |y|= 3.3e+000
Detailed timing (sec)
Pre IPM Post
1.200E-002 1.420E-001 1.600E-002
Max-norms: ||b||=0, ||c|| = 3.999992e+000,
Cholesky |add|=0, |skip| = 1, ||L.L|| = 1.
sol =
yalmiptime: 0.1460
solvertime: 0.1780
info: 'No problems detected (SeDuMi-1.3)'
problem: 0
dimacs: [2.3687e-013 0 0 7.6738e-014 9.6170e-013 1.4390e-012]
pres =
1.9509
>> double(S)
ans =
2.0000 0 1.1955
0 0.6090 1.0000
1.1955 1.0000 4.0000
>> sdisplay(z'*double(S)*z)
2.0000*x1^4+3.0000*x1^2*x2^2+2.0000*x1*x2^3+4.0000*x2^4
>> double(U)
ans =
0.0000 0.0000 -1.4591
0.0000 2.9183 0.0000
-1.4591 0.0000 -0.0000
>> sdisplay(z'*double(U)*z)
-2.9932e-013*x1^2*x2^2
>> sdisplay(clean(z'*double(U)*z,1e-6)) %%% 1e-6 以下の係数は無視
0
>> Q = S + U;
>> double(Q)
ans =
2.0000 0.0000 -0.2636
0.0000 3.5272 1.0000
-0.2636 1.0000 4.0000
>> min(eig(double(Q))) %%% Q の最小固有値
ans =
1.9509
>> sdisplay(z'*double(Q)*z)
2.0000*x1^4+3.0000*x1^2*x2^2+2.0000*x1*x2^3+4.0000*x2^4