|
IBAMR
IBAMR version 0.19.
|
Class CCPoissonPointRelaxationFACOperator is a concrete PoissonFACPreconditionerStrategy for solving elliptic equations of the form \( \mbox{$L u$} = \mbox{$(C I + \nabla \cdot D \nabla) u$} = f \) using a globally second-order accurate cell-centered finite-volume discretization, with support for Robin and periodic boundary conditions. More...
#include <ibtk/CCPoissonPointRelaxationFACOperator.h>

Public Member Functions | |
| CCPoissonPointRelaxationFACOperator (const std::string &object_name, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, std::string default_options_prefix) | |
| Constructor. More... | |
| ~CCPoissonPointRelaxationFACOperator () | |
| Destructor. More... | |
| 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. More... | |
| 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... | |
| const std::string & | getName () const |
| Return the object name. More... | |
| virtual bool | getIsInitialized () const |
| Return whether the operator is initialized. More... | |
| virtual void | setFACPreconditioner (SAMRAI::tbox::ConstPointer< FACPreconditioner > preconditioner) |
| Method to allow the FACPreconditioner object to register itself with the concrete FACPreconditionerStrategy. More... | |
| virtual void | setHomogeneousBc (bool homogeneous_bc) |
| Set whether the solver should use homogeneous boundary conditions. More... | |
| virtual bool | getHomogeneousBc () const |
| Return whether the solver is using homogeneous boundary conditions. More... | |
| virtual void | setSolutionTime (double solution_time) |
| Set the time at which the solution is to be evaluated. More... | |
| virtual double | getSolutionTime () const |
| Get the time at which the solution is being evaluated. More... | |
| virtual void | setTimeInterval (double current_time, double new_time) |
| Set the current time interval. More... | |
| virtual std::pair< double, double > | getTimeInterval () const |
| Get the current time interval. More... | |
| virtual double | getDt () const |
| Get the current time step size. More... | |
| virtual void | allocateScratchData () |
| Allocate scratch data. More... | |
| virtual void | deallocateScratchData () |
| Deallocate scratch data. More... | |
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. More... | |
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 CCPoissonPointRelaxationFACOperator FAC strategy. More... | |
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< SAMRAI::hier::BoxList< NDIM > > > | d_patch_bc_box_overlap |
| std::vector< std::vector< std::map< int, SAMRAI::hier::Box< NDIM > > > > | d_patch_neighbor_overlap |
| std::string | d_data_refine_type = "NONE" |
| bool | d_use_cf_interpolation = true |
| std::string | d_data_coarsen_type = "CUBIC_COARSEN" |
| std::string | d_bdry_extrap_type = "LINEAR" |
| bool | d_use_consistent_type_2_bdry = false |
| 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 a 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. More... | |
| void | deallocateOperatorStateSpecialized (int coarsest_reset_ln, int finest_reset_ln) override |
| Remove implementation-specific hierarchy-dependent data. More... | |
| CCPoissonPointRelaxationFACOperator ()=delete | |
| Default constructor. More... | |
| CCPoissonPointRelaxationFACOperator (const CCPoissonPointRelaxationFACOperator &from)=delete | |
| Copy constructor. More... | |
| CCPoissonPointRelaxationFACOperator & | operator= (const CCPoissonPointRelaxationFACOperator &that)=delete |
| Assignment operator. More... | |
Logging functions. | |
| virtual void | printClassData (std::ostream &stream) |
| Print class data to stream. More... | |
| 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. More... | |
| 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() |
Partial implementation of FACPreconditionerStrategy interface. | |
| void | setToZero (SAMRAI::solv::SAMRAIVectorReal< NDIM, double > &vec, int level_num) override |
| Zero the supplied vector. More... | |
| 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... | |
Methods for executing, caching, and resetting communication | |
| void | xeqScheduleProlongation (int dst_idx, int src_idx, int dst_ln) |
| Execute a refinement schedule for prolonging data. More... | |
| void | xeqScheduleRestriction (int dst_idx, int src_idx, int dst_ln) |
| Execute schedule for restricting solution or residual to the specified level. More... | |
| void | xeqScheduleGhostFillNoCoarse (int dst_idx, int dst_ln) |
| Execute schedule for filling ghosts on the specified level. More... | |
| void | xeqScheduleDataSynch (int dst_idx, int dst_ln) |
| Execute schedule for synchronizing data on the specified level. More... | |
| 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 |
This class provides operators that are used by class FACPreconditioner to solve scalar Poisson-type equations of the form
\[ (C I + \nabla \cdot D \nabla) u = f \]
using a cell-centered, globally second-order accurate finite-volume 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 \( -\nabla^2 u = f \), 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 = "LINEAR_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
solver_type = "PFMG"
num_pre_relax_steps = 0
num_post_relax_steps = 2
}
| IBTK::CCPoissonPointRelaxationFACOperator::CCPoissonPointRelaxationFACOperator | ( | const std::string & | object_name, |
| SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > | input_db, | ||
| std::string | default_options_prefix | ||
| ) |
| IBTK::CCPoissonPointRelaxationFACOperator::~CCPoissonPointRelaxationFACOperator | ( | ) |
|
privatedelete |
|
privatedelete |
| from | The value to copy to this object. |
|
inlinestatic |
|
overridevirtual |
Select from:
"PATCH_GAUSS_SEIDEL" "PROCESSOR_GAUSS_SEIDEL" "RED_BLACK_GAUSS_SEIDEL" Implements IBTK::PoissonFACPreconditionerStrategy.
|
overridevirtual |
Implements IBTK::PoissonFACPreconditionerStrategy.
|
overridevirtual |
| 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.
|
overridevirtual |
| error | error vector |
| residual | residual vector |
| coarsest_ln | coarsest level number |
Implements IBTK::FACPreconditionerStrategy.
|
overridevirtual |
| 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.
|
overrideprotectedvirtual |
Implements IBTK::PoissonFACPreconditionerStrategy.
|
overrideprotectedvirtual |
Implements IBTK::PoissonFACPreconditionerStrategy.
|
privatedelete |
| that | The value to assign to this object. |
|
virtualinherited |
|
virtualinherited |
| bc_coef | Pointer to an object that can set the Robin boundary condition coefficients |
|
virtualinherited |
| bc_coefs | Vector of pointers to objects that can set the Robin boundary condition coefficients |
|
inherited |
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 coarsest_ln = finest_ln = -1.
|
inherited |
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.
|
inherited |
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.
|
inherited |
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.
|
inherited |
|
inherited |
|
overridevirtualinherited |
Implements IBTK::FACPreconditionerStrategy.
|
overridevirtualinherited |
| src | source residual |
| dst | destination residual |
| dst_ln | destination level number |
Implements IBTK::FACPreconditionerStrategy.
Reimplemented in IBTK::VCSCViscousOpPointRelaxationFACOperator.
|
overridevirtualinherited |
| src | source error vector |
| dst | destination error vector |
| dst_ln | destination level number of data transfer |
Implements IBTK::FACPreconditionerStrategy.
|
overridevirtualinherited |
| src | source error vector |
| dst | destination error vector |
| dst_ln | destination level number of data transfer |
Implements IBTK::FACPreconditionerStrategy.
|
overridevirtualinherited |
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:
An unrecoverable error will occur if the rhs vector does not have consistent ghost cell width as that provided in the constructor.
| solution | solution vector u |
| rhs | right hand side vector f |
Reimplemented from IBTK::FACPreconditionerStrategy.
|
overridevirtualinherited |
Remove all hierarchy-dependent data set by initializeOperatorState().
Reimplemented from IBTK::FACPreconditionerStrategy.
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
inherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
virtualinherited |
|
protectedvirtualinherited |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
privateinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
|
protectedinherited |
1.8.17