PipelineManager.hpp

Go to the documentation of this file.

/** \file PipelineManager.hpp
 * \brief Header file for basic interface
 * \ingroup mofem_basic_interface
 *
 * Make basic interface, to speed up problem setup and analysis.
 * See discussion here
 * <a
 * href=https://groups.google.com/d/msg/mofem-group/Vkc00aia4dU/o9RF3ZmPAAAJ>link
 * to google groups</a>
 *
 */
 
#ifndef __BASIC_HPP__
#define __BASIC_HPP__
 
#include "UnknownInterface.hpp"
 
namespace MoFEM {
 
/**
 * \brief PipelineManager interface
 * \ingroup mofem_basic_interface
 */
struct PipelineManager : public UnknownInterface {
 
  /**
   * @brief Query interface for type-safe casting
   *
   * @param type_index Type index of the requested interface
   * @param iface Pointer to interface pointer to be set
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
                                 UnknownInterface **iface) const;
 
  /**
   * @brief Construct PipelineManager
   *
   * @param core Reference to MoFEM core instance
   */
  PipelineManager(const MoFEM::Core &core);
 
  using UserDataOperator = MoFEM::ForcesAndSourcesCore::UserDataOperator;
  using RuleHookFun = MoFEM::ForcesAndSourcesCore::RuleHookFun;
 
  using VolEle = MoFEM::VolumeElementForcesAndSourcesCore;
  using FaceEle = MoFEM::FaceElementForcesAndSourcesCore;
  using EdgeEle = MoFEM::EdgeElementForcesAndSourcesCore;
 
  /**
   * @brief Template struct for dimension-specific finite element types
   *
   * @tparam DIM Problem dimension (1, 2, or 3)
   */
  template <int DIM> struct ElementsAndOpsByDim;
 
  /**
   * @brief Get domain left-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to domain LHS finite element for matrix assembly
   */
  inline boost::shared_ptr<FEMethod> &getDomainLhsFE();
 
  /**
   * @brief Get domain right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to domain RHS finite element for vector assembly
   */
  inline boost::shared_ptr<FEMethod> &getDomainRhsFE();
 
  /**
   * @brief Get boundary left-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to boundary LHS finite element for matrix assembly
   */
  inline boost::shared_ptr<FEMethod> &getBoundaryLhsFE();
 
  /**
   * @brief Get boundary right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to boundary RHS finite element for vector assembly
   */
  inline boost::shared_ptr<FEMethod> &getBoundaryRhsFE();
 
  /**
   * @brief Get domain postprocessing finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to domain postprocessing finite element
   */
  inline boost::shared_ptr<FEMethod> &getDomainPostProcFE();
 
  /**
   * @brief Get boundary postprocessing finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to boundary postprocessing finite element
   */
  inline boost::shared_ptr<FEMethod> &getBoundaryPostProcFE();
 
  /**
   * @brief Get domain evaluation finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to domain evaluation finite element
   */
  inline boost::shared_ptr<FEMethod> &getEvaluationFE();
 
  /**
   * @brief Get skeleton left-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to skeleton LHS finite element for matrix assembly
   */
  inline boost::shared_ptr<FEMethod> &getSkeletonLhsFE();
 
  /**
   * @brief Get skeleton right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to skeleton RHS finite element for vector assembly
   */
  inline boost::shared_ptr<FEMethod> &getSkeletonRhsFE();
 
  /**
   * @brief Get domain explicit right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to domain explicit RHS finite element for IMEX schemes
   */
  inline boost::shared_ptr<FEMethod> &getDomainExplicitRhsFE();
 
  /**
   * @brief Get boundary explicit right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to boundary explicit RHS finite element for IMEX schemes
   */
  inline boost::shared_ptr<FEMethod> &getBoundaryExplicitRhsFE();
 
  /**
   * @brief Get skeleton explicit right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to skeleton explicit RHS finite element for IMEX schemes
   */
  inline boost::shared_ptr<FEMethod> &getSkeletonExplicitRhsFE();
 
  /**
   * @brief Get meshset right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to meshset RHS finite element for vector assembly
   */
  inline boost::shared_ptr<FEMethod> &getMeshsetRhsFE();
 
  /**
   * @brief Get meshset left-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to meshset LHS finite element for matrix assembly
   */
  inline boost::shared_ptr<FEMethod> &getMeshsetLhsFE();
 
  /**
   * @brief Get meshset explicit right-hand side finite element
   *
   * @return boost::shared_ptr<FEMethod>& Reference to meshset explicit RHS finite element for IMEX schemes
   */
  inline boost::shared_ptr<FEMethod> &getMeshsetExplicitRhsFE();
 
  /**
   * @brief Get typed domain left-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed domain LHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastDomainLhsFE();
 
  /**
   * @brief Get typed domain right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed domain RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastDomainRhsFE();
 
  /**
   * @brief Get typed boundary left-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed boundary LHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastBoundaryLhsFE();
 
  /**
   * @brief Get typed boundary right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed boundary RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastBoundaryRhsFE();
 
  /**
   * @brief Get typed domain postprocessing finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed domain postprocessing finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastDomainPostProcFE();
 
  /**
   * @brief Get typed boundary postprocessing finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed boundary postprocessing finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastBoundaryPostProcFE();
 
  /**
   * @brief Get typed domain evaluation finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed domain evaluation finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastEvaluationFE();
 
  /**
   * @brief Get typed skeleton left-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed skeleton LHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastSkeletonLhsFE();
 
  /**
   * @brief Get typed skeleton right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed skeleton RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastSkeletonRhsFE();
 
  /**
   * @brief Get typed domain explicit right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed domain explicit RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastDomainExplicitRhsFE();
 
  /**
   * @brief Get typed boundary explicit right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed boundary explicit RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastBoundaryExplicitRhsFE();
 
  /**
   * @brief Get typed skeleton explicit right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @return auto Shared pointer to typed skeleton explicit RHS finite element
   */
  template <typename T = ForcesAndSourcesCore, int DIM = -1>
  inline auto getCastSkeletonExplicitRhsFE();
 
  /**
   * @brief Get typed meshset right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @return auto Shared pointer to typed meshset RHS finite element
   */
  template <typename T = ForcesAndSourcesCore>
  inline auto getCastMeshsetRhsFE();
 
  /**
   * @brief Get typed meshset left-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @return auto Shared pointer to typed meshset LHS finite element
   */
  template <typename T = ForcesAndSourcesCore>
  inline auto getCastMeshsetLhsFE();
 
  /**
   * @brief Get typed meshset explicit right-hand side finite element
   *
   * @tparam T Target finite element type (default: ForcesAndSourcesCore)
   * @return auto Shared pointer to typed meshset explicit RHS finite element
   */
  template <typename T = ForcesAndSourcesCore>
  inline auto getCastMeshsetExplicitRhsFE();
 
  /**
   * @brief Set integration rule for domain left-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for domain right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for boundary left-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for boundary right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for domain evaluation finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setEvaluationIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for skeleton left-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setSkeletonLhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for skeleton right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setSkeletonRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for domain explicit right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setDomainExplicitRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for boundary explicit right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setBoundaryExplicitRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Set integration rule for skeleton explicit right-hand side finite element
   *
   * @tparam DIM Dimension (-1 for automatic detection from problem)
   * @param rule Integration rule hook function
   * @return MoFEMErrorCode Error code (0 on success)
   */
  template <int DIM = -1>
  inline MoFEMErrorCode setSkeletonExplicitRhsIntegrationRule(RuleHookFun rule);
 
  /**
   * @brief Get the Op Domain Lhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpDomainLhsPipeline();
 
  /**
   * @brief Get the Op Domain Rhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpDomainRhsPipeline();
 
  /**
   * @brief Get the Op Boundary Lhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpBoundaryLhsPipeline();
 
  /**
   * @brief Get the Op Boundary Rhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpBoundaryRhsPipeline();
 
  /**
   * @brief Get the Op Domain PostProc Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpDomainPostProcPipeline();
 
  /**
   * @brief Get the Op Boundary PostProc Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpBoundaryPostProcPipeline();
 
  /**
   * @brief Get the Op Domain Evaluation Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpEvaluationPipeline();
 
  /**
   * @brief Get the Op Skeleton Lhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpSkeletonLhsPipeline();
 
  /**
   * @brief Get the Op Skeleton Rhs Pipeline object
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpSkeletonRhsPipeline();
 
  /**
   * @brief Get the Op Domain Rhs Pipeline object for implicit-explicit G term
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpDomainExplicitRhsPipeline();
 
  /**
   * @brief Get the Op Boundary Rhs Pipeline object for implicit-explicit G term
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpBoundaryExplicitRhsPipeline();
 
  /**
   * @brief Get the Op Skeleton Rhs Pipeline object for implicit-explicit G term
   * @ingroup mofem_basic_interface
   *
   * @tparam -1
   * @return boost::ptr_deque<UserDataOperator>&
   */
  template <int DIM = -1>
  inline boost::ptr_deque<UserDataOperator> &getOpSkeletonExplicitRhsPipeline();
 
  /**
   * @brief Get the Op Meshset Rhs Pipeline object
   * 
   * @return boost::ptr_deque<UserDataOperator>& 
   */
  boost::ptr_deque<UserDataOperator> &getOpMeshsetRhsPipeline();
 
  /**
   * @brief Get the Op Meshset Lhs Pipeline object
   * 
   * @return boost::ptr_deque<UserDataOperator>& 
   */
  boost::ptr_deque<UserDataOperator> &getOpMeshsetLhsPipeline();
  
  /**
   * @brief Get the Op Meshset Explicit Rhs Pipeline object
   * 
   * @return boost::ptr_deque<UserDataOperator>& 
   */
  boost::ptr_deque<UserDataOperator> &getOpMeshsetExplicitRhsPipeline();
 
  /**
   * @brief Iterate finite elements
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode loopFiniteElements(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Iterate postprocessing finite elements
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode loopFiniteElementsPostProc(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Iterate domain postprocessing finite elements only
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode loopDomainPostProc(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Iterate boundary postprocessing finite elements only
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode loopBoundaryPostProc(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Iterate evaluation finite elements
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return MoFEMErrorCode Error code (0 on success)
   */
  MoFEMErrorCode loopFiniteElementsEvaluation(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Create KSP (linear) solver
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<KSP> Smart pointer to PETSc KSP linear solver
   */
  SmartPetscObj<KSP> createKSP(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Create SNES (nonlinear) solver
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<SNES> Smart pointer to PETSc SNES nonlinear solver
   */
  SmartPetscObj<SNES> createSNES(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Enumeration of time solver types
   */
  enum TSType { 
    EX,    ///< Explicit time integration
    IM,    ///< Implicit time integration
    IM2,   ///< Second-order implicit time integration
    IMEX   ///< Implicit-explicit time integration
  };
 
  /**
   * @brief Create TS (time) solver with specified type
   *
   * @param type Type of time solver (EX/IM/IM2/IMEX)
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<TS> Smart pointer to PETSc TS solver
   */
  SmartPetscObj<TS> createTS(const TSType type, SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Create TS (time) explicit solver
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<TS> Smart pointer to PETSc TS explicit time solver
   */
  SmartPetscObj<TS> createTSEX(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @brief Create TS (time) implicit solver
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<TS> Smart pointer to PETSc TS implicit time solver
   */
  SmartPetscObj<TS> createTSIM(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @deprecated  Use version with explicit TS solver type
   */
  DEPRECATED auto createTS(SmartPetscObj<DM> dm = nullptr) {
    return createTSIM(dm);
  }
 
  /**
   * @brief Create TS (time) solver for second order equation in time
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<TS> Smart pointer to PETSc TS second-order time solver
   */
  SmartPetscObj<TS> createTSIM2(SmartPetscObj<DM> dm = nullptr);
 
  /**
   * @deprecated  Change name. Use createTSIM2 instead.
   */
  inline DEPRECATED auto createTS2(SmartPetscObj<DM> dm = nullptr) {
    return createTSIM2(dm);
  }
 
  /**
   * @brief Create TS (time) implicit-explicit solver
   * @ingroup mofem_basic_interface
   *
   * @param dm Optional DM object (default: nullptr)
   * @return SmartPetscObj<TS> Smart pointer to PETSc TS implicit-explicit time solver
   */
  SmartPetscObj<TS> createTSIMEX(SmartPetscObj<DM> dm = nullptr);
 
private:
  MoFEM::Core &cOre; ///< Reference to MoFEM core instance
 
  boost::shared_ptr<FEMethod>
      feDomainRhs; ///< Element to assemble RHS side by integrating domain
  boost::shared_ptr<FEMethod>
      feDomainLhs; ///< Element to assemble LHS side by integrating domain
  boost::shared_ptr<FEMethod>
      feBoundaryRhs; ///< Element to assemble RHS side by integrating boundary
  boost::shared_ptr<FEMethod>
      feDomainPostProc; ///< Element for postprocessing by integrating domain
  boost::shared_ptr<FEMethod>
      feBoundaryPostProc; ///< Element for postprocessing by integrating boundary
  boost::shared_ptr<FEMethod>
      feEvaluation; ///< Element for evaluation by integrating domain
  boost::shared_ptr<FEMethod>
      feBoundaryLhs; ///< Element to assemble LHS side by integrating boundary
  boost::shared_ptr<FEMethod>
      feSkeletonRhs; ///< Element to assemble RHS side by integrating skeleton
  boost::shared_ptr<FEMethod>
      feSkeletonLhs; ///< Element to assemble LHS side by integrating skeleton
 
  boost::shared_ptr<FEMethod>
      feDomainExplicitRhs; ///< Element to assemble explicit Rhs for IMEX solver
  boost::shared_ptr<FEMethod>
      feBoundaryExplicitRhs; ///< Element to assemble explicit Rhs for IMEX
                             ///< solver
  boost::shared_ptr<FEMethod>
      feSkeletonExplicitRhs; ///< Element to assemble explicit Rhs for IMEX
                             ///< solver
 
  boost::shared_ptr<FEMethod> feMeshsetRhs; ///< Element to assemble RHS side by
                                            ///< integrating meshset
  boost::shared_ptr<FEMethod> feMeshsetLhs; ///< Element to assemble LHS side by
                                            ///< integrating meshset
  boost::shared_ptr<FEMethod> feMeshsetExplicitRhs; ///< Element to assemble
                                                    ///< explicit RHS side by
                                                    ///< integrating meshset
 
  /**
   * @brief Create domain finite element pipeline based on dimension
   *
   * @tparam DIM Problem dimension (1, 2, 3, or -1 for automatic detection)
   * @param fe Reference to finite element method shared pointer
   * @return boost::shared_ptr<FEMethod>& Reference to created finite element
   */
  template <int DIM>
  inline boost::shared_ptr<FEMethod> &
  createDomainFEPipeline(boost::shared_ptr<FEMethod> &fe);
 
  /**
   * @brief Create boundary finite element pipeline based on dimension
   *
   * @tparam DIM Problem dimension (1, 2, 3, or -1 for automatic detection)
   * @param fe Reference to finite element method shared pointer
   * @return boost::shared_ptr<FEMethod>& Reference to created finite element
   */
  template <int DIM>
  inline boost::shared_ptr<FEMethod> &
  createBoundaryFEPipeline(boost::shared_ptr<FEMethod> &fe);
 
  /**
   * @brief Create meshset finite element pipeline
   *
   * @param fe Reference to finite element method shared pointer
   * @return boost::shared_ptr<FEMethod>& Reference to created finite element
   */
  inline boost::shared_ptr<FEMethod> &
  createMeshsetFEPipeline(boost::shared_ptr<FEMethod> &fe);
 
private:
  /**
   * @brief Forward declaration for meshset finite element implementation
   */
  struct MeshsetFE;
};
 
template <int DIM>
boost::shared_ptr<FEMethod> &
PipelineManager::createDomainFEPipeline(boost::shared_ptr<FEMethod> &fe) {
  static_assert(DIM == 1 || DIM == 2 || DIM == 3, "not implemented");
  fe = boost::make_shared<FEMethod>();
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createDomainFEPipeline<3>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<VolEle>(cOre);
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createDomainFEPipeline<2>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<FaceEle>(cOre);
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createDomainFEPipeline<1>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<EdgeEle>(cOre);
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createDomainFEPipeline<-1>(boost::shared_ptr<FEMethod> &fe) {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return createDomainFEPipeline<1>(fe);
  case 2:
    return createDomainFEPipeline<2>(fe);
  case 3:
    return createDomainFEPipeline<3>(fe);
  default:
    THROW_MESSAGE("Not implemented");
  }
}
 
template <int DIM>
boost::shared_ptr<FEMethod> &
PipelineManager::createBoundaryFEPipeline(boost::shared_ptr<FEMethod> &fe) {
  static_assert(DIM == 1 || DIM == 2 || DIM == 3, "not implemented");
  fe = boost::make_shared<FEMethod>();
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createBoundaryFEPipeline<3>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<FaceElementForcesAndSourcesCore>(cOre);
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createBoundaryFEPipeline<2>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<EdgeEle>(cOre);
  return fe;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createBoundaryFEPipeline<1>(boost::shared_ptr<FEMethod> &fe) {
  if (!fe)
    fe = boost::make_shared<VertexElementForcesAndSourcesCore>(cOre);
  return fe;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getDomainLhsFE() {
  return feDomainLhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getDomainRhsFE() {
  return feDomainRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getBoundaryLhsFE() {
  return feBoundaryLhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getBoundaryRhsFE() {
  return feBoundaryRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getDomainPostProcFE() {
  return feDomainPostProc;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getBoundaryPostProcFE() {
  return feBoundaryPostProc;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getEvaluationFE() {
  return feEvaluation;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getSkeletonLhsFE() {
  return feSkeletonLhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getSkeletonRhsFE() {
  return feSkeletonRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getDomainExplicitRhsFE() {
  return feDomainExplicitRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getBoundaryExplicitRhsFE() {
  return feBoundaryExplicitRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getSkeletonExplicitRhsFE() {
  return feSkeletonExplicitRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getMeshsetRhsFE() {
  return feMeshsetRhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getMeshsetLhsFE() {
  return feMeshsetLhs;
}
 
boost::shared_ptr<FEMethod> &PipelineManager::getMeshsetExplicitRhsFE() {
  return feMeshsetExplicitRhs;
}
 
template <>
inline boost::shared_ptr<FEMethod> &
PipelineManager::createBoundaryFEPipeline<-1>(boost::shared_ptr<FEMethod> &fe) {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return createBoundaryFEPipeline<1>(fe);
  case 2:
    return createBoundaryFEPipeline<2>(fe);
  case 3:
    return createBoundaryFEPipeline<3>(fe);
  default:
    THROW_MESSAGE("Not implemented");
  }
}
 
template <typename T, int DIM> auto PipelineManager::getCastDomainLhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createDomainFEPipeline<DIM>(feDomainLhs));
}
 
template <typename T, int DIM> auto PipelineManager::getCastDomainRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createDomainFEPipeline<DIM>(feDomainRhs));
}
 
template <typename T, int DIM> auto PipelineManager::getCastBoundaryLhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feBoundaryLhs));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastBoundaryRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feBoundaryRhs));
}
 
template <typename T, int DIM> auto PipelineManager::getCastDomainPostProcFE() {
  return boost::dynamic_pointer_cast<T>(
      createDomainFEPipeline<DIM>(feDomainPostProc));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastBoundaryPostProcFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feBoundaryPostProc));
}
 
template <typename T, int DIM> auto PipelineManager::getCastEvaluationFE() {
  return boost::dynamic_pointer_cast<T>(
      createDomainFEPipeline<DIM>(feEvaluation));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastSkeletonLhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feSkeletonLhs));
}
 
template <typename T, int DIM> auto PipelineManager::getCastSkeletonRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feSkeletonRhs));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastDomainExplicitRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createDomainFEPipeline<DIM>(feDomainExplicitRhs));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastBoundaryExplicitRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feBoundaryExplicitRhs));
}
 
template <typename T, int DIM>
auto PipelineManager::getCastSkeletonExplicitRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createBoundaryFEPipeline<DIM>(feSkeletonExplicitRhs));
}
 
template <typename T> auto PipelineManager::getCastMeshsetRhsFE() {
  return boost::dynamic_pointer_cast<T>(createMeshsetFEPipeline(feMeshsetRhs));
}
 
template <typename T> auto PipelineManager::getCastMeshsetLhsFE() {
  return boost::dynamic_pointer_cast<T>(createMeshsetFEPipeline(feMeshsetLhs));
}
 
template <typename T> auto PipelineManager::getCastMeshsetExplicitRhsFE() {
  return boost::dynamic_pointer_cast<T>(
      createMeshsetFEPipeline(feMeshsetExplicitRhs));
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setDomainLhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createDomainFEPipeline<DIM>(feDomainLhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setDomainLhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBeginHot;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setDomainLhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setDomainLhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setDomainLhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturnHot(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setDomainRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createDomainFEPipeline<DIM>(feDomainRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setDomainRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setDomainRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setDomainRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setDomainRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setBoundaryLhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feBoundaryLhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setBoundaryLhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setBoundaryLhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setBoundaryLhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setBoundaryLhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setBoundaryRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feBoundaryRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setBoundaryRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setBoundaryRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setBoundaryRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setBoundaryRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setEvaluationIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createDomainFEPipeline<DIM>(feEvaluation))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setEvaluationIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setEvaluationIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setEvaluationIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setEvaluationIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setSkeletonLhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feSkeletonLhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setSkeletonLhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setSkeletonLhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setSkeletonLhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setSkeletonLhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setSkeletonRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feSkeletonRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setSkeletonRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setSkeletonRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setSkeletonRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setSkeletonRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setDomainExplicitRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createDomainFEPipeline<DIM>(feDomainExplicitRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode PipelineManager::setDomainExplicitRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setDomainExplicitRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setDomainExplicitRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setDomainExplicitRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setBoundaryExplicitRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feBoundaryExplicitRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode
PipelineManager::setBoundaryExplicitRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setBoundaryExplicitRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setBoundaryExplicitRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setBoundaryExplicitRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
MoFEMErrorCode PipelineManager::setSkeletonExplicitRhsIntegrationRule(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
      createBoundaryFEPipeline<DIM>(feSkeletonExplicitRhs))
      ->getRuleHook = rule;
  MoFEMFunctionReturn(0);
}
 
template <>
inline MoFEMErrorCode
PipelineManager::setSkeletonExplicitRhsIntegrationRule<-1>(
    PipelineManager::RuleHookFun rule) {
  MoFEMFunctionBegin;
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    MoFEMFunctionReturnHot(setSkeletonExplicitRhsIntegrationRule<1>(rule));
  case 2:
    MoFEMFunctionReturnHot(setSkeletonExplicitRhsIntegrationRule<2>(rule));
  case 3:
    MoFEMFunctionReturnHot(setSkeletonExplicitRhsIntegrationRule<3>(rule));
  default:
    THROW_MESSAGE("Not implemented");
  }
  MoFEMFunctionReturn(0);
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainLhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createDomainFEPipeline<DIM>(feDomainLhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainLhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpDomainLhsPipeline<1>();
  case 2:
    return getOpDomainLhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpDomainLhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createDomainFEPipeline<DIM>(feDomainRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpDomainRhsPipeline<1>();
  case 2:
    return getOpDomainRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpDomainRhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryLhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feBoundaryLhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryLhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpBoundaryLhsPipeline<1>();
  case 2:
    return getOpBoundaryLhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpBoundaryLhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feBoundaryRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpBoundaryRhsPipeline<1>();
  case 2:
    return getOpBoundaryRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpBoundaryRhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainPostProcPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createDomainFEPipeline<DIM>(feDomainPostProc))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainPostProcPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpDomainPostProcPipeline<1>();
  case 2:
    return getOpDomainPostProcPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpDomainPostProcPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryPostProcPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feBoundaryPostProc))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryPostProcPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpBoundaryPostProcPipeline<1>();
  case 2:
    return getOpBoundaryPostProcPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpBoundaryPostProcPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpEvaluationPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createDomainFEPipeline<DIM>(feEvaluation))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpEvaluationPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpEvaluationPipeline<1>();
  case 2:
    return getOpEvaluationPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpEvaluationPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonLhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feSkeletonLhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonLhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpSkeletonLhsPipeline<1>();
  case 2:
    return getOpSkeletonLhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpSkeletonLhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feSkeletonRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpSkeletonRhsPipeline<1>();
  case 2:
    return getOpSkeletonRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpSkeletonRhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainExplicitRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createDomainFEPipeline<DIM>(feDomainExplicitRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpDomainExplicitRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpDomainExplicitRhsPipeline<1>();
  case 2:
    return getOpDomainExplicitRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpDomainExplicitRhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryExplicitRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feBoundaryExplicitRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpBoundaryExplicitRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpBoundaryExplicitRhsPipeline<1>();
  case 2:
    return getOpBoundaryExplicitRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpBoundaryExplicitRhsPipeline<3>();
}
 
template <int DIM>
boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonExplicitRhsPipeline() {
  return boost::dynamic_pointer_cast<ForcesAndSourcesCore>(
             createBoundaryFEPipeline<DIM>(feSkeletonExplicitRhs))
      ->getOpPtrVector();
}
 
template <>
inline boost::ptr_deque<PipelineManager::UserDataOperator> &
PipelineManager::getOpSkeletonExplicitRhsPipeline<-1>() {
  switch (cOre.getInterface<Simple>()->getDim()) {
  case 1:
    return getOpSkeletonExplicitRhsPipeline<1>();
  case 2:
    return getOpSkeletonExplicitRhsPipeline<2>();
  case 3:
    break;
  default:
    THROW_MESSAGE("Not implemented");
  }
  return getOpSkeletonExplicitRhsPipeline<3>();
}
 
/**
 * @brief Specialization for 2D problems
 */
template <> struct PipelineManager::ElementsAndOpsByDim<2> {
  using DomainEle = PipelineManager::FaceEle;        ///< 2D domain elements (faces)
  using BoundaryEle = PipelineManager::EdgeEle;      ///< 2D boundary elements (edges)
  using FaceSideEle = FaceElementForcesAndSourcesCoreOnSide;  ///< Face side elements
  using DomainParentEle = FaceElementForcesAndSourcesCoreOnChildParent;  ///< Domain parent-child elements
  using BoundaryParentEle = EdgeElementForcesAndSourcesCoreOnChildParent;  ///< Boundary parent-child elements
 
  using DomianParentEle = DomainParentEle; //< @deprecated Use DomainParentEle instead
};
 
/**
 * @brief Specialization for 3D problems
 */
template <> struct PipelineManager::ElementsAndOpsByDim<3> {
  using DomainEle = PipelineManager::VolEle;         ///< 3D domain elements (volumes)
  using BoundaryEle = PipelineManager::FaceEle;      ///< 3D boundary elements (faces)
  using FaceSideEle = VolumeElementForcesAndSourcesCoreOnSide;  ///< Volume side elements
  // !!! Not yet implemented
  // using DomainParentEle = VolumeElementForcesAndSourcesCoreOnChildParent;
  using BoundaryParentEle = FaceElementForcesAndSourcesCoreOnChildParent;  ///< Boundary parent-child elements
};
 
} // namespace MoFEM
 
#endif // __BASIC_HPP__
 
/**
 * \defgroup mofem_basic_interface PipelineManager interface
 * \brief Implementation of basic interface for rapid problem implementation.
 *
 * \ingroup mofem
 **/