YANE-Framework Tutorial 1.1.0

examples/minprog/src/demo_gaussnewton.cpp

00001 /***************************************************************************
00002  *                                                                         *
00003  * Copyright (C) 2011 by www.nonlinearmpc.com                              *
00004  *                                                                         *
00005  * Authors:                                                                *
00006  *  Juergen Pannek <juergen.pannek@googlemail.com>                         *
00007  *  Michael Schamel <michael.schamel@uni-bayreuth.de>                      *
00008  *                                                                         *
00009  ***************************************************************************/
00010 
00011 #include <yane/utils.h>
00012 #include <yane/minprog.h>
00013 
00014 #include <cmath>
00015 #include <cstring>
00016 #include <fstream>
00017 #include <sstream>
00018 #include <time.h>
00019 #include <stdlib.h>
00020 
00021 using namespace std;
00022 using namespace yane;
00023 
00024 void function ( int * i, int * dim_x, int * dim_a, double * a, bool * achanged,
00025         double * x, double * fx, double * rpar, int * ipar )
00026 {
00027         fx [ 0 ] = a [ 0 ] * sin ( a [ 1 ] * x [ 0 ] ) + a [ 2 ];
00028 }
00029 
00030 void grad ( int * i, int * dim_x, int * dim_a, double * a, bool * achanged,
00031         double * x, double * dfx, double * rpar, int * ipar )
00032 {
00033         dfx [ 0 ] = sin ( a [ 1 ] * x [ 0 ] );
00034         dfx [ 1 ] = a [ 0 ] * x [ 0 ] * cos ( a [ 1 ] * x [ 0 ] );
00035         dfx [ 2 ] = 1;
00036 }
00037 
00038 double dist ( double a )
00039 {
00040         return a + ( double ) ( rand ( ) % 20 - 10 ) / 100.0;
00041 }
00042 
00043 int main ( void )
00044 {
00045         int dim_x, dim_a, dim_y;
00046         double *x, *a, *y;
00047         int *ipar;
00048         double *rpar;
00049         stringstream plot_functions, toplot;
00050         ofstream out;
00051 
00052         // Set dimensions
00053         dim_x = 1;
00054         dim_a = 3;
00055         dim_y = 9;
00056 
00057         // Create data arrays
00058         x = new double [ dim_x * dim_y ];
00059         a = new double [ dim_a ];
00060         y = new double [ dim_y ];
00061 
00062         // Initial data
00063         a [ 0 ] = 1.0;
00064         a [ 1 ] = 1.0;
00065         a [ 2 ] = 0.0;
00066 
00067         srand ( time ( NULL ));
00068 
00069         out.open ( "data.txt" );
00070         for ( int i = 0; i < dim_y; i++ )
00071         {
00072                 x [ i ] = ( double ) i / ( double ) dim_y * 6.0;
00073                 y [ i ] = dist ( a [ 0 ] ) * sin ( dist ( a [ 1 ] ) * x [ i ] + dist (
00074                         a [ 2 ] ) ) + dist ( a [ 3 ] );
00075                 out << x [ i ] << " " << y [ i ] << endl;
00076         }
00077         out.close ( );
00078 
00079         MinProg::GaussNewton * gn = new MinProg::GaussNewton ( dim_x, dim_a, dim_y,
00080                 x, y, function, grad );
00081 
00082         toplot << "plot 'data.txt'";
00083         for ( int i = dim_a; i <= dim_y; i = i + 2 )
00084         {
00085                 gn->calcParams ( a, i, rpar, ipar );
00086                 plot_functions << "f" << i << "(x) =" << a [ 0 ] << " * sin(" << a [ 1 ]
00087                         << " * x) + " << a [ 2 ] << endl;
00088                 toplot << ", f" << i << "(x)";
00089         }
00090 
00091         cout << plot_functions.str ( ) << endl;
00092 
00093         FILE * gnuplotpipe = popen ( "gnuplot", "w" );
00094 
00095         fprintf ( gnuplotpipe, "%s\n", plot_functions.str ( ).c_str ( ) );
00096         fprintf ( gnuplotpipe, "%s\n", toplot.str ( ).c_str ( ) );
00097         fflush ( gnuplotpipe );
00098 
00099         cin.get ( );
00100 
00101         fprintf ( gnuplotpipe, "exit\n" );
00102         pclose ( gnuplotpipe );
00103 
00104         delete gn;
00105         delete [ ] x;
00106         delete [ ] a;
00107         delete [ ] y;
00108 }