00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042 #include "NOX_Solver_LineSearchBased.H"
00043 #include "NOX_Abstract_Vector.H"
00044 #include "NOX_Abstract_Group.H"
00045 #include "NOX_Common.H"
00046 #include "Teuchos_ParameterList.hpp"
00047 #include "NOX_Utils.H"
00048 #include "NOX_GlobalData.H"
00049 #include "NOX_Solver_SolverUtils.H"
00050 #include "NOX_LineSearch_Generic.H"
00051 #include "NOX_LineSearch_Factory.H"
00052 #include "NOX_Direction_Generic.H"
00053 #include "NOX_Direction_Factory.H"
00054
00055 NOX::Solver::LineSearchBased::
00056 LineSearchBased(const Teuchos::RCP<NOX::Abstract::Group>& xGrp,
00057 const Teuchos::RCP<NOX::StatusTest::Generic>& t,
00058 const Teuchos::RCP<Teuchos::ParameterList>& p) :
00059 globalDataPtr(Teuchos::rcp(new NOX::GlobalData(p))),
00060 utilsPtr(globalDataPtr->getUtils()),
00061 solnPtr(xGrp),
00062 oldSolnPtr(xGrp->clone(DeepCopy)),
00063 dirPtr(xGrp->getX().clone(ShapeCopy)),
00064 testPtr(t),
00065 paramsPtr(p),
00066 prePostOperator(utilsPtr, paramsPtr->sublist("Solver Options"))
00067 {
00068 init();
00069 }
00070
00071
00072 void NOX::Solver::LineSearchBased::init()
00073 {
00074
00075 stepSize = 0.0;
00076 nIter = 0;
00077 status = NOX::StatusTest::Unconverged;
00078 checkType = parseStatusTestCheckType(paramsPtr->sublist("Solver Options"));
00079
00080 lineSearchPtr = NOX::LineSearch::
00081 buildLineSearch(globalDataPtr, paramsPtr->sublist("Line Search"));
00082
00083 directionPtr = NOX::Direction::
00084 buildDirection(globalDataPtr, paramsPtr->sublist("Direction"));
00085
00086
00087 if (utilsPtr->isPrintType(NOX::Utils::Parameters))
00088 {
00089 utilsPtr->out() << "\n" << NOX::Utils::fill(72) << "\n";
00090 utilsPtr->out() << "\n-- Parameters Passed to Nonlinear Solver --\n\n";
00091 paramsPtr->print(utilsPtr->out(),5);
00092 }
00093
00094 }
00095
00096 void NOX::Solver::LineSearchBased::
00097 reset(const NOX::Abstract::Vector& initialGuess,
00098 const Teuchos::RCP<NOX::StatusTest::Generic>& t)
00099 {
00100 solnPtr->setX(initialGuess);
00101 testPtr = t;
00102 init();
00103 }
00104
00105 void NOX::Solver::LineSearchBased::
00106 reset(const NOX::Abstract::Vector& initialGuess)
00107 {
00108 solnPtr->setX(initialGuess);
00109 init();
00110 }
00111
00112 NOX::Solver::LineSearchBased::~LineSearchBased()
00113 {
00114
00115 }
00116
00117
00118 NOX::StatusTest::StatusType NOX::Solver::LineSearchBased::getStatus()
00119 {
00120 return status;
00121 }
00122
00123 NOX::StatusTest::StatusType NOX::Solver::LineSearchBased::step()
00124 {
00125 prePostOperator.runPreIterate(*this);
00126
00127
00128 if (nIter == 0) {
00129
00130 NOX::Abstract::Group::ReturnType rtype = solnPtr->computeF();
00131 if (rtype != NOX::Abstract::Group::Ok) {
00132 utilsPtr->out() << "NOX::Solver::LineSearchBased::init - "
00133 << "Unable to compute F" << endl;
00134 throw "NOX Error";
00135 }
00136
00137
00138 status = testPtr->checkStatus(*this, checkType);
00139 if ((status == NOX::StatusTest::Converged) &&
00140 (utilsPtr->isPrintType(NOX::Utils::Warning))) {
00141 utilsPtr->out() << "Warning: NOX::Solver::LineSearchBased::init() - "
00142 << "The solution passed into the solver (either "
00143 << "through constructor or reset method) "
00144 << "is already converged! The solver wil not "
00145 << "attempt to solve this system since status is "
00146 << "flagged as converged." << endl;
00147 }
00148
00149 printUpdate();
00150 }
00151
00152
00153 if (status != NOX::StatusTest::Unconverged) {
00154 prePostOperator.runPostIterate(*this);
00155 printUpdate();
00156 return status;
00157 }
00158
00159
00160 NOX::Abstract::Group& soln = *solnPtr;
00161 NOX::StatusTest::Generic& test = *testPtr;
00162
00163
00164 bool ok;
00165 ok = directionPtr->compute(*dirPtr, soln, *this);
00166 if (!ok)
00167 {
00168 utilsPtr->out() << "NOX::Solver::LineSearchBased::iterate - unable to calculate direction" << endl;
00169 status = NOX::StatusTest::Failed;
00170 prePostOperator.runPostIterate(*this);
00171 printUpdate();
00172 return status;
00173 }
00174
00175
00176 nIter ++;
00177
00178
00179 *oldSolnPtr = *solnPtr;
00180
00181
00182 ok = lineSearchPtr->compute(soln, stepSize, *dirPtr, *this);
00183 if (!ok)
00184 {
00185 if (stepSize == 0.0)
00186 {
00187 utilsPtr->out() << "NOX::Solver::LineSearchBased::iterate - line search failed" << endl;
00188 status = NOX::StatusTest::Failed;
00189 prePostOperator.runPostIterate(*this);
00190 printUpdate();
00191 return status;
00192 }
00193 else if (utilsPtr->isPrintType(NOX::Utils::Warning))
00194 utilsPtr->out() << "NOX::Solver::LineSearchBased::iterate - using recovery step for line search" << endl;
00195 }
00196
00197
00198 NOX::Abstract::Group::ReturnType rtype = soln.computeF();
00199 if (rtype != NOX::Abstract::Group::Ok)
00200 {
00201 utilsPtr->out() << "NOX::Solver::LineSearchBased::iterate - unable to compute F" << endl;
00202 status = NOX::StatusTest::Failed;
00203 prePostOperator.runPostIterate(*this);
00204 printUpdate();
00205 return status;
00206 }
00207
00208
00209 status = test.checkStatus(*this, checkType);
00210
00211 prePostOperator.runPostIterate(*this);
00212
00213 printUpdate();
00214
00215 return status;
00216 }
00217
00218 NOX::StatusTest::StatusType NOX::Solver::LineSearchBased::solve()
00219 {
00220 prePostOperator.runPreSolve(*this);
00221
00222
00223 while (status == NOX::StatusTest::Unconverged)
00224 step();
00225
00226 Teuchos::ParameterList& outputParams = paramsPtr->sublist("Output");
00227 outputParams.set("Nonlinear Iterations", nIter);
00228 outputParams.set("2-Norm of Residual", solnPtr->getNormF());
00229
00230 prePostOperator.runPostSolve(*this);
00231
00232 return status;
00233 }
00234
00235 const NOX::Abstract::Group&
00236 NOX::Solver::LineSearchBased::getSolutionGroup() const
00237 {
00238 return *solnPtr;
00239 }
00240
00241 const NOX::Abstract::Group&
00242 NOX::Solver::LineSearchBased::getPreviousSolutionGroup() const
00243 {
00244 return *oldSolnPtr;
00245 }
00246
00247 int NOX::Solver::LineSearchBased::getNumIterations() const
00248 {
00249 return nIter;
00250 }
00251
00252 const Teuchos::ParameterList&
00253 NOX::Solver::LineSearchBased::getList() const
00254 {
00255 return *paramsPtr;
00256 }
00257
00258
00259 void NOX::Solver::LineSearchBased::printUpdate()
00260 {
00261 double normSoln = 0;
00262 double normStep = 0;
00263
00264
00265 if ((status == NOX::StatusTest::Unconverged) &&
00266 (utilsPtr->isPrintType(NOX::Utils::OuterIterationStatusTest)))
00267 {
00268 utilsPtr->out() << NOX::Utils::fill(72) << "\n";
00269 utilsPtr->out() << "-- Status Test Results --\n";
00270 testPtr->print(utilsPtr->out());
00271 utilsPtr->out() << NOX::Utils::fill(72) << "\n";
00272 }
00273
00274
00275 if (utilsPtr->isPrintType(NOX::Utils::OuterIteration))
00276 {
00277 normSoln = solnPtr->getNormF();
00278 normStep = (nIter > 0) ? dirPtr->norm() : 0;
00279 }
00280
00281
00282 if (utilsPtr->isPrintType(NOX::Utils::OuterIteration))
00283 {
00284 utilsPtr->out() << "\n" << NOX::Utils::fill(72) << "\n";
00285 utilsPtr->out() << "-- Nonlinear Solver Step " << nIter << " -- \n";
00286 utilsPtr->out() << "||F|| = " << utilsPtr->sciformat(normSoln);
00287 utilsPtr->out() << " step = " << utilsPtr->sciformat(stepSize);
00288 utilsPtr->out() << " dx = " << utilsPtr->sciformat(normStep);
00289 if (status == NOX::StatusTest::Converged)
00290 utilsPtr->out() << " (Converged!)";
00291 if (status == NOX::StatusTest::Failed)
00292 utilsPtr->out() << " (Failed!)";
00293 utilsPtr->out() << "\n" << NOX::Utils::fill(72) << "\n" << endl;
00294 }
00295
00296
00297 if ((status != NOX::StatusTest::Unconverged) &&
00298 (utilsPtr->isPrintType(NOX::Utils::OuterIteration)))
00299 {
00300 utilsPtr->out() << NOX::Utils::fill(72) << "\n";
00301 utilsPtr->out() << "-- Final Status Test Results --\n";
00302 testPtr->print(utilsPtr->out());
00303 utilsPtr->out() << NOX::Utils::fill(72) << "\n";
00304 }
00305 }
00306
00307 double NOX::Solver::LineSearchBased::getStepSize() const
00308 {
00309 return stepSize;
00310 }