IBAMR
An adaptive and distributed-memory parallel implementation of the immersed boundary (IB) method
|
Class SCPoissonPointRelaxationFACOperator is a concrete PoissonFACPreconditionerStrategy for solving elliptic equations of the form using a globally second-order accurate side-centered finite-difference discretization, with support for Robin and periodic boundary conditions. More...
#include </home/runner/work/IBAMR/IBAMR/ibtk/include/ibtk/SCPoissonPointRelaxationFACOperator.h>
Public Member Functions | |
SCPoissonPointRelaxationFACOperator (const std::string &object_name, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, std::string default_options_prefix) | |
Constructor. | |
~SCPoissonPointRelaxationFACOperator () | |
Destructor. | |
Functions for configuring the solver. | |
void | setSmootherType (const std::string &smoother_type) override |
Specify the smoother type. More... | |
void | setCoarseSolverType (const std::string &coarse_solver_type) override |
Specify the coarse level solver. | |
Public Member Functions inherited from IBTK::PoissonFACPreconditionerStrategy | |
PoissonFACPreconditionerStrategy (std::string object_name, SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > scratch_var, int ghost_cell_width, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, const std::string &default_options_prefix) | |
Constructor. | |
~PoissonFACPreconditionerStrategy () | |
Destructor. | |
virtual void | setPoissonSpecifications (const SAMRAI::solv::PoissonSpecifications &poisson_spec) |
Set the SAMRAI::solv::PoissonSpecifications object used to specify the coefficients for the scalar-valued or vector-valued Laplace operator. | |
virtual void | setPhysicalBcCoef (SAMRAI::solv::RobinBcCoefStrategy< NDIM > *bc_coef) |
Set the SAMRAI::solv::RobinBcCoefStrategy object used to specify physical boundary conditions. More... | |
virtual void | setPhysicalBcCoefs (const std::vector< SAMRAI::solv::RobinBcCoefStrategy< NDIM > * > &bc_coefs) |
Set the SAMRAI::solv::RobinBcCoefStrategy objects used to specify physical boundary conditions. More... | |
void | setResetLevels (int coarsest_ln, int finest_ln) |
Specify the levels that need to be reset the next time the operator is re-initialized. More... | |
void | setCoarseSolverMaxIterations (int coarse_solver_max_iterations) |
Set the maximum number of iterations for the coarse level solve. More... | |
void | setCoarseSolverAbsoluteTolerance (double coarse_solver_abs_residual_tol) |
Set the absolute residual tolerance for convergence for coarse level solve. More... | |
void | setCoarseSolverRelativeTolerance (double coarse_solver_rel_residual_tol) |
Set the relative residual tolerance for convergence for coarse level solve. More... | |
void | setProlongationMethod (const std::string &prolongation_method) |
Set the name of the prolongation method. | |
void | setRestrictionMethod (const std::string &restriction_method) |
Set the name of the restriction method. | |
void | setToZero (SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &vec, int level_num) override |
Zero the supplied vector. | |
void | restrictResidual (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &src, SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &dst, int dst_ln) override |
Restrict the residual quantity to the specified level from the next finer level. More... | |
void | prolongError (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &src, SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &dst, int dst_ln) override |
Prolong the error quantity to the specified level from the next coarser level. More... | |
void | prolongErrorAndCorrect (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &src, SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &dst, int dst_ln) override |
Prolong the error quantity to the specified level from the next coarser level and apply the correction to the fine-level error. More... | |
void | initializeOperatorState (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &solution, const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &rhs) override |
Compute hierarchy-dependent data. More... | |
void | deallocateOperatorState () override |
Remove all hierarchy-dependent data. More... | |
Public Member Functions inherited from IBTK::FACPreconditionerStrategy | |
FACPreconditionerStrategy (std::string object_name, bool homogeneous_bc=false) | |
Constructor. | |
virtual | ~FACPreconditionerStrategy ()=default |
Empty virtual desctructor. | |
const std::string & | getName () const |
Return the object name. | |
virtual bool | getIsInitialized () const |
Return whether the operator is initialized. | |
virtual void | setFACPreconditioner (SAMRAI::tbox::ConstPointer< FACPreconditioner > preconditioner) |
Method to allow the FACPreconditioner object to register itself with the concrete FACPreconditionerStrategy. | |
virtual void | setHomogeneousBc (bool homogeneous_bc) |
Set whether the solver should use homogeneous boundary conditions. | |
virtual bool | getHomogeneousBc () const |
Return whether the solver is using homogeneous boundary conditions. | |
virtual void | setSolutionTime (double solution_time) |
Set the time at which the solution is to be evaluated. | |
virtual double | getSolutionTime () const |
Get the time at which the solution is being evaluated. | |
virtual void | setTimeInterval (double current_time, double new_time) |
Set the current time interval. | |
virtual std::pair< double, double > | getTimeInterval () const |
Get the current time interval. | |
virtual double | getDt () const |
Get the current time step size. | |
virtual void | allocateScratchData () |
Allocate scratch data. | |
virtual void | deallocateScratchData () |
Deallocate scratch data. | |
virtual void | printClassData (std::ostream &stream) |
Print class data to stream. | |
Static Public Member Functions | |
static SAMRAI::tbox::Pointer< PoissonSolver > | allocate_solver (const std::string &object_name, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, const std::string &default_options_prefix) |
Static function to construct a PoissonFACPreconditioner with a SCPoissonPointRelaxationFACOperator FAC strategy. | |
Implementation of FACPreconditionerStrategy interface. | |
SAMRAI::tbox::Pointer< PoissonSolver > | d_coarse_solver |
SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > | d_coarse_solver_db |
std::vector< std::vector< std::array< SAMRAI::hier::BoxList< NDIM >, NDIM > > > | d_patch_bc_box_overlap |
std::vector< std::vector< std::array< std::map< int, SAMRAI::hier::Box< NDIM > >, NDIM > > > | d_patch_neighbor_overlap |
SAMRAI::tbox::Pointer< StaggeredPhysicalBoundaryHelper > | d_bc_helper |
int | d_mask_idx |
void | smoothError (SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &error, const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &residual, int level_num, int num_sweeps, bool performing_pre_sweeps, bool performing_post_sweeps) override |
Perform a given number of relaxations on the error. More... | |
bool | solveCoarsestLevel (SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &error, const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &residual, int coarsest_ln) override |
Solve the residual equation Ae=r on the coarsest level of the patch hierarchy. More... | |
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 |
Compute composite grid residual on the specified range of levels. More... | |
void | initializeOperatorStateSpecialized (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &solution, const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &rhs, int coarsest_reset_ln, int finest_reset_ln) override |
Compute implementation-specific hierarchy-dependent data. | |
void | deallocateOperatorStateSpecialized (int coarsest_reset_ln, int finest_reset_ln) override |
Remove implementation-specific hierarchy-dependent data. | |
Additional Inherited Members | |
Protected Member Functions inherited from IBTK::PoissonFACPreconditionerStrategy | |
void | xeqScheduleProlongation (int dst_idx, int src_idx, int dst_ln) |
Execute a refinement schedule for prolonging data. | |
void | xeqScheduleRestriction (int dst_idx, int src_idx, int dst_ln) |
Execute schedule for restricting solution or residual to the specified level. | |
void | xeqScheduleGhostFillNoCoarse (int dst_idx, int dst_ln) |
Execute schedule for filling ghosts on the specified level. | |
void | xeqScheduleDataSynch (int dst_idx, int dst_ln) |
Execute schedule for synchronizing data on the specified level. | |
Protected Member Functions inherited from IBTK::FACPreconditionerStrategy | |
virtual SAMRAI::tbox::Pointer< SAMRAI::solv::SAMRAIVectorReal< NDIM, double > > | getLevelSAMRAIVectorReal (const SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &vec, int level_num) const |
Return a SAMRAIVectorReal object that corresponds to the given object but restricted to a single level of the patch hierarchy. | |
Protected Attributes inherited from IBTK::PoissonFACPreconditionerStrategy | |
SAMRAI::tbox::Pointer< SAMRAI::solv::SAMRAIVectorReal< NDIM, double > > | d_solution |
SAMRAI::tbox::Pointer< SAMRAI::solv::SAMRAIVectorReal< NDIM, double > > | d_rhs |
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > | d_hierarchy |
int | d_coarsest_ln = IBTK::invalid_level_number |
int | d_finest_ln = IBTK::invalid_level_number |
std::vector< SAMRAI::tbox::Pointer< SAMRAI::math::HierarchyDataOpsReal< NDIM, double > > > | d_level_data_ops |
std::vector< SAMRAI::tbox::Pointer< IBTK::HierarchyGhostCellInterpolation > > | d_level_bdry_fill_ops |
std::vector< SAMRAI::tbox::Pointer< IBTK::HierarchyMathOps > > | d_level_math_ops |
bool | d_in_initialize_operator_state = false |
int | d_coarsest_reset_ln = IBTK::invalid_level_number |
int | d_finest_reset_ln = IBTK::invalid_level_number |
std::string | d_smoother_type = "DEFAULT" |
std::string | d_prolongation_method = "DEFAULT" |
std::string | d_restriction_method = "DEFAULT" |
std::string | d_coarse_solver_type = "DEFAULT" |
std::string | 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< SAMRAI::hier::VariableContext > | d_context |
int | d_scratch_idx = IBTK::invalid_index |
SAMRAI::solv::PoissonSpecifications | d_poisson_spec |
std::unique_ptr< SAMRAI::solv::RobinBcCoefStrategy< NDIM > > | d_default_bc_coef |
std::vector< SAMRAI::solv::RobinBcCoefStrategy< NDIM > * > | d_bc_coefs |
SAMRAI::hier::IntVector< NDIM > | d_gcw |
SAMRAI::tbox::Pointer< RobinPhysBdryPatchStrategy > | d_bc_op |
SAMRAI::tbox::Pointer< CoarseFineBoundaryRefinePatchStrategy > | d_cf_bdry_op |
SAMRAI::tbox::Pointer< SAMRAI::xfer::VariableFillPattern< NDIM > > | d_op_stencil_fill_pattern |
SAMRAI::tbox::Pointer< SAMRAI::xfer::VariableFillPattern< NDIM > > | d_synch_fill_pattern |
Protected Attributes inherited from IBTK::FACPreconditionerStrategy | |
SAMRAI::tbox::ConstPointer< IBTK::FACPreconditioner > | d_preconditioner |
const std::string | d_object_name |
bool | d_is_initialized = false |
bool | d_homogeneous_bc |
double | d_solution_time = std::numeric_limits<double>::quiet_NaN() |
double | d_current_time = std::numeric_limits<double>::quiet_NaN() |
double | d_new_time = std::numeric_limits<double>::quiet_NaN() |
Class SCPoissonPointRelaxationFACOperator is a concrete PoissonFACPreconditionerStrategy for solving elliptic equations of the form using a globally second-order accurate side-centered finite-difference discretization, with support for Robin and periodic boundary conditions.
This class provides operators that are used by class FACPreconditioner to solve scalar Poisson-type equations of the form
using a side-centered, globally second-order accurate finite-difference discretization, where
Robin boundary conditions may be specified at physical boundaries; see class SAMRAI::solv::RobinBcCoefStrategy.
By default, the class is configured to solve the Poisson problem , subject to homogeneous Dirichlet boundary conditions.
Sample parameters for initialization from database (and their default values):
smoother_type = "PATCH_GAUSS_SEIDEL" // see setSmootherType() prolongation_method = "CONSTANT_REFINE" // see setProlongationMethod() restriction_method = "CONSERVATIVE_COARSEN" // see setRestrictionMethod() coarse_solver_type = "HYPRE_LEVEL_SOLVER" // 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 = 1 // see setCoarseSolverMaxIterations() coarse_solver_db = { ... } // SAMRAI::tbox::Database for initializing coarse level solver
|
overridevirtual |
Compute composite grid residual on the specified range of levels.
residual | residual vector |
solution | solution vector |
rhs | source (right hand side) vector |
coarsest_level_num | coarsest level number |
finest_level_num | finest level number |
Implements IBTK::FACPreconditionerStrategy.
Reimplemented in IBTK::VCSCViscousOpPointRelaxationFACOperator.
|
overridevirtual |
Specify the smoother type.
Select from:
"PATCH_GAUSS_SEIDEL"
"PROCESSOR_GAUSS_SEIDEL"
"RED_BLACK_GAUSS_SEIDEL"
Implements IBTK::PoissonFACPreconditionerStrategy.
|
overridevirtual |
Perform a given number of relaxations on the error.
error | error vector |
residual | residual vector |
level_num | level number |
num_sweeps | number of sweeps to perform |
performing_pre_sweeps | boolean value that is true when pre-smoothing sweeps are being performed |
performing_post_sweeps | boolean value that is true when post-smoothing sweeps are being performed |
Implements IBTK::FACPreconditionerStrategy.
Reimplemented in IBTK::VCSCViscousOpPointRelaxationFACOperator.
|
overridevirtual |
Solve the residual equation Ae=r on the coarsest level of the patch hierarchy.
error | error vector |
residual | residual vector |
coarsest_ln | coarsest level number |
Implements IBTK::FACPreconditionerStrategy.