2009. 3. 11. 16:10
@GSMC/서용덕: Convex Optimization
SeDuMi - Home
optimization package over symmetric cones
SeDuMi = Self-Dual-Minimization/ Optimization over self-dual homogeneous cones
: MatLab Toolbox developed by Jos F. Sturm
1) SeDuMi를 다운로드 받기 위해 홈페이지에 회원 등록
2) 메뉴의 Download에서 적당한 버전 다운로드
3) Matlab을 실행시키고, 메뉴에서 File > Set Path > Add Folder 클릭, SeDuMi가 들어 있는 폴더 지정 후 저장
4) command window에서 'help sedumi' 입력하여 toolbox 추가 확인
optimization package over symmetric cones
SeDuMi = Self-Dual-Minimization/ Optimization over self-dual homogeneous cones
: MatLab Toolbox developed by Jos F. Sturm
1) SeDuMi를 다운로드 받기 위해 홈페이지에 회원 등록
2) 메뉴의 Download에서 적당한 버전 다운로드
3) Matlab을 실행시키고, 메뉴에서 File > Set Path > Add Folder 클릭, SeDuMi가 들어 있는 폴더 지정 후 저장
4) command window에서 'help sedumi' 입력하여 toolbox 추가 확인
>> help sedumi
[x,y,info] = sedumi(A,b,c,K,pars)
SEDUMI Self-Dual-Minimization/ Optimization over self-dual homogeneous
cones.
> X = SEDUMI(A,b,c) yields an optimal solution to the linear program
MINIMIZE c'*x SUCH THAT A*x = b, x >= 0
x is a vector of decision variables.
If size(A,2)==length(b), then it solves the linear program
MINIMIZE c'*x SUCH THAT A'*x = b, x >= 0
> [X,Y,INFO] = SEDUMI(A,b,c) also yields a vector of dual multipliers Y,
and a structure INFO, with the fields INFO.pinf, INFO.dinf and
INFO.numerr.
(1) INFO.pinf=INFO.dinf=0: x is an optimal solution (as above)
and y certifies optimality, viz.\ b'*y = c'*x and c - A'*y >= 0.
Stated otherwise, y is an optimal solution to
MAXIMIZE b'*y SUCH THAT c-A'*y >= 0.
If size(A,2)==length(b), then y solves the linear program
MAXIMIZE b'*y SUCH THAT c-A*y >= 0.
(2) INFO.pinf=1: there cannot be x>=0 with A*x=b, and this is certified
by y, viz. b'*y > 0 and A'*y <= 0. Thus y is a Farkas solution.
(3) INFO.dinf=1: there cannot be y such that c-A'*y >= 0, and this is
certified by x, viz. c'*x <0, A*x = 0, x >= 0. Thus x is a Farkas
solution.
(I) INFO.numerr = 0: desired accuracy achieved (see PARS.eps).
(II) INFO.numerr = 1: numerical problems warning. Results are accurate
merely to the level of PARS.bigeps.
(III) INFO.numerr = 2: complete failure due to numerical problems.
INFO.feasratio is the final value of the feasibility indicator. This
indicator converges to 1 for problems with a complementary solution, and
to -1 for strongly infeasible problems. If feasratio in somewhere in
between, the problem may be nasty (e.g. the optimum is not attained),
if the problem is NOT purely linear (see below). Otherwise, the reason
must lie in numerical problems: try to rescale the problem.
> [X,Y,INFO] = SEDUMI(A,b,0) or SEDUMI(A,b) solves the feasibility problem
FIND x>=0 such that A*x = b
> [X,Y,INFO] = SEDUMI(A,0,c) or SEDUMI(A,c) solves the feasibility problem
FIND y such that A'*y <= c
> [X,Y,INFO] = SEDUMI(A,b,c,K) instead of the constraint "x>=0", this
restricts x to a self-dual homogeneous cone that you describe in the
structure K. Up to 5 fields can be used, called K.f, K.l, K.q, K.r and
K.s, for Free, Linear, Quadratic, Rotated quadratic and Semi-definite.
In addition, there are fields K.xcomplex, K.scomplex and K.ycomplex
for complex-variables.
(1) K.f is the number of FREE, i.e. UNRESTRICTED primal components.
The dual components are restricted to be zero. E.g. if
K.f = 2 then x(1:2) is unrestricted, and z(1:2)=0.
These are ALWAYS the first components in x.
(2) K.l is the number of NONNEGATIVE components. E.g. if K.f=2, K.l=8
then x(3:10) >=0.
(3) K.q lists the dimensions of LORENTZ (quadratic, second-order cone)
constraints. E.g. if K.l=10 and K.q = [3 7] then
x(11) >= norm(x(12:13)),
x(14) >= norm(x(15:20)).
These components ALWAYS immediately follow the K.l nonnegative ones.
If the entries in A and/or c are COMPLEX, then the x-components in
"norm(x(#,#))" take complex-values, whenever that is beneficial.
Use K.ycomplex to impose constraints on the imaginary part of A*x.
(4) K.r lists the dimensions of Rotated LORENTZ
constraints. E.g. if K.l=10, K.q = [3 7] and K.r = [4 6], then
2*x(21)x(22) >= norm(x(23:24))^2,
2*x(25)x(26) >= norm(x(27:30))^2.
These components ALWAYS immediately follow the K.q ones.
Just as for the K.q-variables, the variables in "norm(x(#,#))" are
allowed to be complex, if you provide complex data. Use K.ycomplex
to impose constraints on the imaginary part of A*x.
(5) K.s lists the dimensions of POSITIVE SEMI-DEFINITE (PSD) constraints
E.g. if K.l=10, K.q = [3 7] and K.s = [4 3], then
mat( x(21:36),4 ) is PSD,
mat( x(37:45),3 ) is PSD.
These components are ALWAYS the last entries in x.
(a) K.xcomplex lists the components in f,l,q,r blocks that are allowed
to have nonzero imaginary part in the primal. For f,l blocks, these
(b) K.scomplex lists the PSD blocks that are Hermitian rather than
real symmetric.
(c) Use K.ycomplex to impose constraints on the imaginary part of A*x.
The dual multipliers y have analogous meaning as in the "x>=0" case,
except that instead of "c-A'*y>=0" resp. "-A'*y>=0", one should read that
c-A'*y resp. -A'*y are in the cone that is described by K.l, K.q and K.s.
In the above example, if z = c-A'*y and mat(z(21:36),4) is not symmetric/
Hermitian, then positive semi-definiteness reflects the symmetric/
Hermitian parts, i.e. Z + Z' is PSD.
If the model contains COMPLEX data, then you may provide a list
K.ycomplex, with the following meaning:
y(i) is complex if ismember(i,K.ycomplex)
y(i) is real otherwise
The equality constraints in the primal are then as follows:
A(i,:)*x = b(i) if imag(b(i)) ~= 0 or ismember(i,K.ycomplex)
real(A(i,:)*x) = b(i) otherwise.
Thus, equality constraints on both the real and imaginary part
of A(i,:)*x should be listed in the field K.ycomplex.
You may use EIGK(x,K) and EIGK(c-A'*y,K) to check that x and c-A'*y
are in the cone K.
> [X,Y,INFO] = SEDUMI(A,b,c,K,pars) allows you to override the default
parameter settings, using fields in the structure `pars'.
(1) pars.fid By default, fid=1. If fid=0, then SeDuMi runs quietly,
i.e. no screen output. In general, output is written to a device or
file whose handle is fid. Use fopen to assign a fid to a file.
(2) pars.alg By default, alg=2. If alg=0, then a first-order wide
region algorithm is used, not recommended. If alg=1, then SeDuMi uses
the centering-predictor-corrector algorithm with v-linearization.
If alg=2, then xz-linearization is used in the corrector, similar
to Mehrotra's algorithm. The wide-region centering-predictor-corrector
algorithm was proposed in Chapter 7 of
J.F. Sturm, Primal-Dual Interior Point Approach to Semidefinite Pro-
gramming, TIR 156, Thesis Publishers Amsterdam, 1997.
(3) pars.theta, pars.beta By default, theta=0.25 and beta=0.5. These
are the wide region and neighborhood parameters. Valid choices are
0 < theta <= 1 and 0 < beta < 1. Setting theta=1 restricts the iterates
to follow the central path in an N_2(beta)-neighborhood.
(4) pars.stepdif, pars.w. By default, stepdif = 2 and w=[1 1].
This implements an adaptive heuristic to control ste-differentiation.
You can enable primal/dual step length differentiation by setting stepdif=1 or 0.
If so, it weights the rel. primal, dual and gap residuals as
w(1):w(2):1 in order to find the optimal step differentiation.
(5) pars.eps The desired accuracy. Setting pars.eps=0 lets SeDuMi run
as long as it can make progress. By default eps=1e-8.
(6) pars.bigeps In case the desired accuracy pars.eps cannot be achieved,
the solution is tagged as info.numerr=1 if it is accurate to pars.bigeps,
otherwise it yields info.numerr=2.
(7) pars.maxiter Maximum number of iterations, before termination.
(8) pars.stopat SeDuMi enters debugging mode at the iterations specified in this vector.
(9) pars.cg Various parameters for controling the Preconditioned conjugate
gradient method (CG), which is only used if results from Cholesky are inaccurate.
(a) cg.maxiter Maximum number of CG-iterates (per solve). Theoretically needed
is |add|+2*|skip|, the number of added and skipped pivots in Cholesky.
(Default 49.)
(b) cg.restol Terminates if residual is a "cg.restol" fraction of duality gap.
Should be smaller than 1 in order to make progress (default 5E-3).
(c) cg.refine Number of refinement loops that are allowed. The maximum number
of actual CG-steps will thus be 1+(1+cg.refine)*cg.maxiter. (default 1)
(d) cg.stagtol Terminates if relative function progress less than stagtol (5E-14).
(e) cg.qprec Stores cg-iterates in quadruple precision if qprec=1 (default 0).
(10) pars.chol Various parameters for controling the Cholesky solve.
Subfields of the structure pars.chol are:
(a) chol.canceltol: Rel. tolerance for detecting cancelation during Cholesky (1E-12)
(b) chol.maxu: Adds to diagonal if max(abs(L(:,j))) > chol.maxu otherwise (5E5).
(c) chol.abstol: Skips pivots falling below abstol (1e-20).
(d) chol.maxuden: pivots in dense-column factorization so that these factors
satisfy max(abs(Lk)) <= maxuden (default 5E2).
(11) pars.vplot If this field is 1, then SeDuMi produces a fancy
v-plot, for research purposes. Default: vplot = 0.
(12) pars.errors If this field is 1 then SeDuMi outputs some error
measures as defined in the Seventh DIMACS Challenge. For more details
see the User Guide.
Bug reports can be submitted at http://sedumi.mcmaster.ca.
See also mat, vec, cellK, eyeK, eigK
[x,y,info] = sedumi(A,b,c,K,pars)
SEDUMI Self-Dual-Minimization/ Optimization over self-dual homogeneous
cones.
> X = SEDUMI(A,b,c) yields an optimal solution to the linear program
MINIMIZE c'*x SUCH THAT A*x = b, x >= 0
x is a vector of decision variables.
If size(A,2)==length(b), then it solves the linear program
MINIMIZE c'*x SUCH THAT A'*x = b, x >= 0
> [X,Y,INFO] = SEDUMI(A,b,c) also yields a vector of dual multipliers Y,
and a structure INFO, with the fields INFO.pinf, INFO.dinf and
INFO.numerr.
(1) INFO.pinf=INFO.dinf=0: x is an optimal solution (as above)
and y certifies optimality, viz.\ b'*y = c'*x and c - A'*y >= 0.
Stated otherwise, y is an optimal solution to
MAXIMIZE b'*y SUCH THAT c-A'*y >= 0.
If size(A,2)==length(b), then y solves the linear program
MAXIMIZE b'*y SUCH THAT c-A*y >= 0.
(2) INFO.pinf=1: there cannot be x>=0 with A*x=b, and this is certified
by y, viz. b'*y > 0 and A'*y <= 0. Thus y is a Farkas solution.
(3) INFO.dinf=1: there cannot be y such that c-A'*y >= 0, and this is
certified by x, viz. c'*x <0, A*x = 0, x >= 0. Thus x is a Farkas
solution.
(I) INFO.numerr = 0: desired accuracy achieved (see PARS.eps).
(II) INFO.numerr = 1: numerical problems warning. Results are accurate
merely to the level of PARS.bigeps.
(III) INFO.numerr = 2: complete failure due to numerical problems.
INFO.feasratio is the final value of the feasibility indicator. This
indicator converges to 1 for problems with a complementary solution, and
to -1 for strongly infeasible problems. If feasratio in somewhere in
between, the problem may be nasty (e.g. the optimum is not attained),
if the problem is NOT purely linear (see below). Otherwise, the reason
must lie in numerical problems: try to rescale the problem.
> [X,Y,INFO] = SEDUMI(A,b,0) or SEDUMI(A,b) solves the feasibility problem
FIND x>=0 such that A*x = b
> [X,Y,INFO] = SEDUMI(A,0,c) or SEDUMI(A,c) solves the feasibility problem
FIND y such that A'*y <= c
> [X,Y,INFO] = SEDUMI(A,b,c,K) instead of the constraint "x>=0", this
restricts x to a self-dual homogeneous cone that you describe in the
structure K. Up to 5 fields can be used, called K.f, K.l, K.q, K.r and
K.s, for Free, Linear, Quadratic, Rotated quadratic and Semi-definite.
In addition, there are fields K.xcomplex, K.scomplex and K.ycomplex
for complex-variables.
(1) K.f is the number of FREE, i.e. UNRESTRICTED primal components.
The dual components are restricted to be zero. E.g. if
K.f = 2 then x(1:2) is unrestricted, and z(1:2)=0.
These are ALWAYS the first components in x.
(2) K.l is the number of NONNEGATIVE components. E.g. if K.f=2, K.l=8
then x(3:10) >=0.
(3) K.q lists the dimensions of LORENTZ (quadratic, second-order cone)
constraints. E.g. if K.l=10 and K.q = [3 7] then
x(11) >= norm(x(12:13)),
x(14) >= norm(x(15:20)).
These components ALWAYS immediately follow the K.l nonnegative ones.
If the entries in A and/or c are COMPLEX, then the x-components in
"norm(x(#,#))" take complex-values, whenever that is beneficial.
Use K.ycomplex to impose constraints on the imaginary part of A*x.
(4) K.r lists the dimensions of Rotated LORENTZ
constraints. E.g. if K.l=10, K.q = [3 7] and K.r = [4 6], then
2*x(21)x(22) >= norm(x(23:24))^2,
2*x(25)x(26) >= norm(x(27:30))^2.
These components ALWAYS immediately follow the K.q ones.
Just as for the K.q-variables, the variables in "norm(x(#,#))" are
allowed to be complex, if you provide complex data. Use K.ycomplex
to impose constraints on the imaginary part of A*x.
(5) K.s lists the dimensions of POSITIVE SEMI-DEFINITE (PSD) constraints
E.g. if K.l=10, K.q = [3 7] and K.s = [4 3], then
mat( x(21:36),4 ) is PSD,
mat( x(37:45),3 ) is PSD.
These components are ALWAYS the last entries in x.
(a) K.xcomplex lists the components in f,l,q,r blocks that are allowed
to have nonzero imaginary part in the primal. For f,l blocks, these
(b) K.scomplex lists the PSD blocks that are Hermitian rather than
real symmetric.
(c) Use K.ycomplex to impose constraints on the imaginary part of A*x.
The dual multipliers y have analogous meaning as in the "x>=0" case,
except that instead of "c-A'*y>=0" resp. "-A'*y>=0", one should read that
c-A'*y resp. -A'*y are in the cone that is described by K.l, K.q and K.s.
In the above example, if z = c-A'*y and mat(z(21:36),4) is not symmetric/
Hermitian, then positive semi-definiteness reflects the symmetric/
Hermitian parts, i.e. Z + Z' is PSD.
If the model contains COMPLEX data, then you may provide a list
K.ycomplex, with the following meaning:
y(i) is complex if ismember(i,K.ycomplex)
y(i) is real otherwise
The equality constraints in the primal are then as follows:
A(i,:)*x = b(i) if imag(b(i)) ~= 0 or ismember(i,K.ycomplex)
real(A(i,:)*x) = b(i) otherwise.
Thus, equality constraints on both the real and imaginary part
of A(i,:)*x should be listed in the field K.ycomplex.
You may use EIGK(x,K) and EIGK(c-A'*y,K) to check that x and c-A'*y
are in the cone K.
> [X,Y,INFO] = SEDUMI(A,b,c,K,pars) allows you to override the default
parameter settings, using fields in the structure `pars'.
(1) pars.fid By default, fid=1. If fid=0, then SeDuMi runs quietly,
i.e. no screen output. In general, output is written to a device or
file whose handle is fid. Use fopen to assign a fid to a file.
(2) pars.alg By default, alg=2. If alg=0, then a first-order wide
region algorithm is used, not recommended. If alg=1, then SeDuMi uses
the centering-predictor-corrector algorithm with v-linearization.
If alg=2, then xz-linearization is used in the corrector, similar
to Mehrotra's algorithm. The wide-region centering-predictor-corrector
algorithm was proposed in Chapter 7 of
J.F. Sturm, Primal-Dual Interior Point Approach to Semidefinite Pro-
gramming, TIR 156, Thesis Publishers Amsterdam, 1997.
(3) pars.theta, pars.beta By default, theta=0.25 and beta=0.5. These
are the wide region and neighborhood parameters. Valid choices are
0 < theta <= 1 and 0 < beta < 1. Setting theta=1 restricts the iterates
to follow the central path in an N_2(beta)-neighborhood.
(4) pars.stepdif, pars.w. By default, stepdif = 2 and w=[1 1].
This implements an adaptive heuristic to control ste-differentiation.
You can enable primal/dual step length differentiation by setting stepdif=1 or 0.
If so, it weights the rel. primal, dual and gap residuals as
w(1):w(2):1 in order to find the optimal step differentiation.
(5) pars.eps The desired accuracy. Setting pars.eps=0 lets SeDuMi run
as long as it can make progress. By default eps=1e-8.
(6) pars.bigeps In case the desired accuracy pars.eps cannot be achieved,
the solution is tagged as info.numerr=1 if it is accurate to pars.bigeps,
otherwise it yields info.numerr=2.
(7) pars.maxiter Maximum number of iterations, before termination.
(8) pars.stopat SeDuMi enters debugging mode at the iterations specified in this vector.
(9) pars.cg Various parameters for controling the Preconditioned conjugate
gradient method (CG), which is only used if results from Cholesky are inaccurate.
(a) cg.maxiter Maximum number of CG-iterates (per solve). Theoretically needed
is |add|+2*|skip|, the number of added and skipped pivots in Cholesky.
(Default 49.)
(b) cg.restol Terminates if residual is a "cg.restol" fraction of duality gap.
Should be smaller than 1 in order to make progress (default 5E-3).
(c) cg.refine Number of refinement loops that are allowed. The maximum number
of actual CG-steps will thus be 1+(1+cg.refine)*cg.maxiter. (default 1)
(d) cg.stagtol Terminates if relative function progress less than stagtol (5E-14).
(e) cg.qprec Stores cg-iterates in quadruple precision if qprec=1 (default 0).
(10) pars.chol Various parameters for controling the Cholesky solve.
Subfields of the structure pars.chol are:
(a) chol.canceltol: Rel. tolerance for detecting cancelation during Cholesky (1E-12)
(b) chol.maxu: Adds to diagonal if max(abs(L(:,j))) > chol.maxu otherwise (5E5).
(c) chol.abstol: Skips pivots falling below abstol (1e-20).
(d) chol.maxuden: pivots in dense-column factorization so that these factors
satisfy max(abs(Lk)) <= maxuden (default 5E2).
(11) pars.vplot If this field is 1, then SeDuMi produces a fancy
v-plot, for research purposes. Default: vplot = 0.
(12) pars.errors If this field is 1 then SeDuMi outputs some error
measures as defined in the Seventh DIMACS Challenge. For more details
see the User Guide.
Bug reports can be submitted at http://sedumi.mcmaster.ca.
See also mat, vec, cellK, eyeK, eigK
'@GSMC > 서용덕: Convex Optimization' 카테고리의 다른 글
Homework #1 (0) | 2009.03.11 |
---|---|
[Kolman & Beck] Chapter 1 Introduction to Linear Programming (0) | 2009.03.11 |
Kolman & Beck <Elementary linear programming with applications> 2nd ed. (0) | 2009.03.09 |
[Boyd & Vandenberghe] Chapter 4 Convex Optimization Problems (0) | 2009.03.04 |
[Boyd & Vandenberghe] Chapter 3 Convex Functions (0) | 2009.02.23 |