text: // ---------------------------------------------------------------------
  //
  // Copyright (c) 2014 - 2019 by the IBAMR developers
  // All rights reserved.
  //
  // This file is part of IBAMR.
  //
  // IBAMR is free software and is distributed under the 3-clause BSD
  // license. The full text of the license can be found in the file
  // COPYRIGHT at the top level directory of IBAMR.
  //
  // ---------------------------------------------------------------------
  
  #ifndef included_IBAMR_StaggeredStokesFACPreconditionerStrategy
  #define included_IBAMR_StaggeredStokesFACPreconditionerStrategy
  
  /////////////////////////////// INCLUDES /////////////////////////////////////
  
  #include "ibamr/StaggeredStokesPhysicalBoundaryHelper.h"
  
  #include "ibtk/CartCellRobinPhysBdryOp.h"
  #include "ibtk/CartSideRobinPhysBdryOp.h"
  #include "ibtk/CoarseFineBoundaryRefinePatchStrategy.h"
  #include "ibtk/FACPreconditionerStrategy.h"
  #include "ibtk/ibtk_utilities.h"
  
  #include "CoarsenAlgorithm.h"
  #include "CoarsenOperator.h"
  #include "IntVector.h"
  #include "PatchHierarchy.h"
  #include "PoissonSpecifications.h"
  #include "RefineAlgorithm.h"
  #include "RefineOperator.h"
  #include "RefinePatchStrategy.h"
  #include "SAMRAIVectorReal.h"
  #include "VariableContext.h"
  #include "VariableFillPattern.h"
  #include "tbox/Database.h"
  #include "tbox/Pointer.h"
  
  #include <string>
  #include <utility>
  #include <vector>
  
  namespace IBAMR
  {
  class StaggeredStokesSolver;
  } // namespace IBAMR
  namespace IBTK
  {
  class HierarchyGhostCellInterpolation;
  class HierarchyMathOps;
  } // namespace IBTK
  namespace SAMRAI
  {
  namespace solv
  {
  template <int DIM>
  class RobinBcCoefStrategy;
  } // namespace solv
  namespace xfer
  {
  template <int DIM>
  class CoarsenSchedule;
  template <int DIM>
  class RefineSchedule;
  } // namespace xfer
  } // namespace SAMRAI
  
  /////////////////////////////// CLASS DEFINITION /////////////////////////////
  
  namespace IBAMR
  {
  /*!
   * \brief Class StaggeredStokesFACPreconditionerStrategy is an abstract
   * FACPreconditionerStrategy implementing many of the operations required by
   * smoothers for staggered-grid (MAC) discretizations of the incompressible
   * Stokes equations and related problems.
   *
   * Sample parameters for initialization from database (and their default
   * values): \verbatim
  
   smoother_type = "ADDITIVE"                     // see setSmootherType()
   U_prolongation_method = "CONSTANT_REFINE"      // see setProlongationMethods()
   P_prolongation_method = "LINEAR_REFINE"        // see setProlongationMethods()
   U_restriction_method = "CONSERVATIVE_COARSEN"  // see setRestrictionMethods()
   P_restriction_method = "CONSERVATIVE_COARSEN"  // see setRestrictionMethods()
   coarse_solver_type = "LEVEL_SMOOTHER"          // see setCoarseSolverType()
   coarse_solver_rel_residual_tol = 1.0e-5        // see setCoarseSolverRelativeTolerance()
   coarse_solver_abs_residual_tol = 1.0e-50       // see setCoarseSolverAbsoluteTolerance()
   coarse_solver_max_iterations = 10              // see setCoarseSolverMaxIterations()
   coarse_solver_db = { ... }                     // SAMRAI::tbox::Database for initializing
   coarse
   level solver
   \endverbatim
  */
  class StaggeredStokesFACPreconditionerStrategy : public IBTK::FACPreconditionerStrategy
  {
  public:
      /*!
       * \brief Constructor.
       */
      StaggeredStokesFACPreconditionerStrategy(const std::string& object_name,
                                               int ghost_cell_width,
                                               SAMRAI::tbox::Pointer<SAMRAI::tbox::Database> input_db,
                                               const std::string& default_options_prefix);
  
      /*!
       * \brief Destructor.
       */
      ~StaggeredStokesFACPreconditionerStrategy();
  
      /*!
       * \brief Set the PoissonSpecifications object used to specify the
       * coefficients for the momentum equation in the incompressible Stokes
       * operator.
       */
      virtual void setVelocityPoissonSpecifications(const SAMRAI::solv::PoissonSpecifications& U_problem_coefs);
  
      /*!
       * \brief Set if velocity and pressure have nullspace.
       */
      virtual void setComponentsHaveNullspace(const bool has_velocity_nullspace, const bool has_pressure_nullspace);
  
      /*!
       * \brief Set the SAMRAI::solv::RobinBcCoefStrategy objects used to specify
       * physical boundary conditions.
       *
       * \note Any of the elements of \a U_bc_coefs may be NULL.  In this case,
       * homogeneous Dirichlet boundary conditions are employed for that data
       * depth.  \a P_bc_coef may also be NULL; in that case, homogeneous Neumann
       * boundary conditions are employed for the pressure.
       *
       * \param U_bc_coefs  IBTK::Vector of pointers to objects that can set the Robin boundary
       *condition coefficients for the velocity
       * \param P_bc_coef   Pointer to object that can set the Robin boundary condition
       *coefficients
       *for the pressure
       */
      virtual void setPhysicalBcCoefs(const std::vector<SAMRAI::solv::RobinBcCoefStrategy<NDIM>*>& U_bc_coefs,
                                      SAMRAI::solv::RobinBcCoefStrategy<NDIM>* P_bc_coef);
  
      /*!
       * \brief Set the StokesSpecifications object and timestep size used to specify
       * the coefficients for the time-dependent incompressible Stokes operator.
       */
      virtual void setPhysicalBoundaryHelper(SAMRAI::tbox::Pointer<StaggeredStokesPhysicalBoundaryHelper> bc_helper);
  
      /*!
       * \name Functions for configuring the solver.
       */
      //\{
  
      /*!
       * \brief Specify the levels that need to be reset the next time the
       * operator is re-initialized.
       *
       * When the operator is initialized, then only the specified range of levels
       * are reset in the operator state the next time that the operator is
       * initialized.  If the operator is not initialized, this method has no
       * effect.
       *
       * To ensure the range of levels that is reset includes all levels in the
       * patch hierarchy, use \a coarsest_ln = \a finest_ln = \p -1.
       *
       * \note This function is used to save some unnecessary computations when
       * the hierarchy is regridded.  The range of levels specified must include
       * all levels which need to be reset by
       * SAMRAI::mesh::StandardTagAndInitStrategy::resetHierarchyConfiguration().
       * Any data residing outside of this range of levels will not be reset.
       * This \b is \b not what you want to have happen if, for instance, the
       * Poisson specifications changes.
       */
      void setResetLevels(int coarsest_ln, int finest_ln);
  
      /*!
       * \brief Specify the smoother type.
       */
      void setSmootherType(const std::string& smoother_type);
  
      /*!
       * \brief Specify the coarse level solver.
       */
      void setCoarseSolverType(const std::string& coarse_solver_type);
  
      /*!
       * \brief Set the maximum number of iterations for the coarse level solve.
       *
       * If the coarse level solver uses a maximum number of iterations parameter,
       * the specified value is used.  If the coarse level solver does not use
       * such a stopping parameter, implementations are free to ignore this value.
       */
      void setCoarseSolverMaxIterations(int coarse_solver_max_iterations);
  
      /*!
       * \brief Set the absolute residual tolerance for convergence for coarse
       * level solve.
       *
       * If the coarse level solver uses a absolute convergence tolerance
       * parameter, the specified value is used.  If the coarse level solver does
       * not use such a stopping parameter, implementations are free to ignore
       * this value.
       */
      void setCoarseSolverAbsoluteTolerance(double coarse_solver_abs_residual_tol);
  
      /*!
       * \brief Set the relative residual tolerance for convergence for coarse
       * level solve.
       *
       * If the coarse level solver uses a relative convergence tolerance
       * parameter, the specified value is used.  If the coarse level solver does
       * not use such a stopping parameter, implementations are free to ignore
       * this value.
       */
      void setCoarseSolverRelativeTolerance(double coarse_solver_rel_residual_tol);
  
      /*!
       * \brief Set the prolongation methods.
       */
      void setProlongationMethods(const std::string& U_prolongation_method, const std::string& P_prolongation_method);
  
      /*!
       * \brief Set the restriction methods.
       */
      void setRestrictionMethods(const std::string& U_restriction_method, const std::string& P_restriction_method);
  
      //\}
  
      /*!
       * \name Partial implementation of FACPreconditionerStrategy interface.
       */
      //\{
  
      /*!
       * \brief Zero the supplied vector.
       */
      void setToZero(SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& vec, int level_num) override;
  
      /*!
       * \brief Restrict the residual quantity to the specified level from the
       * next finer level.
       *
       * \param src source residual
       * \param dst destination residual
       * \param dst_ln destination level number
       */
      void restrictResidual(const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& src,
                            SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& dst,
                            int dst_ln) override;
  
      /*!
       * \brief Prolong the error quantity to the specified level from the next
       * coarser level.
       *
       * \param src source error vector
       * \param dst destination error vector
       * \param dst_ln destination level number of data transfer
       */
      void prolongError(const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& src,
                        SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& dst,
                        int dst_ln) override;
  
      /*!
       * \brief Prolong the error quantity to the specified level from the next
       * coarser level and apply the correction to the fine-level error.
       *
       * \param src source error vector
       * \param dst destination error vector
       * \param dst_ln destination level number of data transfer
       */
      void prolongErrorAndCorrect(const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& src,
                                  SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& dst,
                                  int dst_ln) override;
  
      /*!
       * \brief Solve the residual equation Ae=r on the coarsest level of the
       * patch hierarchy.
       *
       * \param error error vector
       * \param residual residual vector
       * \param coarsest_ln coarsest level number
       */
      bool solveCoarsestLevel(SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& error,
                              const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& residual,
                              int coarsest_ln) override;
  
      /*!
       * \brief Compute the composite-grid residual on the specified range of
       * levels of the patch hierarchy.
       */
      void computeResidual(SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& residual,
                           const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& solution,
                           const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& rhs,
                           int coarsest_level_num,
                           int finest_level_num) override;
  
      /*!
       * \brief Compute hierarchy-dependent data.
       *
       * Note that although the vector arguments given to other methods in this
       * class may not necessarily be the same as those given to this method,
       * there will be similarities, including:
       *
       * - hierarchy configuration (hierarchy pointer and level range)
       * - number, type and alignment of vector component data
       * - ghost cell width of data in the solution (or solution-like) vector
       *
       * \param solution solution vector u
       * \param rhs right hand side vector f
       */
      void initializeOperatorState(const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& solution,
                                   const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& rhs) override;
  
      /*!
       * \brief Remove all hierarchy-dependent data.
       *
       * Remove all hierarchy-dependent data set by initializeOperatorState().
       *
       * \see initializeOperatorState
       */
      void deallocateOperatorState() override;
  
      //\}
  
  protected:
      /*!
       * \brief Compute implementation-specific hierarchy-dependent data.
       */
      virtual void initializeOperatorStateSpecialized(const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& solution,
                                                      const SAMRAI::solv::SAMRAIVectorReal<NDIM, double>& rhs,
                                                      int coarsest_reset_ln,
                                                      int finest_reset_ln) = 0;
  
      /*!
       * \brief Remove implementation-specific hierarchy-dependent data.
       */
      virtual void deallocateOperatorStateSpecialized(int coarsest_reset_ln, int finest_reset_ln) = 0;
  
      /*!
       * \name Methods for executing, caching, and resetting communication
       * schedules.
       */
      //\{
  
      /*!
       * \brief Execute a refinement schedule for prolonging data.
       */
      void xeqScheduleProlongation(const std::pair<int, int>& dst_idxs, const std::pair<int, int>& src_idxs, int dst_ln);
  
      /*!
       * \brief Execute schedule for restricting solution or residual to the
       * specified level.
       */
      void xeqScheduleRestriction(const std::pair<int, int>& dst_idxs, const std::pair<int, int>& src_idxs, int dst_ln);
  
      /*!
       * \brief Execute schedule for filling ghosts on the specified level.
       */
      void xeqScheduleGhostFillNoCoarse(const std::pair<int, int>& dst_idxs, int dst_ln);
  
      /*!
       * \brief Execute schedule for synchronizing data on the specified level.
       */
      void xeqScheduleDataSynch(int dst_idx, int dst_ln);
  
      //\}
  
      /*
       * Problem specification.
       */
      SAMRAI::solv::PoissonSpecifications d_U_problem_coefs;
      SAMRAI::solv::RobinBcCoefStrategy<NDIM>* d_default_U_bc_coef;
      std::vector<SAMRAI::solv::RobinBcCoefStrategy<NDIM>*> d_U_bc_coefs;
      SAMRAI::solv::RobinBcCoefStrategy<NDIM>* d_default_P_bc_coef;
      SAMRAI::solv::RobinBcCoefStrategy<NDIM>* d_P_bc_coef;
  
      /*
       * Boundary condition helper object.
       */
      SAMRAI::tbox::Pointer<StaggeredStokesPhysicalBoundaryHelper> d_bc_helper;
  
      /*
       * Ghost cell width.
       */
      const int d_gcw;
  
      /*!
       * \name Hierarchy-dependent objects.
       */
      //\{
  
      /*
       * Solution and rhs vectors.
       */
      SAMRAI::tbox::Pointer<SAMRAI::solv::SAMRAIVectorReal<NDIM, double> > d_solution, d_rhs;
  
      /*
       * Reference patch hierarchy and range of levels involved in the solve.
       *
       * This variable is non-null between the initializeOperatorState() and
       * deallocateOperatorState() calls.  It is not truly needed, because the
       * hierarchy is obtainable through variables in most function argument
       * lists.  We use it to enforce working on one hierarchy at a time.
       */
      SAMRAI::tbox::Pointer<SAMRAI::hier::PatchHierarchy<NDIM> > d_hierarchy;
      int d_coarsest_ln = IBTK::invalid_level_number, d_finest_ln = IBTK::invalid_level_number;
  
      /*
       * Level operators, used to compute composite-grid residuals.
       */
      std::vector<SAMRAI::tbox::Pointer<IBTK::HierarchyGhostCellInterpolation> > d_level_bdry_fill_ops;
      std::vector<SAMRAI::tbox::Pointer<IBTK::HierarchyMathOps> > d_level_math_ops;
  
      /*
       * Range of levels to be reset the next time the operator is initialized.
       */
      bool d_in_initialize_operator_state = false;
      int d_coarsest_reset_ln = IBTK::invalid_level_number, d_finest_reset_ln = IBTK::invalid_level_number;
  
      //\}
  
      /*!
       * \name Solver configuration variables.
       */
      //\{
  
      /*
       * The kind of smoothing to perform.
       */
      std::string d_smoother_type = "ADDITIVE";
  
      /*
       * The names of the refinement operators used to prolong the coarse grid
       * correction.
       */
      std::string d_U_prolongation_method = "CONSTANT_REFINE", d_P_prolongation_method = "LINEAR_REFINE";
  
      /*
       * The names of the coarsening operators used to restrict the fine grid
       * error or residual.
       */
      std::string d_U_restriction_method = "CONSERVATIVE_COARSEN", d_P_restriction_method = "CONSERVATIVE_COARSEN";
  
      /*
       * Coarse level solver parameters.
       */
      bool d_coarse_solver_init_subclass = false;
      std::string d_coarse_solver_type = "LEVEL_SMOOTHER", d_coarse_solver_default_options_prefix;
      double d_coarse_solver_rel_residual_tol = 1.0e-5;
      double d_coarse_solver_abs_residual_tol = 1.0e-50;
      int d_coarse_solver_max_iterations = 10;
      SAMRAI::tbox::Pointer<IBAMR::StaggeredStokesSolver> d_coarse_solver;
      SAMRAI::tbox::Pointer<SAMRAI::tbox::Database> d_coarse_solver_db;
  
      /*
       * Nullspace info.
       */
      bool d_has_velocity_nullspace = false, d_has_pressure_nullspace = false;
  
      //\}
  
      /*!
       * \name Internal context and scratch data.
       */
      //\{
  
      /*
       * Variable context for internally maintained hierarchy data.
       */
      SAMRAI::tbox::Pointer<SAMRAI::hier::VariableContext> d_context;
  
      /*
       * Patch descriptor indices for scratch data.
       */
      int d_side_scratch_idx = IBTK::invalid_index, d_cell_scratch_idx = IBTK::invalid_index;
  
      //\}
  
      /*!
       * \name Various refine and coarsen objects.
       */
      //\{
  
      /*
       * Physical boundary operators.
       */
      SAMRAI::tbox::Pointer<IBTK::CartSideRobinPhysBdryOp> d_U_bc_op;
      SAMRAI::tbox::Pointer<IBTK::CartCellRobinPhysBdryOp> d_P_bc_op;
  
      /*
       * Coarse-fine interface interpolation objects.
       */
      SAMRAI::tbox::Pointer<IBTK::CoarseFineBoundaryRefinePatchStrategy> d_U_cf_bdry_op, d_P_cf_bdry_op;
  
      /*
       * Variable fill pattern object.
       */
      SAMRAI::tbox::Pointer<SAMRAI::xfer::VariableFillPattern<NDIM> > d_U_op_stencil_fill_pattern,
          d_P_op_stencil_fill_pattern, d_U_synch_fill_pattern;
  
      //\}
  
  private:
      /*!
       * \brief Default constructor.
       *
       * \note This constructor is not implemented and should not be used.
       */
      StaggeredStokesFACPreconditionerStrategy() = delete;
  
      /*!
       * \brief Copy constructor.
       *
       * \note This constructor is not implemented and should not be used.
       *
       * \param from The value to copy to this object.
       */
      StaggeredStokesFACPreconditionerStrategy(const StaggeredStokesFACPreconditionerStrategy& from) = delete;
  
      /*!
       * \brief Assignment operator.
       *
       * \note This operator is not implemented and should not be used.
       *
       * \param that The value to assign to this object.
       *
       * \return A reference to this object.
       */
      StaggeredStokesFACPreconditionerStrategy& operator=(const StaggeredStokesFACPreconditionerStrategy& that) = delete;
  
      /*
       * Combined U & P physical boundary operator.
       */
      SAMRAI::xfer::RefinePatchStrategy<NDIM>* d_U_P_bc_op;
  
      /*
       * Error prolongation (refinement) operator.
       */
      SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineOperator<NDIM> > d_U_prolongation_refine_operator,
          d_P_prolongation_refine_operator;
      SAMRAI::tbox::Pointer<SAMRAI::xfer::RefinePatchStrategy<NDIM> > d_prolongation_refine_patch_strategy;
      SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineAlgorithm<NDIM> > d_prolongation_refine_algorithm;
      std::vector<SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineSchedule<NDIM> > > d_prolongation_refine_schedules;
  
      /*
       * Residual restriction (coarsening) operator.
       */
      SAMRAI::tbox::Pointer<SAMRAI::xfer::CoarsenOperator<NDIM> > d_U_restriction_coarsen_operator,
          d_P_restriction_coarsen_operator;
      SAMRAI::tbox::Pointer<SAMRAI::xfer::CoarsenAlgorithm<NDIM> > d_restriction_coarsen_algorithm;
      std::vector<SAMRAI::tbox::Pointer<SAMRAI::xfer::CoarsenSchedule<NDIM> > > d_restriction_coarsen_schedules;
  
      /*
       * Refine operator for side and cell data from same level.
       */
      SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineAlgorithm<NDIM> > d_ghostfill_nocoarse_refine_algorithm;
      std::vector<SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineSchedule<NDIM> > > d_ghostfill_nocoarse_refine_schedules;
  
      /*
       * Operator for side data synchronization on same level.
       */
      SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineAlgorithm<NDIM> > d_synch_refine_algorithm;
      std::vector<SAMRAI::tbox::Pointer<SAMRAI::xfer::RefineSchedule<NDIM> > > d_synch_refine_schedules;
  };
  } // namespace IBAMR
  
  //////////////////////////////////////////////////////////////////////////////
  
  #endif //#ifndef included_IBAMR_StaggeredStokesFACPreconditionerStrategy
  