Thyra Version of the Day
Thyra_DampenedNewtonNonlinearSolver.hpp
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //    Thyra: Interfaces and Support for Abstract Numerical Algorithms
00005 //                 Copyright (2004) 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 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov) 
00038 // 
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 #ifndef THYRA_DAMPENED_NEWTON_NONLINEAR_SOLVER_HPP
00043 #define THYRA_DAMPENED_NEWTON_NONLINEAR_SOLVER_HPP
00044 
00045 #include "Thyra_NonlinearSolverBase.hpp"
00046 #include "Thyra_ModelEvaluatorHelpers.hpp"
00047 #include "Thyra_TestingTools.hpp"
00048 #include "Teuchos_StandardMemberCompositionMacros.hpp"
00049 #include "Teuchos_StandardCompositionMacros.hpp"
00050 #include "Teuchos_VerboseObject.hpp"
00051 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
00052 #include "Teuchos_StandardParameterEntryValidators.hpp"
00053 #include "Teuchos_as.hpp"
00054 
00055 
00056 namespace Thyra {
00057 
00058 
00079 template <class Scalar>
00080 class DampenedNewtonNonlinearSolver : public NonlinearSolverBase<Scalar> {
00081 public:
00082 
00084   typedef Teuchos::ScalarTraits<Scalar> ST;
00086   typedef typename ST::magnitudeType ScalarMag;
00088   typedef Teuchos::ScalarTraits<ScalarMag> SMT;
00089 
00091   STANDARD_MEMBER_COMPOSITION_MEMBERS( ScalarMag, defaultTol );
00092 
00094   STANDARD_MEMBER_COMPOSITION_MEMBERS( int, defaultMaxNewtonIterations );
00095 
00097   STANDARD_MEMBER_COMPOSITION_MEMBERS( bool, useDampenedLineSearch  );
00098   
00100   STANDARD_MEMBER_COMPOSITION_MEMBERS( Scalar, armijoConstant );
00101   
00103   STANDARD_MEMBER_COMPOSITION_MEMBERS( int, maxLineSearchIterations );
00104 
00106   DampenedNewtonNonlinearSolver(
00107     const ScalarMag defaultTol = 1e-2
00108     ,const int defaultMaxNewtonIterations = 1000
00109     ,const bool useDampenedLineSearch = true
00110     ,const Scalar armijoConstant = 1e-4
00111     ,const int maxLineSearchIterations = 20
00112     );
00113 
00115   static RCP<const Teuchos::ParameterList>
00116   getValidSolveCriteriaExtraParameters();
00117 
00120 
00122   void setParameterList(RCP<Teuchos::ParameterList> const& paramList);
00124   RCP<Teuchos::ParameterList> getNonconstParameterList();
00126   RCP<Teuchos::ParameterList> unsetParameterList();
00128   RCP<const Teuchos::ParameterList> getParameterList() const;
00130   RCP<const Teuchos::ParameterList> getValidParameters() const;
00131 
00133 
00136 
00138   void setModel(
00139     const RCP<const ModelEvaluator<Scalar> > &model
00140     );
00142   RCP<const ModelEvaluator<Scalar> > getModel() const;
00144   SolveStatus<Scalar> solve(
00145     VectorBase<Scalar> *x,
00146     const SolveCriteria<Scalar> *solveCriteria,
00147     VectorBase<Scalar> *delta
00148     );
00150   RCP<const VectorBase<Scalar> > get_current_x() const;
00152   bool is_W_current() const;
00154   RCP<LinearOpWithSolveBase<Scalar> > get_nonconst_W(const bool forceUpToDate);
00156   RCP<const LinearOpWithSolveBase<Scalar> > get_W() const;
00158   void set_W_is_current(bool W_is_current);
00159 
00161 
00162 private:
00163 
00164   RCP<Teuchos::ParameterList> paramList_;
00165   RCP<const ModelEvaluator<Scalar> > model_;
00166   RCP<LinearOpWithSolveBase<Scalar> > J_;
00167   RCP<VectorBase<Scalar> > current_x_;
00168   bool J_is_current_;
00169 
00170 };
00171 
00172 // ////////////////////////
00173 // Defintions
00174 
00175 template <class Scalar>
00176 DampenedNewtonNonlinearSolver<Scalar>::DampenedNewtonNonlinearSolver(
00177   const ScalarMag my_defaultTol
00178   ,const int my_defaultMaxNewtonIterations
00179   ,const bool my_useDampenedLineSearch
00180   ,const Scalar my_armijoConstant
00181   ,const int my_maxLineSearchIterations
00182   )
00183   :defaultTol_(my_defaultTol)
00184   ,defaultMaxNewtonIterations_(my_defaultMaxNewtonIterations)
00185   ,useDampenedLineSearch_(my_useDampenedLineSearch)
00186   ,armijoConstant_(my_armijoConstant)
00187   ,maxLineSearchIterations_(my_maxLineSearchIterations)
00188   ,J_is_current_(false)
00189 {}
00190 
00191 template <class Scalar>
00192 RCP<const Teuchos::ParameterList>
00193 DampenedNewtonNonlinearSolver<Scalar>::getValidSolveCriteriaExtraParameters()
00194 {
00195   static RCP<const Teuchos::ParameterList> validSolveCriteriaExtraParameters;
00196   if(!validSolveCriteriaExtraParameters.get()) {
00197     RCP<Teuchos::ParameterList>
00198       paramList = Teuchos::rcp(new Teuchos::ParameterList);
00199     paramList->set("Max Iters",int(1000));
00200     validSolveCriteriaExtraParameters = paramList;
00201   }
00202   return validSolveCriteriaExtraParameters;
00203 }
00204 
00205 // Overridden from Teuchos::ParameterListAcceptor
00206 
00207 template<class Scalar>
00208 void DampenedNewtonNonlinearSolver<Scalar>::setParameterList(
00209   RCP<Teuchos::ParameterList> const& paramList
00210   )
00211 {
00212   using Teuchos::get;
00213   TEUCHOS_TEST_FOR_EXCEPT(is_null(paramList));
00214   paramList->validateParametersAndSetDefaults(*getValidParameters(),0);
00215   paramList_ = paramList;
00216   TEUCHOS_TEST_FOR_EXCEPT("ToDo: Implement!");
00217   Teuchos::readVerboseObjectSublist(&*paramList_,this);
00218 #ifdef TEUCHOS_DEBUG
00219   paramList_->validateParameters(*getValidParameters(),0);
00220 #endif // TEUCHOS_DEBUG
00221 }
00222 
00223 template<class Scalar>
00224 RCP<Teuchos::ParameterList>
00225 DampenedNewtonNonlinearSolver<Scalar>::getNonconstParameterList()
00226 {
00227   return paramList_;
00228 }
00229 
00230 template<class Scalar>
00231 RCP<Teuchos::ParameterList>
00232 DampenedNewtonNonlinearSolver<Scalar>::unsetParameterList()
00233 {
00234   RCP<Teuchos::ParameterList> _paramList = paramList_;
00235   paramList_ = Teuchos::null;
00236   return _paramList;
00237 }
00238 
00239 template<class Scalar>
00240 RCP<const Teuchos::ParameterList>
00241 DampenedNewtonNonlinearSolver<Scalar>::getParameterList() const
00242 {
00243   return paramList_;
00244 }
00245 
00246 template<class Scalar>
00247 RCP<const Teuchos::ParameterList>
00248 DampenedNewtonNonlinearSolver<Scalar>::getValidParameters() const
00249 {
00250   using Teuchos::setDoubleParameter; using Teuchos::setIntParameter;
00251   static RCP<const Teuchos::ParameterList> validPL;
00252   if (is_null(validPL)) {
00253     RCP<Teuchos::ParameterList>
00254       pl = Teuchos::parameterList();
00255     TEUCHOS_TEST_FOR_EXCEPT("ToDo: Implement!");
00256     Teuchos::setupVerboseObjectSublist(&*pl);
00257     validPL = pl;
00258   }
00259   return validPL;
00260 }
00261 
00262 // Overridden from NonlinearSolverBase
00263 
00264 template <class Scalar>
00265 void DampenedNewtonNonlinearSolver<Scalar>::setModel(
00266   const RCP<const ModelEvaluator<Scalar> > &model
00267   )
00268 {
00269   TEUCHOS_TEST_FOR_EXCEPT(model.get()==NULL);
00270   model_ = model;
00271   J_ = Teuchos::null;
00272   current_x_ = Teuchos::null;
00273   J_is_current_ = false;
00274 }
00275 
00276 template <class Scalar>
00277 RCP<const ModelEvaluator<Scalar> >
00278 DampenedNewtonNonlinearSolver<Scalar>::getModel() const
00279 {
00280   return model_;
00281 }
00282 
00283 template <class Scalar>
00284 SolveStatus<Scalar>
00285 DampenedNewtonNonlinearSolver<Scalar>::solve(
00286   VectorBase<Scalar> *x_inout
00287   ,const SolveCriteria<Scalar> *solveCriteria
00288   ,VectorBase<Scalar> *delta
00289   ) 
00290 {
00291 
00292   using std::endl;
00293   using Teuchos::as;
00294 
00295   // Validate input
00296 #ifdef TEUCHOS_DEBUG
00297   TEUCHOS_TEST_FOR_EXCEPT(0==x_inout);
00298   THYRA_ASSERT_VEC_SPACES(
00299     "DampenedNewtonNonlinearSolver<Scalar>::solve(...)",
00300     *x_inout->space(), *model_->get_x_space() );
00301 #endif
00302 
00303   // Get the output stream and verbosity level
00304   const RCP<Teuchos::FancyOStream> out = this->getOStream();
00305   const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
00306   const bool showNewtonIters = (verbLevel==Teuchos::VERB_LOW);
00307   const bool showLineSearchIters = (as<int>(verbLevel) >= as<int>(Teuchos::VERB_MEDIUM));
00308   const bool showNewtonDetails = (as<int>(verbLevel) >= as<int>(Teuchos::VERB_HIGH));
00309   const bool dumpAll = (as<int>(verbLevel) == as<int>(Teuchos::VERB_EXTREME)); 
00310   TEUCHOS_OSTAB;
00311   if(out.get() && showNewtonIters) {
00312     *out << "\nBeginning dampended Newton solve of model = " << model_->description() << "\n\n";
00313     if (!useDampenedLineSearch())
00314       *out << "\nDoing undampened newton ...\n\n";
00315   }
00316 
00317   // Initialize storage for algorithm
00318   if(!J_.get()) J_ = model_->create_W();
00319   RCP<VectorBase<Scalar> > f = createMember(model_->get_f_space());
00320   RCP<VectorBase<Scalar> > x = Teuchos::rcp(x_inout,false);
00321   RCP<VectorBase<Scalar> > dx = createMember(model_->get_x_space());
00322   RCP<VectorBase<Scalar> > x_new = createMember(model_->get_x_space());
00323   RCP<VectorBase<Scalar> > ee = createMember(model_->get_x_space());
00324   V_S(ee.ptr(),ST::zero());
00325 
00326   // Get convergence criteria
00327   ScalarMag tol = this->defaultTol();
00328   int maxIters = this->defaultMaxNewtonIterations();
00329   if(solveCriteria && !solveCriteria->solveMeasureType.useDefault()) {
00330     TEUCHOS_TEST_FOR_EXCEPTION(
00331       !solveCriteria->solveMeasureType(SOLVE_MEASURE_NORM_RESIDUAL,SOLVE_MEASURE_NORM_RHS), CatastrophicSolveFailure
00332       ,"DampenedNewtonNonlinearSolver<Scalar>::solve(...): Error, can only support resudual-based"
00333       " convergence criteria!");
00334     tol = solveCriteria->requestedTol;
00335     if(solveCriteria->extraParameters.get()) {
00336       solveCriteria->extraParameters->validateParameters(*getValidSolveCriteriaExtraParameters());
00337       maxIters = solveCriteria->extraParameters->get("Max Iters",int(maxIters));
00338     }
00339   }
00340 
00341   if(out.get() && showNewtonDetails)
00342     *out << "\nCompute the initial starting point ...\n";
00343 
00344   eval_f_W( *model_, *x, &*f, &*J_ );
00345   if(out.get() && dumpAll) {
00346     *out << "\nInitial starting point:\n";
00347     *out << "\nx =\n" << *x;
00348     *out << "\nf =\n" << *f;
00349     *out << "\nJ =\n" << *J_;
00350   }
00351 
00352   // Peform the Newton iterations
00353   int newtonIter, num_residual_evals = 1;
00354   SolveStatus<Scalar> solveStatus;
00355   solveStatus.solveStatus = SOLVE_STATUS_UNCONVERGED;
00356 
00357   for( newtonIter = 1; newtonIter <= maxIters; ++newtonIter ) {
00358 
00359     if(out.get() && showNewtonDetails) *out << "\n*** newtonIter = " << newtonIter << endl;
00360 
00361     // Check convergence
00362     if(out.get() && showNewtonDetails) *out << "\nChecking for convergence ... : ";
00363     const Scalar phi = scalarProd(*f,*f), sqrt_phi = ST::squareroot(phi); // merit function: phi(f) = <f,f>
00364     solveStatus.achievedTol = sqrt_phi;
00365     const bool isConverged = sqrt_phi <= tol;
00366     if(out.get() && showNewtonDetails) *out
00367       << "sqrt(phi) = sqrt(<f,f>) = ||f|| = " << sqrt_phi << ( isConverged ? " <= " : " > " ) << "tol = " << tol << endl;
00368     if(out.get() && showNewtonIters) *out
00369       << "newton_iter="<<newtonIter<<": Check convergence: ||f|| = "
00370       << sqrt_phi << ( isConverged ? " <= " : " > " ) << "tol = " << tol << ( isConverged ? ", Converged!!!" : "" ) << endl;
00371     if(isConverged) {
00372       if(x_inout != x.get()) assign( ptr(x_inout), *x ); // Assign the solution if we have to
00373       if(out.get() && showNewtonDetails) {
00374         *out << "\nWe have converged :-)\n"
00375              << "\n||x||inf = " << norm_inf(*x) << endl;
00376         if(dumpAll) *out << "\nx =\n" << *x;
00377         *out << "\nExiting SimpleNewtonSolver::solve(...)\n";
00378       }
00379       std::ostringstream oss;
00380       oss << "Converged! ||f|| = " << sqrt_phi << ", num_newton_iters="<<newtonIter<<", num_residual_evals="<<num_residual_evals<<".";
00381       solveStatus.solveStatus = SOLVE_STATUS_CONVERGED;
00382       solveStatus.message = oss.str();
00383       break;
00384     }
00385     if(out.get() && showNewtonDetails) *out << "\nWe have to keep going :-(\n";
00386 
00387     // Compute the Jacobian if we have not already
00388     if(newtonIter > 1) {
00389       if(out.get() && showNewtonDetails) *out << "\nComputing the Jacobian J_ at current point ...\n";
00390       eval_f_W<Scalar>( *model_, *x, NULL, &*J_ );
00391       if(out.get() && dumpAll) *out << "\nJ =\n" << *J_;
00392     }
00393 
00394     // Compute the newton step: dx = -inv(J)*f
00395     if(out.get() && showNewtonDetails) *out << "\nComputing the Newton step: dx = - inv(J)*f ...\n";
00396     if(out.get() && showNewtonIters) *out << "newton_iter="<<newtonIter<<": Computing Newton step ...\n";
00397     assign( dx.ptr(), ST::zero() ); // Initial guess for the linear solve
00398     J_->solve(NOTRANS,*f,dx.ptr()); // Solve: J*dx = f
00399     Vt_S( dx.ptr(), Scalar(-ST::one()) ); // dx *= -1.0
00400     Vp_V( ee.ptr(), *dx); // ee += dx
00401     if(out.get() && showNewtonDetails) *out << "\n||dx||inf = " << norm_inf(*dx) << endl;
00402     if(out.get() && dumpAll) *out << "\ndy =\n" << *dx;
00403 
00404     // Perform backtracking armijo line search
00405     if(out.get() && showNewtonDetails) *out << "\nStarting backtracking line search iterations ...\n";
00406     if(out.get() && showNewtonIters) *out << "newton_iter="<<newtonIter<<": Starting backtracking line search ...\n";
00407     const Scalar Dphi = -2.0*phi; // D(phi(x+alpha*dx))/D(alpha) at alpha=0.0 => -2.0*<f,c>: where dx = -inv(J)*f
00408     Scalar alpha = 1.0; // Try a full step initially since it will eventually be accepted near solution
00409     int lineSearchIter;
00410     ++num_residual_evals;
00411     for( lineSearchIter = 1; lineSearchIter <= maxLineSearchIterations(); ++lineSearchIter, ++num_residual_evals ) {
00412       TEUCHOS_OSTAB_DIFF(lineSearchIter);
00413       if(out.get() && showNewtonDetails) *out << "\n*** lineSearchIter = " << lineSearchIter << endl;
00414       // x_new = x + alpha*dx
00415       assign( x_new.ptr(), *x ); Vp_StV( x_new.ptr(), alpha, *dx );
00416       if(out.get() && showNewtonDetails) *out << "\n||x_new||inf = " << norm_inf(*x_new) << endl;
00417       if(out.get() && dumpAll) *out << "\nx_new =\n" << *x_new;
00418       // Compute the residual at the updated point
00419       eval_f(*model_,*x_new,&*f);
00420       if(out.get() && dumpAll) *out << "\nf_new =\n" << *f;
00421       const Scalar phi_new = scalarProd(*f,*f), phi_frac = phi + alpha * armijoConstant() * Dphi;
00422       if(out.get() && showNewtonDetails) *out << "\nphi_new = <f_new,f_new> = " << phi_new << endl;
00423       if( Teuchos::ScalarTraits<Scalar>::isnaninf(phi_new) ) {
00424         if(out.get() && showNewtonDetails) *out << "\nphi_new is not a valid number, backtracking (alpha = 0.1*alpha) ...\n";
00425         alpha *= 0.1;
00426         continue;
00427       }
00428       const bool acceptPoint = (phi_new <= phi_frac);
00429       if(out.get() && showNewtonDetails) *out
00430         << "\nphi_new = " << phi_new << ( acceptPoint ? " <= " : " > " )
00431         << "phi + alpha * eta * Dphi = " << phi << " + " << alpha << " * " << armijoConstant() << " * " << Dphi
00432         << " = " << phi_frac << endl;
00433       if(out.get() && (showLineSearchIters || (showNewtonIters && acceptPoint))) *out
00434         << "newton_iter="<<newtonIter<<", ls_iter="<<lineSearchIter<<" : "
00435         << "phi(alpha="<<alpha<<") = "<<phi_new<<(acceptPoint ? " <=" : " >")<<" armijo_cord = " << phi_frac << endl;
00436       if (out.get() && showNewtonDetails && !useDampenedLineSearch())
00437         *out << "\nUndamped newton, always accpeting the point!\n";
00438       if( acceptPoint || !useDampenedLineSearch() ) {
00439         if(out.get() && showNewtonDetails) *out << "\nAccepting the current step with step length alpha = " << alpha << "!\n";
00440         break;
00441       }
00442       if(out.get() && showNewtonDetails) *out << "\nBacktracking (alpha = 0.5*alpha) ...\n";
00443       alpha *= 0.5;
00444     }
00445 
00446     // Check for line search failure
00447     if( lineSearchIter > maxLineSearchIterations() ) {
00448       std::ostringstream oss;
00449       oss
00450         << "lineSearchIter = " << lineSearchIter << " > maxLineSearchIterations = " << maxLineSearchIterations()
00451         << ": Linear search failure! Algorithm terminated!";
00452       solveStatus.message = oss.str();
00453       if(out.get() && (showNewtonIters || showNewtonDetails)) *out << endl << oss.str() << endl;
00454       goto exit;
00455     }
00456 
00457     // Take the Newton step
00458     std::swap<RCP<VectorBase<Scalar> > >( x_new, x ); // Now x is current point!
00459 
00460   }
00461 
00462 exit:
00463 
00464   if(out.get() && showNewtonIters) *out
00465     << "\n[Final] newton_iters="<<newtonIter<<", num_residual_evals="<<num_residual_evals<<"\n";
00466 
00467   if(newtonIter > maxIters) {
00468     std::ostringstream oss;
00469     oss
00470       << "newton_iter = " << newtonIter << " > maxIters = " << maxIters
00471       << ": Newton algorithm terminated!";
00472     solveStatus.message = oss.str();
00473     if( out.get() && (showNewtonIters || showNewtonDetails)) *out << endl << oss.str() << endl;
00474   }
00475 
00476   if(x_inout != x.get()) assign( ptr(x_inout), *x ); // Assign the final point
00477   if(delta != NULL) assign( ptr(delta), *ee );
00478   current_x_ = x_inout->clone_v(); // Remember the final point
00479   J_is_current_ = newtonIter==1; // J is only current with x if initial point was converged!
00480 
00481   if(out.get() && showNewtonDetails) *out
00482     << "\n*** Ending dampended Newton solve." << endl; 
00483 
00484   return solveStatus;
00485 
00486 }
00487 
00488 template <class Scalar>
00489 RCP<const VectorBase<Scalar> >
00490 DampenedNewtonNonlinearSolver<Scalar>::get_current_x() const
00491 {
00492   return current_x_;
00493 }
00494 
00495 template <class Scalar>
00496 bool DampenedNewtonNonlinearSolver<Scalar>::is_W_current() const
00497 {
00498   return J_is_current_;
00499 }
00500 
00501 template <class Scalar>
00502 RCP<LinearOpWithSolveBase<Scalar> >
00503 DampenedNewtonNonlinearSolver<Scalar>::get_nonconst_W(const bool forceUpToDate)
00504 {
00505   if (forceUpToDate) {
00506     TEUCHOS_TEST_FOR_EXCEPT(forceUpToDate);
00507   }
00508   return J_;
00509 }
00510 
00511 template <class Scalar>
00512 RCP<const LinearOpWithSolveBase<Scalar> >
00513 DampenedNewtonNonlinearSolver<Scalar>::get_W() const
00514 {
00515   return J_;
00516 }
00517 
00518 template <class Scalar>
00519 void DampenedNewtonNonlinearSolver<Scalar>::set_W_is_current(bool W_is_current)
00520 {
00521   J_is_current_ = W_is_current;
00522 }
00523 
00524 
00525 } // namespace Thyra
00526 
00527 
00528 #endif // THYRA_DAMPENED_NEWTON_NONLINEAR_SOLVER_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines