IBAMR
An adaptive and distributed-memory parallel implementation of the immersed boundary (IB) method
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
IBAMR::IBHierarchyIntegrator Class Reference

Class IBHierarchyIntegrator provides an abstract interface for a time integrator for various versions of the immersed boundary method on an AMR grid hierarchy, along with basic data management for variables defined on that hierarchy. More...

#include </home/runner/work/IBAMR/IBAMR/include/ibamr/IBHierarchyIntegrator.h>

Inheritance diagram for IBAMR::IBHierarchyIntegrator:
Inheritance graph
[legend]

Classes

class  IBEulerianForceFunction
 A class to communicate the Eulerian body force computed by class IBHierarchyIntegrator to the incompressible Navier-Stokes solver. More...
 
class  IBEulerianSourceFunction
 A class to communicate the Eulerian fluid source-sink distribution computed by class IBHierarchyIntegrator to the incompressible Navier-Stokes solver. More...
 

Public Member Functions

 ~IBHierarchyIntegrator ()=default
 
TimeSteppingType getTimeSteppingType () const
 
SAMRAI::tbox::Pointer< IBStrategygetIBStrategy () const
 
void registerBodyForceFunction (SAMRAI::tbox::Pointer< IBTK::CartGridFunction > F_fcn)
 
virtual void registerLoadBalancer (SAMRAI::tbox::Pointer< SAMRAI::mesh::LoadBalancer< NDIM > > load_balancer) override
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > getVelocityVariable () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > getPressureVariable () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > getBodyForceVariable () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > getFluidSourceVariable () const
 
IBTK::RobinPhysBdryPatchStrategygetVelocityPhysBdryOp () const
 
void preprocessIntegrateHierarchy (double current_time, double new_time, int num_cycles=1) override
 
void postprocessIntegrateHierarchy (double current_time, double new_time, bool skip_synchronize_new_state_data, int num_cycles=1) override
 
void initializeHierarchyIntegrator (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > > gridding_alg) override
 
void initializePatchHierarchy (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > > gridding_alg) override
 
- Public Member Functions inherited from IBTK::HierarchyIntegrator
 HierarchyIntegrator (std::string object_name, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, bool register_for_restart)
 
 ~HierarchyIntegrator ()
 
const std::stringgetName () const
 
virtual void advanceHierarchy (double dt)
 
double getMinimumTimeStepSize ()
 
double getMaximumTimeStepSize ()
 
void synchronizeHierarchyData (VariableContextType ctx_type)
 
void resetTimeDependentHierarchyData (double new_time)
 
void resetIntegratorToPreadvanceState ()
 
virtual void regridHierarchy ()
 
bool atRegridPoint () const
 
double getIntegratorTime () const
 
double getStartTime () const
 
double getEndTime () const
 
int getIntegratorStep () const
 
int getMaxIntegratorSteps () const
 
bool stepsRemaining () const
 
void updateWorkloadEstimates ()
 
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > getPatchHierarchy () const
 
SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > > getGriddingAlgorithm () const
 
int getWorkloadDataIndex () const
 
void registerVisItDataWriter (SAMRAI::tbox::Pointer< SAMRAI::appu::VisItDataWriter< NDIM > > visit_writer)
 
SAMRAI::tbox::Pointer< SAMRAI::appu::VisItDataWriter< NDIM > > getVisItDataWriter () const
 
void setupPlotData ()
 
virtual int getNumberOfCycles () const
 
virtual int getCurrentCycleNumber () const
 
virtual double getCurrentTimeStepSize () const
 
void integrateHierarchy (double current_time, double new_time, int cycle_num=0)
 
void skipCycle (double current_time, double new_time, int cycle_num=0)
 
void registerPreprocessIntegrateHierarchyCallback (PreprocessIntegrateHierarchyCallbackFcnPtr callback, void *ctx=nullptr)
 
void registerIntegrateHierarchyCallback (IntegrateHierarchyCallbackFcnPtr callback, void *ctx=nullptr)
 
void registerPostprocessIntegrateHierarchyCallback (PostprocessIntegrateHierarchyCallbackFcnPtr callback, void *ctx=nullptr)
 
void registerApplyGradientDetectorCallback (ApplyGradientDetectorCallbackFcnPtr callback, void *ctx=nullptr)
 
void registerRegridHierarchyCallback (RegridHierarchyCallbackFcnPtr, void *ctx=nullptr)
 
void initializeCompositeHierarchyData (double init_data_time, bool initial_time)
 
void initializeLevelData (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double init_data_time, bool can_be_refined, bool initial_time, SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchLevel< NDIM > > old_level=SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchLevel< NDIM > >(NULL), bool allocate_data=true) override
 
void resetHierarchyConfiguration (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int coarsest_level, int finest_level) override
 
void applyGradientDetector (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double error_data_time, int tag_index, bool initial_time, bool uses_richardson_extrapolation_too) override
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextgetContext (VariableContextType ctx_type) const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextgetCurrentContext () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextgetNewContext () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextgetScratchContext () const
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextgetPlotContext () const
 
bool isAllocatedPatchData (int data_idx, int coarsest_ln=invalid_level_number, int finest_ln=invalid_level_number) const
 
void allocatePatchData (int data_idx, double data_time, int coarsest_ln=invalid_level_number, int finest_ln=invalid_level_number) const
 
void deallocatePatchData (int data_idx, int coarsest_ln=invalid_level_number, int finest_ln=invalid_level_number) const
 
SAMRAI::tbox::Pointer< HierarchyMathOpsgetHierarchyMathOps () const
 
void registerVariable (int &current_idx, int &new_idx, int &scratch_idx, SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > variable, const SAMRAI::hier::IntVector< NDIM > &scratch_ghosts=SAMRAI::hier::IntVector< NDIM >(0), const std::string &coarsen_name="NO_COARSEN", const std::string &refine_name="NO_REFINE", SAMRAI::tbox::Pointer< CartGridFunction > init_fcn=SAMRAI::tbox::Pointer< CartGridFunction >(NULL), const bool register_for_restart=true)
 
void registerVariable (int &idx, SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > variable, const SAMRAI::hier::IntVector< NDIM > &ghosts=SAMRAI::hier::IntVector< NDIM >(0), SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContext > ctx=SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContext >(NULL), const bool register_for_restart=true)
 
void putToDatabase (SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > db) override
 
- Public Member Functions inherited from SAMRAI::mesh::StandardTagAndInitStrategy< NDIM >
virtual double getLevelDt (const tbox::Pointer< hier::BasePatchLevel< NDIM > > level, const double dt_time, const bool initial_time)
 
virtual double advanceLevel (const tbox::Pointer< hier::BasePatchLevel< NDIM > > level, const tbox::Pointer< hier::BasePatchHierarchy< NDIM > > hierarchy, const double current_time, const double new_time, const bool first_step, const bool last_step, const bool regrid_advance=false)
 
virtual void resetTimeDependentData (const tbox::Pointer< hier::BasePatchLevel< NDIM > > level, const double new_time, const bool can_be_refined)
 
virtual void resetDataToPreadvanceState (const tbox::Pointer< hier::BasePatchLevel< NDIM > > level)
 
virtual void applyRichardsonExtrapolation (const tbox::Pointer< hier::PatchLevel< NDIM > > level, const double error_data_time, const int tag_index, const double deltat, const int error_coarsen_ratio, const bool initial_time, const bool uses_gradient_detector_too)
 
virtual void coarsenDataForRichardsonExtrapolation (const tbox::Pointer< hier::PatchHierarchy< NDIM > > hierarchy, const int level_number, const tbox::Pointer< hier::PatchLevel< NDIM > > coarser_level, const double coarsen_data_time, const bool before_advance)
 

Protected Member Functions

void regridHierarchyBeginSpecialized () override
 
void regridHierarchyEndSpecialized () override
 
 IBHierarchyIntegrator (const std::string &object_name, SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db, SAMRAI::tbox::Pointer< IBStrategy > ib_method_ops, SAMRAI::tbox::Pointer< INSHierarchyIntegrator > ins_hier_integrator, bool register_for_restart=true)
 
bool atRegridPointSpecialized () const override
 
void initializeLevelDataSpecialized (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double init_data_time, bool can_be_refined, bool initial_time, SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchLevel< NDIM > > old_level, bool allocate_data) override
 
void resetHierarchyConfigurationSpecialized (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int coarsest_level, int finest_level) override
 
void applyGradientDetectorSpecialized (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double error_data_time, int tag_index, bool initial_time, bool uses_richardson_extrapolation_too) override
 
void putToDatabaseSpecialized (SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > db) override
 
virtual void addWorkloadEstimate (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, const int workload_data_idx) override
 
- Protected Member Functions inherited from IBTK::HierarchyIntegrator
virtual void integrateHierarchySpecialized (double current_time, double new_time, int cycle_num=0)=0
 
virtual double getMinimumTimeStepSizeSpecialized ()
 
virtual double getMaximumTimeStepSizeSpecialized ()
 
virtual void synchronizeHierarchyDataSpecialized (VariableContextType ctx_type)
 
virtual void resetTimeDependentHierarchyDataSpecialized (double new_time)
 
virtual void resetIntegratorToPreadvanceStateSpecialized ()
 
virtual void setupPlotDataSpecialized ()
 
virtual void initializeCompositeHierarchyDataSpecialized (double init_data_time, bool initial_time)
 
virtual void executePreprocessIntegrateHierarchyCallbackFcns (double current_time, double new_time, int num_cycles)
 
virtual void executeIntegrateHierarchyCallbackFcns (double current_time, double new_time, int cycle_num)
 
virtual void executePostprocessIntegrateHierarchyCallbackFcns (double current_time, double new_time, bool skip_synchronize_new_state_data, int num_cycles)
 
virtual void executeApplyGradientDetectorCallbackFcns (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double error_data_time, int tag_index, bool initial_time, bool uses_richardson_extrapolation_too)
 
virtual void executeRegridHierarchyCallbackFcns (SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, double data_time, bool initial_time)
 
void registerGhostfillRefineAlgorithm (const std::string &name, SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > ghostfill_alg, std::unique_ptr< SAMRAI::xfer::RefinePatchStrategy< NDIM > > ghostfill_patch_strategy=nullptr)
 
void registerProlongRefineAlgorithm (const std::string &name, SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > prolong_alg, std::unique_ptr< SAMRAI::xfer::RefinePatchStrategy< NDIM > > prolong_patch_strategy=nullptr)
 
void registerCoarsenAlgorithm (const std::string &name, SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenAlgorithm< NDIM > > coarsen_alg, std::unique_ptr< SAMRAI::xfer::CoarsenPatchStrategy< NDIM > > coarsen_patch_strategy=nullptr)
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > getGhostfillRefineAlgorithm (const std::string &name) const
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > getProlongRefineAlgorithm (const std::string &name) const
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenAlgorithm< NDIM > > getCoarsenAlgorithm (const std::string &name) const
 
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > & getGhostfillRefineSchedules (const std::string &name) const
 
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > & getProlongRefineSchedules (const std::string &name) const
 
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenSchedule< NDIM > > > & getCoarsenSchedules (const std::string &name) const
 
void registerChildHierarchyIntegrator (HierarchyIntegrator *child_integrator)
 
void registerParentHierarchyIntegrator (HierarchyIntegrator *parent_integrator)
 
SAMRAI::tbox::Pointer< HierarchyMathOpsbuildHierarchyMathOps (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy)
 
void setupTagBuffer (SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > > gridding_alg)
 
bool regriddingHierarchy () const
 
bool atRegridTimeStep () const
 

Protected Attributes

bool d_integrator_is_initialized = false
 
TimeSteppingType d_time_stepping_type = MIDPOINT_RULE
 
bool d_use_structure_predictor
 
bool d_error_on_dt_change = true
 
bool d_warn_on_dt_change = false
 
SAMRAI::tbox::Pointer< INSHierarchyIntegratord_ins_hier_integrator
 
double d_regrid_fluid_cfl_interval = -1.0
 
double d_regrid_structure_cfl_interval = -1.0
 
double d_regrid_fluid_cfl_estimate = 0.0
 
double d_regrid_structure_cfl_estimate = 0.0
 
SAMRAI::tbox::Pointer< IBStrategyd_ib_method_ops
 
SAMRAI::tbox::Pointer< SAMRAI::math::HierarchyDataOpsReal< NDIM, double > > d_hier_velocity_data_ops
 
SAMRAI::tbox::Pointer< SAMRAI::math::HierarchyDataOpsReal< NDIM, double > > d_hier_pressure_data_ops
 
SAMRAI::tbox::Pointer< SAMRAI::math::HierarchyCellDataOpsReal< NDIM, double > > d_hier_cc_data_ops
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > d_u_var
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > d_p_var
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > d_f_var
 
SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > d_q_var
 
int d_u_idx = IBTK::invalid_index
 
int d_p_idx = IBTK::invalid_index
 
int d_f_idx = IBTK::invalid_index
 
int d_f_current_idx = IBTK::invalid_index
 
int d_q_idx = IBTK::invalid_index
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_ib_context
 
SAMRAI::hier::ComponentSelector d_ib_data
 
IBTK::RobinPhysBdryPatchStrategyd_u_phys_bdry_op
 
IBTK::RobinPhysBdryPatchStrategyd_p_phys_bdry_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > d_u_ghostfill_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > d_f_prolong_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > d_p_ghostfill_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineAlgorithm< NDIM > > d_q_prolong_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineOperator< NDIM > > d_u_ghostfill_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineOperator< NDIM > > d_f_prolong_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineOperator< NDIM > > d_p_ghostfill_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineOperator< NDIM > > d_q_prolong_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenAlgorithm< NDIM > > d_u_coarsen_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenAlgorithm< NDIM > > d_p_coarsen_alg
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenOperator< NDIM > > d_u_coarsen_op
 
SAMRAI::tbox::Pointer< SAMRAI::xfer::CoarsenOperator< NDIM > > d_p_coarsen_op
 
SAMRAI::tbox::Pointer< IBTK::CartGridFunctiond_body_force_fcn
 
- Protected Attributes inherited from IBTK::HierarchyIntegrator
std::string d_object_name
 
bool d_registered_for_restart = false
 
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > d_hierarchy
 
SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > > d_gridding_alg
 
SAMRAI::tbox::Pointer< SAMRAI::mesh::LoadBalancer< NDIM > > d_load_balancer
 
SAMRAI::tbox::Pointer< SAMRAI::pdat::CellVariable< NDIM, double > > d_workload_var
 
int d_workload_idx = IBTK::invalid_index
 
bool d_hierarchy_is_initialized = false
 
bool d_may_need_to_reset_hierarchy_configuration = false
 
HierarchyIntegratord_parent_integrator = nullptr
 
std::set< HierarchyIntegrator * > d_child_integrators
 
SAMRAI::tbox::Pointer< SAMRAI::appu::VisItDataWriter< NDIM > > d_visit_writer
 
double d_integrator_time = std::numeric_limits<double>::quiet_NaN()
 
double d_start_time = 0.0
 
double d_end_time = std::numeric_limits<double>::max()
 
double d_dt_init = std::numeric_limits<double>::max()
 
double d_dt_min = 0.0
 
double d_dt_max = std::numeric_limits<double>::max()
 
double d_dt_growth_factor = 2.0
 
int d_integrator_step = 0
 
int d_max_integrator_steps = std::numeric_limits<int>::max()
 
std::deque< double > d_dt_previous
 
int d_num_cycles = 1
 
int d_current_num_cycles = -1
 
int d_current_cycle_num = -1
 
double d_current_dt = std::numeric_limits<double>::quiet_NaN()
 
int d_regrid_interval = 1
 
RegridMode d_regrid_mode = STANDARD
 
bool d_enable_logging = false
 
bool d_enable_logging_solver_iterations = false
 
std::string d_bdry_extrap_type = "LINEAR"
 
SAMRAI::tbox::Array< int > d_tag_buffer = { 0 }
 
SAMRAI::tbox::Pointer< HierarchyMathOpsd_hier_math_ops
 
bool d_manage_hier_math_ops = true
 
std::list< SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > > d_state_variables
 
std::list< SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > > d_scratch_variables
 
std::list< SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > > d_copy_scratch_to_current_fast
 
std::list< SAMRAI::tbox::Pointer< SAMRAI::hier::Variable< NDIM > > > d_copy_scratch_to_current_slow
 
SAMRAI::hier::ComponentSelector d_current_data
 
SAMRAI::hier::ComponentSelector d_new_data
 
SAMRAI::hier::ComponentSelector d_scratch_data
 
SAMRAI::hier::ComponentSelector d_plot_data
 
std::map< SAMRAI::hier::Variable< NDIM > *, SAMRAI::tbox::Pointer< CartGridFunction > > d_state_var_init_fcns
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_current_context
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_new_context
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_scratch_context
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_plot_context
 
SAMRAI::hier::ComponentSelector d_fill_after_regrid_bc_idxs
 
SAMRAI::xfer::RefineAlgorithm< NDIM > d_fill_after_regrid_prolong_alg
 
std::unique_ptr< SAMRAI::xfer::RefinePatchStrategy< NDIM > > d_fill_after_regrid_phys_bdry_bc_op
 
std::vector< PreprocessIntegrateHierarchyCallbackFcnPtrd_preprocess_integrate_hierarchy_callbacks
 
std::vector< void * > d_preprocess_integrate_hierarchy_callback_ctxs
 
std::vector< IntegrateHierarchyCallbackFcnPtrd_integrate_hierarchy_callbacks
 
std::vector< void * > d_integrate_hierarchy_callback_ctxs
 
std::vector< PostprocessIntegrateHierarchyCallbackFcnPtrd_postprocess_integrate_hierarchy_callbacks
 
std::vector< void * > d_postprocess_integrate_hierarchy_callback_ctxs
 
std::vector< ApplyGradientDetectorCallbackFcnPtrd_apply_gradient_detector_callbacks
 
std::vector< void * > d_apply_gradient_detector_callback_ctxs
 
std::vector< RegridHierarchyCallbackFcnPtrd_regrid_hierarchy_callbacks
 
std::vector< void * > d_regrid_hierarchy_callback_ctxs
 

Friends

class IBStrategy
 
class IBEulerianForceFunction
 
class IBEulerianSourceFunction
 

Additional Inherited Members

- Public Types inherited from IBTK::HierarchyIntegrator
using PreprocessIntegrateHierarchyCallbackFcnPtr = void(*)(double current_time, double new_time, int num_cycles, void *ctx)
 
using IntegrateHierarchyCallbackFcnPtr = void(*)(double current_time, double new_time, int cycle_num, void *ctx)
 
using PostprocessIntegrateHierarchyCallbackFcnPtr = void(*)(double current_time, double new_time, bool skip_synchronize_new_state_data, int num_cycles, void *ctx)
 
using ApplyGradientDetectorCallbackFcnPtr = void(*)(SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, int level_number, double error_data_time, int tag_index, bool initial_time, bool uses_richardson_extrapolation_too, void *ctx)
 
using RegridHierarchyCallbackFcnPtr = void(*)(SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > > hierarchy, double data_time, bool initial_time, void *ctx)
 
- Static Protected Attributes inherited from IBTK::HierarchyIntegrator
static const std::string SYNCH_CURRENT_DATA_ALG = "SYNCH_CURRENT_DATA"
 
static const std::string SYNCH_NEW_DATA_ALG = "SYNCH_NEW_DATA"
 

Detailed Description

Class IBHierarchyIntegrator provides an abstract interface for a time integrator for various versions of the immersed boundary method on an AMR grid hierarchy, along with basic data management for variables defined on that hierarchy.

Options Controlling Regridding

Most IBAMR applications involve structure meshes defined on the finest level of the patch hierarchy. These structures move: in particular, they can potentially move off the finest grid level, causing the interaction routines to no longer work.

This class offers two different strategies for calculating how much the structure has moved (which it then uses to specify that a regrid is required):

  1. Estimation based on the fluid: The structure is assumed to move at the maximum velocity of the fluid at each time step. This class will regrid once a single fluid point has moved at least regrid_fluid_cfl_interval cell widths since the last regrid.

  2. Estimation based on the structure: The structure's displacement is calculated directly from its position vector and that value is used to determine the number of cell widths the structure has moved. In this case we will regrid once a point on the structure has moved at least regrid_structure_cfl_interval cell widths since the last regrid.

Both regrid_fluid_cfl_interval and regrid_structure_cfl_interval can be specified in the input database. For backwards compatibility the value regrid_cfl_interval is equivalent to regrid_fluid_cfl_interval. At the present time regrid_structure_cfl_interval is not implemented for all IBStrategy classes.

Alternatively, one can request that the solver (regardless of any computed displacement or velocity) regrid every time a fixed number of timesteps have occurred by specifying regrid_interval in the input database. If either regrid_structure_cfl_interval or regrid_fluid_cfl_interval are provided in the input database then regrid_interval is ignored.

Constructor & Destructor Documentation

◆ ~IBHierarchyIntegrator()

IBAMR::IBHierarchyIntegrator::~IBHierarchyIntegrator ( )
default

The destructor for class IBHierarchyIntegrator unregisters the integrator object with the restart manager when the object is so registered.

◆ IBHierarchyIntegrator()

IBAMR::IBHierarchyIntegrator::IBHierarchyIntegrator ( const std::string object_name,
SAMRAI::tbox::Pointer< SAMRAI::tbox::Database input_db,
SAMRAI::tbox::Pointer< IBStrategy ib_method_ops,
SAMRAI::tbox::Pointer< INSHierarchyIntegrator ins_hier_integrator,
bool  register_for_restart = true 
)
protected

The constructor for class IBHierarchyIntegrator sets some default values, reads in configuration information from input and restart databases, and registers the integrator object with the restart manager when requested.

Member Function Documentation

◆ addWorkloadEstimate()

void IBAMR::IBHierarchyIntegrator::addWorkloadEstimate ( SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
const int  workload_data_idx 
)
overrideprotectedvirtual

Add the work contributions (excluding the background grid) for the current hierarchy into the variable with index workload_data_idx. The only direct workload contribution of this hierarchy manager is usually the work done by the IBStrategy object.

Reimplemented from IBTK::HierarchyIntegrator.

◆ applyGradientDetectorSpecialized()

void IBAMR::IBHierarchyIntegrator::applyGradientDetectorSpecialized ( SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > >  hierarchy,
int  level_number,
double  error_data_time,
int  tag_index,
bool  initial_time,
bool  uses_richardson_extrapolation_too 
)
overrideprotectedvirtual

Set integer tags to "one" in cells where refinement of the given level should occur according to the magnitude of the fluid vorticity.

Reimplemented from IBTK::HierarchyIntegrator.

◆ atRegridPointSpecialized()

bool IBAMR::IBHierarchyIntegrator::atRegridPointSpecialized ( ) const
overrideprotectedvirtual

Function to determine whether regridding should occur at the current time step.

Reimplemented from IBTK::HierarchyIntegrator.

◆ getBodyForceVariable()

Pointer< Variable< NDIM > > IBAMR::IBHierarchyIntegrator::getBodyForceVariable ( ) const

Return a pointer to the body force variable.

◆ getFluidSourceVariable()

Pointer< Variable< NDIM > > IBAMR::IBHierarchyIntegrator::getFluidSourceVariable ( ) const

Return a pointer to the source strength variable.

◆ getIBStrategy()

Pointer< IBStrategy > IBAMR::IBHierarchyIntegrator::getIBStrategy ( ) const

Return a pointer to the IBStrategy object registered with this integrator.

◆ getPressureVariable()

Pointer< Variable< NDIM > > IBAMR::IBHierarchyIntegrator::getPressureVariable ( ) const

Return a pointer to the fluid pressure state variable.

◆ getTimeSteppingType()

TimeSteppingType IBAMR::IBHierarchyIntegrator::getTimeSteppingType ( ) const

Return the time stepping scheme.

◆ getVelocityPhysBdryOp()

IBTK::RobinPhysBdryPatchStrategy * IBAMR::IBHierarchyIntegrator::getVelocityPhysBdryOp ( ) const

Return a pointer to the velocity physical boundary conditions

◆ getVelocityVariable()

Pointer< Variable< NDIM > > IBAMR::IBHierarchyIntegrator::getVelocityVariable ( ) const

Return a pointer to the fluid velocity variable.

◆ initializeHierarchyIntegrator()

void IBAMR::IBHierarchyIntegrator::initializeHierarchyIntegrator ( SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > >  gridding_alg 
)
overridevirtual

Initialize the variables, basic communications algorithms, solvers, and other data structures used by this time integrator object.

This method is called automatically by initializePatchHierarchy() prior to the construction of the patch hierarchy. It is also possible for users to make an explicit call to initializeHierarchyIntegrator() prior to calling initializePatchHierarchy().

Implements IBTK::HierarchyIntegrator.

Reimplemented in IBAMR::IBInterpolantHierarchyIntegrator, and IBAMR::IBImplicitStaggeredHierarchyIntegrator.

◆ initializeLevelDataSpecialized()

void IBAMR::IBHierarchyIntegrator::initializeLevelDataSpecialized ( SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > >  hierarchy,
int  level_number,
double  init_data_time,
bool  can_be_refined,
bool  initial_time,
SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchLevel< NDIM > >  old_level,
bool  allocate_data 
)
overrideprotectedvirtual

Initialize data on a new level after it is inserted into an AMR patch hierarchy by the gridding algorithm.

Reimplemented from IBTK::HierarchyIntegrator.

◆ initializePatchHierarchy()

void IBAMR::IBHierarchyIntegrator::initializePatchHierarchy ( SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
SAMRAI::tbox::Pointer< SAMRAI::mesh::GriddingAlgorithm< NDIM > >  gridding_alg 
)
overridevirtual

Initialize the AMR patch hierarchy and data defined on the hierarchy at the start of a computation. If the computation is begun from a restart file, the patch hierarchy and patch data are read from the hierarchy database. Otherwise, the patch hierarchy and patch data are initialized by the gridding algorithm associated with the integrator object.

The implementation of this function assumes that the hierarchy exists upon entry to the function, but that it contains no patch levels. On return from this function, the state of the integrator object will be such that it is possible to step through time via the advanceHierarchy() function.

Reimplemented from IBTK::HierarchyIntegrator.

◆ postprocessIntegrateHierarchy()

void IBAMR::IBHierarchyIntegrator::postprocessIntegrateHierarchy ( double  current_time,
double  new_time,
bool  skip_synchronize_new_state_data,
int  num_cycles = 1 
)
overridevirtual

◆ preprocessIntegrateHierarchy()

void IBAMR::IBHierarchyIntegrator::preprocessIntegrateHierarchy ( double  current_time,
double  new_time,
int  num_cycles = 1 
)
overridevirtual

Basic functions to prepare to advance data from current_time to new_time.

A default implementation is provided that sets the current values of num_cycles and the time step size and checks to see if the time step size has changed.

Reimplemented from IBTK::HierarchyIntegrator.

Reimplemented in IBAMR::IBInterpolantHierarchyIntegrator, and IBAMR::IBImplicitStaggeredHierarchyIntegrator.

◆ putToDatabaseSpecialized()

void IBAMR::IBHierarchyIntegrator::putToDatabaseSpecialized ( SAMRAI::tbox::Pointer< SAMRAI::tbox::Database db)
overrideprotectedvirtual

Write out specialized object state to the given database.

Reimplemented from IBTK::HierarchyIntegrator.

Reimplemented in IBAMR::IBInterpolantHierarchyIntegrator, and IBAMR::IBImplicitStaggeredHierarchyIntegrator.

◆ registerBodyForceFunction()

void IBAMR::IBHierarchyIntegrator::registerBodyForceFunction ( SAMRAI::tbox::Pointer< IBTK::CartGridFunction F_fcn)

Supply a body force (optional).

◆ registerLoadBalancer()

void IBAMR::IBHierarchyIntegrator::registerLoadBalancer ( SAMRAI::tbox::Pointer< SAMRAI::mesh::LoadBalancer< NDIM > >  load_balancer)
overridevirtual

Register a load balancer for non-uniform load balancing.

Reimplemented from IBTK::HierarchyIntegrator.

◆ regridHierarchyBeginSpecialized()

void IBAMR::IBHierarchyIntegrator::regridHierarchyBeginSpecialized ( )
overrideprotectedvirtual

Perform necessary data movement, workload estimation, and logging prior to regridding.

Reimplemented from IBTK::HierarchyIntegrator.

◆ regridHierarchyEndSpecialized()

void IBAMR::IBHierarchyIntegrator::regridHierarchyEndSpecialized ( )
overrideprotectedvirtual

Perform necessary data movement and logging after regridding.

Reimplemented from IBTK::HierarchyIntegrator.

◆ resetHierarchyConfigurationSpecialized()

void IBAMR::IBHierarchyIntegrator::resetHierarchyConfigurationSpecialized ( SAMRAI::tbox::Pointer< SAMRAI::hier::BasePatchHierarchy< NDIM > >  hierarchy,
int  coarsest_level,
int  finest_level 
)
overrideprotectedvirtual

Reset cached hierarchy dependent data.

Reimplemented from IBTK::HierarchyIntegrator.

Member Data Documentation

◆ d_error_on_dt_change

bool IBAMR::IBHierarchyIntegrator::d_error_on_dt_change = true
protected

Flags to determine whether warnings or error messages should be emitted when time step size changes are encountered.

◆ d_ib_context

SAMRAI::tbox::Pointer<SAMRAI::hier::VariableContext> IBAMR::IBHierarchyIntegrator::d_ib_context
protected

Context containing all patch data indices relevant to IB operations.

◆ d_ib_data

SAMRAI::hier::ComponentSelector IBAMR::IBHierarchyIntegrator::d_ib_data
protected

ComponentSelector corresponding to d_ib_context. Also contains patch data indices for relevant cloned indices (which, as they are clones, cannot be placed in the Context).

◆ d_regrid_fluid_cfl_estimate

double IBAMR::IBHierarchyIntegrator::d_regrid_fluid_cfl_estimate = 0.0
protected

Estimation on the maximum fraction of fluid cells the structure has moved based on the maximum fluid velocity.

◆ d_regrid_fluid_cfl_interval

double IBAMR::IBHierarchyIntegrator::d_regrid_fluid_cfl_interval = -1.0
protected

The regrid CFL interval indicates the number of meshwidths a particle may move in any coordinate direction between invocations of the regridding process.

Note
Currently, when the CFL-based regrid interval is specified, it is always used instead of the fixed-step regrid interval.

◆ d_regrid_structure_cfl_estimate

double IBAMR::IBHierarchyIntegrator::d_regrid_structure_cfl_estimate = 0.0
protected

Estimation on the maximum fraction of fluid cells the structure has moved based on the infinity norm of the structure's displacement.

◆ d_time_stepping_type

TimeSteppingType IBAMR::IBHierarchyIntegrator::d_time_stepping_type = MIDPOINT_RULE
protected

Enum indicating the time integration employed for the IB equations.

◆ d_use_structure_predictor

bool IBAMR::IBHierarchyIntegrator::d_use_structure_predictor
protected

Flag indicating whether to use an explicit predictor for the structure configuration in the time stepping scheme.


The documentation for this class was generated from the following files: