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

Class IBStandardForceGen is a concrete IBLagrangianForceStrategy that is intended to be used in conjunction with curvilinear mesh data generated by class IBStandardInitializer. More...

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

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

Public Member Functions

 IBStandardForceGen (SAMRAI::tbox::Pointer< SAMRAI::tbox::Database > input_db=SAMRAI::tbox::Pointer< SAMRAI::tbox::Database >())
 Default constructor.
 
void registerSpringForceFunction (int force_fcn_index, const SpringForceFcnPtr spring_force_fcn_ptr, const SpringForceDerivFcnPtr spring_force_deriv_fcn_ptr=nullptr)
 Register a spring force specification function with the force generator. More...
 
void setUniformBodyForce (IBTK::Vector F, int structure_id, int level_number)
 Set a uniform body force that is applied on each point in the structure with the given structure_id.
 
void initializeLevelData (SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, int level_number, double init_data_time, bool initial_time, IBTK::LDataManager *l_data_manager) override
 Setup the data needed to compute the forces on the specified level of the patch hierarchy.
 
void computeLagrangianForce (SAMRAI::tbox::Pointer< IBTK::LData > F_data, SAMRAI::tbox::Pointer< IBTK::LData > X_data, SAMRAI::tbox::Pointer< IBTK::LData > U_data, SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, int level_number, double data_time, IBTK::LDataManager *l_data_manager) override
 Compute the force generated by the Lagrangian structure on the specified level of the patch hierarchy. More...
 
void computeLagrangianForceJacobianNonzeroStructure (std::vector< int > &d_nnz, std::vector< int > &o_nnz, SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, int level_number, IBTK::LDataManager *l_data_manager) override
 Compute the non-zero structure of the force Jacobian matrix. More...
 
void computeLagrangianForceJacobian (Mat &J_mat, MatAssemblyType assembly_type, double X_coef, SAMRAI::tbox::Pointer< IBTK::LData > X_data, double U_coef, SAMRAI::tbox::Pointer< IBTK::LData > U_data, SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, int level_number, double data_time, IBTK::LDataManager *l_data_manager) override
 Compute the Jacobian of the force with respect to the present structure configuration. More...
 
double computeLagrangianEnergy (SAMRAI::tbox::Pointer< IBTK::LData > X_data, SAMRAI::tbox::Pointer< IBTK::LData > U_data, SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > > hierarchy, int level_number, double data_time, IBTK::LDataManager *l_data_manager) override
 Compute the potential energy with respect to the present structure configuration and velocity.
 
- Public Member Functions inherited from IBAMR::IBLagrangianForceStrategy
 IBLagrangianForceStrategy ()=default
 Default constructor.
 
virtual ~IBLagrangianForceStrategy ()=default
 Virtual destructor.
 
virtual void setTimeInterval (double current_time, double new_time)
 Set the current and new times for the present timestep. More...
 

Detailed Description

Class IBStandardForceGen is a concrete IBLagrangianForceStrategy that is intended to be used in conjunction with curvilinear mesh data generated by class IBStandardInitializer.

Class IBStandardForceGen provides support for linear and nonlinear spring forces, linear beam forces, and target point penalty forces.

Note
By default, function default_linear_spring_force() is associated with spring force_fcn_idx 0. This is the default spring force function index used by class IBStandardInitializer in cases in which a force function index is not specified in a spring input file. Users may override this default force function with any function that implements the interface required by registerSpringForceFunction(). Users may also specify additional force functions that may be associated with arbitrary integer indices.

Member Function Documentation

◆ computeLagrangianForce()

void IBAMR::IBStandardForceGen::computeLagrangianForce ( SAMRAI::tbox::Pointer< IBTK::LData F_data,
SAMRAI::tbox::Pointer< IBTK::LData X_data,
SAMRAI::tbox::Pointer< IBTK::LData U_data,
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
int  level_number,
double  data_time,
IBTK::LDataManager l_data_manager 
)
overridevirtual

Compute the force generated by the Lagrangian structure on the specified level of the patch hierarchy.

Note
Nodal forces computed by this method are added to the force vector.

Reimplemented from IBAMR::IBLagrangianForceStrategy.

◆ computeLagrangianForceJacobian()

void IBAMR::IBStandardForceGen::computeLagrangianForceJacobian ( Mat &  J_mat,
MatAssemblyType  assembly_type,
double  X_coef,
SAMRAI::tbox::Pointer< IBTK::LData X_data,
double  U_coef,
SAMRAI::tbox::Pointer< IBTK::LData U_data,
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
int  level_number,
double  data_time,
IBTK::LDataManager l_data_manager 
)
overridevirtual

Compute the Jacobian of the force with respect to the present structure configuration.

Note
The elements of the Jacobian should be "accumulated" in the provided matrix J.

Reimplemented from IBAMR::IBLagrangianForceStrategy.

◆ computeLagrangianForceJacobianNonzeroStructure()

void IBAMR::IBStandardForceGen::computeLagrangianForceJacobianNonzeroStructure ( std::vector< int > &  d_nnz,
std::vector< int > &  o_nnz,
SAMRAI::tbox::Pointer< SAMRAI::hier::PatchHierarchy< NDIM > >  hierarchy,
int  level_number,
IBTK::LDataManager l_data_manager 
)
overridevirtual

Compute the non-zero structure of the force Jacobian matrix.

Note
Elements indices must be global PETSc indices.

Reimplemented from IBAMR::IBLagrangianForceStrategy.

◆ registerSpringForceFunction()

void IBAMR::IBStandardForceGen::registerSpringForceFunction ( int  force_fcn_index,
const SpringForceFcnPtr  spring_force_fcn_ptr,
const SpringForceDerivFcnPtr  spring_force_deriv_fcn_ptr = nullptr 
)

Register a spring force specification function with the force generator.

These functions are employed to compute the force generated by a particular spring for the specified displacement, spring constant, rest length, and Lagrangian index.

Note
By default, function default_linear_spring_force() is associated with force_fcn_idx 0.

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