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 #include "LOCA_Stepper.H"
00042 #include "NOX_StatusTest_Generic.H"
00043
00044
00045 #include "NOX_Utils.H"
00046 #include "NOX_Solver_Factory.H"
00047 #include "LOCA_ErrorCheck.H"
00048 #include "LOCA_GlobalData.H"
00049 #include "LOCA_Factory.H"
00050 #include "LOCA_Parameter_Vector.H"
00051 #include "LOCA_Parameter_SublistParser.H"
00052 #include "LOCA_MultiPredictor_AbstractStrategy.H"
00053 #include "LOCA_MultiContinuation_AbstractStrategy.H"
00054 #include "LOCA_MultiContinuation_AbstractGroup.H"
00055 #include "LOCA_MultiContinuation_ExtendedGroup.H"
00056 #include "LOCA_MultiContinuation_ExtendedVector.H"
00057 #include "LOCA_StepSize_AbstractStrategy.H"
00058 #include "LOCA_Eigensolver_AbstractStrategy.H"
00059 #include "LOCA_SaveEigenData_AbstractStrategy.H"
00060 #include "LOCA_MultiContinuation_ConstrainedGroup.H"
00061 #include "LOCA_MultiContinuation_ConstraintInterface.H"
00062
00063 LOCA::Stepper::Stepper(
00064 const Teuchos::RCP<LOCA::GlobalData>& global_data,
00065 const Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup>& initialGuess,
00066 const Teuchos::RCP<NOX::StatusTest::Generic>& t,
00067 const Teuchos::RCP<Teuchos::ParameterList>& p) :
00068 LOCA::Abstract::Iterator(),
00069 globalData(),
00070 parsedParams(),
00071 predictor(),
00072 curGroupPtr(),
00073 prevGroupPtr(),
00074 eigensolver(),
00075 saveEigenData(),
00076 bifGroupPtr(),
00077 statusTestPtr(),
00078 paramListPtr(),
00079 stepperList(),
00080 solverPtr(),
00081 curPredictorPtr(),
00082 prevPredictorPtr(),
00083 stepSizeStrategyPtr(),
00084 conParamName(),
00085 conParamIDs(1),
00086 startValue(0.0),
00087 maxValue(0.0),
00088 minValue(0.0),
00089 stepSize(0.0),
00090 maxNonlinearSteps(15),
00091 targetValue(0.0),
00092 isTargetStep(false),
00093 doTangentFactorScaling(false),
00094 tangentFactor(1.0),
00095 minTangentFactor(0.1),
00096 tangentFactorExponent(1.0),
00097 calcEigenvalues(false),
00098 return_failed_on_max_steps(true)
00099 {
00100 reset(global_data, initialGuess, t, p);
00101 }
00102
00103 LOCA::Stepper::~Stepper()
00104 {
00105 }
00106
00107 bool
00108 LOCA::Stepper::reset(
00109 const Teuchos::RCP<LOCA::GlobalData>& global_data,
00110 const Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup>& initialGuess,
00111 const Teuchos::RCP<NOX::StatusTest::Generic>& t,
00112 const Teuchos::RCP<Teuchos::ParameterList>& p)
00113 {
00114 globalData = global_data;
00115 paramListPtr = p;
00116 statusTestPtr = t;
00117
00118
00119 parsedParams = Teuchos::rcp(new LOCA::Parameter::SublistParser(globalData));
00120 parsedParams->parseSublists(paramListPtr);
00121
00122
00123 stepperList = parsedParams->getSublist("Stepper");
00124
00125
00126 LOCA::Abstract::Iterator::resetIterator(*stepperList);
00127
00128
00129 Teuchos::RCP<Teuchos::ParameterList> predictorParams =
00130 parsedParams->getSublist("Predictor");
00131 predictor = globalData->locaFactory->createPredictorStrategy(
00132 parsedParams,
00133 predictorParams);
00134
00135
00136 Teuchos::RCP<Teuchos::ParameterList> eigenParams =
00137 parsedParams->getSublist("Eigensolver");
00138 eigensolver = globalData->locaFactory->createEigensolverStrategy(
00139 parsedParams,
00140 eigenParams);
00141
00142
00143 saveEigenData = globalData->locaFactory->createSaveEigenDataStrategy(
00144 parsedParams,
00145 eigenParams);
00146
00147
00148 Teuchos::RCP<Teuchos::ParameterList> stepsizeParams =
00149 parsedParams->getSublist("Step Size");
00150 stepSizeStrategyPtr = globalData->locaFactory->createStepSizeStrategy(
00151 parsedParams,
00152 stepsizeParams);
00153
00154
00155 if (stepperList->isParameter("Initial Value"))
00156 startValue = stepperList->get("Initial Value", 0.0);
00157 else {
00158 globalData->locaErrorCheck->throwError(
00159 "LOCA::Stepper::reset()",
00160 "\"Initial Value\" of continuation parameter is not set!");
00161 }
00162
00163
00164 if (stepperList->isParameter("Continuation Parameter")) {
00165 conParamName = stepperList->get("Continuation Parameter","None");
00166 initialGuess->setParam(conParamName,startValue);
00167 }
00168 else {
00169 globalData->locaErrorCheck->throwError(
00170 "LOCA::Stepper::reset()",
00171 "\"Continuation Parameter\" name is not set!");
00172 }
00173
00174
00175 const LOCA::ParameterVector& pv = initialGuess->getParams();
00176 conParamIDs[0] = pv.getIndex(conParamName);
00177
00178
00179 if (stepperList->isParameter("Max Value"))
00180 maxValue = stepperList->get("Max Value", 0.0);
00181 else {
00182 globalData->locaErrorCheck->throwError(
00183 "LOCA::Stepper::reset()",
00184 "\"Maximum Value\" of continuation parameter is not set!");
00185 }
00186 if (stepperList->isParameter("Min Value"))
00187 minValue = stepperList->get("Min Value", 0.0);
00188 else {
00189 globalData->locaErrorCheck->throwError(
00190 "LOCA::Stepper::reset()",
00191 "\"Minimum Value\" of continuation parameter is not set!");
00192 }
00193
00194
00195
00196 stepSize = stepSizeStrategyPtr->getStartStepSize();
00197 maxNonlinearSteps =
00198 stepperList->get("Max Nonlinear Iterations", 15);
00199
00200 targetValue = 0.0;
00201 isTargetStep = false;
00202 tangentFactor = 1.0;
00203 doTangentFactorScaling =
00204 stepperList->get("Enable Tangent Factor Step Size Scaling",
00205 false);
00206 minTangentFactor = stepperList->get("Min Tangent Factor",0.1);
00207 tangentFactorExponent =
00208 stepperList->get("Tangent Factor Exponent",1.0);
00209 calcEigenvalues = stepperList->get("Compute Eigenvalues",false);
00210 return_failed_on_max_steps =
00211 stepperList->get("Return Failed on Reaching Max Steps", true);
00212
00213
00214
00215 Teuchos::RCP<Teuchos::ParameterList> firstStepperParams =
00216 Teuchos::rcp(new Teuchos::ParameterList(*stepperList));
00217 firstStepperParams->set("Continuation Method", "Natural");
00218
00219
00220 Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup> constraintsGrp
00221 = buildConstrainedGroup(initialGuess);
00222
00223
00224 Teuchos::RCP<Teuchos::ParameterList> bifurcationParams =
00225 parsedParams->getSublist("Bifurcation");
00226 bifGroupPtr = globalData->locaFactory->createBifurcationStrategy(
00227 parsedParams,
00228 bifurcationParams,
00229 constraintsGrp);
00230
00231
00232 curGroupPtr = globalData->locaFactory->createContinuationStrategy(
00233 parsedParams,
00234 firstStepperParams,
00235 bifGroupPtr, predictor,
00236 conParamIDs);
00237
00238
00239 curGroupPtr->setStepSize(0.0);
00240
00241
00242 curGroupPtr->setPrevX(curGroupPtr->getX());
00243
00244
00245 solverPtr = NOX::Solver::buildSolver(curGroupPtr, statusTestPtr,
00246 parsedParams->getSublist("NOX"));
00247
00248 printInitializationInfo();
00249
00250 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperParameters))
00251 paramListPtr->print(globalData->locaUtils->out());
00252
00253 return true;
00254 }
00255
00256 LOCA::Abstract::Iterator::IteratorStatus
00257 LOCA::Stepper::start() {
00258 NOX::StatusTest::StatusType solverStatus;
00259 string callingFunction = "LOCA::Stepper::start()";
00260
00261
00262 curGroupPtr->preProcessContinuationStep(LOCA::Abstract::Iterator::Successful);
00263
00264 printStartStep();
00265
00266
00267 solverStatus = solverPtr->solve();
00268
00269
00270 if (solverStatus == NOX::StatusTest::Converged)
00271 curGroupPtr->postProcessContinuationStep(LOCA::Abstract::Iterator::Successful);
00272 else
00273 curGroupPtr->postProcessContinuationStep(LOCA::Abstract::Iterator::Unsuccessful);
00274
00275
00276 const LOCA::MultiContinuation::ExtendedGroup& constSolnGrp =
00277 dynamic_cast<const LOCA::MultiContinuation::ExtendedGroup&>(
00278 solverPtr->getSolutionGroup());
00279 Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup> underlyingGroup
00280 = Teuchos::rcp_const_cast<LOCA::MultiContinuation::AbstractGroup>(constSolnGrp.getUnderlyingGroup());
00281
00282
00283 curGroupPtr = globalData->locaFactory->createContinuationStrategy(
00284 parsedParams,
00285 stepperList,
00286 underlyingGroup,
00287 predictor,
00288 conParamIDs);
00289
00290
00291 if (solverStatus == NOX::StatusTest::Failed)
00292 printEndStep(LOCA::Abstract::Iterator::Unsuccessful);
00293 else
00294 printEndStep(LOCA::Abstract::Iterator::Successful);
00295
00296
00297 curGroupPtr->setStepSize(stepSize);
00298
00299 prevGroupPtr =
00300 Teuchos::rcp_dynamic_cast<LOCA::MultiContinuation::AbstractStrategy>(
00301 curGroupPtr->clone());
00302
00303
00304
00305 if (solverStatus != NOX::StatusTest::Converged)
00306 return LOCA::Abstract::Iterator::Failed;
00307
00308
00309 curGroupPtr->printSolution();
00310
00311
00312 if (calcEigenvalues) {
00313 Teuchos::RCP< std::vector<double> > evals_r;
00314 Teuchos::RCP< std::vector<double> > evals_i;
00315 Teuchos::RCP< NOX::Abstract::MultiVector > evecs_r;
00316 Teuchos::RCP< NOX::Abstract::MultiVector > evecs_i;
00317 eigensolver->computeEigenvalues(
00318 *curGroupPtr->getBaseLevelUnderlyingGroup(),
00319 evals_r, evals_i, evecs_r, evecs_i);
00320
00321 saveEigenData->save(evals_r, evals_i, evecs_r, evecs_i);
00322 }
00323
00324
00325 NOX::Abstract::Group::ReturnType predictorStatus =
00326 curGroupPtr->computePredictor();
00327 globalData->locaErrorCheck->checkReturnType(predictorStatus,
00328 callingFunction);
00329 curPredictorPtr =
00330 Teuchos::rcp_dynamic_cast<LOCA::MultiContinuation::ExtendedVector>(curGroupPtr->getPredictorTangent()[0].clone(NOX::DeepCopy));
00331 prevPredictorPtr =
00332 Teuchos::rcp_dynamic_cast<LOCA::MultiContinuation::ExtendedVector>(curGroupPtr->getPredictorTangent()[0].clone(NOX::ShapeCopy));
00333
00334
00335 solverPtr = NOX::Solver::buildSolver(curGroupPtr, statusTestPtr,
00336 parsedParams->getSublist("NOX"));
00337
00338 return LOCA::Abstract::Iterator::NotFinished;
00339 }
00340
00341 LOCA::Abstract::Iterator::IteratorStatus
00342 LOCA::Stepper::finish(LOCA::Abstract::Iterator::IteratorStatus itStatus)
00343 {
00344 string callingFunction = "LOCA::Stepper::finish()";
00345
00346
00347
00348
00349
00350
00351
00352
00353 curGroupPtr->copy(solverPtr->getSolutionGroup());
00354
00355
00356 if (itStatus == LOCA::Abstract::Iterator::Failed)
00357 return itStatus;
00358
00359 bool do_target = stepperList->get("Hit Continuation Bound", true);
00360 if (!do_target)
00361 return LOCA::Abstract::Iterator::Finished;
00362
00363
00364 double value = curGroupPtr->getContinuationParameter();
00365
00366 if (fabs(value-targetValue) > 1.0e-15*(1.0 + fabs(targetValue))) {
00367
00368 isTargetStep = true;
00369
00370
00371 prevGroupPtr->copy(*curGroupPtr);
00372
00373
00374 Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup> underlyingGrp
00375 = curGroupPtr->getUnderlyingGroup();
00376
00377
00378 Teuchos::RCP<Teuchos::ParameterList> lastStepPredictorParams =
00379 parsedParams->getSublist("Last Step Predictor");
00380
00381 lastStepPredictorParams->get("Method", "Constant");
00382 predictor = globalData->locaFactory->createPredictorStrategy(
00383 parsedParams,
00384 lastStepPredictorParams);
00385
00386
00387
00388 Teuchos::RCP<Teuchos::ParameterList> lastStepperParams =
00389 Teuchos::rcp(new Teuchos::ParameterList(*stepperList));
00390 lastStepperParams->set("Continuation Method", "Natural");
00391
00392
00393 curGroupPtr = globalData->locaFactory->createContinuationStrategy(
00394 parsedParams,
00395 lastStepperParams,
00396 underlyingGrp,
00397 predictor,
00398 conParamIDs);
00399
00400
00401 stepSize = targetValue - value;
00402 curGroupPtr->setStepSize(stepSize);
00403
00404
00405 NOX::Abstract::Group::ReturnType predictorStatus =
00406 curGroupPtr->computePredictor();
00407 globalData->locaErrorCheck->checkReturnType(predictorStatus,
00408 callingFunction);
00409 *curPredictorPtr = curGroupPtr->getPredictorTangent()[0];
00410
00411
00412 curGroupPtr->setPrevX(curGroupPtr->getX());
00413
00414
00415 curGroupPtr->computeX(*curGroupPtr, *curPredictorPtr, stepSize);
00416
00417
00418 curGroupPtr->preProcessContinuationStep(LOCA::Abstract::Iterator::Successful);
00419
00420 printStartStep();
00421
00422
00423 solverPtr = NOX::Solver::buildSolver(curGroupPtr, statusTestPtr,
00424 parsedParams->getSublist("NOX"));
00425
00426
00427 NOX::StatusTest::StatusType solverStatus = solverPtr->solve();
00428
00429
00430 if (solverStatus == NOX::StatusTest::Converged)
00431 curGroupPtr->postProcessContinuationStep(LOCA::Abstract::Iterator::Successful);
00432 else
00433 curGroupPtr->postProcessContinuationStep(LOCA::Abstract::Iterator::Unsuccessful);
00434
00435
00436 curGroupPtr->copy(solverPtr->getSolutionGroup());
00437
00438 if (solverStatus != NOX::StatusTest::Converged) {
00439 printEndStep(LOCA::Abstract::Iterator::Unsuccessful);
00440 return LOCA::Abstract::Iterator::Failed;
00441 }
00442
00443 printEndStep(LOCA::Abstract::Iterator::Successful);
00444
00445 curGroupPtr->printSolution();
00446
00447 }
00448
00449 return LOCA::Abstract::Iterator::Finished;
00450 }
00451
00452 LOCA::Abstract::Iterator::StepStatus
00453 LOCA::Stepper::preprocess(LOCA::Abstract::Iterator::StepStatus stepStatus)
00454 {
00455 if (stepStatus == LOCA::Abstract::Iterator::Unsuccessful) {
00456
00457
00458 curGroupPtr->copy(*prevGroupPtr);
00459 }
00460 else {
00461
00462
00463 prevGroupPtr->copy(*curGroupPtr);
00464 }
00465
00466
00467 stepStatus = computeStepSize(stepStatus, stepSize);
00468
00469
00470 curGroupPtr->setStepSize(stepSize);
00471
00472
00473 curGroupPtr->setPrevX(prevGroupPtr->getX());
00474
00475
00476 curGroupPtr->computeX(*prevGroupPtr, *curPredictorPtr, stepSize);
00477
00478
00479 curGroupPtr->preProcessContinuationStep(stepStatus);
00480
00481
00482 solverPtr = NOX::Solver::buildSolver(curGroupPtr, statusTestPtr,
00483 parsedParams->getSublist("NOX"));
00484
00485 return stepStatus;
00486 }
00487
00488 LOCA::Abstract::Iterator::StepStatus
00489 LOCA::Stepper::compute(LOCA::Abstract::Iterator::StepStatus stepStatus)
00490 {
00491 NOX::StatusTest::StatusType solverStatus;
00492
00493
00494 printStartStep();
00495
00496
00497 solverStatus = solverPtr->solve();
00498
00499
00500 if (solverStatus == NOX::StatusTest::Failed) {
00501 printEndStep(LOCA::Abstract::Iterator::Unsuccessful);
00502 return LOCA::Abstract::Iterator::Unsuccessful;
00503 }
00504
00505
00506 curGroupPtr->copy(solverPtr->getSolutionGroup());
00507
00508
00509 printEndStep(LOCA::Abstract::Iterator::Successful);
00510
00511 return LOCA::Abstract::Iterator::Successful;
00512 }
00513
00514 LOCA::Abstract::Iterator::StepStatus
00515 LOCA::Stepper::postprocess(LOCA::Abstract::Iterator::StepStatus stepStatus)
00516 {
00517 string callingFunction = "LOCA::Stepper::postprocess()";
00518
00519
00520 curGroupPtr->postProcessContinuationStep(stepStatus);
00521
00522 if (stepStatus == LOCA::Abstract::Iterator::Unsuccessful)
00523 return stepStatus;
00524
00525 *prevPredictorPtr = *curPredictorPtr;
00526
00527 NOX::Abstract::Group::ReturnType predictorStatus =
00528 curGroupPtr->computePredictor();
00529 globalData->locaErrorCheck->checkReturnType(predictorStatus,
00530 callingFunction);
00531 *curPredictorPtr = curGroupPtr->getPredictorTangent()[0];
00532
00533 if (doTangentFactorScaling && (getStepNumber() > 1)) {
00534 tangentFactor = curGroupPtr->computeScaledDotProduct(*curPredictorPtr,
00535 *prevPredictorPtr) /
00536 sqrt(curGroupPtr->computeScaledDotProduct(*curPredictorPtr,
00537 *curPredictorPtr) *
00538 curGroupPtr->computeScaledDotProduct(*prevPredictorPtr,
00539 *prevPredictorPtr));
00540
00541 if (tangentFactor < minTangentFactor) {
00542 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperDetails)) {
00543 globalData->locaUtils->out()
00544 << "\n\tTangent factor scaling: Failing step! Tangent factor "
00545 << "less than" << std::endl << "\t\tspecified bound: "
00546 << globalData->locaUtils->sciformat(tangentFactor) << " < "
00547 << globalData->locaUtils->sciformat(minTangentFactor) << std::endl;
00548 }
00549 return LOCA::Abstract::Iterator::Unsuccessful;
00550 }
00551 }
00552
00553
00554 curGroupPtr->printSolution();
00555
00556
00557 if (calcEigenvalues) {
00558 Teuchos::RCP< std::vector<double> > evals_r;
00559 Teuchos::RCP< std::vector<double> > evals_i;
00560 Teuchos::RCP< NOX::Abstract::MultiVector > evecs_r;
00561 Teuchos::RCP< NOX::Abstract::MultiVector > evecs_i;
00562 eigensolver->computeEigenvalues(
00563 *curGroupPtr->getBaseLevelUnderlyingGroup(),
00564 evals_r, evals_i, evecs_r, evecs_i);
00565
00566 saveEigenData->save(evals_r, evals_i, evecs_r, evecs_i);
00567 }
00568
00569 return stepStatus;
00570 }
00571
00572 LOCA::Abstract::Iterator::IteratorStatus
00573 LOCA::Stepper::stop(LOCA::Abstract::Iterator::StepStatus stepStatus)
00574 {
00575
00576
00577 if (LOCA::Abstract::Iterator::numTotalSteps
00578 >= LOCA::Abstract::Iterator::maxSteps) {
00579 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00580 globalData->locaUtils->out()
00581 << "\n\tContinuation run stopping: reached maximum number of steps "
00582 << LOCA::Abstract::Iterator::maxSteps << std::endl;
00583 }
00584 if (return_failed_on_max_steps)
00585 return LOCA::Abstract::Iterator::Failed;
00586 else
00587 return LOCA::Abstract::Iterator::Finished;
00588
00589 }
00590
00591 if (stepStatus == LOCA::Abstract::Iterator::Successful) {
00592
00593 double value = curGroupPtr->getContinuationParameter();
00594 double paramStep = value - prevGroupPtr->getContinuationParameter();
00595
00596
00597 if ( value >= maxValue*(1.0 - 1.0e-15) && paramStep > 0 ) {
00598 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00599 globalData->locaUtils->out()
00600 << "\n\tContinuation run stopping: parameter reached bound of "
00601 << globalData->locaUtils->sciformat(maxValue) << std::endl;
00602 }
00603 targetValue = maxValue;
00604 return LOCA::Abstract::Iterator::Finished;
00605 }
00606 if ( value <= minValue*(1.0 + 1.0e-15) && paramStep < 0 ) {
00607 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00608 globalData->locaUtils->out()
00609 << "\n\tContinuation run stopping: parameter reached bound of "
00610 << globalData->locaUtils->sciformat(minValue) << std::endl;
00611 }
00612 targetValue = minValue;
00613 return LOCA::Abstract::Iterator::Finished;
00614 }
00615
00616
00617 if (isLastIteration()) {
00618
00619
00620 if (withinThreshold()) {
00621 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00622 globalData->locaUtils->out()
00623 << "\n\tContinuation run stopping: parameter stepped to bound"
00624 << std::endl;
00625 }
00626 return LOCA::Abstract::Iterator::Finished;
00627 }
00628 else
00629 return LOCA::Abstract::Iterator::NotFinished;
00630 }
00631 }
00632 else if (isLastIteration())
00633 return LOCA::Abstract::Iterator::NotFinished;
00634
00635 return LOCA::Abstract::Iterator::NotFinished;
00636 }
00637
00638 Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup>
00639 LOCA::Stepper::buildConstrainedGroup(
00640 const Teuchos::RCP<LOCA::MultiContinuation::AbstractGroup>& grp)
00641 {
00642
00643 Teuchos::RCP<Teuchos::ParameterList> constraintsList =
00644 parsedParams->getSublist("Constraints");
00645
00646
00647 if (!constraintsList->isParameter("Constraint Object"))
00648 return grp;
00649
00650 string methodName = "LOCA::Stepper::buildConstrainedGroup()";
00651
00652 Teuchos::RCP<LOCA::MultiContinuation::ConstraintInterface> constraints;
00653 Teuchos::RCP< vector<string> > constraintParamNames;
00654
00655
00656 if ((*constraintsList).INVALID_TEMPLATE_QUALIFIER
00657 isType< Teuchos::RCP<LOCA::MultiContinuation::ConstraintInterface> >("Constraint Object"))
00658 constraints = (*constraintsList).INVALID_TEMPLATE_QUALIFIER
00659 get< Teuchos::RCP<LOCA::MultiContinuation::ConstraintInterface> >("Constraint Object");
00660 else
00661 globalData->locaErrorCheck->throwError(methodName,
00662 "\"Constraint Object\" parameter is not of type Teuchos::RCP<LOCA::MultiContinuation::ConstraintInterface>!");
00663
00664
00665 if ((*constraintsList).INVALID_TEMPLATE_QUALIFIER
00666 isType< Teuchos::RCP< vector<string> > > ("Constraint Parameter Names"))
00667 constraintParamNames = (*constraintsList).INVALID_TEMPLATE_QUALIFIER
00668 get< Teuchos::RCP< vector<string> > > ("Constraint Parameter Names");
00669 else
00670 globalData->locaErrorCheck->throwError(methodName,
00671 "\"Constraint Parameter Names\" parameter is not of type Teuchos::RCP< vector<string> >!");
00672
00673
00674 vector<int> constraintParamIDs(constraintParamNames->size());
00675 const LOCA::ParameterVector& pvec = grp->getParams();
00676 for (unsigned int i=0; i<constraintParamIDs.size(); i++)
00677 constraintParamIDs[i] = pvec.getIndex((*constraintParamNames)[i]);
00678
00679
00680 return
00681 Teuchos::rcp(new LOCA::MultiContinuation::ConstrainedGroup(
00682 globalData,
00683 parsedParams,
00684 constraintsList,
00685 grp,
00686 constraints,
00687 constraintParamIDs));
00688 }
00689
00690 LOCA::Abstract::Iterator::StepStatus
00691 LOCA::Stepper::computeStepSize(LOCA::Abstract::Iterator::StepStatus stepStatus,
00692 double& stepSz)
00693 {
00694 NOX::Abstract::Group::ReturnType res =
00695 stepSizeStrategyPtr->computeStepSize(*curGroupPtr, *curPredictorPtr,
00696 *solverPtr,
00697 stepStatus, *this, stepSz);
00698
00699 if (res == NOX::Abstract::Group::Failed)
00700 return LOCA::Abstract::Iterator::Provisional;
00701
00702 if (doTangentFactorScaling) {
00703 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperDetails)) {
00704 globalData->locaUtils->out()
00705 << "\n\tTangent factor scaling: Rescaling step size by "
00706 << globalData->locaUtils->sciformat(pow(fabs(tangentFactor),
00707 tangentFactorExponent))
00708 << std::endl;
00709 }
00710
00711 stepSz *= pow(fabs(tangentFactor), tangentFactorExponent);
00712 }
00713
00714
00715 double prevValue = curGroupPtr->getContinuationParameter();
00716 double dpds = curPredictorPtr->getScalar(0);
00717 if ( (prevValue+stepSz*dpds > maxValue*(1.0 - 1.0e-15)) ) {
00718 stepSz = (maxValue - prevValue)/dpds;
00719 targetValue = maxValue;
00720 setLastIteration();
00721 }
00722 if ( (prevValue+stepSz*dpds < minValue*(1.0 + 1.0e-15)) ) {
00723 stepSz = (minValue - prevValue)/dpds;
00724 targetValue = minValue;
00725 setLastIteration();
00726 }
00727
00728 return LOCA::Abstract::Iterator::Successful;
00729 }
00730
00731 Teuchos::RCP<const LOCA::MultiContinuation::AbstractGroup>
00732 LOCA::Stepper::getSolutionGroup() const
00733 {
00734 return curGroupPtr->getBaseLevelUnderlyingGroup();
00735 }
00736
00737 Teuchos::RCP<const LOCA::MultiContinuation::AbstractGroup>
00738 LOCA::Stepper::getBifurcationGroup() const
00739 {
00740 return curGroupPtr->getUnderlyingGroup();
00741 }
00742
00743 Teuchos::RCP<const Teuchos::ParameterList>
00744 LOCA::Stepper::getList() const
00745 {
00746 return paramListPtr;
00747 }
00748
00749 Teuchos::RCP<const NOX::Solver::Generic>
00750 LOCA::Stepper::getSolver() const
00751 {
00752 if (solverPtr.get() == NULL) {
00753 globalData->locaErrorCheck->throwError(
00754 "LOCA::Stepper::getSolver()",
00755 "Solver has not been constructed yet!");
00756 }
00757
00758 return solverPtr;
00759 }
00760
00761 void
00762 LOCA::Stepper::printInitializationInfo()
00763 {
00764 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00765 globalData->locaUtils->out()
00766 << std::endl
00767 << globalData->locaUtils->fill(72, '~') << std::endl;
00768
00769 globalData->locaUtils->out()
00770 << "Beginning Continuation Run \n"
00771 << "Stepper Method: "
00772 << stepperList->get("Continuation Method", "Arc Length")
00773 << "\n"
00774 << "Initial Parameter Value = "
00775 << globalData->locaUtils->sciformat(startValue)
00776 << "\n"
00777 << "Maximum Parameter Value = "
00778 << globalData->locaUtils->sciformat(maxValue) << "\n"
00779 << "Minimum Parameter Value = "
00780 << globalData->locaUtils->sciformat(minValue) << "\n"
00781 << "Maximum Number of Continuation Steps = "
00782 << LOCA::Abstract::Iterator::maxSteps
00783 << std::endl;
00784
00785 globalData->locaUtils->out()
00786 << globalData->locaUtils->fill(72, '~') << std::endl << std::endl;
00787 }
00788 }
00789
00790 void
00791 LOCA::Stepper::printStartStep()
00792 {
00793 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00794 globalData->locaUtils->out()
00795 << std::endl << globalData->locaUtils->fill(72, '~') << std::endl;
00796
00797 globalData->locaUtils->out()
00798 << "Start of Continuation Step " << stepNumber << " : ";
00799 if (stepNumber==0) {
00800 globalData->locaUtils->out()
00801 << "Attempting to converge initial guess at initial parameter "
00802 << "values." << std::endl;
00803 }
00804 else if (isTargetStep) {
00805 globalData->locaUtils->out()
00806 << "Attempting to hit final target value "
00807 << globalData->locaUtils->sciformat(targetValue) << std::endl;
00808 }
00809 else {
00810 globalData->locaUtils->out()
00811 << "Parameter: " << conParamName
00812 << " = "
00813 << globalData->locaUtils->sciformat(curGroupPtr->getContinuationParameter())
00814 << " from "
00815 << globalData->locaUtils->sciformat(prevGroupPtr->getContinuationParameter())
00816 << std::endl;
00817 globalData->locaUtils->out()
00818 << "Continuation Method: "
00819 << stepperList->get("Continuation Method", "Arc Length")
00820 << std::endl;
00821 globalData->locaUtils->out()
00822 << "Current step size = "
00823 << globalData->locaUtils->sciformat(stepSize) << " "
00824 << "Previous step size = "
00825 << globalData->locaUtils->sciformat(stepSizeStrategyPtr->getPrevStepSize())
00826 << std::endl;
00827 }
00828 globalData->locaUtils->out()
00829 << globalData->locaUtils->fill(72, '~') << std::endl << std::endl;
00830 }
00831 }
00832
00833 void
00834 LOCA::Stepper::printEndStep(LOCA::Abstract::Iterator::StepStatus stepStatus)
00835 {
00836 if (stepStatus == LOCA::Abstract::Iterator::Successful) {
00837
00838 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00839 globalData->locaUtils->out()
00840 << std::endl << globalData->locaUtils->fill(72, '~') << std::endl;
00841 globalData->locaUtils->out()
00842 << "End of Continuation Step " << stepNumber << " : "
00843 << "Parameter: " << conParamName << " = "
00844 << globalData->locaUtils->sciformat(curGroupPtr->getContinuationParameter());
00845 if (stepNumber != 0)
00846 globalData->locaUtils->out()
00847 << " from "
00848 << globalData->locaUtils->sciformat(prevGroupPtr->getContinuationParameter());
00849 globalData->locaUtils->out()
00850 << std::endl << "--> Step Converged in "
00851 << solverPtr->getNumIterations()
00852 <<" Nonlinear Solver Iterations!\n";
00853 globalData->locaUtils->out()
00854 << globalData->locaUtils->fill(72, '~') << std::endl << std::endl;
00855 }
00856 }
00857 else {
00858 if (globalData->locaUtils->isPrintType(NOX::Utils::StepperIteration)) {
00859
00860
00861 globalData->locaUtils->out()
00862 << std::endl << globalData->locaUtils->fill(72, '~') << std::endl;
00863 globalData->locaUtils->out()
00864 << "Continuation Step Number " << stepNumber
00865 << " experienced a convergence failure in\n"
00866 << "the nonlinear solver after "<< solverPtr->getNumIterations()
00867 <<" Iterations\n";
00868 globalData->locaUtils->out()
00869 << "Value of continuation parameter at failed step = "
00870 << globalData->locaUtils->sciformat(curGroupPtr->getContinuationParameter());
00871 if (stepNumber != 0)
00872 globalData->locaUtils->out()
00873 << " from "
00874 << globalData->locaUtils->sciformat(prevGroupPtr->getContinuationParameter());
00875 globalData->locaUtils->out()
00876 << std::endl << globalData->locaUtils->fill(72, '~') << endl;
00877 }
00878 }
00879 }
00880
00881 void
00882 LOCA::Stepper::printEndInfo()
00883 {
00884
00885 }
00886
00887 bool
00888 LOCA::Stepper::withinThreshold()
00889 {
00890 Teuchos::RCP<Teuchos::ParameterList> stepSizeList =
00891 parsedParams->getSublist("Step Size");
00892 double relt = stepperList->get("Relative Stopping Threshold", 0.9);
00893 double initialStep = stepSizeList->get("Initial Step Size", 1.0);
00894 double conParam = curGroupPtr->getContinuationParameter();
00895
00896 return (fabs(conParam-targetValue) < relt*fabs(initialStep));
00897 }