MoochoPack_ThyraModelEvaluatorSolver.cpp

00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
00005 //                  Copyright (2003) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // 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 this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #include "MoochoPack_ThyraModelEvaluatorSolver.hpp"
00030 #include "NLPInterfacePack_NLPDirectThyraModelEvaluator.hpp"
00031 #include "NLPInterfacePack_NLPFirstOrderThyraModelEvaluator.hpp"
00032 #include "Thyra_DefaultFiniteDifferenceModelEvaluator.hpp"
00033 #include "Thyra_DefaultStateEliminationModelEvaluator.hpp"
00034 #include "Thyra_DefaultEvaluationLoggerModelEvaluator.hpp"
00035 #include "Thyra_SpmdMultiVectorFileIO.hpp"
00036 #include "Thyra_DampenedNewtonNonlinearSolver.hpp"
00037 #include "Thyra_VectorStdOps.hpp"
00038 
00039 namespace {
00040 
00041 typedef AbstractLinAlgPack::value_type  Scalar;
00042 
00043 Teuchos::RefCountPtr<Thyra::VectorBase<Scalar> >
00044 readVectorFromFile(
00045   const std::string                                                   &fileNameBase
00046   ,const Teuchos::RefCountPtr<const Thyra::VectorSpaceBase<Scalar> >  &vs
00047   ,const double                                                       scaleBy = 1.0
00048   )
00049 {
00050   Thyra::SpmdMultiVectorFileIO<Scalar> fileIO;
00051   Teuchos::RefCountPtr<Thyra::VectorBase<Scalar> >
00052     vec = fileIO.readVectorFromFile(fileNameBase,vs);
00053   Thyra::Vt_S(&*vec,scaleBy);
00054   return vec;
00055 }
00056 
00057 void writeVectorToFile(
00058   const Thyra::VectorBase<Scalar>    &vec
00059   ,const std::string                 &fileNameBase
00060   )
00061 {
00062   Thyra::SpmdMultiVectorFileIO<Scalar> fileIO;
00063   fileIO.writeToFile(vec,fileNameBase);
00064 }
00065 
00066 const int numFDMethodTypes = 7;
00067 
00068 Thyra::DirectionalFiniteDiffCalculatorTypes::EFDMethodType fdMethodTypeValues[numFDMethodTypes] =
00069 {
00070   Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_ONE
00071   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_TWO
00072   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_TWO_CENTRAL
00073   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_TWO_AUTO
00074   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_FOUR
00075   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_FOUR_CENTRAL
00076   ,Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_FOUR_AUTO
00077 };
00078 
00079 const char* fdMethodTypeNames[numFDMethodTypes] =
00080 {
00081   "order-one"
00082   ,"order-two"
00083   ,"order-two-central"
00084   ,"order-two-auto"
00085   ,"order-four"
00086   ,"order-four-central"
00087   ,"order-four-auto"
00088 };
00089 
00090 } // namespace
00091 
00092 
00093 namespace MoochoPack {
00094 
00095 ThyraModelEvaluatorSolver::ThyraModelEvaluatorSolver()
00096   :insertStateElimCommandLineOptions_(false)
00097   ,insertFiniteDiffCommandLineOptions_(false)
00098   ,do_sim_(false)
00099   ,use_direct_(false)
00100   ,use_black_box_(false)
00101   ,use_finite_diff_(false)
00102   ,fd_method_type_(Thyra::DirectionalFiniteDiffCalculatorTypes::FD_ORDER_ONE)
00103   ,fd_step_len_(1e-4)
00104   ,fwd_newton_tol_(-1.0)
00105   ,fwd_newton_max_iters_(20)
00106   ,stateGuessFileBase_("")
00107   ,scaleStateGuess_(1.0)
00108   ,paramGuessFileBase_("")
00109   ,scaleParamGuess_(1.0)
00110   ,stateSoluFileBase_("")
00111   ,paramSoluFileBase_("")
00112 {}
00113 
00114 void ThyraModelEvaluatorSolver::setupCLP(
00115   Teuchos::CommandLineProcessor *clp
00116   )
00117 {
00118   solver_.setup_commandline_processor(clp);
00119   clp->setOption( "do-sim", "do-opt",  &do_sim_, "Flag for if only the square constraints are solved" );
00120   clp->setOption( "use-direct", "use-first-order",  &use_direct_, "Flag for if we use the NLPDirect or NLPFirstOrderInfo implementation." );
00121   if(insertStateElimCommandLineOptions()) {
00122     clp->setOption( "black-box", "all-at-once",  &use_black_box_, "Flag for if we do black-box NAND or all-at-once SAND." );
00123     clp->setOption( "fwd-newton-tol", &fwd_newton_tol_, "Nonlinear residual tolerance for the forward Newton solve." );
00124     clp->setOption( "fwd-newton-max-iters", &fwd_newton_max_iters_, "Max number of iterations for the forward Newton solve." );
00125   }
00126   if(insertFiniteDiffCommandLineOptions()) {
00127     clp->setOption( "use-finite-diff", "no-use-finite-diff",  &use_finite_diff_, "Flag for if we are using finite differences or not." );
00128     clp->setOption(
00129       "finite-diff-type", &fd_method_type_
00130       ,numFDMethodTypes,fdMethodTypeValues,fdMethodTypeNames
00131       ,"The type of finite differences used [---use-finite-diff]"
00132       );
00133     clp->setOption( "fd-step-len", &fd_step_len_, "The finite-difference step length ( < 0 determine automatically )." );
00134   }
00135   clp->setOption( "x-guess-file", &stateGuessFileBase_, "Base file name to read the guess of the state x from." );
00136   clp->setOption( "scale-x-guess", &scaleStateGuess_, "Amount to scale the guess for x read in by --x-guess-file." );
00137   clp->setOption( "p-guess-file", &paramGuessFileBase_, "Base file name to read the guess of the parameters p from." );
00138   clp->setOption( "scale-p-guess", &scaleParamGuess_, "Amount to scale the guess for p read in by --p-guess-file." );
00139   clp->setOption( "x-solu-file", &stateSoluFileBase_, "Base file name to write the state solution x to." );
00140   clp->setOption( "p-solu-file", &paramSoluFileBase_, "Base file name to write the parameter solution p to." );
00141 }
00142 
00143 MoochoSolver& ThyraModelEvaluatorSolver::getSolver()
00144 {
00145   return solver_;
00146 }
00147 
00148 const MoochoSolver& ThyraModelEvaluatorSolver::getSolver() const
00149 {
00150   return solver_;
00151 }
00152 
00153 void ThyraModelEvaluatorSolver::setDoSim( const bool doSim )
00154 {
00155   do_sim_ = doSim;
00156 }
00157 
00158 bool ThyraModelEvaluatorSolver::getDoSim() const
00159 {
00160   return do_sim_;
00161 }
00162 
00163 void ThyraModelEvaluatorSolver::setModel(
00164   const Teuchos::RefCountPtr<Thyra::ModelEvaluator<value_type> > &model
00165   ,const int                                                     p_idx
00166   ,const int                                                     g_idx
00167   )
00168 {
00169 
00170   using Teuchos::rcp;
00171   using Teuchos::RefCountPtr;
00172   using NLPInterfacePack::NLP;
00173   using NLPInterfacePack::NLPDirectThyraModelEvaluator;
00174   using NLPInterfacePack::NLPFirstOrderThyraModelEvaluator;
00175 
00176   origModel_ = model;
00177   p_idx_ = p_idx;
00178   g_idx_ = g_idx;
00179 
00180   const int procRank = Teuchos::GlobalMPISession::getRank();
00181   const int numProcs = Teuchos::GlobalMPISession::getNProc();
00182 
00183   //
00184   // Wrap the orginal model
00185   //
00186 
00187   Teuchos::RefCountPtr<std::ostream>
00188     modelEvalLogOut;
00189   if(procRank==0)
00190     modelEvalLogOut = rcp(new std::ofstream("ModelEvaluationLog.out"));
00191   else
00192     modelEvalLogOut = rcp(new Teuchos::oblackholestream());
00193   Teuchos::RefCountPtr<Thyra::DefaultEvaluationLoggerModelEvaluator<Scalar> >
00194     loggerThyraModel
00195     = rcp(
00196       new Thyra::DefaultEvaluationLoggerModelEvaluator<Scalar>(
00197         origModel_,modelEvalLogOut
00198         )
00199       );
00200   nominalModel_
00201     = rcp(
00202       new Thyra::DefaultNominalBoundsOverrideModelEvaluator<Scalar>(loggerThyraModel,Teuchos::null)
00203       );
00204   finalPointModel_
00205     = rcp(
00206       new Thyra::DefaultFinalPointCaptureModelEvaluator<value_type>(nominalModel_)
00207       );
00208   outerModel_ = finalPointModel_;
00209 
00210   //
00211   // Create the NLP
00212   //
00213     
00214   Teuchos::RefCountPtr<NLP> nlp;
00215 
00216   if(do_sim_) {
00217     nlp = rcp(new NLPFirstOrderThyraModelEvaluator(outerModel_,-1,-1));
00218   }
00219   else {
00220     // Setup finite difference object
00221     RefCountPtr<Thyra::DirectionalFiniteDiffCalculator<Scalar> > direcFiniteDiffCalculator;
00222     if(use_finite_diff_) {
00223       direcFiniteDiffCalculator = rcp(new Thyra::DirectionalFiniteDiffCalculator<Scalar>());
00224       direcFiniteDiffCalculator->fd_method_type(fd_method_type_);
00225       direcFiniteDiffCalculator->fd_step_size(fd_step_len_);
00226     }
00227     if( use_black_box_ ) {
00228       // Create a Thyra::NonlinearSolverBase object to solve and eliminate the
00229       // state variables and the state equations
00230       Teuchos::RefCountPtr<Thyra::DampenedNewtonNonlinearSolver<Scalar> >
00231         stateSolver = rcp(new Thyra::DampenedNewtonNonlinearSolver<Scalar>()); // ToDo: Replace with MOOCHO!
00232       stateSolver->defaultTol(fwd_newton_tol_);
00233       stateSolver->defaultMaxNewtonIterations(fwd_newton_max_iters_);
00234       // Create the reduced  Thyra::ModelEvaluator object for p -> g_hat(p)
00235       Teuchos::RefCountPtr<Thyra::DefaultStateEliminationModelEvaluator<Scalar> >
00236         reducedThyraModel = rcp(new Thyra::DefaultStateEliminationModelEvaluator<Scalar>(outerModel_,stateSolver));
00237       Teuchos::RefCountPtr<Thyra::ModelEvaluator<Scalar> >
00238         finalReducedThyraModel;
00239       if(use_finite_diff_) {
00240         // Create the finite-difference wrapped Thyra::ModelEvaluator object
00241         Teuchos::RefCountPtr<Thyra::DefaultFiniteDifferenceModelEvaluator<Scalar> >
00242           fdReducedThyraModel = rcp(
00243             new Thyra::DefaultFiniteDifferenceModelEvaluator<Scalar>(
00244               reducedThyraModel,direcFiniteDiffCalculator)
00245             );
00246         finalReducedThyraModel = fdReducedThyraModel;
00247       }
00248       else {
00249         finalReducedThyraModel = reducedThyraModel;
00250       }
00251       // Wrap the reduced NAND Thyra::ModelEvaluator object in an NLP object
00252       nlp = rcp(new NLPFirstOrderThyraModelEvaluator(finalReducedThyraModel,p_idx,g_idx));
00253     }
00254     else {
00255       if(use_direct_) {
00256         Teuchos::RefCountPtr<NLPDirectThyraModelEvaluator>
00257           nlpDirect = rcp(new NLPDirectThyraModelEvaluator(outerModel_,p_idx,g_idx));
00258         if(use_finite_diff_) {
00259           nlpDirect->set_direcFiniteDiffCalculator(direcFiniteDiffCalculator);
00260         }
00261         nlp = nlpDirect;
00262       }
00263       else {
00264         nlp = rcp(new NLPFirstOrderThyraModelEvaluator(outerModel_,p_idx,g_idx));
00265       }
00266     }
00267   }
00268     
00269   // Set the NLP
00270   solver_.set_nlp(nlp);
00271 
00272 }
00273 
00274 void ThyraModelEvaluatorSolver::readInitialGuess(
00275   std::ostream *out_arg
00276   )
00277 {
00278   using Teuchos::OSTab;
00279   Teuchos::RefCountPtr<Teuchos::FancyOStream>
00280     out = Teuchos::getFancyOStream(Teuchos::rcp(out_arg,false));
00281   Teuchos::RefCountPtr<Thyra::ModelEvaluatorBase::InArgs<value_type> >
00282     initialGuess = Teuchos::rcp(new Thyra::ModelEvaluatorBase::InArgs<value_type>(origModel_->createInArgs()));
00283   if(stateGuessFileBase_ != "") {
00284     if(out.get())
00285       *out << "\nReading the guess of the state \'x\' from the file(s) with base name \""<<stateGuessFileBase_<<"\" ...\n";
00286     initialGuess->set_x(readVectorFromFile(stateGuessFileBase_,origModel_->get_x_space(),scaleStateGuess_));
00287   }
00288   else {
00289     if(out.get())
00290       *out << "\nNot reading the guess of the state \'x\' from a file!\n";
00291   }
00292   if(paramGuessFileBase_ != "") {
00293     if(out.get())
00294       *out << "\nReading the guess of the parameters \'p\' from the file(s) with base name \""<<paramGuessFileBase_<<"\" ...\n";
00295     initialGuess->set_p(p_idx_,readVectorFromFile(paramGuessFileBase_,origModel_->get_p_space(p_idx_),scaleParamGuess_));
00296   }
00297   else {
00298     if(out.get())
00299       *out << "\nNot reading the guess of the parameters \'p\' from a file!\n";
00300   }
00301   nominalModel_->setNominalValues(initialGuess);
00302 }
00303 
00304 void ThyraModelEvaluatorSolver::setInitialGuess(
00305   const Teuchos::RefCountPtr<const Thyra::ModelEvaluatorBase::InArgs<value_type> > &initialGuess
00306   )
00307 {
00308   nominalModel_->setNominalValues(initialGuess);
00309 }
00310 
00311 void ThyraModelEvaluatorSolver::setInitialGuess(
00312   const Thyra::ModelEvaluatorBase::InArgs<value_type> &initialGuess
00313   )
00314 {
00315   nominalModel_->setNominalValues(
00316     Teuchos::rcp(new Thyra::ModelEvaluatorBase::InArgs<value_type>(initialGuess))
00317     );
00318 }
00319   
00320 MoochoSolver::ESolutionStatus ThyraModelEvaluatorSolver::solve()
00321 {
00322   return solver_.solve_nlp();
00323 }
00324 
00325 void ThyraModelEvaluatorSolver::writeFinalSolution(
00326   std::ostream *out_arg
00327   ) const
00328 {
00329   using Teuchos::OSTab;
00330   Teuchos::RefCountPtr<Teuchos::FancyOStream>
00331     out = Teuchos::getFancyOStream(Teuchos::rcp(out_arg,false));
00332   if( stateSoluFileBase_ != "" && finalPointModel_->getFinalPoint().get_x().get() ) {
00333     if(out.get())
00334       *out << "\nWriting the state solution \'x\' to the file(s) with base name \""<<stateSoluFileBase_<<"\" ...\n";
00335     writeVectorToFile(*finalPointModel_->getFinalPoint().get_x(),stateSoluFileBase_);
00336   }
00337   if( paramSoluFileBase_ != "" ) {
00338     if(out.get())
00339       *out << "\nWriting the parameter solution \'p\' to the file(s) with base name \""<<paramSoluFileBase_<<"\" ...\n";
00340     writeVectorToFile(*finalPointModel_->getFinalPoint().get_p(p_idx_),paramSoluFileBase_);
00341   }
00342 }
00343 
00344 const Thyra::ModelEvaluatorBase::InArgs<value_type>&
00345 ThyraModelEvaluatorSolver::getFinalPoint() const
00346 {
00347   return finalPointModel_->getFinalPoint();
00348 }
00349 
00350 } // namespace MoochoPack

Generated on Thu Sep 18 12:34:37 2008 for MOOCHO/Thyra Adapter Software by doxygen 1.3.9.1