Rythmos - Transient Integration for Differential Equations Version of the Day
diagonalTransientMain.cpp
00001 //@HEADER
00002 
00003 // ***********************************************************************
00004 //
00005 //                     Rythmos Package
00006 //                 Copyright (2006) Sandia Corporation
00007 //
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 //
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
00026 //
00027 // ***********************************************************************
00028 //@HEADER
00029 
00030 #include "EpetraExt_DiagonalTransientModel.hpp"
00031 #include "Rythmos_BackwardEulerStepper.hpp"
00032 #include "Rythmos_ImplicitBDFStepper.hpp"
00033 #include "Rythmos_ImplicitRKStepper.hpp"
00034 #include "Rythmos_ForwardSensitivityStepper.hpp"
00035 #include "Rythmos_TimeStepNonlinearSolver.hpp"
00036 #include "Rythmos_DefaultIntegrator.hpp"
00037 #include "Rythmos_SimpleIntegrationControlStrategy.hpp"
00038 #include "Rythmos_StepperAsModelEvaluator.hpp"
00039 #include "Rythmos_RKButcherTableauBuilder.hpp"
00040 #include "Stratimikos_DefaultLinearSolverBuilder.hpp"
00041 #include "Thyra_EpetraModelEvaluator.hpp"
00042 #include "Thyra_DirectionalFiniteDiffCalculator.hpp"
00043 #include "Thyra_TestingTools.hpp"
00044 #include "Teuchos_StandardCatchMacros.hpp"
00045 #include "Teuchos_GlobalMPISession.hpp"
00046 #include "Teuchos_DefaultComm.hpp"
00047 #include "Teuchos_VerboseObject.hpp"
00048 #include "Teuchos_XMLParameterListHelpers.hpp"
00049 #include "Teuchos_CommandLineProcessor.hpp"
00050 #include "Teuchos_VerbosityLevelCommandLineProcessorHelpers.hpp"
00051 #include "Teuchos_as.hpp"
00052 
00053 #ifdef HAVE_MPI
00054 #  include "Epetra_MpiComm.h"
00055 #else
00056 #  include "Epetra_SerialComm.h"
00057 #endif // HAVE_MPI
00058 
00059 
00060 namespace {
00061 
00062 
00063 const std::string TimeStepNonlinearSolver_name = "TimeStepNonlinearSolver";
00064 
00065 const std::string Stratimikos_name = "Stratimikos";
00066 
00067 const std::string DiagonalTransientModel_name = "DiagonalTransientModel";
00068 
00069 const std::string RythmosStepper_name = "Rythmos Stepper";
00070 
00071 const std::string RythmosIntegrator_name = "Rythmos Integrator";
00072 
00073 const std::string FdCalc_name = "FD Calc";
00074 
00075 
00076 Teuchos::RCP<const Teuchos::ParameterList>
00077 getValidParameters()
00078 {
00079   using Teuchos::RCP; using Teuchos::ParameterList;
00080   static RCP<const ParameterList> validPL;
00081   if (is_null(validPL)) {
00082     RCP<ParameterList> pl = Teuchos::parameterList();
00083     pl->sublist(TimeStepNonlinearSolver_name);
00084     pl->sublist(Stratimikos_name);
00085     pl->sublist(DiagonalTransientModel_name);
00086     pl->sublist(RythmosStepper_name);
00087     pl->sublist(RythmosIntegrator_name);
00088     pl->sublist(FdCalc_name);
00089     validPL = pl;
00090   }
00091   return validPL;
00092 }
00093 
00094 
00095 } // namespace
00096 
00097 
00098 int main(int argc, char *argv[])
00099 {
00100 
00101   using std::endl;
00102   typedef double Scalar;
00103   typedef double ScalarMag;
00104   using Teuchos::describe;
00105   using Teuchos::RCP;
00106   using Teuchos::rcp;
00107   using Teuchos::rcp_implicit_cast;
00108   using Teuchos::rcp_dynamic_cast;
00109   using Teuchos::as;
00110   using Teuchos::ParameterList;
00111   using Teuchos::CommandLineProcessor;
00112   typedef Teuchos::ParameterList::PrintOptions PLPrintOptions;
00113   typedef Thyra::ModelEvaluatorBase MEB;
00114   typedef Thyra::DefaultMultiVectorProductVectorSpace<Scalar> DMVPVS;
00115   using Thyra::productVectorBase;
00116 
00117   bool result, success = true;
00118 
00119   Teuchos::GlobalMPISession mpiSession(&argc,&argv);
00120 
00121   RCP<Epetra_Comm> epetra_comm;
00122 #ifdef HAVE_MPI
00123   epetra_comm = rcp( new Epetra_MpiComm(MPI_COMM_WORLD) );
00124 #else
00125   epetra_comm = rcp( new Epetra_SerialComm );
00126 #endif // HAVE_MPI
00127 
00128   RCP<Teuchos::FancyOStream>
00129     out = Teuchos::VerboseObjectBase::getDefaultOStream();
00130 
00131   try {
00132 
00133     //
00134     // Read commandline options
00135     //
00136 
00137     CommandLineProcessor clp;
00138     clp.throwExceptions(false);
00139     clp.addOutputSetupOptions(true);
00140 
00141     std::string paramsFileName = "";
00142     clp.setOption( "params-file", &paramsFileName,
00143       "File name for XML parameters" );
00144 
00145     std::string extraParamsString = "";
00146     clp.setOption( "extra-params", &extraParamsString,
00147       "Extra XML parameters" );
00148 
00149     std::string extraParamsFile = "";
00150     clp.setOption( "extra-params-file", &extraParamsFile, "File containing extra parameters in XML format.");
00151 
00152     double maxStateError = 1e-6;
00153     clp.setOption( "max-state-error", &maxStateError,
00154       "The maximum allowed error in the integrated state in relation to the exact state solution" );
00155 
00156     double finalTime = 1e-3;
00157     clp.setOption( "final-time", &finalTime,
00158       "Final integration time (initial time is 0.0)" );
00159 
00160     int numTimeSteps = 10;
00161     clp.setOption( "num-time-steps", &numTimeSteps,
00162       "Number of (fixed) time steps.  If <= 0.0, then variable time steps are taken" );
00163 
00164     bool useBDF = false;
00165     clp.setOption( "use-BDF", "use-BE", &useBDF,
00166       "Use BDF or Backward Euler (BE)" );
00167 
00168     bool useIRK = false;
00169     clp.setOption( "use-IRK", "use-other", &useIRK,
00170       "Use IRK or something" );
00171 
00172     bool doFwdSensSolve = false;
00173     clp.setOption( "fwd-sens-solve", "state-solve", &doFwdSensSolve,
00174       "Do the forward sensitivity solve or just the state solve" );
00175 
00176     bool doFwdSensErrorControl = false;
00177     clp.setOption( "fwd-sens-err-cntrl", "no-fwd-sens-err-cntrl", &doFwdSensErrorControl,
00178       "Do error control on the forward sensitivity solve or not" );
00179 
00180     double maxRestateError = 0.0;
00181     clp.setOption( "max-restate-error", &maxRestateError,
00182       "The maximum allowed error between the state integrated by itself verses integrated along with DxDp" );
00183 
00184     double maxSensError = 1e-4;
00185     clp.setOption( "max-sens-error", &maxSensError,
00186       "The maximum allowed error in the integrated sensitivity in relation to"
00187       " the finite-difference sensitivity" );
00188 
00189     Teuchos::EVerbosityLevel verbLevel = Teuchos::VERB_DEFAULT;
00190     setVerbosityLevelOption( "verb-level", &verbLevel,
00191       "Top-level verbosity level.  By default, this gets deincremented as you go deeper into numerical objects.",
00192       &clp );
00193 
00194     bool testExactSensitivity = false;
00195     clp.setOption( "test-exact-sens", "no-test-exact-sens", &testExactSensitivity,
00196       "Test the exact sensitivity with finite differences or not." );
00197 
00198     bool dumpFinalSolutions = false;
00199     clp.setOption(
00200       "dump-final-solutions", "no-dump-final-solutions", &dumpFinalSolutions,
00201       "Determine if the final solutions are dumpped or not." );
00202 
00203     CommandLineProcessor::EParseCommandLineReturn parse_return = clp.parse(argc,argv);
00204     if( parse_return != CommandLineProcessor::PARSE_SUCCESSFUL ) return parse_return;
00205 
00206     if ( Teuchos::VERB_DEFAULT == verbLevel )
00207       verbLevel = Teuchos::VERB_LOW;
00208 
00209     const Teuchos::EVerbosityLevel
00210       solnVerbLevel = ( dumpFinalSolutions ? Teuchos::VERB_EXTREME : verbLevel );
00211 
00212     //
00213     // Get the base parameter list that all other parameter lists will be read
00214     // from.
00215     //
00216 
00217     RCP<ParameterList>
00218       paramList = Teuchos::parameterList();
00219     if (paramsFileName.length())
00220       updateParametersFromXmlFile( paramsFileName, paramList.ptr() );
00221     if(extraParamsFile.length())
00222       Teuchos::updateParametersFromXmlFile( "./"+extraParamsFile, paramList.ptr() );
00223     if (extraParamsString.length())
00224       updateParametersFromXmlString( extraParamsString, paramList.ptr() );
00225 
00226     if (testExactSensitivity) {
00227       paramList->sublist(DiagonalTransientModel_name).set("Exact Solution as Response",true);
00228     }
00229 
00230     paramList->validateParameters(*getValidParameters(),0); // Only validate top level lists!
00231 
00232     //
00233     // Create the Stratimikos linear solver factory.
00234     //
00235     // This is the linear solve strategy that will be used to solve for the
00236     // linear system with the W.
00237     //
00238 
00239     Stratimikos::DefaultLinearSolverBuilder linearSolverBuilder;
00240     linearSolverBuilder.setParameterList(sublist(paramList,Stratimikos_name));
00241     RCP<Thyra::LinearOpWithSolveFactoryBase<Scalar> >
00242       W_factory = createLinearSolveStrategy(linearSolverBuilder);
00243 
00244     //
00245     // Create the underlying EpetraExt::ModelEvaluator
00246     //
00247 
00248     RCP<EpetraExt::DiagonalTransientModel>
00249       epetraStateModel = EpetraExt::diagonalTransientModel(
00250         epetra_comm,
00251         sublist(paramList,DiagonalTransientModel_name)
00252         );
00253 
00254     *out <<"\nepetraStateModel valid options:\n";
00255     epetraStateModel->getValidParameters()->print(
00256       *out, PLPrintOptions().indent(2).showTypes(true).showDoc(true)
00257       );
00258 
00259     //
00260     // Create the Thyra-wrapped ModelEvaluator
00261     //
00262 
00263     RCP<Thyra::ModelEvaluator<double> >
00264       stateModel = epetraModelEvaluator(epetraStateModel,W_factory);
00265 
00266     *out << "\nParameter names = " << *stateModel->get_p_names(0) << "\n";
00267 
00268     //
00269     // Create the Rythmos stateStepper
00270     //
00271 
00272     RCP<Rythmos::TimeStepNonlinearSolver<double> >
00273       nonlinearSolver = Rythmos::timeStepNonlinearSolver<double>();
00274     RCP<ParameterList>
00275       nonlinearSolverPL = sublist(paramList,TimeStepNonlinearSolver_name);
00276     nonlinearSolverPL->get("Default Tol",1e-3*maxStateError); // Set default if not set
00277     nonlinearSolver->setParameterList(nonlinearSolverPL);
00278 
00279     RCP<Rythmos::StepperBase<Scalar> > stateStepper;
00280 
00281     if (useBDF) {
00282       stateStepper = rcp(
00283         new Rythmos::ImplicitBDFStepper<double>(
00284           stateModel, nonlinearSolver
00285           )
00286         );
00287     }
00288     else if (useIRK) {
00289       // We need a separate LOWSFB object for the IRK stepper
00290       RCP<Thyra::LinearOpWithSolveFactoryBase<Scalar> >
00291         irk_W_factory = createLinearSolveStrategy(linearSolverBuilder);
00292       RCP<Rythmos::RKButcherTableauBase<double> > irkbt = Rythmos::createRKBT<double>("Backward Euler");
00293       stateStepper = Rythmos::implicitRKStepper<double>(
00294         stateModel, nonlinearSolver, irk_W_factory, irkbt
00295         );
00296     }
00297     else {
00298       stateStepper = rcp(
00299         new Rythmos::BackwardEulerStepper<double>(
00300           stateModel, nonlinearSolver
00301           )
00302         );
00303     }
00304 
00305     *out <<"\nstateStepper:\n" << describe(*stateStepper,verbLevel);
00306     *out <<"\nstateStepper valid options:\n";
00307     stateStepper->getValidParameters()->print(
00308       *out, PLPrintOptions().indent(2).showTypes(true).showDoc(true)
00309       );
00310 
00311     stateStepper->setParameterList(sublist(paramList,RythmosStepper_name));
00312 
00313     //
00314     // Setup finite difference objects that will be used for tests
00315     //
00316 
00317     Thyra::DirectionalFiniteDiffCalculator<Scalar> fdCalc;
00318     fdCalc.setParameterList(sublist(paramList,FdCalc_name));
00319     fdCalc.setOStream(out);
00320     fdCalc.setVerbLevel(verbLevel);
00321 
00322     //
00323     // Use a StepperAsModelEvaluator to integrate the state
00324     //
00325 
00326     const MEB::InArgs<Scalar>
00327       state_ic = stateModel->getNominalValues();
00328     *out << "\nstate_ic:\n" << describe(state_ic,verbLevel);
00329 
00330     RCP<Rythmos::IntegratorBase<Scalar> > integrator;
00331     {
00332       RCP<ParameterList>
00333         integratorPL = sublist(paramList,RythmosIntegrator_name);
00334       integratorPL->set( "Take Variable Steps", as<bool>(numTimeSteps < 0) );
00335       integratorPL->set( "Fixed dt", as<double>((finalTime - state_ic.get_t())/numTimeSteps) );
00336       RCP<Rythmos::IntegratorBase<Scalar> >
00337         defaultIntegrator = Rythmos::controlledDefaultIntegrator<Scalar>(
00338           Rythmos::simpleIntegrationControlStrategy<Scalar>(integratorPL)
00339           );
00340       integrator = defaultIntegrator;
00341     }
00342 
00343     RCP<Rythmos::StepperAsModelEvaluator<Scalar> >
00344       stateIntegratorAsModel = Rythmos::stepperAsModelEvaluator(
00345         stateStepper, integrator, state_ic
00346         );
00347     stateIntegratorAsModel->setVerbLevel(verbLevel);
00348 
00349     *out << "\nUse the StepperAsModelEvaluator to integrate state x(p,finalTime) ... \n";
00350 
00351     RCP<Thyra::VectorBase<Scalar> > x_final;
00352 
00353     {
00354 
00355       Teuchos::OSTab tab(out);
00356 
00357       x_final = createMember(stateIntegratorAsModel->get_g_space(0));
00358 
00359       eval_g(
00360         *stateIntegratorAsModel,
00361         0, *state_ic.get_p(0),
00362         finalTime,
00363         0, &*x_final
00364         );
00365 
00366       *out
00367         << "\nx_final = x(p,finalTime) evaluated using stateIntegratorAsModel:\n"
00368         << describe(*x_final,solnVerbLevel);
00369 
00370     }
00371 
00372     //
00373     // Test the integrated state against the exact analytical state solution
00374     //
00375 
00376     RCP<const Thyra::VectorBase<Scalar> >
00377       exact_x_final = create_Vector(
00378         epetraStateModel->getExactSolution(finalTime),
00379         stateModel->get_x_space()
00380         );
00381 
00382     result = Thyra::testRelNormDiffErr(
00383       "exact_x_final", *exact_x_final, "x_final", *x_final,
00384       "maxStateError", maxStateError, "warningTol", 1.0, // Don't warn
00385       &*out, solnVerbLevel
00386       );
00387     if (!result) success = false;
00388 
00389     //
00390     // Solve and test the forward sensitivity computation
00391     //
00392 
00393     if (doFwdSensSolve) {
00394 
00395       //
00396       // Create the forward sensitivity stepper
00397       //
00398 
00399       RCP<Rythmos::ForwardSensitivityStepper<Scalar> > stateAndSensStepper =
00400         Rythmos::forwardSensitivityStepper<Scalar>();
00401       if (doFwdSensErrorControl) {
00402         stateAndSensStepper->initializeDecoupledSteppers(
00403           stateModel, 0, stateModel->getNominalValues(),
00404           stateStepper, nonlinearSolver,
00405           integrator->cloneIntegrator(), finalTime
00406           );
00407       }
00408       else {
00409         stateAndSensStepper->initializeSyncedSteppers(
00410           stateModel, 0, stateModel->getNominalValues(),
00411           stateStepper, nonlinearSolver
00412           );
00413         // The above call will result in stateStepper and nonlinearSolver being
00414         // cloned.  This helps to ensure consistency between the state and
00415         // sensitivity computations!
00416       }
00417 
00418       //
00419       // Set the initial condition for the state and forward sensitivities
00420       //
00421 
00422       RCP<Thyra::VectorBase<Scalar> > s_bar_init
00423         = createMember(stateAndSensStepper->getFwdSensModel()->get_x_space());
00424       assign( s_bar_init.ptr(), 0.0 );
00425       RCP<Thyra::VectorBase<Scalar> > s_bar_dot_init
00426         = createMember(stateAndSensStepper->getFwdSensModel()->get_x_space());
00427       assign( s_bar_dot_init.ptr(), 0.0 );
00428       // Above, I believe that these are the correct initial conditions for
00429       // s_bar and s_bar_dot given how the EpetraExt::DiagonalTransientModel
00430       // is currently implemented!
00431 
00432       RCP<const Rythmos::StateAndForwardSensitivityModelEvaluator<Scalar> >
00433         stateAndSensModel = stateAndSensStepper->getStateAndFwdSensModel();
00434 
00435       MEB::InArgs<Scalar>
00436         state_and_sens_ic = stateAndSensStepper->getModel()->createInArgs();
00437 
00438       // Copy time, parameters etc.
00439       state_and_sens_ic.setArgs(state_ic);
00440       // Set initial condition for x_bar = [ x; s_bar ]
00441       state_and_sens_ic.set_x(
00442         stateAndSensModel->create_x_bar_vec(state_ic.get_x(),s_bar_init)
00443         );
00444       // Set initial condition for x_bar_dot = [ x_dot; s_bar_dot ]
00445       state_and_sens_ic.set_x_dot(
00446         stateAndSensModel->create_x_bar_vec(state_ic.get_x_dot(),s_bar_dot_init)
00447         );
00448 
00449       *out << "\nstate_and_sens_ic:\n" << describe(state_and_sens_ic,verbLevel);
00450 
00451       stateAndSensStepper->setInitialCondition(state_and_sens_ic);
00452 
00453       //
00454       // Use a StepperAsModelEvaluator to integrate the state+sens
00455       //
00456 
00457       RCP<Rythmos::StepperAsModelEvaluator<Scalar> >
00458         stateAndSensIntegratorAsModel = Rythmos::stepperAsModelEvaluator(
00459           rcp_implicit_cast<Rythmos::StepperBase<Scalar> >(stateAndSensStepper),
00460           integrator, state_and_sens_ic
00461           );
00462       stateAndSensIntegratorAsModel->setVerbLevel(verbLevel);
00463 
00464       *out << "\nUse the StepperAsModelEvaluator to integrate state + sens x_bar(p,finalTime) ... \n";
00465 
00466       RCP<Thyra::VectorBase<Scalar> > x_bar_final;
00467 
00468       {
00469 
00470         Teuchos::OSTab tab(out);
00471 
00472         x_bar_final = createMember(stateAndSensIntegratorAsModel->get_g_space(0));
00473 
00474         eval_g(
00475           *stateAndSensIntegratorAsModel,
00476           0, *state_ic.get_p(0),
00477           finalTime,
00478           0, &*x_bar_final
00479           );
00480 
00481         *out
00482           << "\nx_bar_final = x_bar(p,finalTime) evaluated using stateAndSensIntegratorAsModel:\n"
00483           << describe(*x_bar_final,solnVerbLevel);
00484 
00485       }
00486 
00487       //
00488       // Test that the state computed above is same as computed initially!
00489       //
00490 
00491       *out << "\nChecking that x(p,finalTime) computed as part of x_bar above is the same ...\n";
00492 
00493       {
00494 
00495         Teuchos::OSTab tab(out);
00496 
00497         RCP<const Thyra::VectorBase<Scalar> >
00498           x_in_x_bar_final = productVectorBase<Scalar>(x_bar_final)->getVectorBlock(0);
00499 
00500         result = Thyra::testRelNormDiffErr<Scalar>(
00501           "x_final", *x_final,
00502           "x_in_x_bar_final", *x_in_x_bar_final,
00503           "maxRestateError", maxRestateError,
00504           "warningTol", 1.0, // Don't warn
00505           &*out, solnVerbLevel
00506           );
00507         if (!result) success = false;
00508 
00509       }
00510 
00511       //
00512       // Compute DxDp using finite differences
00513       //
00514 
00515       *out << "\nApproximating DxDp(p,t) using directional finite differences of integrator for x(p,t) ...\n";
00516 
00517       RCP<Thyra::MultiVectorBase<Scalar> > DxDp_fd_final;
00518 
00519       {
00520 
00521         Teuchos::OSTab tab(out);
00522 
00523 
00524         MEB::InArgs<Scalar>
00525           fdBasePoint = stateIntegratorAsModel->createInArgs();
00526 
00527         fdBasePoint.set_t(finalTime);
00528         fdBasePoint.set_p(0,stateModel->getNominalValues().get_p(0));
00529 
00530         DxDp_fd_final = createMembers(
00531           stateIntegratorAsModel->get_g_space(0),
00532           stateIntegratorAsModel->get_p_space(0)->dim()
00533           );
00534 
00535         typedef Thyra::DirectionalFiniteDiffCalculatorTypes::SelectedDerivatives
00536           SelectedDerivatives;
00537 
00538         MEB::OutArgs<Scalar> fdOutArgs =
00539           fdCalc.createOutArgs(
00540             *stateIntegratorAsModel,
00541             SelectedDerivatives().supports(MEB::OUT_ARG_DgDp,0,0)
00542             );
00543         fdOutArgs.set_DgDp(0,0,DxDp_fd_final);
00544 
00545         // Silence the model evaluators that are called.  The fdCal object
00546         // will show all of the inputs and outputs for each call.
00547         stateStepper->setVerbLevel(Teuchos::VERB_NONE);
00548         stateIntegratorAsModel->setVerbLevel(Teuchos::VERB_NONE);
00549 
00550         fdCalc.calcDerivatives(
00551           *stateIntegratorAsModel, fdBasePoint,
00552           stateIntegratorAsModel->createOutArgs(), // Don't bother with function value
00553           fdOutArgs
00554           );
00555 
00556         *out
00557           << "\nFinite difference DxDp_fd_final = DxDp(p,finalTime): "
00558           << describe(*DxDp_fd_final,solnVerbLevel);
00559 
00560       }
00561 
00562       //
00563       // Test that the integrated sens and the F.D. sens are similar
00564       //
00565 
00566       *out << "\nChecking that integrated DxDp(p,finalTime) and finite-diff DxDp(p,finalTime) are similar ...\n";
00567 
00568       {
00569 
00570         Teuchos::OSTab tab(out);
00571 
00572         RCP<const Thyra::VectorBase<Scalar> >
00573           DxDp_vec_final = Thyra::productVectorBase<Scalar>(x_bar_final)->getVectorBlock(1);
00574 
00575         RCP<const Thyra::VectorBase<Scalar> >
00576           DxDp_fd_vec_final = Thyra::multiVectorProductVector(
00577             rcp_dynamic_cast<const Thyra::DefaultMultiVectorProductVectorSpace<Scalar> >(
00578               DxDp_vec_final->range()
00579               ),
00580             DxDp_fd_final
00581             );
00582 
00583         result = Thyra::testRelNormDiffErr(
00584           "DxDp_vec_final", *DxDp_vec_final,
00585           "DxDp_fd_vec_final", *DxDp_fd_vec_final,
00586           "maxSensError", maxSensError,
00587           "warningTol", 1.0, // Don't warn
00588           &*out, solnVerbLevel
00589           );
00590         if (!result) success = false;
00591 
00592       }
00593 
00594     }
00595 
00596   }
00597   TEUCHOS_STANDARD_CATCH_STATEMENTS(true,*out,success);
00598 
00599   if(success)
00600     *out << "\nEnd Result: TEST PASSED" << endl;
00601   else
00602     *out << "\nEnd Result: TEST FAILED" << endl;
00603 
00604   return ( success ? 0 : 1 );
00605 
00606 } // end main() [Doxygen looks for this!]
00607 
 All Classes Functions Variables Typedefs Friends