Index: mttroot/mtt/bin/mtt ================================================================== --- mttroot/mtt/bin/mtt +++ mttroot/mtt/bin/mtt @@ -15,10 +15,13 @@ ############################################################### ## Version control history ############################################################### ## $Header$ ## $Log$ +## Revision 1.322 2001/08/01 04:06:07 geraint +## Added -i dassl for -cc and -oct. +## ## Revision 1.321 2001/07/27 23:43:34 geraint ## Added -cc to usage options (required for use with xmtt). ## ## Revision 1.320 2001/07/27 23:38:38 geraint ## Changed comment below (# SUMMARY) - fixes xmtt. @@ -2354,12 +2357,13 @@ ar -cr \$@ \$^ $1_ode2odes_dassl.oct : $1_ode.oct $1_odeo.oct mtt_dassl.oct @echo > /dev/null mtt_Solver.cc: mtt_Solver.hh -$1_ode2odes_${algebraic_solver}.cc: mtt_Solver.cc mtt_${algebraic_solver}.hh mtt_${algebraic_solver}.cc -$1_ode2odes_${algebraic_solver}.o: mtt_Solver.o mtt_${algebraic_solver}.o +mtt_AlgebraicSolver.cc: mtt_AlgebraicSolver.hh +$1_ode2odes_${algebraic_solver}.cc: mtt_Solver.cc mtt_AlgebraicSolver.cc mtt_${algebraic_solver}.hh mtt_${algebraic_solver}.cc +$1_ode2odes_${algebraic_solver}.o: mtt_Solver.o mtt_AlgebraicSolver.o mtt_${algebraic_solver}.o @echo "Creating $1_ode2odes_${algebraic_solver}.o" ar -cr \$@ \$^ #SUMMARY numpar numerical parameter declaration (m) $1_numpar.m: $1_numpar.txt $1_sympars.txt @@ -2701,15 +2705,16 @@ ${MTT_CXX} ${MTT_CXXFLAGS} ${MTT_CXXINCS} -c $1_ode2odes.cc -DSTANDALONE $1_ode2odes.oct: $1_ode2odes.cc $1_ode2odes_common.oct $1_ode2odes_${integration_method}.oct $1_ode2odes_${algebraic_solver}.o touch $1_ode2odes.m echo Creating $1_ode2odes.oct - $MKOCTFILE $1_ode2odes.cc mtt_${algebraic_solver}.cc mtt_Solver.cc + $MKOCTFILE $1_ode2odes.cc mtt_${algebraic_solver}.cc mtt_Solver.cc mtt_AlgebraicSolver.cc $1_ode2odes.cc: $1_def.r $1_sympars.txt\ $1_ode2odes_common.m $1_ode2odes_common.cc\ - $1_ode2odes_${integration_method}.m $1_ode2odes_${integration_method}.cc mtt_Solver.cc mtt_${algebraic_solver}.cc mtt_${algebraic_solver}.hh + $1_ode2odes_${integration_method}.m $1_ode2odes_${integration_method}.cc\ + mtt_Solver.cc mtt_AlgebraicSolver.cc mtt_${algebraic_solver}.cc mtt_${algebraic_solver}.hh touch $1_ode2odes.m make_ode2odes $1 cc $integration_method $algebraic_solver #Conversion of m to p to c #SUMMARY ode ordinary differential equations (c) Index: mttroot/mtt/bin/trans/make_ode2odes ================================================================== --- mttroot/mtt/bin/trans/make_ode2odes +++ mttroot/mtt/bin/trans/make_ode2odes @@ -7,10 +7,13 @@ ############################################################### ## Version control history ############################################################### ## $Id$ ## $Log$ +## Revision 1.63 2001/08/07 04:39:24 geraint +## Consolidated dassl and residual functions. +## ## Revision 1.62 2001/08/01 22:14:32 geraint ## Bug fix for dassl. ## ## Revision 1.61 2001/08/01 04:06:07 geraint ## Added -i dassl for -cc and -oct. @@ -535,11 +538,11 @@ mtt_input (ColumnVector &x, ColumnVector &y, const double &t, ColumnVector &par) { - static ${algebraic_solver} ae(mtt_ae,MTTNPAR,MTTNU,MTTNX,MTTNY,MTTNYZ); + static MTT::${algebraic_solver} ae(MTTNPAR,MTTNU,MTTNX,MTTNY,MTTNYZ); static ColumnVector u (MTTNU); #ifdef STANDALONE u = F${sys}_input (x, y, t, par); #else ADDED mttroot/mtt/lib/cc/mtt_AlgebraicSolver.cc Index: mttroot/mtt/lib/cc/mtt_AlgebraicSolver.cc ================================================================== --- /dev/null +++ mttroot/mtt/lib/cc/mtt_AlgebraicSolver.cc @@ -0,0 +1,33 @@ + +#include "mtt_AlgebraicSolver.hh" + +ColumnVector +MTT::AlgebraicSolver::solve (const ColumnVector &x, + const ColumnVector &u, + const double &t, + const ColumnVector &par) +{ + if (_nyz > 0) + { + _x = x; + _uui.insert(u,0); + _t = t; + _par = par; + _ui = ColumnVector(_nyz,1.0); + Solve(); + _uui.insert(_ui,_nu); + } + else + { + _uui = u; + } + return _uui; +} + +ColumnVector +MTT::AlgebraicSolver::eval (const ColumnVector &ui) +{ + if (_nyz > 0) + _uui.insert(ui,_nu); + return mtt_ae(_x,_uui,_t,_par); +} ADDED mttroot/mtt/lib/cc/mtt_AlgebraicSolver.hh Index: mttroot/mtt/lib/cc/mtt_AlgebraicSolver.hh ================================================================== --- /dev/null +++ mttroot/mtt/lib/cc/mtt_AlgebraicSolver.hh @@ -0,0 +1,49 @@ + +#ifndef MTT_ALGEBRAICSOLVER +#define MTT_ALGEBRAICSOLVER + + +#include "mtt_Solver.hh" + + +namespace MTT +{ + class AlgebraicSolver : public MTT::Solver + { + public: + + AlgebraicSolver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) + : MTT::Solver (npar,nu,nx,ny,nyz) + {;} + + ColumnVector + solve (const ColumnVector &x, + const ColumnVector &u, + const double &t, + const ColumnVector &par); + + ColumnVector + eval (const ColumnVector &ui); + + virtual ~AlgebraicSolver (void) {}; + + protected: + + virtual void + Solve (void) = 0; + }; +} + + +extern ColumnVector +mtt_ae(ColumnVector &x, + ColumnVector &u, + const double &t, + ColumnVector &par); + + +#endif // MTT_ALGEBRAICSOLVER Index: mttroot/mtt/lib/cc/mtt_HJ_Solver.hh ================================================================== --- mttroot/mtt/lib/cc/mtt_HJ_Solver.hh +++ mttroot/mtt/lib/cc/mtt_HJ_Solver.hh @@ -1,53 +1,65 @@ -#include "mtt_Solver.hh" - -class HJ_Solver : public Solver { - - // http://www.netlib.org/opt/hooke.c - // Hooke and Jeeves solution - -public: - - HJ_Solver (sys_ae ae, - const int npar, - const int nu, - const int nx, - const int ny, - const int nyz) - : Solver (ae,npar,nu,nx,ny,nyz) - { static_ptr = this; VARS = nyz; }; - - static double - f (double tryUi[], int nyz); - - ~HJ_Solver (void) {}; - -protected: - - void - Solve (void); - - double - best_nearby (double delta[], - double point[], - double prevbest, - int nvars); - - int - hooke (int nvars, // MTTNYZ - double startpt[], // user's initial guess - double endpt[], // result - double rho = 0.05, // geometric shrink factor - double epsilon = 1e-3, // end value stepsize - int itermax = 5000); // max # iterations - -private: - - int VARS; - -public: - - static HJ_Solver *static_ptr; - -}; - +#ifndef MTT_HJSOLVER +#define MTT_HJSOLVER + + +#include "mtt_AlgebraicSolver.hh" + + +namespace MTT +{ + class HJ_Solver : public MTT::AlgebraicSolver + { + // http://www.netlib.org/opt/hooke.c + // Hooke and Jeeves solution + + public: + + HJ_Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) + : MTT::AlgebraicSolver (npar,nu,nx,ny,nyz) + { + static_ptr = this; + VARS = nyz; + } + + static double + f (double tryUi[], int nyz); + + ~HJ_Solver (void) {}; + + protected: + + void + Solve (void); + + double + best_nearby (double delta[], + double point[], + double prevbest, + int nvars); + + int + hooke (int nvars, // MTTNYZ + double startpt[], // user's initial guess + double endpt[], // result + double rho = 0.05, // geometric shrink factor + double epsilon = 1e-3, // end value stepsize + int itermax = 5000); // max # iterations + + private: + + int VARS; + + public: + + static HJ_Solver *static_ptr; + + }; +} + + +#endif // MTT_HJSOLVER Index: mttroot/mtt/lib/cc/mtt_Hybrd_Solver.cc ================================================================== --- mttroot/mtt/lib/cc/mtt_Hybrd_Solver.cc +++ mttroot/mtt/lib/cc/mtt_Hybrd_Solver.cc @@ -1,22 +1,31 @@ #include "mtt_Hybrd_Solver.hh" // http://www.netlib.org/minpack/hybrd.f // used by Octave's fsolve - -Hybrd_Solver *Hybrd_Solver::static_ptr; + +MTT::Hybrd_Solver *MTT::Hybrd_Solver::static_ptr; + +MTT::Hybrd_Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) +{ + static_ptr = this; +} ColumnVector -Hybrd_Solver::f_hybrd (const ColumnVector &tryUi) +MTT::Hybrd_Solver::f_hybrd (const ColumnVector &tryUi) { - Hybrd_Solver::static_ptr->_yz = Hybrd_Solver::static_ptr->eval(tryUi); - return Hybrd_Solver::static_ptr->_yz; + MTT::Hybrd_Solver::static_ptr->_yz = MTT::Hybrd_Solver::static_ptr->eval(tryUi); + return MTT::Hybrd_Solver::static_ptr->_yz; } void -Hybrd_Solver::Solve (void) +MTT::Hybrd_Solver::Solve (void) { int info; static int input_errors; static int user_errors; static int convergences; @@ -62,6 +71,5 @@ << " limit (" << limit_errors << ")" << " unknown (" << unknown_errors << ")" << " (max error = " << std::abs(eval(_ui).max()) << ")" << std::endl; } } - Index: mttroot/mtt/lib/cc/mtt_Hybrd_Solver.hh ================================================================== --- mttroot/mtt/lib/cc/mtt_Hybrd_Solver.hh +++ mttroot/mtt/lib/cc/mtt_Hybrd_Solver.hh @@ -1,38 +1,43 @@ -#include "mtt_Solver.hh" -#include - -class Hybrd_Solver : public Solver { - - // http://www.netlib.org/minpack/hybrd.f - // used by Octave's fsolve - -public: - - Hybrd_Solver (sys_ae ae, - const int npar, - const int nu, - const int nx, - const int ny, - const int nyz) - : Solver (ae,npar,nu,nx,ny,nyz) - { - static_ptr = this; - } - - static ColumnVector - f_hybrd (const ColumnVector &tryUi); - - ~Hybrd_Solver (void) {}; - -protected: - - void - Solve (void); - -public: - - static Hybrd_Solver *static_ptr; - -}; - +#ifndef MTT_HYBRDSOLVER +#define MTT_HYBRDSOLVER + + +#include +#include "mtt_AlgebraicSolver.hh" + + +namespace MTT +{ + class Hybrd_Solver : public MTT::AlgebraicSolver + { + // http://www.netlib.org/minpack/hybrd.f + // used by Octave's fsolve + + public: + + Hybrd_Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) + : MTT::AlgebraicSolver (npar,nu,nx,ny,nyz); + + static ColumnVector + f_hybrd (const ColumnVector &tryUi); + + ~Hybrd_Solver (void) {}; + + protected: + + void + Solve (void); + + public: + + static Hybrd_Solver *static_ptr; + }; +} + + +#endif // MTT_HYBRDSOLVER Index: mttroot/mtt/lib/cc/mtt_Reduce_Solver.cc ================================================================== --- mttroot/mtt/lib/cc/mtt_Reduce_Solver.cc +++ mttroot/mtt/lib/cc/mtt_Reduce_Solver.cc @@ -1,20 +1,17 @@ #include "mtt_Reduce_Solver.hh" - void -Reduce_Solver::Solve (void) +MTT::Reduce_Solver::Solve (void) { - // std::cerr << "Error:" - // << " Symbolic solution of equations failed during model build" << std::endl - // << " Try using one of the other algebraic solution methods" << std::endl; + ; } ColumnVector -Reduce_Solver::solve (const ColumnVector &x, - const ColumnVector &u, - const double &t, - const ColumnVector &par) +MTT::Reduce_Solver::solve (const ColumnVector &x, + const ColumnVector &u, + const double &t, + const ColumnVector &par) { return u; } Index: mttroot/mtt/lib/cc/mtt_Reduce_Solver.hh ================================================================== --- mttroot/mtt/lib/cc/mtt_Reduce_Solver.hh +++ mttroot/mtt/lib/cc/mtt_Reduce_Solver.hh @@ -1,33 +1,41 @@ -#include "mtt_Solver.hh" - -class Reduce_Solver : public Solver { - - // Dummy class - // This will not be used unless the Reduce solver has failed earlier - // in the model build process - -public: - - Reduce_Solver (sys_ae ae, - const int npar, - const int nu, - const int nx, - const int ny, - const int nyz) - : Solver (ae,npar,nu,nx,ny,nyz) - { ; }; - - void - Solve (void); - - ColumnVector - solve (const ColumnVector &x, - const ColumnVector &u, - const double &t, - const ColumnVector &par); - - ~Reduce_Solver (void) {}; - -}; - +#ifndef MTT_REDUCESOLVER +#define MTT_REDUCESOLVER + + +#include "mtt_AlgebraicSolver.hh" + + +namespace MTT +{ + class Reduce_Solver : public MTT::AlgebraicSolver + { + // Dummy class + // This will not be used unless the Reduce solver has failed earlier + // in the model build process + + public: + + Reduce_Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) + : AlgebraicSolver (npar,nu,nx,ny,nyz) + {;} + + void + Solve (void); + + ColumnVector + solve (const ColumnVector &x, + const ColumnVector &u, + const double &t, + const ColumnVector &par); + + ~Reduce_Solver (void) {}; + }; +} + + +#endif // MTT_REDUCESOLVER Index: mttroot/mtt/lib/cc/mtt_Solver.cc ================================================================== --- mttroot/mtt/lib/cc/mtt_Solver.cc +++ mttroot/mtt/lib/cc/mtt_Solver.cc @@ -1,42 +1,16 @@ #include "mtt_Solver.hh" -Solver::Solver (sys_ae ae, - const int npar, - const int nu, - const int nx, - const int ny, - const int nyz) -{ - _ae = ae; +MTT::Solver::Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz) +{ _np = npar; _nu = nu; _nx = nx; _ny = ny; _nyz = nyz; - _uui = ColumnVector (_nu+_nyz); }; - -ColumnVector -Solver::solve (const ColumnVector &x, - const ColumnVector &u, - const double &t, - const ColumnVector &par) -{ - _x = x; - _uui.insert(u,0); - _t = t; - _par = par; - _ui = ColumnVector(_nyz,1.0); - Solve (); - _uui.insert(_ui,_nu); - return _uui; -} - -ColumnVector -Solver::eval (const ColumnVector &ui) -{ - _uui.insert(ui,_nu); - return _ae (_x, _uui, _t, _par); -} Index: mttroot/mtt/lib/cc/mtt_Solver.hh ================================================================== --- mttroot/mtt/lib/cc/mtt_Solver.hh +++ mttroot/mtt/lib/cc/mtt_Solver.hh @@ -1,63 +1,47 @@ #ifndef MTT_SOLVER #define MTT_SOLVER + #include #include #include #include -class Solver { - -protected: - typedef ColumnVector (*sys_ae) // pointer to F${sys}_ae function - (ColumnVector &,ColumnVector &,const double &t,ColumnVector &); - -public: - - Solver (sys_ae ae, - const int npar, - const int nu, - const int nx, - const int ny, - const int nyz); - - ColumnVector - solve (const ColumnVector &x, - const ColumnVector &u, - const double &t, - const ColumnVector &par); - - ColumnVector - eval (const ColumnVector &ui); - - virtual ~Solver (void) {}; - -protected: - - virtual void - Solve (void) = 0; - -protected: - - ColumnVector _x; - ColumnVector _uui; - double _t; - ColumnVector _par; - - ColumnVector _ui; - ColumnVector _yz; - - int _nu; - int _np; - int _nx; - int _ny; - int _nyz; - - sys_ae _ae; - -}; - + +namespace MTT +{ + class Solver + { + public: + + Solver (const int npar, + const int nu, + const int nx, + const int ny, + const int nyz); + + virtual ~Solver (void) {}; + + protected: + + ColumnVector _x; + ColumnVector _uui; + double _t; + ColumnVector _par; + + ColumnVector _ui; + ColumnVector _yz; + + int _nu; + int _np; + int _nx; + int _ny; + int _nyz; + }; +} + + #endif // MTT_SOLVER