YANE-Framework Tutorial 1.1.0

examples/minprog/src/demo_nlpsolver.cpp

00001 /***************************************************************************
00002  *                                                                         *
00003  * Copyright (C) 2011 by www.nonlinearmpc.com                              *
00004  *                                                                         *
00005  * Authors:                                                                *
00006  *  Michael Bodenschatz <michael.bodenschatz@uni-bayreuth.de>             *
00007  *  Juergen Pannek <juergen.pannek@googlemail.com>                         *
00008  *                                                                         *
00009  ***************************************************************************/
00010 
00011 #include <iostream>
00012 #include <iomanip>
00013 #include <cstdio>
00014 #include <cstdlib>
00015 #include <string>
00016 
00017 #include <yane.h>
00018 
00019 #include <yaneipopt.h>
00020 #include <sqpnagc3.h>
00021 #include <sqpfortran.h>
00022 
00023 using namespace std;
00024 using namespace yane::Utils;
00025 
00026 using namespace yane::YaneIpopt;
00027 using namespace yane::SqpNagC3;
00028 using namespace yane::SqpFortran;
00029 
00030 void objectivFunction ( const double *x, double *fx, void *params );
00031 void restrictionFunction ( const double *x, double *restx, void *params );
00032 void boxConstraints ( int dimension, double * lb, double * ub );
00033 
00034 int main ( void )
00035 {
00036         Exception::enableDebugMessage ( true );
00037         RTClock timer;
00038         yane::MinProg::NLP *sqpf, *sqpnagc, *ipopt, **minimizer;
00039 
00040         string * names = new string [ 3 ];
00041 
00042         double duration;
00043         int n = 4; // Dimension of the problem
00044         int nrest = 2; // Number of constrains (both equality and inequality)
00045         int nrest_eq = 1; // Number of equality constraints
00046         double optval = 0.0; // Optimal solution
00047         double *x0 = new double [ n ]; // Startvalue
00048         double *x = new double [ n ]; // Optimization variable
00049         double *lb = new double [ n ]; // Lower Bound
00050         double *ub = new double [ n ]; // Upper Bound
00051         double *lambda = new double [ nrest ]; // Array for Lagrangian multipliers of the nonlinear constraints
00052         double *rest_jac = new double [ n * nrest ]; // Array for the Jacobian matrix of the constraints
00053 
00054         // Initialization of the startvalue and the boxed constraints
00055         x0 [ 0 ] = 1.0;
00056         x0 [ 1 ] = 5.0;
00057         x0 [ 2 ] = 5.0;
00058         x0 [ 3 ] = 1.0;
00059         boxConstraints ( n, lb, ub );
00060 
00061         // Creation and initalization of the optimization routine objects
00062         names [ 0 ] = "SqpFortran";
00063         sqpf = new SqpFortran ( n, objectivFunction, restrictionFunction, nrest,
00064                 nrest_eq );
00065         sqpf->setAccuracy ( 1.0e-8 );
00066         sqpf->setMaxIterations ( 1000 );
00067         sqpf->setAbortTime ( 0.002 );
00068 
00069         names [ 1 ] = "SqpNagC";
00070         sqpnagc = new SqpNagC ( n, objectivFunction, restrictionFunction, nrest,
00071                 nrest_eq );
00072         sqpnagc->setAccuracy ( 1.0e-8 );
00073         sqpnagc->setMaxIterations ( 1000 );
00074 
00075         names [ 2 ] = "YaneIpopt";
00076         ipopt = new YaneIpopt ( n, objectivFunction, restrictionFunction, nrest,
00077                 nrest_eq );
00078         ipopt->setAccuracy ( 1.0e-8 );
00079         ipopt->setMaxIterations ( 1000 );
00080 
00081         // Auxiliary MinProg-Array
00082         minimizer = new yane::MinProg::NLP* [ 3 ];
00083         minimizer [ 0 ] = sqpf;
00084         minimizer [ 1 ] = sqpnagc;
00085         minimizer [ 2 ] = ipopt;
00086 
00087         for ( int k = 0; k < 3; k++ )
00088         {
00089                 for ( int i = 0; i < n; i++ )
00090                         x [ i ] = x0 [ i ];
00091                 timer.reset ( );
00092 
00093                 try
00094                 {
00095                         minimizer [ k ]->calcMin ( x, & optval, lb, ub );
00096                 }
00097                 catch ( yane::MinProg::SolverWarning e )
00098                 {
00099                         cout << e.what ( ) << endl;
00100                 }
00101                 duration = timer.elapsedSeconds ( ) * 1000.0;
00102 
00103                 try
00104                 {
00105                         minimizer [ k ]->getLagrangeParameters ( lambda );
00106                 }
00107                 catch ( yane::MinProg::MinProgException e )
00108                 {
00109                         cout << e.what ( ) << endl;
00110                 }
00111 
00112                 try
00113                 {
00114                         minimizer [ k ]->getRestrictionJacobi ( rest_jac );
00115                 }
00116                 catch ( yane::Utils::Exception e )
00117                 {
00118                         cout << e.what ( ) << endl;
00119                 }
00120 
00121                 objectivFunction ( x, & optval, ( void* ) 0 ); // calculation of the optimal value
00122 
00123                 // Print solution and other information
00124                 cout << names [ k ] << endl;
00125                 cout << "Calculation Time: " << fixed << setw ( 6 )
00126                         << setprecision ( 4 ) << duration << " ms" << endl;
00127                 cout << "Optimalvalue: " << setw ( 12 ) << setprecision ( 7 ) << optval
00128                         << endl;
00129                 cout << "Optimalsolution:" << endl;
00130                 for ( int i = 0; i < n; i++ )
00131                         cout << setw ( 12 ) << setprecision ( 7 ) << x [ i ];
00132                 cout << endl;
00133 
00134                 cout << "Lagrangian Multipliers of the constraints: " << endl;
00135                 for ( int i = 0; i < nrest; i++ )
00136                         cout << setw ( 12 ) << setprecision ( 7 ) << lambda [ i ];
00137                 cout << endl;
00138 
00139                 cout << "Jacobian matrix of the constraints:" << endl;
00140                 switch ( minimizer [ k ]->memoryModel ( ) )
00141                 {
00142                         case yane::MinProg::BYLINE :
00143                                 for ( int i = 0; i < n; i++ )
00144                                 {
00145                                         for ( int j = 0; j < nrest; j++ )
00146                                                 cout << setw ( 12 ) << setprecision ( 7 )
00147                                                         << rest_jac [ i + j * n ];
00148                                         cout << endl;
00149                                 }
00150                                 break;
00151                         case yane::MinProg::BYCOLUMN :
00152                                 for ( int i = 0; i < n; i++ )
00153                                 {
00154                                         for ( int j = 0; j < nrest; j++ )
00155                                                 cout << setw ( 12 ) << setprecision ( 7 )
00156                                                         << rest_jac [ i * nrest + j ];
00157                                         cout << endl;
00158                                 }
00159                                 break;
00160                 }
00161                 cout << endl << endl;
00162         }
00163 
00164         for ( int i = 0; i < 3; i++ )
00165                 delete minimizer [ i ];
00166         delete [ ] minimizer;
00167         delete [ ] names;
00168         delete [ ] x0;
00169         delete [ ] x;
00170         delete [ ] lb;
00171         delete [ ] ub;
00172 
00173         return 0;
00174 }
00175 
00176 // Function which contains the objective function
00177 void objectivFunction ( const double *x, double *fx, void *params )
00178 {
00179         fx [ 0 ] = x [ 0 ] * x [ 3 ] * ( x [ 0 ] + x [ 1 ] + x [ 2 ] ) + x [ 2 ];
00180 }
00181 
00182 // Function which contains the constraints (equality constraints listet first)
00183 void restrictionFunction ( const double *x, double *restx, void *params )
00184 {
00185         restx [ 0 ] = 40.0 - ( x [ 0 ] * x [ 0 ] + x [ 1 ] * x [ 1 ] + x [ 2 ]
00186                 * x [ 2 ] + x [ 3 ] * x [ 3 ] ); // = 0.0
00187         restx [ 1 ] = - 25.0 + x [ 0 ] * x [ 1 ] * x [ 2 ] * x [ 3 ]; // >= 0.0
00188         // restx[2] = -4.0 + x[2];              // => x[2] >= 4.0;
00189 }
00190 
00191 // Function which contains the box constraints of the optimization variable
00192 void boxConstraints ( int dimension, double * lb, double * ub )
00193 {
00194         for ( int i = 0; i < dimension; i++ )
00195         {
00196                 lb [ i ] = 1.0;
00197                 ub [ i ] = 5.0;
00198         }
00199 }