|
IBAMR
IBAMR version 0.19.
|
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>

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 |
| SystemDofMapCache * | getDofMapCache (const std::string &system_name) |
| SystemDofMapCache * | getDofMapCache (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 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 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< FEData > | d_fe_data |
| std::shared_ptr< FEProjector > | d_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 | |
| 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< SAMRAIDataCache > | d_eulerian_data_cache |
| SAMRAI::tbox::Pointer< SAMRAI::hier::VariableContext > | d_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 InterpSpec & | getDefaultInterpSpec () const |
| const SpreadSpec & | getDefaultSpreadSpec () 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... | |
| FEDataManager & | operator= (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 () |
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:
node_outside_permit: Permit nodes to be outside the finest patch level. node_outside_warn: Permit nodes to be outside the finest patch level, but log a warning whenever this is detected. 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
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.
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.
Alias FEData::SystemDofMapCache for backwards compatibility.
|
protected |
|
protected |
|
privatedelete |
|
privatedelete |
| from | The value to copy to this object. |
|
inline |
|
inline |
|
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.
|
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.
| void IBTK::FEDataManager::setEquationSystems | ( | libMesh::EquationSystems * | equation_systems, |
| int | level_number | ||
| ) |
| libMesh::EquationSystems* IBTK::FEDataManager::getEquationSystems | ( | ) | const |
| SystemDofMapCache* IBTK::FEDataManager::getDofMapCache | ( | const std::string & | system_name | ) |
| SystemDofMapCache* IBTK::FEDataManager::getDofMapCache | ( | unsigned int | system_num | ) |
| void IBTK::FEDataManager::setLoggingEnabled | ( | bool | enable_logging = true | ) |
| bool IBTK::FEDataManager::getLoggingEnabled | ( | ) | const |
| 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.
| SAMRAI::tbox::Pointer<SAMRAI::hier::PatchHierarchy<NDIM> > IBTK::FEDataManager::getPatchHierarchy | ( | ) | const |
| int IBTK::FEDataManager::getCoarsestPatchLevelNumber | ( | ) | const |
Get the coarsest patch level number on which elements are assigned.
| int IBTK::FEDataManager::getFinestPatchLevelNumber | ( | ) | const |
Get the finest patch level number on which elements are assigned.
| const SAMRAI::hier::IntVector<NDIM>& IBTK::FEDataManager::getGhostCellWidth | ( | ) | const |
| const InterpSpec& IBTK::FEDataManager::getDefaultInterpSpec | ( | ) | const |
| const SpreadSpec& IBTK::FEDataManager::getDefaultSpreadSpec | ( | ) | const |
| const std::vector<std::vector<libMesh::Elem*> >& IBTK::FEDataManager::getActivePatchElementMap | ( | ) | const |
| const std::vector<std::vector<libMesh::Node*> >& IBTK::FEDataManager::getActivePatchNodeMap | ( | ) | const |
| void IBTK::FEDataManager::reinitElementMappings | ( | ) |
| libMesh::NumericVector<double>* IBTK::FEDataManager::getSolutionVector | ( | const std::string & | system_name | ) | const |
| libMesh::NumericVector<double>* IBTK::FEDataManager::buildGhostedSolutionVector | ( | const std::string & | system_name, |
| bool | localize_data = true |
||
| ) |
localize_data is true.| std::unique_ptr<libMesh::PetscVector<double> > IBTK::FEDataManager::buildIBGhostedVector | ( | const std::string & | system_name | ) |
| libMesh::NumericVector<double>* IBTK::FEDataManager::getCoordsVector | ( | ) | const |
| libMesh::NumericVector<double>* IBTK::FEDataManager::buildGhostedCoordsVector | ( | bool | localize_data = true | ) |
| std::shared_ptr<FEData>& IBTK::FEDataManager::getFEData | ( | ) |
| const std::shared_ptr<FEData>& IBTK::FEDataManager::getFEData | ( | ) | const |
| void IBTK::FEDataManager::spread | ( | int | f_data_idx, |
| libMesh::NumericVector< double > & | F, | ||
| libMesh::NumericVector< double > & | X, | ||
| const std::string & | system_name | ||
| ) |
| 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 | ||
| ) |
| [in] | f_data_idx | SAMRAI index into the SAMRAI::hier::PatchHierarchy object owned by this class. The spread force values will be added into this index. |
| [in] | F | Finite element solution vector containing the field which will be spread onto the Eulerian grid. |
| [in] | X | Finite element solution vector containing the current position of the mesh. |
| [in] | system_name | name 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.
| 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 |
||
| ) |
| 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 |
||
| ) |
| 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 |
||
| ) |
| 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::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), |
||
| double | fill_data_time = 0.0, |
||
| bool | close_F = true, |
||
| bool | close_X = true |
||
| ) |
| [in] | f_data_idx | Index of the variable being projected. |
| [in] | IB-ghosted | position 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] | F | Vector 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);
IBTK_CHKERRQ(ierr);
ierr = VecGhostUpdateEnd(F.vec(), INSERT_VALUES, SCATTER_FORWARD);
IBTK_CHKERRQ(ierr);
|
| [in] | system_name | Name of the libMesh system corresponding to the vector F. |
| [in] | f_refine_scheds | Refinement schedules to process before actually running this function. |
| [in] | fill_data_time | Time at which the data in f_data_idx was filled. |
| [in] | Whether | or not to close F after assembly. |
| [in] | Whether | or not to close X before assembly. |
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.| 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::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), |
||
| double | fill_data_time = 0.0, |
||
| bool | close_F = true, |
||
| bool | close_X = true |
||
| ) |
| 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::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), |
||
| double | fill_data_time = 0.0, |
||
| bool | close_X = true |
||
| ) |
| 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::Pointer< SAMRAI::xfer::RefineSchedule< NDIM > > >(), |
||
| double | fill_data_time = 0.0, |
||
| bool | close_X = true |
||
| ) |
| 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 |
||
| ) |
| libMesh::NumericVector<double>* IBTK::FEDataManager::buildDiagonalL2MassMatrix | ( | const std::string & | system_name | ) |
| libMesh::PetscVector<double>* IBTK::FEDataManager::buildIBGhostedDiagonalL2MassMatrix | ( | const std::string & | system_name | ) |
| 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 |
||
| ) |
|
static |
Update the quadrature rule for the current element. If the provided qrule is already configured appropriately, it is not modified.
|
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.
|
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.
| 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 |
||
| ) |
| 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.
|
overridevirtual |
Write out object state to the given database.
When assertion checking is active, database pointer must be non-null.
Implements SAMRAI::tbox::Serializable.
|
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.
|
privatedelete |
| that | The value to assign to this object. |
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.
|
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:
level_number - the level number in the patch hierarchy on which we are identifying intersections. 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) 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.
|
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.
|
private |
Get the patch level on which an element lives.
|
private |
Collect all ghost DOF indices for the specified collection of elements.
|
private |
Reinitialize IB ghosted DOF data structures for the specified system.
|
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:
|
protected |
FEData object that contains the libMesh data structures.
|
protected |
FEProjector object that handles L2 projection functionality.
|
protected |
IB ghosted diagonal mass matrix representations.
| std::string& IBTK::FEDataManager::COORDINATES_SYSTEM_NAME |
|
static |
|
static |
|
static |
|
static |
|
static |
|
static |
|
static |
|
private |
Store the association between subdomain ids and patch levels.
|
staticprivate |
Static data members used to control access to and destruction of singleton data manager instance.
|
staticprivate |
|
staticprivate |
|
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.
|
private |
|
private |
Whether or not to log data to the screen: see FEDataManager::setLoggingEnabled() and FEDataManager::getLoggingEnabled().
|
private |
Grid hierarchy information.
|
private |
Maximum possible level number in the patch hierarchy.
|
private |
Cached Eulerian data to reduce the number of allocations/deallocations.
|
private |
SAMRAI::hier::VariableContext object used for data management.
|
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.
|
private |
|
private |
The default parameters used during workload calculations.
|
private |
The default kernel functions and quadrature rule used to mediate Lagrangian-Eulerian interaction.
|
private |
|
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.
|
private |
SAMRAI::hier::IntVector object which determines the required ghost cell width of this class.
|
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.
|
private |
Data to manage mappings between mesh elements and grid patches.
|
private |
|
private |
|
private |
|
private |
Ghost vectors for the various equation systems.
|
private |
Exemplar relevant IB-ghosted vectors for the various equation systems. These vectors are cloned for fast initialization in buildIBGhostedVector.
1.8.17