IBAMR  IBAMR version 0.19.
Classes | Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Attributes | List of all members
IBTK::FEDataManager Class Reference

Class FEDataManager coordinates data required for Lagrangian-Eulerian interaction between a Lagrangian finite element (FE) mesh. In particular, the FEData member object stores the necessary finite element data while this class stores additional data dependent on the Eulerian grid. More...

#include <ibtk/FEDataManager.h>

Inheritance diagram for IBTK::FEDataManager:
Inheritance graph
[legend]

Classes

struct  InterpSpec
 Struct InterpSpec encapsulates data needed to specify the manner in which Eulerian-to-Lagrangian interpolation is performed when using an FE structural discretization. More...
 
struct  SpreadSpec
 Struct SpreadSpec encapsulates data needed to specify the manner in which Lagrangian-to-Eulerian spreading is performed when using an FE structural discretization. More...
 
struct  WorkloadSpec
 Struct WorkloadSpec encapsulates the parameters used to calculate workload estimates (i.e., the input to the load balancing algorithm). More...
 

Public Types

using SystemDofMapCache = FEData::SystemDofMapCache
 

Public Member Functions

const std::string & getCurrentCoordinatesSystemName () const
 The name of the equation system which stores the spatial position data. The actual string is stored by FEData. More...
 
void setCurrentCoordinatesSystemName (const std::string &coordinates_system_name) const
 Set name of the equation system which stores the spatial position data. The actual string is stored by FEData. More...
 
void setEquationSystems (libMesh::EquationSystems *equation_systems, int level_number)
 Set the equations systems object that is associated with the FEData object. Currently, each set of equation systems must be assigned to a particular level of the AMR grid. More...
 
libMesh::EquationSystems * getEquationSystems () const
 
SystemDofMapCachegetDofMapCache (const std::string &system_name)
 
SystemDofMapCachegetDofMapCache (unsigned int system_num)
 
void setLoggingEnabled (bool enable_logging=true)
 Enable or disable logging. More...
 
bool getLoggingEnabled () const
 Determine whether logging is enabled or disabled. More...
 

Static Public Member Functions

static FEDataManagergetManager (std::shared_ptr< FEData > fe_data, const std::string &name, const SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > &input_db, const int max_levels, const InterpSpec &default_interp_spec, const SpreadSpec &default_spread_spec, const WorkloadSpec &default_workload_spec, const SAMRAI::hier::IntVector< NDIM > &min_ghost_width=SAMRAI::hier::IntVector< NDIM >(0), std::shared_ptr< SAMRAIDataCache > eulerian_data_cache=nullptr, bool register_for_restart=true)
 
static void freeAllManagers ()
 

Public Attributes

std::string & COORDINATES_SYSTEM_NAME
 The name of the equation system which stores the spatial position data. The actual string is stored by FEData. More...
 

Static Public Attributes

static const libMesh::boundary_id_type ZERO_DISPLACEMENT_X_BDRY_ID
 The libMesh boundary IDs to use for specifying essential boundary conditions. More...
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_Y_BDRY_ID
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_Z_BDRY_ID
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_XY_BDRY_ID
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_XZ_BDRY_ID
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_YZ_BDRY_ID
 
static const libMesh::boundary_id_type ZERO_DISPLACEMENT_XYZ_BDRY_ID
 

Protected Attributes

std::shared_ptr< FEDatad_fe_data
 
std::shared_ptr< FEProjectord_fe_projector
 
std::map< std::string, std::unique_ptr< libMesh::PetscVector< double > > > d_L2_proj_matrix_diag_ghost
 

Methods to set and get the patch hierarchy and range of patch

levels associated with this manager class.

SubdomainToPatchLevelTranslation d_level_lookup
 
std::string d_object_name
 
bool d_registered_for_restart
 
bool d_enable_logging = false
 
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > d_hierarchy
 
int d_max_level_number = IBTK::invalid_level_number
 
std::shared_ptr< SAMRAIDataCached_eulerian_data_cache
 
SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContextd_context
 
SAMRAI::tbox::Pointer< SAMRAI::pdat::CellVariable< NDIM, double > > d_qp_count_var
 
int d_qp_count_idx
 
const WorkloadSpec d_default_workload_spec
 
const InterpSpec d_default_interp_spec
 
const SpreadSpec d_default_spread_spec
 
NodeOutsidePatchCheckType d_node_patch_check = NODE_OUTSIDE_ERROR
 
const SAMRAI::hier::IntVector< NDIM > d_ghost_width
 
const SAMRAI::hier::IntVector< NDIM > d_associated_elem_ghost_width = SAMRAI::hier::IntVector<NDIM>(1)
 
std::vector< std::vector< std::vector< libMesh::Elem * > > > d_active_patch_elem_map
 
std::vector< std::vector< std::vector< libMesh::Node * > > > d_active_patch_node_map
 
std::map< std::string, std::vector< unsigned int > > d_active_patch_ghost_dofs
 
std::vector< libMesh::Elem * > d_active_elems
 
std::map< std::string, std::unique_ptr< libMesh::NumericVector< double > > > d_system_ghost_vec
 
std::map< std::string, std::unique_ptr< libMesh::PetscVector< double > > > d_system_ib_ghost_vec
 
static std::map< std::string, FEDataManager * > s_data_manager_instances
 
static bool s_registered_callback
 
static unsigned char s_shutdown_priority
 
void setPatchHierarchy (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy)
 Reset patch hierarchy over which operations occur. More...
 
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > getPatchHierarchy () const
 Get the patch hierarchy used by this object. More...
 
int getCoarsestPatchLevelNumber () const
 
int getFinestPatchLevelNumber () const
 
const SAMRAI::hier::IntVector< NDIM > & getGhostCellWidth () const
 
const InterpSpecgetDefaultInterpSpec () const
 
const SpreadSpecgetDefaultSpreadSpec () const
 
const std::vector< std::vector< libMesh::Elem * > > & getActivePatchElementMap () const
 
const std::vector< std::vector< libMesh::Node * > > & getActivePatchNodeMap () const
 
void reinitElementMappings ()
 Reinitialize the mappings from elements to Cartesian grid patches. More...
 
libMesh::NumericVector< double > * getSolutionVector (const std::string &system_name) const
 
libMesh::NumericVector< double > * buildGhostedSolutionVector (const std::string &system_name, bool localize_data=true)
 
std::unique_ptr< libMesh::PetscVector< double > > buildIBGhostedVector (const std::string &system_name)
 
libMesh::NumericVector< double > * getCoordsVector () const
 
libMesh::NumericVector< double > * buildGhostedCoordsVector (bool localize_data=true)
 
std::shared_ptr< FEData > & getFEData ()
 
const std::shared_ptr< FEData > & getFEData () const
 
void spread (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name)
 Spread a density from the FE mesh to the Cartesian grid using the default spreading spec. More...
 
void spread (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const SpreadSpec &spread_spec)
 Spread a density from the FE mesh to the Cartesian grid. More...
 
void spread (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, RobinPhysBdryPatchStrategy *f_phys_bdry_op, double fill_data_time, bool close_F=true, bool close_X=true)
 Spread a density from the FE mesh to the Cartesian grid using the default spreading spec. This is a convenience overload of spread where f_phys_bdry_op is used to correctly accumulate force values spread outside the physical domain into it. More...
 
void spread (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const SpreadSpec &spread_spec, RobinPhysBdryPatchStrategy *f_phys_bdry_op, double fill_data_time, bool close_F=true, bool close_X=true)
 Spread a density from the FE mesh to the Cartesian grid using a specified spreading spec. This is a convenience overload of spread where f_phys_bdry_op is used to correctly accumulate force values spread outside the physical domain into it. More...
 
void prolongData (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, bool is_density=true, bool accumulate_on_grid=true, bool close_F=true, bool close_X=true)
 Prolong a value or a density from the FE mesh to the Cartesian grid. More...
 
void interpWeighted (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &f_refine_scheds=std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), double fill_data_time=0.0, bool close_F=true, bool close_X=true)
 Set up the right-hand side F of an L2 projection problem where Eulerian data given by f_data_idx is projected onto the finite element space given by system_name. More...
 
void interpWeighted (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const InterpSpec &interp_spec, const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &f_refine_scheds=std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), double fill_data_time=0.0, bool close_F=true, bool close_X=true)
 Interpolate a value from the Cartesian grid to the FE mesh using a specified interpolation spec. This interpolation function does NOT do an L2-projection of the interpolated quantity. It does however weighs/filters the interpolated quantity at the quadrature points to the nodes. Here, the basis functions of the deformational field is used as the filter. More...
 
void interp (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &f_refine_scheds=std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), double fill_data_time=0.0, bool close_X=true)
 Interpolate a value from the Cartesian grid to the FE mesh using the default interpolation spec. More...
 
void interp (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, const InterpSpec &interp_spec, const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &f_refine_scheds=std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), double fill_data_time=0.0, bool close_X=true)
 Interpolate a value from the Cartesian grid to the FE mesh using a specified interpolation spec. More...
 
void restrictData (int f_data_idx, libMesh::NumericVector< double > &F, libMesh::NumericVector< double > &X, const std::string &system_name, bool use_consistent_mass_matrix=true, bool close_X=true)
 Restrict a value from the Cartesian grid to the FE mesh. More...
 
libMesh::NumericVector< double > * buildDiagonalL2MassMatrix (const std::string &system_name)
 
libMesh::PetscVector< double > * buildIBGhostedDiagonalL2MassMatrix (const std::string &system_name)
 
bool computeL2Projection (libMesh::NumericVector< double > &U, libMesh::NumericVector< double > &F, const std::string &system_name, bool consistent_mass_matrix=true, bool close_U=true, bool close_F=true, double tol=1.0e-6, unsigned int max_its=100)
 Set U to be the L2 projection of F. More...
 
void addWorkloadEstimate (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, const int workload_data_idx, const int coarsest_ln=invalid_level_number, const int finest_ln=invalid_level_number)
 Update the cell workload estimate by adding a value (see the main documentation of this class for information on how this is computed) to the d_workload_idx cell variable. More...
 
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)
 
void putToDatabase (SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > db) override
 
static bool updateQuadratureRule (std::unique_ptr< libMesh::QBase > &qrule, libMesh::QuadratureType quad_type, libMesh::Order quad_order, bool use_adaptive_quadrature, double point_density, bool allow_rules_with_negative_weights, const libMesh::Elem *elem, const boost::multi_array< double, 2 > &X_node, double dx_min)
 
static bool updateInterpQuadratureRule (std::unique_ptr< libMesh::QBase > &qrule, const InterpSpec &spec, const libMesh::Elem *elem, const boost::multi_array< double, 2 > &X_node, double dx_min)
 
static bool updateSpreadQuadratureRule (std::unique_ptr< libMesh::QBase > &qrule, const SpreadSpec &spec, const libMesh::Elem *elem, const boost::multi_array< double, 2 > &X_node, double dx_min)
 
static void zeroExteriorValues (const SAMRAI::geom::CartesianPatchGeometry< NDIM > &patch_geom, const std::vector< double > &X_qp, std::vector< double > &F_qp, int n_vars)
 Zero the values corresponding to points on the given patch (described by patch_geometry) that might be duplicated on another patch. More...
 
 FEDataManager (std::string object_name, const SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > &input_db, const int max_levels, InterpSpec default_interp_spec, SpreadSpec default_spread_spec, WorkloadSpec default_workload_spec, SAMRAI::hier::IntVector< NDIM > ghost_width, std::shared_ptr< SAMRAIDataCache > eulerian_data_cache, std::shared_ptr< FEData > fe_data, bool register_for_restart=true)
 Constructor, where the FEData object owned by this class may be co-owned by other objects. More...
 
 ~FEDataManager ()
 The FEDataManager destructor cleans up any allocated data objects. More...
 
 FEDataManager ()=delete
 Default constructor. More...
 
 FEDataManager (const FEDataManager &from)=delete
 Copy constructor. More...
 
FEDataManageroperator= (const FEDataManager &that)=delete
 Assignment operator. More...
 
void updateQuadPointCountData (int coarsest_ln, int finest_ln)
 
void collectActivePatchElements (std::vector< std::vector< libMesh::Elem * > > &active_patch_elems, int level_number, int coarsest_elem_ln, int finest_elem_ln)
 
void collectActivePatchNodes (std::vector< std::vector< libMesh::Node * > > &active_patch_nodes, const std::vector< std::vector< libMesh::Elem * > > &active_patch_elems)
 
int getPatchLevel (const libMesh::Elem *elem) const
 
void collectGhostDOFIndices (std::vector< unsigned int > &ghost_dofs, const std::vector< libMesh::Elem * > &active_elems, const std::string &system_name)
 
void reinitializeIBGhostedDOFs (const std::string &system_name)
 
virtual void getFromRestart ()
 

Detailed Description

Parameters read from the input database

node_outside_patch_check: parameter controlling how this class responds to mesh nodes outside the finest patch level. In all cases, for backwards compatibility, nodes that are outside the computational domain are permitted and are ignored by this check. Possible values are:

  1. node_outside_permit: Permit nodes to be outside the finest patch level.
  2. node_outside_warn: Permit nodes to be outside the finest patch level, but log a warning whenever this is detected.
  3. node_outside_error: Abort the program when nodes are detected outside the finest patch level.

The default value is node_outside_error.

subdomain_ids_on_levels: a database correlating libMesh subdomain IDs to patch levels. A possible value for this is

subdomain_ids_on_levels
{
level_1 = 4
level_3 = 10, 12, 14
level_5 = 1000, 1003, 1006, 1009
}

This particular input will associate all elements with subdomain id 4 with patch level 1, all elements with subdomain ids 10, 12, or 14 with patch level 3, etc. All unspecified subdomain ids will be associated with the finest patch level. All inputs in this database for levels finer than the finest level are ignored (e.g., if the maximum patch level number is 4, then the values given in the example for level 5 ultimately end up on level 4). This feature is experimental: at the current time it is known that it produces some artifacts at the coarse-fine interface, but that these generally don't effect the overall solution quality.

Parameters effecting workload estimate calculations

FEDataManager can estimate the amount of work done in IBFE calculations (such as FEDataManager::spread). Since most calculations use a variable number of quadrature points on each libMesh element this estimate can vary quite a bit over different Eulerian cells corresponding to a single mesh. The current implementation estimates the workload on each cell of the background Eulerian grid by applying a background value representing the work on the Eulerian cell itself and a weight times the number of quadrature points on that cell. These values are set at the time of object construction through the FEDataManager::WorkloadSpec object, which contains reasonable defaults.

Note
Multiple FEDataManager objects may be instantiated simultaneously.

Member Typedef Documentation

◆ SystemDofMapCache

Alias FEData::SystemDofMapCache for backwards compatibility.

Constructor & Destructor Documentation

◆ FEDataManager() [1/3]

IBTK::FEDataManager::FEDataManager ( std::string  object_name,
const SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > &  input_db,
const int  max_levels,
InterpSpec  default_interp_spec,
SpreadSpec  default_spread_spec,
WorkloadSpec  default_workload_spec,
SAMRAI::hier::IntVector< NDIM >  ghost_width,
std::shared_ptr< SAMRAIDataCache eulerian_data_cache,
std::shared_ptr< FEData fe_data,
bool  register_for_restart = true 
)
protected

◆ ~FEDataManager()

IBTK::FEDataManager::~FEDataManager ( )
protected

◆ FEDataManager() [2/3]

IBTK::FEDataManager::FEDataManager ( )
privatedelete
Note
This constructor is not implemented and should not be used.

◆ FEDataManager() [3/3]

IBTK::FEDataManager::FEDataManager ( const FEDataManager from)
privatedelete
Note
This constructor is not implemented and should not be used.
Parameters
fromThe value to copy to this object.

Member Function Documentation

◆ getCurrentCoordinatesSystemName()

const std::string& IBTK::FEDataManager::getCurrentCoordinatesSystemName ( ) const
inline
Note
The default value for this string is "coordinates system".

◆ setCurrentCoordinatesSystemName()

void IBTK::FEDataManager::setCurrentCoordinatesSystemName ( const std::string &  coordinates_system_name) const
inline
Note
The default value for this string is "coordinates system".

◆ getManager()

static FEDataManager* IBTK::FEDataManager::getManager ( std::shared_ptr< FEData fe_data,
const std::string &  name,
const SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > &  input_db,
const int  max_levels,
const InterpSpec default_interp_spec,
const SpreadSpec default_spread_spec,
const WorkloadSpec default_workload_spec,
const SAMRAI::hier::IntVector< NDIM > &  min_ghost_width = SAMRAI::hier::IntVector< NDIM >(0),
std::shared_ptr< SAMRAIDataCache eulerian_data_cache = nullptr,
bool  register_for_restart = true 
)
static

Return a pointer to the instance of the Lagrangian data manager corresponding to the specified name. Access to FEDataManager objects is mediated by the getManager() function.

Note
When a manager is accessed for the first time, the FEDataManager::freeAllManagers() static method is registered with the SAMRAI::tbox::ShutdownRegistry class. Consequently, all allocated managers are freed at program completion. Thus, users of this class do not explicitly allocate or deallocate the FEDataManager instances.
Returns
A pointer to the data manager instance.

◆ freeAllManagers()

static void IBTK::FEDataManager::freeAllManagers ( )
static

Deallocate all of the FEDataManager instances.

It is not necessary to call this function at program termination since it is automatically called by the SAMRAI::tbox::ShutdownRegistry class.

◆ setEquationSystems()

void IBTK::FEDataManager::setEquationSystems ( libMesh::EquationSystems *  equation_systems,
int  level_number 
)
Deprecated:
This function is deprecated since the FEData constructor now requires an EquationSystems argument.

◆ getEquationSystems()

libMesh::EquationSystems* IBTK::FEDataManager::getEquationSystems ( ) const
Returns
A pointer to the equations systems object that is associated with the FEData object.

◆ getDofMapCache() [1/2]

SystemDofMapCache* IBTK::FEDataManager::getDofMapCache ( const std::string &  system_name)
Returns
The DofMapCache for a specified system.

◆ getDofMapCache() [2/2]

SystemDofMapCache* IBTK::FEDataManager::getDofMapCache ( unsigned int  system_num)
Returns
The DofMapCache for a specified system.

◆ setLoggingEnabled()

void IBTK::FEDataManager::setLoggingEnabled ( bool  enable_logging = true)
Note
This is usually set by the IBFEMethod which owns the current FEDataManager, which reads the relevant boolean from the database.

◆ getLoggingEnabled()

bool IBTK::FEDataManager::getLoggingEnabled ( ) const

◆ setPatchHierarchy()

void IBTK::FEDataManager::setPatchHierarchy ( SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy)

The patch hierarchy must be fully set up (i.e., contain all the levels it is expected to have) at the point this function is called. If you need to tag cells for refinement to create the initial hierarchy then use applyGradientDetector, which does not use the stored patch hierarchy.

◆ getPatchHierarchy()

SAMRAI::tbox::Pointer<SAMRAI::hier::PatchHierarchy<NDIM> > IBTK::FEDataManager::getPatchHierarchy ( ) const

◆ getCoarsestPatchLevelNumber()

int IBTK::FEDataManager::getCoarsestPatchLevelNumber ( ) const

Get the coarsest patch level number on which elements are assigned.

◆ getFinestPatchLevelNumber()

int IBTK::FEDataManager::getFinestPatchLevelNumber ( ) const

Get the finest patch level number on which elements are assigned.

◆ getGhostCellWidth()

const SAMRAI::hier::IntVector<NDIM>& IBTK::FEDataManager::getGhostCellWidth ( ) const
Returns
The ghost cell width used for quantities that are to be interpolated from the Cartesian grid to the FE mesh.

◆ getDefaultInterpSpec()

const InterpSpec& IBTK::FEDataManager::getDefaultInterpSpec ( ) const
Returns
The specifications of the scheme used for interpolating from the Cartesian grid to the FE mesh.

◆ getDefaultSpreadSpec()

const SpreadSpec& IBTK::FEDataManager::getDefaultSpreadSpec ( ) const
Returns
The specifications of the scheme used for spreading densities from the FE mesh to the Cartesian grid

◆ getActivePatchElementMap()

const std::vector<std::vector<libMesh::Elem*> >& IBTK::FEDataManager::getActivePatchElementMap ( ) const
Returns
A const reference to the map from local patch number to local active elements.

◆ getActivePatchNodeMap()

const std::vector<std::vector<libMesh::Node*> >& IBTK::FEDataManager::getActivePatchNodeMap ( ) const
Returns
A const reference to the map from local patch number to local active nodes.
Note
The local active nodes are the nodes of the local active elements.

◆ reinitElementMappings()

void IBTK::FEDataManager::reinitElementMappings ( )

◆ getSolutionVector()

libMesh::NumericVector<double>* IBTK::FEDataManager::getSolutionVector ( const std::string &  system_name) const
Returns
A pointer to the unghosted solution vector associated with the specified system.

◆ buildGhostedSolutionVector()

libMesh::NumericVector<double>* IBTK::FEDataManager::buildGhostedSolutionVector ( const std::string &  system_name,
bool  localize_data = true 
)
Returns
A pointer to the ghosted solution vector associated with the specified system. The vector contains positions for values in the relevant IB ghost region which are populated if localize_data is true.
Note
The vector returned by pointer is owned by this class (i.e., no copying is done).
Deprecated:
Use buildIBGhostedVector() instead which clones a vector with the same ghost region.

◆ buildIBGhostedVector()

std::unique_ptr<libMesh::PetscVector<double> > IBTK::FEDataManager::buildIBGhostedVector ( const std::string &  system_name)
Returns
A pointer to a vector, with ghost entries corresponding to relevant IB data, associated with the specified system.

◆ getCoordsVector()

libMesh::NumericVector<double>* IBTK::FEDataManager::getCoordsVector ( ) const
Returns
A pointer to the unghosted coordinates (nodal position) vector.

◆ buildGhostedCoordsVector()

libMesh::NumericVector<double>* IBTK::FEDataManager::buildGhostedCoordsVector ( bool  localize_data = true)
Returns
A pointer to the ghosted coordinates (nodal position) vector.
Deprecated:
Use buildIBGhostedVector() instead.

◆ getFEData() [1/2]

std::shared_ptr<FEData>& IBTK::FEDataManager::getFEData ( )
Returns
The shared pointer to the object managing the Lagrangian data.

◆ getFEData() [2/2]

const std::shared_ptr<FEData>& IBTK::FEDataManager::getFEData ( ) const
Returns
The shared pointer to the object managing the Lagrangian data.

◆ spread() [1/4]

void IBTK::FEDataManager::spread ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name 
)
Note
This function may spread forces from points near the physical boundary into ghost cells outside the physical domain. To account for these forces one should usually call IBTK::RobinPhysBdryPatchStrategy::accumulateFromPhysicalBoundaryData() after this function to move said force values into the physical domain.

◆ spread() [2/4]

void IBTK::FEDataManager::spread ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const SpreadSpec spread_spec 
)
Parameters
[in]f_data_idxSAMRAI index into the SAMRAI::hier::PatchHierarchy object owned by this class. The spread force values will be added into this index.
[in]FFinite element solution vector containing the field which will be spread onto the Eulerian grid.
[in]XFinite element solution vector containing the current position of the mesh.
[in]system_namename of the system corresponding to F.

Both X and F should contain ghost values corresponding to the IB partitioning of the Lagrangian data, i.e., vectors returned from buildIBGhostedVector.

Note
This function may spread forces from points near the physical boundary into ghost cells outside the physical domain. To account for these forces one should usually call IBTK::RobinPhysBdryPatchStrategy::accumulateFromPhysicalBoundaryData() after this function to move said force values into the physical domain.

◆ spread() [3/4]

void IBTK::FEDataManager::spread ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
RobinPhysBdryPatchStrategy f_phys_bdry_op,
double  fill_data_time,
bool  close_F = true,
bool  close_X = true 
)

◆ spread() [4/4]

void IBTK::FEDataManager::spread ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const SpreadSpec spread_spec,
RobinPhysBdryPatchStrategy f_phys_bdry_op,
double  fill_data_time,
bool  close_F = true,
bool  close_X = true 
)

◆ prolongData()

void IBTK::FEDataManager::prolongData ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
bool  is_density = true,
bool  accumulate_on_grid = true,
bool  close_F = true,
bool  close_X = true 
)

◆ interpWeighted() [1/2]

void IBTK::FEDataManager::interpWeighted ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &  f_refine_scheds = std::vector< SAMRAI::tbox::PointerSAMRAI::xfer::RefineSchedule< NDIM > > >(),
double  fill_data_time = 0.0,
bool  close_F = true,
bool  close_X = true 
)
Parameters
[in]f_data_idxIndex of the variable being projected.
[in]IB-ghostedposition vector containing the current position of all dofs whose basis functions have support (in the deformed configuration) on the locally owned set of patches.
[out]FVector into which the L2-projection RHS vector will be assembled. If this vector is ghosted then off-processor entries are assembled into the ghost value region of the vector and callers should complete the assembly process with, e.g.,
ierr = VecGhostUpdateBegin(F.vec(), ADD_VALUES, SCATTER_REVERSE);
ierr = VecGhostUpdateEnd(F.vec(), INSERT_VALUES, SCATTER_FORWARD);
otherwise values are added into the vector directly (which, for a libMesh::PetscVector, will involve the use of the PETSc VecCache object to track off-processor entries). It is strongly recommended, for performance reasons, that one assemble into a ghosted vector instead of the approach with VecCache.
[in]system_nameName of the libMesh system corresponding to the vector F.
[in]f_refine_schedsRefinement schedules to process before actually running this function.
[in]fill_data_timeTime at which the data in f_data_idx was filled.
[in]Whetheror not to close F after assembly.
[in]Whetheror not to close X before assembly.
Note
We recommend against using the last two booleans: it is usually better to do any vector communication before calling this function.
This function is poorly named: it actually sets up a dual-space vector F which is the right-hand side of an L2 projection problem. Callers will still need to solve the resulting linear system. The result, therefore, is projected, not interpolated.
The vector set up by this function is set up without applying any constraints (e.g., hanging node or periodicity constraints) during assembly. This is because we do not store the constraints on the IB data partitioning, only on libMesh's own data partitioning. For more information on how to assemble RHS vectors in this case, see the documentation of the function apply_transposed_constraint_matrix().

◆ interpWeighted() [2/2]

void IBTK::FEDataManager::interpWeighted ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const InterpSpec interp_spec,
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &  f_refine_scheds = std::vector< SAMRAI::tbox::PointerSAMRAI::xfer::RefineSchedule< NDIM > > >(),
double  fill_data_time = 0.0,
bool  close_F = true,
bool  close_X = true 
)

◆ interp() [1/2]

void IBTK::FEDataManager::interp ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &  f_refine_scheds = std::vector< SAMRAI::tbox::PointerSAMRAI::xfer::RefineSchedule< NDIM > > >(),
double  fill_data_time = 0.0,
bool  close_X = true 
)

◆ interp() [2/2]

void IBTK::FEDataManager::interp ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
const InterpSpec interp_spec,
const std::vector< SAMRAI::tbox::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > > &  f_refine_scheds = std::vector< SAMRAI::tbox::PointerSAMRAI::xfer::RefineSchedule< NDIM > > >(),
double  fill_data_time = 0.0,
bool  close_X = true 
)

◆ restrictData()

void IBTK::FEDataManager::restrictData ( int  f_data_idx,
libMesh::NumericVector< double > &  F,
libMesh::NumericVector< double > &  X,
const std::string &  system_name,
bool  use_consistent_mass_matrix = true,
bool  close_X = true 
)

◆ buildDiagonalL2MassMatrix()

libMesh::NumericVector<double>* IBTK::FEDataManager::buildDiagonalL2MassMatrix ( const std::string &  system_name)
Returns
Pointer to vector representation of diagonal L2 mass matrix.

◆ buildIBGhostedDiagonalL2MassMatrix()

libMesh::PetscVector<double>* IBTK::FEDataManager::buildIBGhostedDiagonalL2MassMatrix ( const std::string &  system_name)
Returns
Pointer to IB ghosted vector representation of diagonal L2 mass matrix.

◆ computeL2Projection()

bool IBTK::FEDataManager::computeL2Projection ( libMesh::NumericVector< double > &  U,
libMesh::NumericVector< double > &  F,
const std::string &  system_name,
bool  consistent_mass_matrix = true,
bool  close_U = true,
bool  close_F = true,
double  tol = 1.0e-6,
unsigned int  max_its = 100 
)

◆ updateQuadratureRule()

static bool IBTK::FEDataManager::updateQuadratureRule ( std::unique_ptr< libMesh::QBase > &  qrule,
libMesh::QuadratureType  quad_type,
libMesh::Order  quad_order,
bool  use_adaptive_quadrature,
double  point_density,
bool  allow_rules_with_negative_weights,
const libMesh::Elem *  elem,
const boost::multi_array< double, 2 > &  X_node,
double  dx_min 
)
static

Update the quadrature rule for the current element. If the provided qrule is already configured appropriately, it is not modified.

Returns
true if the quadrature rule is updated or otherwise requires reinitialization (e.g. because the element type or p_level changed); false otherwise.

◆ updateInterpQuadratureRule()

static bool IBTK::FEDataManager::updateInterpQuadratureRule ( std::unique_ptr< libMesh::QBase > &  qrule,
const InterpSpec spec,
const libMesh::Elem *  elem,
const boost::multi_array< double, 2 > &  X_node,
double  dx_min 
)
static

Update the quadrature rule for the current element used by the Lagrangian-Eulerian interaction scheme. If the provided qrule is already configured appropriately, it is not modified.

Returns
true if the quadrature rule is updated or otherwise requires reinitialization (e.g. because the element type or p_level changed); false otherwise.

◆ updateSpreadQuadratureRule()

static bool IBTK::FEDataManager::updateSpreadQuadratureRule ( std::unique_ptr< libMesh::QBase > &  qrule,
const SpreadSpec spec,
const libMesh::Elem *  elem,
const boost::multi_array< double, 2 > &  X_node,
double  dx_min 
)
static

Update the quadrature rule for the current element used by the Lagrangian-Eulerian interaction scheme. If the provided qrule is already configured appropriately, it is not modified.

Returns
true if the quadrature rule is updated or otherwise requires reinitialization (e.g. because the element type or p_level changed); false otherwise.

◆ addWorkloadEstimate()

void IBTK::FEDataManager::addWorkloadEstimate ( SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
const int  workload_data_idx,
const int  coarsest_ln = invalid_level_number,
const int  finest_ln = invalid_level_number 
)

◆ applyGradientDetector()

void IBTK::FEDataManager::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 
)

Set integer tags to "one" in cells where refinement of the given level should occur due to the presence of Lagrangian data. The double time argument is the regrid time. The integer "tag_index" argument is the patch descriptor index of the cell centered integer tag array on each patch in the hierarchy. The boolean argument initial_time indicates whether the level is being subject to refinement at the initial simulation time. If it is false, then the error estimation process is being invoked at some later time after the AMR hierarchy was initially constructed. The boolean argument uses_richardson_extrapolation_too is true when Richardson extrapolation error estimation is used in addition to the gradient detector, and false otherwise. This argument helps the user to manage multiple regridding criteria.

When assertion checking is active, an unrecoverable exception will result if the hierarchy pointer is null or the level number does not match any existing level in the hierarchy.

Note
This function is analogous to SAMRAI::mesh::StandardTagAndInitStrategy::applyGradientDetector() and is only meant to be called from IBAMR::IBFEMethod::applyGradientDetector().

◆ putToDatabase()

void IBTK::FEDataManager::putToDatabase ( SAMRAI::tbox::Pointer< SAMRAI::tbox::Database db)
overridevirtual

Write out object state to the given database.

When assertion checking is active, database pointer must be non-null.

Implements SAMRAI::tbox::Serializable.

◆ zeroExteriorValues()

static void IBTK::FEDataManager::zeroExteriorValues ( const SAMRAI::geom::CartesianPatchGeometry< NDIM > &  patch_geom,
const std::vector< double > &  X_qp,
std::vector< double > &  F_qp,
int  n_vars 
)
static

In the process of computing either forces on the Lagrangian mesh to spread to the Eulerian grid or interpolating values from the Eulerian grid to use in the Lagrangian structure we may, under extraordinary circumstances, double-count a point that lies on the boundary between two patches. Note that due to the CFL condition assumed by this class and the width of the tag buffer used by IBAMR::IBFEMethod we can assume that no IB point will ever lie on the boundary of the finest grid level: hence, if a point lies on the boundary of a patch, it must also lie on the boundary of a neighboring patch.

Hence, to establish uniqueness, we zero data that lies on the 'upper' face (corresponding to the unit normal vector pointing towards infinity) since we are guaranteed, by the previous assumptions, that that data must also lie on the lower face of a neighboring patch.

◆ operator=()

FEDataManager& IBTK::FEDataManager::operator= ( const FEDataManager that)
privatedelete
Note
This operator is not implemented and should not be used.
Parameters
thatThe value to assign to this object.
Returns
A reference to this object.

◆ updateQuadPointCountData()

void IBTK::FEDataManager::updateQuadPointCountData ( int  coarsest_ln,
int  finest_ln 
)
private

Compute the quadrature point counts in each cell of the level in which the FE mesh is embedded. Also zeros out node count data for other levels within the specified range of level numbers.

◆ collectActivePatchElements()

void IBTK::FEDataManager::collectActivePatchElements ( std::vector< std::vector< libMesh::Elem * > > &  active_patch_elems,
int  level_number,
int  coarsest_elem_ln,
int  finest_elem_ln 
)
private

Collect all of the active elements which are located within a local Cartesian grid patch grown by a ghost width of 1 (like IBTK::LEInteractor::getMinimumGhostWidth(), we assume that IB points are allowed to move no more than one cell width between regridding operations).

The parameters refer to the levels of different objects:

  1. level_number - the level number in the patch hierarchy on which we are identifying intersections.
  2. coarsest_elem_ln - The minimum level number of elements we should consider (see the main documentation of this class for an explanation on how elements are assigned to particular levels)
  3. finest_elem_ln - The maximum level number of elements we should consider.

All three parameters are necessary because we use this function both to tag cells for refinement (i.e., we want to refine cells containing elements on levels higher than the present level) and to do IB calculations (where all three numbers will be the same).

In this method, the determination as to whether an element is local or not is based on the position of the bounding box of the element.

◆ collectActivePatchNodes()

void IBTK::FEDataManager::collectActivePatchNodes ( std::vector< std::vector< libMesh::Node * > > &  active_patch_nodes,
const std::vector< std::vector< libMesh::Elem * > > &  active_patch_elems 
)
private

Collect all of the nodes of the active elements that are located within a local Cartesian grid patch grown by the specified ghost cell width.

◆ getPatchLevel()

int IBTK::FEDataManager::getPatchLevel ( const libMesh::Elem *  elem) const
private

Get the patch level on which an element lives.

◆ collectGhostDOFIndices()

void IBTK::FEDataManager::collectGhostDOFIndices ( std::vector< unsigned int > &  ghost_dofs,
const std::vector< libMesh::Elem * > &  active_elems,
const std::string &  system_name 
)
private

Collect all ghost DOF indices for the specified collection of elements.

◆ reinitializeIBGhostedDOFs()

void IBTK::FEDataManager::reinitializeIBGhostedDOFs ( const std::string &  system_name)
private

Reinitialize IB ghosted DOF data structures for the specified system.

◆ getFromRestart()

virtual void IBTK::FEDataManager::getFromRestart ( )
privatevirtual

Read object state from the restart file and initialize class data members. The database from which the restart data is read is determined by the object_name specified in the constructor.

Unrecoverable Errors:

  • The database corresponding to object_name is not found in the restart file.
  • The class version number and restart version number do not match.

Member Data Documentation

◆ d_fe_data

std::shared_ptr<FEData> IBTK::FEDataManager::d_fe_data
protected

FEData object that contains the libMesh data structures.

Note
multiple FEDataManager objects may use the same FEData object, usually combined with different hierarchies.

◆ d_fe_projector

std::shared_ptr<FEProjector> IBTK::FEDataManager::d_fe_projector
protected

FEProjector object that handles L2 projection functionality.

◆ d_L2_proj_matrix_diag_ghost

std::map<std::string, std::unique_ptr<libMesh::PetscVector<double> > > IBTK::FEDataManager::d_L2_proj_matrix_diag_ghost
protected

IB ghosted diagonal mass matrix representations.

◆ COORDINATES_SYSTEM_NAME

std::string& IBTK::FEDataManager::COORDINATES_SYSTEM_NAME
Note
The default value for this string is "coordinates system".

◆ ZERO_DISPLACEMENT_X_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_X_BDRY_ID
static

◆ ZERO_DISPLACEMENT_Y_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_Y_BDRY_ID
static

◆ ZERO_DISPLACEMENT_Z_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_Z_BDRY_ID
static

◆ ZERO_DISPLACEMENT_XY_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_XY_BDRY_ID
static

◆ ZERO_DISPLACEMENT_XZ_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_XZ_BDRY_ID
static

◆ ZERO_DISPLACEMENT_YZ_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_YZ_BDRY_ID
static

◆ ZERO_DISPLACEMENT_XYZ_BDRY_ID

const libMesh::boundary_id_type IBTK::FEDataManager::ZERO_DISPLACEMENT_XYZ_BDRY_ID
static

◆ d_level_lookup

SubdomainToPatchLevelTranslation IBTK::FEDataManager::d_level_lookup
private

Store the association between subdomain ids and patch levels.

◆ s_data_manager_instances

std::map<std::string, FEDataManager*> IBTK::FEDataManager::s_data_manager_instances
staticprivate

Static data members used to control access to and destruction of singleton data manager instance.

◆ s_registered_callback

bool IBTK::FEDataManager::s_registered_callback
staticprivate

◆ s_shutdown_priority

unsigned char IBTK::FEDataManager::s_shutdown_priority
staticprivate

◆ d_object_name

std::string IBTK::FEDataManager::d_object_name
private

The object name is used as a handle to databases stored in restart files and for error reporting purposes. The boolean is used to control restart file writing operations.

◆ d_registered_for_restart

bool IBTK::FEDataManager::d_registered_for_restart
private

◆ d_enable_logging

bool IBTK::FEDataManager::d_enable_logging = false
private

Whether or not to log data to the screen: see FEDataManager::setLoggingEnabled() and FEDataManager::getLoggingEnabled().

Note
This is usually set by IBFEMethod, which reads the relevant boolean from the database.

◆ d_hierarchy

SAMRAI::tbox::Pointer<SAMRAI::hier::PatchHierarchy<NDIM> > IBTK::FEDataManager::d_hierarchy
private

Grid hierarchy information.

◆ d_max_level_number

int IBTK::FEDataManager::d_max_level_number = IBTK::invalid_level_number
private

Maximum possible level number in the patch hierarchy.

◆ d_eulerian_data_cache

std::shared_ptr<SAMRAIDataCache> IBTK::FEDataManager::d_eulerian_data_cache
private

Cached Eulerian data to reduce the number of allocations/deallocations.

◆ d_context

SAMRAI::tbox::Pointer<SAMRAI::hier::VariableContext> IBTK::FEDataManager::d_context
private

SAMRAI::hier::VariableContext object used for data management.

◆ d_qp_count_var

SAMRAI::tbox::Pointer<SAMRAI::pdat::CellVariable<NDIM, double> > IBTK::FEDataManager::d_qp_count_var
private

SAMRAI::hier::Variable pointer and patch data descriptor indices for the cell variable used to keep track of the count of the quadrature points in each cell.

◆ d_qp_count_idx

int IBTK::FEDataManager::d_qp_count_idx
private

◆ d_default_workload_spec

const WorkloadSpec IBTK::FEDataManager::d_default_workload_spec
private

The default parameters used during workload calculations.

◆ d_default_interp_spec

const InterpSpec IBTK::FEDataManager::d_default_interp_spec
private

The default kernel functions and quadrature rule used to mediate Lagrangian-Eulerian interaction.

◆ d_default_spread_spec

const SpreadSpec IBTK::FEDataManager::d_default_spread_spec
private

◆ d_node_patch_check

NodeOutsidePatchCheckType IBTK::FEDataManager::d_node_patch_check = NODE_OUTSIDE_ERROR
private

after reassociating patches with elements a node may still lie outside all patches on the finest level in unusual circumstances (like when the parent integrator class does not regrid sufficiently frequently and has more than one patch level). This enum controls what we do when this problem is detected.

◆ d_ghost_width

const SAMRAI::hier::IntVector<NDIM> IBTK::FEDataManager::d_ghost_width
private

SAMRAI::hier::IntVector object which determines the required ghost cell width of this class.

◆ d_associated_elem_ghost_width

const SAMRAI::hier::IntVector<NDIM> IBTK::FEDataManager::d_associated_elem_ghost_width = SAMRAI::hier::IntVector<NDIM>(1)
private

SAMRAI::hier::IntVector object which determines how many ghost cells we should enlarge a patch by when associating an element with a patch. An element is associated with a patch when its bounding box (defined as the bounding box of both its nodes and quadrature points) intersects the bounding box (including ghost cells) of that patch.

Note
At the present time this is always 1, which matches the assumption made by IBTK::LEInteractor::getMinimumGhostWidth().

◆ d_active_patch_elem_map

std::vector<std::vector<std::vector<libMesh::Elem*> > > IBTK::FEDataManager::d_active_patch_elem_map
private

Data to manage mappings between mesh elements and grid patches.

◆ d_active_patch_node_map

std::vector<std::vector<std::vector<libMesh::Node*> > > IBTK::FEDataManager::d_active_patch_node_map
private

◆ d_active_patch_ghost_dofs

std::map<std::string, std::vector<unsigned int> > IBTK::FEDataManager::d_active_patch_ghost_dofs
private

◆ d_active_elems

std::vector<libMesh::Elem*> IBTK::FEDataManager::d_active_elems
private

◆ d_system_ghost_vec

std::map<std::string, std::unique_ptr<libMesh::NumericVector<double> > > IBTK::FEDataManager::d_system_ghost_vec
private

Ghost vectors for the various equation systems.

◆ d_system_ib_ghost_vec

std::map<std::string, std::unique_ptr<libMesh::PetscVector<double> > > IBTK::FEDataManager::d_system_ib_ghost_vec
private

Exemplar relevant IB-ghosted vectors for the various equation systems. These vectors are cloned for fast initialization in buildIBGhostedVector.


The documentation for this class was generated from the following file:
IBTK_CHKERRQ
#define IBTK_CHKERRQ(ierr)
Throw an error exception from within any C++ source code.
Definition: IBTK_CHKERRQ.h:39