PENNON User's Guide (Version 0.9)
Contents
1. 6 Examples 6 1 NLP SDP example Consider the following simple NLP SDP example in matrix variable X S 3 min 2 Hy 8 subject to TrX 1 X gt 0 where 22 1 1 0 H gt 11 19 1 15 0 Al We will treat the matrix variable as a sparse tri diagonal matrix 6 1 1 AMPL interface nlpsdp mod var 1l 5 default OF param mie oer minimize O0J sumir amn ao en aa 2 subject to IRES Sees See m ee data param h ie gag 2 elol 3 a 2 ile nlpsdp sdp Nr of sdp blocks 1 t Ni of non lin sde blocks 0 t Nee of lin sde blocks il Block sizes 3 lower eigenvalue bounds o upper eigenvalue bounds 1 0E38 Constraint types 0 nonzeroes per block 5 nonzero structure per block 1 0 0 21 1 0 1 dl in 1 al 1 2 al 2 2 6 1 2 MATLAB interface Eem function EX a x Reena Sale ee ia eles e E x reshape x length x Le x D B diem fonction maz 1nd val af z By gt 2 ae A aa a Ge x reshape x length x 1 fing length x ind 1002 val 2 x xh Ea fonction nnz Cow Col vall hit ss nnz length x row l nnz col lnnzi yal 2 x ones nnz g m function aqua x gx X x t5 1 0 dg m function nnz ind vall dq i x nnz 3 ind iS val dele iis hg m function nhz row Col val h2. the value of the primal variable at the computed optimum status exit information see below iresults a 4x1 matrix with elements as described below dresults a 5x1 matrix with elements as described below IRESULTS meaning iresults 1 number of outer iterations iresults 2 number of inner iterations iresults 3 number of linesearch steps iresults 4 ellapsed time in seconds DRESULTS meaning dresults 1 primal objective dresults 2 relative precision at Zopt dresults 3 feasibility at topt dresults 4 complementary slackness at Zopt dresults 5 gradient of augmented lagrangian at Zopt STATUS meaning status 0 converged optimal solution status 1 converged suboptimal solution gradient large status 2 converged solution primal infeasible status 3 aborted no progress problem may be primal infeasible status 4 aborted primal unbounded or initial multipliers too small status 5 aborted iteration limit exceeded status 6 aborted line search failure status 7 aborted aborted cholesky solver failed status 8 aborted wrong parameters status 9 aborted resource limit status 10 aborted internal error please contact PENOPT Gbr contact penopt com status 11 aborted error in user s memory allocation status 12 aborted error in user supplied routines 20
3. PENNON User s Guide Version 0 9 Michal Ko vara Michael Sting www penopt com June 11 2008 Contents 1 Installation 2 1h APA 4 0 5 Shes Bek OS eee Soe ee ee Bes 2 1 2 COP AS Wire ty A Ga Ate teen ot wh te 3 2 The problem 3 3 The algorithm 3 4 AMPL interface 5 4i Matix variables MAMP Lenne Scat a e eS Ad 5 42 Preparing AMPL inputdata ae e eee e Se GS oe 6 42 1 Thesdpilen e moe a Oe ae a 6 43 Redundant constraints lt e p e a la E 8 Aas Running PENNON s e oee SOME a e a Se a1 a a eee as 8 ARS Pro srant opugns n ay aeaaaee wish ee ati 2S Sy eee nse eS Pe 8 5 MATLAB interface 12 5 1 Calling PENNONM from MATLAB s sece arer eere eaa n 12 Sab Usen provided MncHons eos ee 12 52 Theperanput structure mi MATLAB o tose to a aes 15 59 The PENNONM tunctioncall o s rea eee a 19 6 Examples 21 EE NEPS DP example oin sont a r dak ang a a Haniel 21 Gosh AMEE interlace e manm da es aaa a o 21 612 MATLAB mterace aaora a a e dae 22 6 2 Correlation matrix with the constrained condition number 25 6 3 Truss topology optimization a e e eaen ae o e e a a e E 28 Copyright 2002 2008 Kocvara amp Stingl PENOPT GbR 1 Installation 1 1 Unpacking UNIX versions The distribution is packed in file pennon tar gz Put this file in an arbitrary directory After uncompressing the file to pennon tar by command gunzip pennon tar gz the files are extracted by command tar xvf pennon tar Win32 version The
4. e Y S Y S are the matrix variables k symmetric matrices of dimen sions p X Pi Pk X Pk e we denote Y Yj Yp f gi and h are C functions from R x SP x x S to R e A and A are the lower and upper bounds respectively on the eigenvalues of Y ld 3 The algorithm To simplify the presentation of the algorithm we only consider inequality constraints For the treatment of the equality constraints see 2 The algorithm is based on a choice of penalty barrier functions y4 R R that penalize the inequality constraints and Pp S S penalizing the matrix inequalities These functions satisfy a number of properties see 2 that guarantee that for any p gt 0 i 1 myg we have gla 0 pei 6 L 00 and Z30 lt Op Z X00 ZeES This means that for any p gt 0 problem NLP SDP has the same solution as the following augmented problem min f x Y subject to pypy g 2 Y p lt 0 i 1 m p A J Yi lt 0 i 1 k NLP SDP4 p Y AL 0 The Lagrangian of NLP SDP can be viewed as a generalized augmented Lagrangian of NLP SDP F a u U U p f Y wipieg aia Y pi 1 k k Un At Y Us PW AD 1 4 1 here u R s and U U are Lagrange multipliers associated with the standard and the matrix inequality constraints respectively The algorithm combines ideas of the exterior penalty
5. output to file 0 no 1 yes check density of the Hessian 0 automatic 1 dense ignore initial solutions O do not ignore 1 do ignore equilibrate linear system matrix 0 no 1 yes maximum number of outer iterations restriction factor u of multiplier update maximum number of iterations in the inner loop Newton or Trust region method 9 1 0E 1 1 0e0 1 5 0e 2 1 0e0 zil 100 0 5 100 nwtmode nwtstopcrit objno ordering ia penalty penup penupmode pes pinit paved precision PECAR pmecond SDE file Cmi Urnie Mindo WIinitnS umin usebarrier usesdpbarrier PENNON AMPL options cont linear system solver O Cholesky method 1 conjugate gradient method 2 conjugate gradient method with approximate Hessian calculation 1 conjugate gradient method to dual system stopping criterium for the inner loop Movie a L VL a 2 lt a luf uf llo 2 VERA la a VEA le objective number in the AMPL mod file ordering for MA57 output level 0 silent mode 1 brief output 2 full output penalty function 0 logarithmic barrier quadratic penalty 1 reciprocal barrier quadratic penalty penalty update penalty update is performed 0 adaptively _ L after each outer iteration minimal penalty initial penalty pivot tolerance for MA27 57 required final precision required final precision of the KKT conditions precondition
6. 0 X 1 1000 1 8333 1 1000 0 0 1 1000 2 0333 The corresponding MATLAB files for this example are included in the directory mat lab nlpsdp 24 6 2 Correlation matrix with the constrained condition number We consider the problem of finding a nearest correlation matrix n mn 2 Xe H 9 subject to gt li he alle X 0 This problem is based on a practical application see 4 Assume that the original correla tion matrix is 1 044 0 20 0 81 0 46 Wa 1 0 87 0 38 0 81 Horig 0 20 0 87 1 0 17 0 65 osi 038 ir 1 i7 0 46 0 81 0 65 0 37 1 When we solve problem 9 with H Horig the solution will be as expected the original matrix Horig We now add a new asset class that means we add one row and column to the original matrix The new data is based on a different frequency than the original part of the matrix which means that the new matrix is no longer positive definite 1 0 44 0 20 0 81 0 46 0 05 0 44 1 0 87 0 38 0 81 0 58 H 0 20 87 1 0 17 0 65 0 56 ext 0 81 0 38 0 17 1 03 0 15 0 46 0 81 0 65 0 37 1 0 08 0 05 0 58 0 56 0 15 0 08 1 Let us find the nearest correlation matrix to Hex by solving 9 We obtain the following result for the presentation of results we will use matlab output in short precision x 1 0000 0 4420 O 2000 0 8096 0 4585 0 0513 0 4420 1 0000 0 8704 OS a Ont IIS 0 3549 0 2000 0 8
7. For the sample problem NLP1 the structure can be as follows ini ivy UE SO m De pen nvars n pen nlin 1 pen neonstT 1 pen nsdp 1 pen blks 3 pen nnz_gradient n pen nnz_hessian n pen dby Infinity lt onesi n 1 pen ubv Infinity xones n 1 pen lbc 0 pen ubc 0 pen lbmv 0 pen ubmv Infinity pen mtype 0 pen mnzs 5 pen mrow 07071172 pen meol Oe 1 I2 2 pen xinit 1L 0 1071 Pen E A pen my_f_gradient df pen my_f_hessian hf pen my_g g 17 pen my_g_gradient dg pen my_g_hessian hg pen ioptions pen doptions A sample implementation is included in the files nlp m f m df m hf m g m dg m 100 100 2 0001001 00 0 1 0 1 Ql ES TORO OO a Ob 5 ES 5 OI lei i QE 6 ELA 1 0e 7 0 05 10 10 1 017 and hg m in directory mat lab IOPTIONS name value meaning default ioptions 1 maxit maximum numbers of outer iterations 100 ioptions 2 nwtiters maximum number of iterations in inner loop 100 ioptions 3 outlev output level 2 0 no output 1 only options are displayed 2 brief output 3 full output ioptions 4 hessianmode check density of hessian 0 0 automatic 1 dense ioptions 5 autoscale automatic scaling 0 0 on 1 off ioptions 6 convex convex problem 0 0 no 1 yes HOSE wens 9 eqltymode the way to
8. Program options The options are summarized in Table 1 Recommendations e Whenever you know that the problem is convex use convex 1 e When you have problems with convergence of the algorithm try to nonconvex problems if an initial guess of the solution is available increase decrease uinit e g uinit 10000 switch to Trust Region algorithm by ncmode 1 decrease alpha e g alpha le 3 change stopping criterion for inner loop by setting nwt stopcrt 1 PENNON AMPL options decrease pinit e g pinit 0 01 This should be particulary helpful for option meaning default alpha alphaupd autoini autoscale egtolmia cgtolup cmaxnzs Convex eqltymode filerep hessianmode ignoreinit KKTscale maxit mu nwtiters stopping parameter a for the Newton Trust region method in the inner loop update of a automatic initialization of multipliers 0 off 1 nonlinear nonconvex mode 2 lp qp mode automatic scaling 0 on 1 off minimum tolerance for the conjugate gradient algorithm update of tolerance for the conjugate gradient algorithm tuning parameter for Hessian assembling in nwtmode put gt 0 to switch it on convex problem 0 generally nonconvex 1 convex initialization of equality multipliers 0 two inequalities symmetric 1 two inequalities unsymmetric 2 augmented lagrangian 3 direct 4 direct only nonlinear equalities
9. Ut U p enri 1 AF Y ds T 1f Y where e is by default 107 parameter precision Optionally the user can choose a stopping critrion based on the KKT error lt 4 AMPL interface AMPL is a comfortable modelling language for optimization problems For a description of AMPL we refer to 3 or www amp1 com 4 1 Matrix variables in AMPL AMPL does not support matrix variables However the format of our problem 2 allows us to use them e within an AMPL script file lt name gt mod matrix variables are treated as vectors using the function svec S R defined by 411 Q12 aim a22 ES 02m FP svec a11 Q12 022 01m 42m Gram sym Umm Example Assume we have a matrix variable X S T T2 Ta XS m m Ta Tg XG and a constraint 0 0 1 Tr XA 3 wthA 0 1 0 1 0 0 The matrix variable is treated as a vector svec X 1 2 6 and the above constraint is equivalent to the following constraint gts 2r4 3 The corresponding lines in the lt name gt mod file would be var IO subject to C RS 22 4 537 e The order of the variables should be 1 nonlinear matrix the matrix or its elements are invovled in a nonlinear expres sion 2 standard real variables 3 linear matrix variables Example Consider a problem with constraints XPX T yiAi y242 0 in variables X S and y R The second constraint should be re written
10. and interior barrier methods with the Augmented Lagrangian method Algorithm 3 1 Let xt Y and alu be given Let p gt 0 i L My and P gt Q For k 1 2 repeat till a stopping criterium is reached i Find x 1 and Y such that V F x Y u ye U p lt K A a a 14 Nm CS orf Y UF i 1 k C2 ie E ya e iii lt 1 2 0 PA The approximate unconstrained minimization in Step i is performed by the Newton method with line search or by a variant of the Trust Region method for details see 2 The min imization is optionally stopped when either X k MBA ET E or k Ver att Y uk 0 0 p lla lt a lug zos gs a pF Ile or k k Aa A ee oy ere a a lle with optional parameter a by default 1071 The multipliers calculated in Step ii are restricted in order to satisfy T 1 AA u i p with some positive 4 lt 1 by default y 0 3 Similar estimates are satisfied for the matrix multipliers The update of the penalty parameter p in Step iii is performed in two different ways Either the penalty vector is updated by some constant factor dependent on the initial penalty parameter 7 or the penalty vector is updated only in case the progress of the method becomes too slow In either case the penalty update is stopped if Peps by default 1079 is reached Algorithm 3 1 is stopped when both of the inequalities holds et YE Fat Y uk
11. distribution is packed in file pennon zip Put this file in an arbitrary directory and extract the files by PKZIP In both cases the directory PENNONO 9 containing the following files and subdirectories will be created LICENSE bin lib matlab matlab cond matlab corr matlab fmo matlab nlpsdp matlab truss file containing the PENNON license agreement directory containing the files pennon0 9 exe the binary executable with AMPL interface nipsdp mod a model file of a sample problem in AMPL format nipsdp sdp an sdp file of a sample problem in AMPL format nipsdp nl an nl file created by AMPL from nlpsdp mod cond mod a model file of a sample problem in AMPL format cond sdp an sdp file of a sample problem in AMPL format cond nl an nl file created by AMPL from cond mod corr mod a model file of a sample problem in AMPL format corr sdp an sdp file of a sample problem in AMPL format corr nl an nl file created by AMPL from corr mod fmo mod a model file of a sample problem in AMPL format fmo dat a data file of a sample problem in AMPL format fmo sdp an sdp file of a sample problem in AMPL format fmo nl an nl file created by AMPL from fmo mod directory containing the PENNON libraries directory containing the files pennonm c the MATLAB interface file penoutm c MATLAB version of penout c make_pennonm m M file containing MEX link command nip m f m df m hf m g m dg m hg m M
12. files defining a sample problem in PEN format nipsdp bfgs cond directories containing examples from the last section directory containing the files for example cond directory containing the files for example corr directory containing the files for example fmo directory containing the files for example nlpsdp directory containing the files for example truss 1 2 Compilation Requirements For successful compilation and linkage depending on the operating system and the pro gram to be created the following software packages have to be installed UNIX versions e MATLAB version 5 0 or later including MEX compiler package and gcc compiler package MATLAB dynamic link library pennonm x Win32 version e MATLAB version 5 0 or later including MEX compiler package and VISUAL C version 6 0 or later MATLAB dynamic link library pennonm x To build a MATLAB dynamic link library pennonm x Start MATLAB go to directory mat 1ab and invoke link command by make_pennonm In case the user wants to use his her own LAPACK BLAS or ATLAS implementations the M file in directory mat lab has to be modified appropriately 2 The problem We solve optimization problems with nonlinear objective subject to nonlinear inequality and equality constraints and semidefinite bound constraints 2 Yi ae ESP Ha Y subjectto g 1 Y lt 0 li Ley 0 O NLP SDP A aad alleen ae e x R is the vector variable
13. treat equality constraints 3 0 as two inequalities unsymmetric initialization 1 as two inequalities symmetric initialization 2 handled by standard augmented lagrangian 3 direct handling all equalities ioptions 8 ignoreinit ignore initial solutions 0 0 do not ignore 1 ignore ioptions 9 cholmode cholesky system mode 0 0 Solve directly 1 Solve augmented system Z Split into two systems ioptions 10 nwtstopcrit stopping criterion for the inner loop 2 0 IVL aF l2 lt a 1 VL DIla lt a gt luf ula 2 WL ya lt a VE 2 la ioptions 11 penalty penalty function 0 0 logarithmic barrier quadratic penalty 1 reciprocal barrier quadratic penalty Tope Lons 12 nwtmode mode of solving the Newton system 0 0 cholesky standard 1 cg with exact hessian 2 cg with appr hessian 3 cg with user provided Hessian vector routine 18 ioptions 13 prec preconditioner for the cg method 0 no precond 1 diagonal precond 2 bfgs precond 3 appr inverse precond 4 sgs precond ioptions 14 cmaxnzs tuning parameter for Hessian assembling in nwt 1 mode 1 3 1 off gt 0 on options 15 autoini automatic initialization of multipliers 0 off 1 nonlinear nonconvex mode 2 Ip qp mode ioptions 16 penup penalty parameter update is performed 0 adaptively 1 after each outer iteration ioptions 17 usebarrier box constraint mode 0 no sp
14. 704 1 0000 0 1699 0 6497 0 3397 0 8096 043714 O 1699 1 0000 0 3766 0 1445 0 4505 Cielo 08197 U 2000 1 OOOO 0 0608 0 0513 O 5949 D oot O 1445 0 0608 1 0000 with eigenvalues eigen 0 0000 o isi Oe 2120 UMSS levls2 oly ay As we can see one eigenvalue of the nearest correlation matrix is zero This is highly undesirable from the application point of view To avoid this we can add lower and upper bounds on the matrix variable i e constraints Ne eae Vi This would be reflected in the following lines in the nlpsdp m matlab code from page 22 25 pen lbmv lambda_min pen ubmv lambda_max However the application requires a more general approach when we only want to bound the condition number of the nearest correlation matrix This can be guaranteed by introducing a pair of new variables y z R and adding the following set of constraints to 9 Xz 10 AS 11 y lt kz 12 where is the required condition number Notice that the above constraints do not fit into our required NLP SDP problem structure see page 3 The standard way to transorm the condition constraints into our standard structure is to introduce two new slack matrix variables say P and Q and replace 10 and 11 by X 2I P 0 X yl Q 0 P gt 0 Q 0 We may not like the additional matrix variables as they increase the problem size consid erably There are two ways how avoid using them First we may use the automatic sl
15. a tr gammal pen xinit 0 1 ones m 1 zeros nnn 1 pen my ft E truss pen my Gradient di truss 7 pen my_f_hessian hf_truss pen my q g truss pen my_g_gradient dg_truss pen my_g_hessian hg_truss pen topetrons 100 100 2 0001 0 0 1 0 0 0 1 0 121 pen doptioas 1 0E 2 1050 1 0n 0 1 0E 2 5 06 11 SS 06 10h 6 T 0E 12 Ae 0 05 1 0 1 0 1 0 w1 w2 pennonm pen We further show the files defining the objective function and constraints E ECUSS function fx f Esuessix global par m par m fx sum x 1 m GUervss function lox g cruss ir X global par A m par m n par n nl par nl ifa lt nit yenila gx 0 for j 1 m gx gx A j itl x 9 end 29 gx ox alan else gx so em end References 1 W Achtziger and M Ko vara Structural Topology Optimization with Eigenvalues SIAM J Optimization 18 4 1129 1164 2007 2 M Ko vara and M Stingl PENNON a code for convex nonlinear and semidefinite programming Optimization Methods and Software 8 3 317 333 2003 3 R Fourer D M Gay and B W Kernighan AMPL a modelling language for mathematical programming Scientific Press 1993 4 R Werner and K Sch ttle Calibration or corellation matrices SDP or not SDP Submitted 30
16. ack removal option mtype 4 The above con straints then should be reformulated as X zI P 0 X yl Q 0 P 0 Q 0 see Remark 5 2 on page 15 All data and m files are prepared with the matrices present in the formulation then we set pen mtype 0 4 4 end PENNON will automatically remove the slacks from the formulation such that the actual calculations will only use variables y z and X The corersponding m files can be found in the directory mat Llab cond_slack Second we can rewrite constraints 10 11 as fad aa 13 assuming that y kz and using the transormation of the variable X X X Now of course we also have to change the other constraints and the objective function the new problem of finding the nearest correlation matrix with a boudned condition number reads as follows min 3 es H 14 ZX i j 1 subject to 2 X4 1 i 1 n Pax 26 26 The new problem now has the NLP SDP problem structure required by PENNLP see page 3 When solving the problem by PENNLP with 10 we get the solution after 11 outer and 37 inner iterations The solution is z 0 2866 and Xtilde 3 4886 12 3170 05 77380 2 4761 49102 0 2456 Sal SEO 3 4986 DAS SOs 20926 ATS O 7 T80 AS 3 4886 0 5392 129269 e bas 2A ToN il 1005 AMS 3 4886 SUSO MESAS 4202 20926 LS BENS Se SS 3 4886 0 2 010 SUIZOS AS IES ASS 0 2008 3 4886 After the back substitution X 1X we get the nearest correlation m
17. atrix X 1 0000 SNS O 2230 0 7098 SAAT 0 0 704 O oT 1 0000 076990 0 2155 073999 0 4218 D 2230 On 6930 1 0000 O 1546 OV S023 0 291 4 0 7098 O 3155 SAS 1 0000 03857 0 1294 0 4272 0 5998 0 5923 0 3857 1 0000 0 0576 0 0704 0 4218 04914 0 1294 O 0376 1 0000 with eigenvalues eigenvals 0 2866 0 2866 EAS OMGER TGO 2 8664 and the condition number equal to 10 Below we show the corresponding AMPL files cond mod and cond sdp These can be found in directory bin cond mod param hil 21 set Indi within fl 21 var ade 24 default 0 var z 3 minimises Obj sumi am D21 zsz a subject to bit in eee Se exe li LO oe pJ Z Z lt 10000 Dae in indit zx il 1 data param h 1 SOO e eS OC A SOA a a OO Ol amp Usa 2 O17 O 00 MORAG 12 eg 13 O65 TA O 3 AO e005 U oo cy SO 0 0e 2 N0 set indi 1 3 6 10 15 21 cond sdp 27 Nr of sdp blocks lt Ne of non lin sdp blocks dE Nr of lin sdp blocks 0 Block sizes 6 lower eigenvalue bounds 1 upper eigenvalue bounds HON Constraint types 0 nonzeroes per block 21 nonzero structure per block Important As mentioned earlier AMPL may reorder the variables It is thus important to add redundant nonlinear constraints that include all variables In this case we added constraints pfi im ls ar Ss 00005 DI z z lt 10000 The corresponding MATLAB files for this exampl
18. e are included in the directory mat lab cond 6 3 Truss topology optimization The single load truss topology optimization problem can be formulated as a linear semidef inite program see e g 1 mint as i l subject to 14 f A e 7 1 y saad Uae TE Here A S i 1 m are given symmetric matrices and f R y R given data To cast the problem into our canonical NLP SDP form we introduce a slack variable Sest ia Dt i 1 subject to g_ 2 t4 0 _ Pr q TA 0 0 Ss 0 ti 0 E 1 N Below we will present part of the MATLAB interface Again the complete m files can be found in directory mat lab truss We have two variables in the problem the matrix S only involved in linear expressions and the vector t In MATLAB interface vector variables should preceede linear matrix variables We thus introduce a joint variable x Re n 2 2 ses The relevant part after the definition of the data matrices of the file truss m reads as nan mlr2 ni 2 nl is the size of matrices AI Infinity 1 0538 pen nvars m nnn pen nlin 0 Pen aconster ANN pen nsdp 1 pen blks nl 1 pen lbmy 0 pen ubmv Infinity pen mtype 0 pen nnz_gradient m 1 pen nnz_hessian m pen lbv 0 0 x ones m 1 Infinity ones nnn 1 pen ubv Infinity ones m 1 Infinity xones nnn 1 pen ibe z res mltlieni Z 1 tt gammal pen ubc zeros inmi
19. e elements the numbers are zero based mcol element column indices per block not required if integer varies all matrices are dense for each sparse block with s non zero elements give s numbers of column indices of the elements the numbers are zero based nae initial guess for the solution double nvars my E actual name of the file f m char my f gradient actual name of the file df m char my_f_hessian actual name of the file hf m char my_g actual name of the file g m char my_g_gradient actual name of the file dg m char my_g_hessian actual name of the file hg m char ioptions integer valued options integer 18 doptions real valued options double 14 16 e order the matrix variables as X P Q or X Q P inlbmv ubmv mtype etc If only of of the marices P Q is to be removed it has to be the last one For instance if only Q is to be removed as a slack variable the order must be X P Q and the array mtype will be mtype 0 0 4 e constraints 4 and 5 must be reformulated such that P and Q have positive signs i e X 2I P 0 6 X yI Q 0 7 PO Q lt 0 e constraint 7 must further be reformulated such that the slack Q is positive definite not negative definite as above The final version of the constraints is then X zI P 0 X yl Q 0 Pie 0 Q 0 Remark 5 3 For a dense block there is no need to give data mrow mcol Hence if all the blocks are dense arrays nsdp mrow mcol can be omitted
20. ecial treatment 1 use strict barrier function 2 use strict modified barrier function ioptions 18 dercheck derivative check 0 no derivative check 1 check gradients 2 check hessians DOPTIONS name value meaning default doptions 1 precision required final precision 1 0e 7 doptions 2 uinit initial multilplier scaling factor 1 0 doptions 3 pinit initial penalty 1 0 doptions 4 alpha stopping parameter alpha for the Newton Trust re 0 01 gion method in the inner loop doptions 5 mu restriction factor of multiplier update 0 5 doptions 6 penup penalty update 0 1 doptions 7 peps minimal penalty 1 0e 8 doptions 8 umin minimal multiplier 1 0e 12 doptions 9 preckkt precision of the KKT conditions 1 0e 1 doptions 10 cgtolmin minimum tolerance of the conjugate gradient algo 5 0e 2 rithm doptions 11 cgtolup update of tolerance of the conjugate gradient algo 1 0e0 rithm doptions 12 uinitbox initial multiplier box constraints 1 0e0 doptions 13 uinitnc initial multiplier nonlinear constraints 1 0e0 doptions 14 uinitsdp initial multiplier sdp constraints 1 0e0 5 3 The PENNONM function call In MATLAB PENNONM is called with the following arguments f x u status iresults dresults pennonm pen where pen the input structure described in the next section 19 f the value of the objective function at the computed optimum X the value of the dual variable at the computed optimum u
21. er type 0 no preconditioner 1 diagonal 2 L BFGS 3 approximate inverse 4 symmetric Gauss Seidel name of the SDP input file timing destination 0 no 1 stdout 2 stderr 3 both initial multiplier scaling factor initial multiplier scaling factor for box constraints initial multiplier scaling factor for nonlinear constraints minimal multiplier Use mod barrier approach for boxes 0 no 1 barrier 2 strict modified barrier Use barrier approach for SDP variables 0 no 1 yes 10 1 0E 7 1 0E0 1 0E 2 1 0E 7 1 0E 5 1 0 1 0 1 0E 10 version wantsol PENNON AMPL options cont report PENNON version 0 yes a tepon without AMPL Sum of O do not write sol file write sol file print primal variable print dual variable do not print solution message on 11 5 MATLAB interface 5 1 Calling PENNONM from MATLAB 5 1 1 User provided functions The user is required to provide six MATLAB functions The names of the functions can be chosen by the user here we use the following names f evaluates the objective function df evaluates the gradient of objective function hf evaluates the Hessian of objective function g evaluates the constraints dg evaluates the gradient of constraints hg evaluates the Hessian of constraints Similarly as in the AMPL interface matrix variables are treated as vectors using the functio
22. gs bigs 4 x naz 0 row 0 col 0 val 0 nlpsdp m 22 n 5 infinity 1 038 pen nvars n pen nlin 1 pen nconstr 1 pen nsdp 1 pen blks 3 pen nnz_gradient n pen nnz_hessian n pen lby Inftinity sones n 1 pen ubv Infinity x ones n 1 pen be 1015 pen tube 015 pen lbmv 0 pen bmy Infinity pen mtype 0 pen mnzs 5 pen mrow 07 0 1 172 pen meol Oy le T 22l pen indir LP Or 0 Li Pen my ES pen my_f_gradient df pen my_f hessian hf pen my_g g pen my_g_gradient dg pen my_g_hessian hg pen ioptions pen doptions wl w2 pennonm pen Below is an output of the command gt gt niipsdp Variables 5 NEU Constraints 0 B constralntes 1 Number of bounds O 1 0 u 0 Dense Problem Number of nonlinear constraints Number of linear constraints Number of Variables bounded below bounded above Number of Matrix variables degrees of freedom bounded below bounded above Number of Nonlinear Equalities Number of Nonlinear Inequalities 23 ineq A OA Oe OE LOS 2 1 OE0 I ESO a Ob e SUE Ela IESO 1 0E 12 DS 0 05 O O SISSI G ee Ca Number of Linear Equalities 1 Number of Linear Inequalities 0 KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK RS PENNON 0 9 ES Max Min Cin Mult 1 000000 7 1000000 maximal penalty 1 000000 penalt
23. infeasible at some iterations 1 lower strict the lower bounds will always be strictly feasible during the iretative process will be treated by an inner barrier 2 upper strict same as above now for the upper bounds 3 lower upper strict both constraints will always remain strictly feasible 4 slack the block is a slack variable non zeroes per block k numbers if the block is dense input the number p 1 p 2 where p is the block size if the block is known to be sparse give the actual number of its nonzero ele ments all other elements of this block will be ignored by the code nonzero structure per block not required if all matrices are dense for each sparse block with s non zero elements give s lines one for each non zero element where each line looks as follows block_number row_index column_index The numbering of the row column indices is zero based For dense block there is no need to give these data Hence is all the blocks are dense this part of the sdp file will be void Below is an example of the sdp file with three blocks the first one dense and the other two sparse diagonal EFE aE TE FE HH aE TE FE EE EE EH EE EEE HEHEHE HEH sdp example 1 PENNON 0 9 EFE HE TE FE HH a TE FE EE EEE EERE EEE EEE HEH Ne of sdp blocks 3 t Nr of non lin sdp blocks 0 Nra of lin sdp blocks 3 Block sizes 3 2 2 lower eigenvalue bounds Oe Ole e up
24. input nnz variable storing number of non zeros of Vg output ind nnz x 1 matrix storing non zero structure of V g output val nnz x 1 matrix storing non zero values of Vg output Description 1 Compute Vg xit 2 Assign non zero structure to ind and the corresponding values to val Note e Vector x should not be modified by the user e Linear constraints should be specified after nonlinear constraints e Non zero indices should be zero based e Non zero structure should be constant hg function nnz Low Col val hgs bigs i x mag 105 row 0 col 0 val 0 14 Arguments a variable storing constraint number input x vector of length n storing current iterate x input nnz variable storing number of non zeros of V g output row nnz x 1 matrix non zero row indices of V g output col nnz x 1 matrix non zero column indices of Vg output val nnz x 1 matrix storing non zero values of Vg output Description 1 Compute V7gi it 2 Assign non zero structure to row and col and the corresponding values to val Note Vector x should not be modified by the user Linear constraints should be specified after nonlinear constraints Non zero indices should be zero based e Non zero structure should be constant Values should be assigned to lower triangular part of V f xj only 5 2 The pen input structure in MATLAB The user must create a MATLAB structure array with field
25. n svec S R 2 defined by 411 12 Aim 0492 a2m T svec a a11 0412 022 Alm 42m Gmm sym Umm Important The order of the variables should be 1 standard real variables 2 matrix variables The parameter nvar contains the number of all variables For instance if we have two matrix mariables one dense of size 2 i e 3 unknowns and one sparse with 7 nonzeros in the triangular part and 4 standard real variables then nvar 3 7 4 14 The specification of the user defined functions will be explained using the sample problem NLP SDP2 E fonction fx f x in 2 24 SL le 1 9 Si 2 1 oye x reshape x length x x xn e Arguments x vector of length n storing current iterate x input fx variable storing f x output Description e Compute f 2 and store it in fx Note e Vector x should not be modified by the user 12 at function mz ana vall df x he 1222 li 19 2 1165 x reshape x length x 1 nnz length x ind 1 nnz val 2 xh Arguments X vector of length n storing current iterate x input nnz variable storing number of non zeros of V f output ind nnz x 1 matrix storing non zero structure of V f output val nnz x 1 matrix storing non zero values of V f output Description 1 Compute V f ait 2 Assign non zero structure to ind and the corresponding values to val Note e Vecto
26. per eigenvalue bounds 1 0E28 IES Cle Constraint types O standard l lower strict 2 upper strict 3 lower upper strict 4 slack not implemented yet nonzeroes per block 6 2 2 nonzero structure per block Wu N NY SoS So 1 O 4 3 Redundant constraints Important AMPL may reorder the variables which would lead to a collaps of our treat ment of the matrix variables To avoid the reordering all variables should be involved in nonlinear constraints If this is not the case it is important to add redundant nonlinear constraints that include all variables We recommend to add constraints of the type rea i in le ns xlil 211 1000000 This kind of redundant constraints will not change the sparsity structure of the Hessian and will not influence the efficiency of the code 4 4 Running PENNON PENNON is called in the standard AMPL style i e either by a sequence like gt model example mod gt data example dat gt option presolve 0 gt gt options solver pennon options pennon_options sdpfile example sdp outlev 2 for instance solve within the AMPL environment or from the command line by emp OS aio emo e decano cate gt pennon example nl sdpfile example sdp outlev 2 It is necessary to surpress AMPL preprocessing either by the command option presolve 0 within AMPL or using option P when running AMPL from the command line Sample files are included in directory bin 4 5
27. r x should not be modified by the user e Non zero indices should be zero based e Non zero structure should be constant E function nnz row col vall mt ss nnz lemoth 2 row noz col I nnz val 2 x ones nnz Arguments x vector of length n storing current iterate x input nnz variable storing number of non zeros of V f output row nnz x 1 matrix non zero row indices of V f output col nnz x 1 matrix nnz non zero column indices of V f output val nnz x 1 matrix nnz storing non zero values of V f output Description 1 Compute V f xi 2 Assign non zero structure to row and col and the corresponding values to val Note e Vector x should not be modified by the user e Non zero indices should be zero based e Non zero structure should be constant Values should be assigned to lower triangular part of V f xj only 13 g function lox iy xX q a aloe Arguments a variable storing constraint number input x vector of length n storing current iterate x input gx variable storing f x output Description Compute g x and store it in gx Note e Vector x should not be modified by the user e Linear constraints should be specified after nonlinear constraints dg function maz ind vall dq nnz 3 ind 325 yal pizi Arguments i variable storing constraint number input x vector of length n storing current iterate xj
28. s described below in Table 2 The structure also specifies the names of the user defined functions Remark 5 1 Arrays 1bv and ubv define lower and upper bounds on variables in case of matrix variables these are bounds on the elements of the matrices Arrays 1bmv and ubmv on the other hand define spectral bounds on the matrices i e bounds on the smallest and largest eigenvalues Remark 5 2 When using the automatic removal of slack variables array mt ype it is important to follow these points e the slack variables should be tha last ones in the list of variables e the equality constraint that defines the slack variable should be formulated in such a way that the slack variable has a positive sign e the constraints on the slack variables are always of the type S 0 Example Consider a problem in variables y z R X SP with constraints X eal 2 Xxyl 3 To transorm the constraints into our standard structure we need to introduce two slack matrix variables P and Q and replace 2 and 3 by X 2I P 0 4 X yI Q 0 5 P gt 0 o o If we want to use the automatic slack removal option mt ype 4 we have to 15 Table 2 The pen structure name description type length nvars number of variables integer number nconstr number of constraints including linear integer number da number of linear constraints integer n
29. umber nsdp number of matrix variables integer number blks row dimensions of matrix variables integer nsdp nnz_gradient maximal number of non zero entries in user spec integer number ified gradients nnz_hessian maximal number of non zero entries in user spec integer number ified Hessians lbv lower bounds on variables double nvars ubv upper bounds on variables double nvars be lower bounds on constraints double nconstr ube upper bounds on constraints double nconstr lbmv lower bounds on matrix variables min eigenval double nsdp ues ubmv upper bounds on matrix variables max eigenval double nsdp ues mtype type of matrix constraints double nsdp 0 standard constraints may be infeasible at some iterations 1 lower strict the lower bounds will always be strictly feasible during the iterative process 2 upper strict same as above now for the upper bounds 3 lower upper strict both constraints will always remain strictly feasible 4 slack the block is a slack variable mnzs non zeros per block if the block is dense input integer nsdp the number p 1 p 2 where p is the block size if the block is known to be sparse give the actual number of its nonzero elements all other elements of this block will be ignored by the code mrow element row indices per block not required if all integer varies matrices are dense for each sparse block with s non zero elements give s numbers of row indices of th
30. using a slack variable S S XPX I y1A1 y2A2 S S 0 In this example X is a nonlinear matrix variable y a standard real variable and S a linear matrix variable So the definition of these variables in the lt name gt mod file must follow the order var x 1 6 vectorized upper triangle of X var yil 2 var s l 3 wectorized upper triangle of e data needed to identify the matrix variables their number and size are included in a file lt name gt sdp that is read by the PENNON 4 2 Preparing AMPL input data 4 2 1 The sdp file File lt name gt sdp includes the following data The order of the data is compulsory 1 Number of blocks Give the total number of blocks in the block matrix variable say k 2 Number of nonlinear blocks The number of blocks involved in nonlinear constraints 3 Number of linear blocks The number of blocks involved in linear constraints 4 Block sizes k numbers the sizes of the single blocks If a block is a p x p matrix its size is p 5 Lower eigenvalue bounds k numbers the lower bounds on the eigenvalues in the single blocks numbers A in problem NLP SDP 6 Upper eigenvalue bounds k numbers the upper bounds on the eigenvalues in the single blocks numbers A in problem NLP SDP Constraint types k numbers identifying the treatment of the lower upper bound constraints in the al gorithm O standard treatment of the constraints the constraints may be
31. y update 0 562341 KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK Ae obj UG Cy MaE feas pmin Nwt Fact KKEKKKKKKKKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK AAA AAA AAA ARA AAA AAA A Ol 5 88000e 000 4 7e 000 1 2e 001 3 0e 000 DEF ODO 0 0 2 S21 se 002 2 8e 000 1 7e 004 1 7e 001 1 0e 000 3 AS SO Ge 0072 6 0e 001 7 3e 003 2 1e 002 y 0e 00 4 1 il 3 1 64866e 002 Il lS 001 1 8e 003 1 0e 002 2 5620017 il dl A 1 36399e 002 1 9e 002 1 2e 004 1 8e 003 1 3e 001 i il ar 335825002 2 5e 003 8 1e 007 1 8e 004 6 3e 002 dl il 6 1 33345e 002 2 1e 004 4 8e 010 8 9e 006 3 le 002 1 il 7 1 33334e 002 1 3e 005 5 9e 014 2 2e 007 6e 002 1 il KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK Objective 1 33333627497024 04e 002 Relative Precision 12990484 US 763813 7E 7003 Gradient Augm Lagrangian 5 9044538720747749E 014 Complementary Slackness 1 2950484187658137E 005 Feasibility 2 200 121623 1608824058 007 Feasibility LMI 0 0000000000000000E 000 Outer Iterations Inner Iterations 3 Linesearch steps 9 Start time sun Jan 16 19 05 23 2008 End time sun Jan 06 19 05 24 2008 Process time O bh O min 1 sec KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK wl Or OM sa w2 2 1333 ale 1000 1 8323 1 1000 220383 gt gt The optimal matrix is thus 2 1333 1 1000 0
Download Pdf Manuals
Related Search