00001 /* 00002 * This file is part of ACADO Toolkit. 00003 * 00004 * ACADO Toolkit -- A Toolkit for Automatic Control and Dynamic Optimization. 00005 * Copyright (C) 2008-2014 by Boris Houska, Hans Joachim Ferreau, 00006 * Milan Vukov, Rien Quirynen, KU Leuven. 00007 * Developed within the Optimization in Engineering Center (OPTEC) 00008 * under supervision of Moritz Diehl. All rights reserved. 00009 * 00010 * ACADO Toolkit is free software; you can redistribute it and/or 00011 * modify it under the terms of the GNU Lesser General Public 00012 * License as published by the Free Software Foundation; either 00013 * version 3 of the License, or (at your option) any later version. 00014 * 00015 * ACADO Toolkit is distributed in the hope that it will be useful, 00016 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00017 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00018 * Lesser General Public License for more details. 00019 * 00020 * You should have received a copy of the GNU Lesser General Public 00021 * License along with ACADO Toolkit; if not, write to the Free Software 00022 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 00023 * 00024 */ 00025 00026 00027 00035 #include <acado_toolkit.hpp> 00036 #include <acado_gnuplot.hpp> 00037 00038 00039 int main( ){ 00040 00041 USING_NAMESPACE_ACADO 00042 00043 00044 DifferentialState phi, omega; // the states of the pendulum 00045 Parameter l, alpha ; // its length and the friction 00046 const double g = 9.81 ; // the gravitational constant 00047 DifferentialEquation f ; // the model equations 00048 Function h ; // the measurement function 00049 00050 VariablesGrid measurements; // read the measurements 00051 measurements = readFromFile( "data.txt" ); // from a file. 00052 00053 // -------------------------------------- 00054 OCP ocp(measurements.getTimePoints()); // construct an OCP 00055 h << phi ; // the state phi is measured 00056 ocp.minimizeLSQ( h, measurements ) ; // fit h to the data 00057 00058 f << dot(phi ) == omega ; // a symbolic implementation 00059 f << dot(omega) == -(g/l) *sin(phi ) // of the model 00060 - alpha*omega ; // equations 00061 00062 ocp.subjectTo( f ); // solve OCP s.t. the model, 00063 ocp.subjectTo( 0.0 <= alpha <= 4.0 ); // the bounds on alpha 00064 ocp.subjectTo( 0.0 <= l <= 2.0 ); // and the bounds on l. 00065 // -------------------------------------- 00066 00067 GnuplotWindow window; 00068 window.addSubplot( phi , "The angle phi", "time [s]", "angle [rad]" ); 00069 window.addSubplot( omega, "The angular velocity dphi" ); 00070 window.addData( 0, measurements(0) ); 00071 00072 ParameterEstimationAlgorithm algorithm(ocp); // the parameter estimation 00073 algorithm << window; 00074 algorithm.solve(); // algorithm solves the problem. 00075 00076 00077 return 0; 00078 } 00079 00080 00081