Tempus_WrapStaggeredFSAModelEvaluator.hpp

Go to the documentation of this file.

// @HEADER
// ****************************************************************************
//                Tempus: Copyright (2017) Sandia Corporation
//
// Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
// ****************************************************************************
// @HEADER

#ifndef Tempus_WrapStaggeredFSAModelEvaluator_hpp
#define Tempus_WrapStaggeredFSAModelEvaluator_hpp

#include "Tempus_config.hpp"
#include "Tempus_SensitivityModelEvaluatorBase.hpp"
#include "Tempus_StaggeredForwardSensitivityModelEvaluator.hpp"
#include "Tempus_WrapperModelEvaluatorPairIMEX_Basic.hpp"
#include "Tempus_WrapperModelEvaluatorPairIMEX_StaggeredFSA.hpp"
#include "Tempus_WrapperModelEvaluatorPairPartIMEX_StaggeredFSA.hpp"

namespace Tempus {

template <typename Scalar>
Teuchos::RCP<SensitivityModelEvaluatorBase<Scalar> >
wrapStaggeredFSAModelEvaluator(
    const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model,
    const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >&
        sens_residual_model,
    const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& sens_solve_model,
    const bool is_pseudotransient,
    const Teuchos::RCP<const Teuchos::ParameterList>& pList = Teuchos::null)
{
  using Teuchos::RCP;
  using Teuchos::rcp_dynamic_cast;

  RCP<SensitivityModelEvaluatorBase<Scalar> > wrapped_model;

  // Test if model is an IMEX pair
  RCP<const WrapperModelEvaluatorPairIMEX_Basic<Scalar> > modelPairIMEX =
      rcp_dynamic_cast<const WrapperModelEvaluatorPairIMEX_Basic<Scalar> >(
          model);
  RCP<const WrapperModelEvaluatorPairPartIMEX_Basic<Scalar> >
      modelPairPartIMEX = rcp_dynamic_cast<
          const WrapperModelEvaluatorPairPartIMEX_Basic<Scalar> >(model);

  // It isn't clear how to handle separate sensitivity and solve model
  // evaluators in the IMEX case, so punt for now.  Would they also be IMEX MEs?
  if ((modelPairIMEX != Teuchos::null || modelPairPartIMEX != Teuchos::null) &&
      (model.ptr() != sens_residual_model.ptr()) &&
      (model.ptr() != sens_solve_model.ptr()))
    TEUCHOS_TEST_FOR_EXCEPTION(
        true, std::logic_error,
        "Unique model evaluators for state, sensitivity residual, and "
        "sensitivity solve is not supported for IMEX");

  if (modelPairIMEX != Teuchos::null) {
    wrapped_model = rcp(new WrapperModelEvaluatorPairIMEX_StaggeredFSA<Scalar>(
        modelPairIMEX, is_pseudotransient, pList));
  }
  else if (modelPairPartIMEX != Teuchos::null) {
    wrapped_model =
        rcp(new WrapperModelEvaluatorPairPartIMEX_StaggeredFSA<Scalar>(
            modelPairPartIMEX, is_pseudotransient, pList));
  }
  else {
    wrapped_model = rcp(new StaggeredForwardSensitivityModelEvaluator<Scalar>(
        model, sens_residual_model, sens_solve_model, is_pseudotransient,
        pList));
  }

  return wrapped_model;
}

template <typename Scalar>
Teuchos::RCP<SensitivityModelEvaluatorBase<Scalar> >
wrapStaggeredFSAModelEvaluator(
    const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& model,
    const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& sens_residual_model,
    const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& sens_solve_model,
    const bool is_pseudotransient,
    const Teuchos::RCP<const Teuchos::ParameterList>& pList = Teuchos::null)
{
  Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> > cmodel = model;
  Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> > csens_residual_model =
      sens_residual_model;
  Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> > csens_solve_model =
      sens_solve_model;
  return wrapStaggeredFSAModelEvaluator(cmodel, csens_residual_model,
                                        csens_solve_model, is_pseudotransient,
                                        pList);
}

}  // namespace Tempus

#endif