IBAMR
An adaptive and distributed-memory parallel implementation of the immersed boundary (IB) method
|
Class CIBStrategy is a lightweight abstract strategy class which provides support for constraint based IB methods for rigid bodies. More...
#include </home/runner/work/IBAMR/IBAMR/include/ibamr/CIBStrategy.h>
Public Member Functions | |
CIBStrategy (const unsigned int parts) | |
Constructor of the class. | |
virtual | ~CIBStrategy () |
Destructor of the class. | |
virtual void | setConstraintForce (Vec L, const double data_time, const double scale=1.0)=0 |
Prepare the implementation class for sprading constraint force. In particular, set the constraint Lagrangian force in the internal data structure of the class. More... | |
virtual void | getConstraintForce (Vec *L, const double data_time)=0 |
Get the constraint rigid body force at the specified time within the current time interval. Generally, the implementation class maintains and stores the constraint force. This routine is called by constraint solver to update the contraint force after the (converged) solution is obtained. More... | |
virtual void | getFreeRigidVelocities (Vec *U, const double data_time) |
Get the free rigid velocities (DOFs) at the specified time within the current time interval. This routine is called by constraint solver to update the free rigid DOFs after the (converged) solution is obtained. More... | |
virtual void | getNetExternalForceTorque (Vec *F, const double data_time) |
Get net external force and torque at the specified time within the current time interval. This routine is called by constraint solver to form the appropriate RHS. More... | |
virtual void | subtractMeanConstraintForce (Vec L, int f_data_idx, const double scale=1.0)=0 |
Subtract the mean of constraint force from the background Eulerian grid. This is required for certain cases like periodic steady Stokes. More... | |
virtual void | setInterpolatedVelocityVector (Vec V, const double data_time) |
Prepare the implementation class for getting the interpolated fluid velocity on the Lagrangian vector V . More... | |
virtual void | getInterpolatedVelocity (Vec V, const double data_time, const double scale=1.0)=0 |
Get the interpolated velocity from the Eulerian grid at the specified time. More... | |
virtual void | computeMobilityRegularization (Vec D, Vec L, const double scale=1.0)=0 |
Compute regularization vector for the mobility problem. More... | |
unsigned int | getNumberOfRigidStructures () const |
Get number of rigid structures registered with this class. | |
virtual unsigned int | getNumberOfNodes (const unsigned int part) const =0 |
Get number of nodes associated with the particular structure. | |
void | setInitialCenterOfMass (const unsigned int part, const Eigen::Vector3d &XCOM_0) |
Set the initial center of mass location for the structures. More... | |
void | setSolveRigidBodyVelocity (const unsigned int part, const IBTK::FreeRigidDOFVector &solve_rigid_dofs) |
Set what rigid DOFs need to be solved for this particular structure. More... | |
const IBTK::FreeRigidDOFVector & | getSolveRigidBodyVelocity (const unsigned int part, int &num_free_dofs) const |
Query what rigid DOFs need to be solved for. | |
virtual void | setRigidBodyVelocity (const unsigned int part, const IBTK::RigidDOFVector &U, Vec V)=0 |
Set the rigid body velocity at the nodal/marker points contained in the Vec V. More... | |
virtual void | setRigidBodyVelocity (const unsigned int part, Vec U, Vec V) |
Set the rigid body velocity at the nodal/marker points contained in the Vec V. More... | |
virtual void | setRigidBodyVelocity (Vec U, Vec V, const bool only_free_dofs, const bool only_imposed_dofs, const bool all_dofs=false) |
Set the rigid body velocity at the nodal/marker points contained in the Vec V. More... | |
virtual void | computeNetRigidGeneralizedForce (const unsigned int part, Vec L, IBTK::RigidDOFVector &F)=0 |
Compute total force and torque on the structure. More... | |
virtual void | computeNetRigidGeneralizedForce (const unsigned int part, Vec L, Vec F) |
Compute total force and torque on the structure. More... | |
virtual void | computeNetRigidGeneralizedForce (Vec L, Vec F, const bool only_free_dofs, const bool only_imposed_dofs, const bool all_dofs=false) |
Compute total force and torque on the structure. More... | |
const IBTK::RigidDOFVector & | getNetRigidGeneralizedForce (const unsigned int part) |
Get total torque and force on the structure at new_time within the current time interval. More... | |
void | updateFreeDOFsMapping () |
Update the mapping of free DOFs for all structures if they are collected in a global vector. | |
void | updateNewRigidBodyVelocity (const unsigned int part, const IBTK::RigidDOFVector &U) |
Update the rigid body velocity obtained from the constraint Stokes solver for free-moving case. | |
void | updateNewRigidBodyVelocity (const unsigned int part, Vec U) |
Update the rigid body velocity obtained from the constraint Stokes solver for free-moving case. | |
void | updateNewRigidBodyVelocity (Vec U, const bool only_free_dofs, const bool only_imposed_dofs, const bool all_dofs=false) |
Update the rigid body velocity obtained from the constraint Stokes solver for free-moving case. | |
virtual void | copyVecToArray (Vec b, double *array, const std::vector< unsigned > &struct_ids, const int data_depth, const int array_rank) |
Copy data from distributed PETSc Vec for specified stucture indices to an array defined on a single processor. A default empty implementation is provided. More... | |
virtual void | copyFreeDOFsVecToArray (Vec b, double *array, const std::vector< unsigned > &struct_ids, const int array_rank) |
Copy data from distributed PETSc Vec for specified stucture indices to an array defined on a single processor. A default implementation is provided. More... | |
virtual void | copyArrayToVec (Vec b, double *array, const std::vector< unsigned > &struct_ids, const int data_depth, const int array_rank) |
Copy data from array defined on a single processor for specified stucture indices to distributed PETScVec. A default empty implementation is provided. More... | |
virtual void | copyFreeDOFsArrayToVec (Vec b, double *array, const std::vector< unsigned > &struct_ids, const int array_rank) |
Copy data from array defined on a single processor for specified stucture indices to distributed PETScVec. A default implementation is provided. More... | |
void | getCurrentRigidBodyVelocity (const unsigned int part, IBTK::RigidDOFVector &U) |
Get the rigid body translational velocity at the beginning of the timestep. | |
void | getNewRigidBodyVelocity (const unsigned int part, IBTK::RigidDOFVector &U) |
Get the rigid body translational velocity at the end of the timestep. | |
const Eigen::Vector3d & | getCurrentBodyCenterOfMass (const unsigned int part) |
Get body center of mass at the current time step. | |
const Eigen::Vector3d & | getNewBodyCenterOfMass (const unsigned int part) |
Get body center of mass at half time step. | |
virtual void | constructMobilityMatrix (const std::string &mat_name, MobilityMatrixType mat_type, Mat &mobility_mat, const std::vector< unsigned > &prototype_struct_ids, const double *grid_dx, const double *domain_extents, const bool initial_time, double rho, double mu, const std::pair< double, double > &scale, double f_periodic_corr, const int managing_rank) |
Construct dense mobility matrix for the prototypical structures identified by their indices. More... | |
virtual void | constructGeometricMatrix (const std::string &mat_name, Mat &geometric_mat, const std::vector< unsigned > &prototype_struct_ids, const bool initial_time, const int managing_rank) |
Construct a geometric matrix for the prototypical structures identified by their indices. A geometric matrix maps center of mass rigid body velocity to nodal velocities. Geometric matrix is generally used with a dense mobility matrices to construct an associated body-mobility matrix algebrically. More... | |
virtual void | rotateArray (double *array, const std::vector< unsigned > &struct_ids, const bool use_transpose, const int managing_rank, const int depth) |
Rotate vector using rotation matrix to/from the reference frame of the structures (which is at the initial time of the simulation). More... | |
Static Public Member Functions | |
static void | vecToRDV (Vec U, IBTK::RigidDOFVector &Ur) |
Set the DOFs from PETSc Vec U to RigidDOFVector Ur . | |
static void | rdvToVec (const IBTK::RigidDOFVector &Ur, Vec &U) |
Set the DOFs from RigidDOFVector Ur to PETSc Vec U . | |
static void | eigenToRDV (const Eigen::Vector3d &U, const Eigen::Vector3d &W, IBTK::RigidDOFVector &UW) |
Set the DOFs from Eigen::Vector3d U and W to RigidDOFVector UW . | |
static void | rdvToEigen (const IBTK::RigidDOFVector &UW, Eigen::Vector3d &U, Eigen::Vector3d &W) |
Set the DOFs from RigidDOFVector UW to Eigen::Vector3d U and W . | |
Protected Member Functions | |
void | setRotationMatrix (const IBTK::EigenAlignedVector< Eigen::Vector3d > &rot_vel, const IBTK::EigenAlignedVector< Eigen::Quaterniond > &q_old, IBTK::EigenAlignedVector< Eigen::Quaterniond > &q_new, IBTK::EigenAlignedVector< Eigen::Matrix3d > &rot_mat, const double dt) |
Fill the rotation matrix. More... | |
Protected Attributes | |
unsigned int | d_num_rigid_parts |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_center_of_mass_initial |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_center_of_mass_current |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_center_of_mass_half |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_center_of_mass_new |
std::vector< bool > | d_compute_center_of_mass_initial |
IBTK::EigenAlignedVector< Eigen::Quaterniond > | d_quaternion_current |
IBTK::EigenAlignedVector< Eigen::Quaterniond > | d_quaternion_half |
IBTK::EigenAlignedVector< Eigen::Quaterniond > | d_quaternion_new |
std::vector< IBTK::FRDV > | d_solve_rigid_vel |
Vec | d_U |
Vec | d_F |
std::vector< std::pair< int, int > > | d_free_dofs_map |
bool | d_free_dofs_map_updated |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_trans_vel_current |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_trans_vel_half |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_trans_vel_new |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_rot_vel_current |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_rot_vel_half |
IBTK::EigenAlignedVector< Eigen::Vector3d > | d_rot_vel_new |
std::vector< IBTK::RigidDOFVector > | d_net_rigid_generalized_force |
Class CIBStrategy is a lightweight abstract strategy class which provides support for constraint based IB methods for rigid bodies.
|
pure virtual |
Compute regularization vector for the mobility problem.
D | Vector containing the regularization for the mobility problem. |
L | Vector from which regularization is to be computed. |
Implemented in IBAMR::CIBMethod.
|
pure virtual |
Compute total force and torque on the structure.
L | The Lagrange multiplier vector. |
F | Vector RDV storing the net generalized force. |
F | RDV storing the net generalized force. |
Implemented in IBAMR::CIBMethod.
|
virtual |
Compute total force and torque on the structure.
part | The structure index. |
L | The Lagrange multiplier vector. |
F | Vec storing the net generalized force. |
|
virtual |
Compute total force and torque on the structure.
L | The Lagrange multiplier vector. |
F | Vec storing the net generalized force. |
only_free_dofs | Boolean indicating if the net generalized force and torque is to be computed only for free dofs of all bodies. |
only_imposed_dofs | Boolean indicating if the net generalized force and torque is to be computed for imposed dofs of all bodies. |
all_dofs | Boolean indicating if the net generalized force and torque is to be computed for all dofs of all bodies. |
|
virtual |
Construct a geometric matrix for the prototypical structures identified by their indices. A geometric matrix maps center of mass rigid body velocity to nodal velocities. Geometric matrix is generally used with a dense mobility matrices to construct an associated body-mobility matrix algebrically.
mat_name | Matrix handle. |
geometric_mat | Dense sequential geometric matrix. The matrix is stored in column-major(FORTRAN) order. |
prototype_struct_ids | Indices of the structures as registered with |
initial_time | Boolean to indicate if the corresponding geometric matrix is to be generated for the initial position of the structures. |
managing_rank | Rank of the processor managing this dense matrix. |
Reimplemented in IBAMR::CIBMethod.
|
virtual |
Construct dense mobility matrix for the prototypical structures identified by their indices.
mat_name | Matrix handle. |
mat_type | Mobility matrix type, e.g., RPY, EMPIRICAL, etc. |
mobility_mat | Dense sequential mobility matrix. The matrix is stored in column-major(FORTRAN) order. |
prototype_struct_ids | Indices of the structures as registered with |
grid_dx | NDIM vector of grid spacing of structure level. |
domain_extents | NDIM vector of domain length. |
initial_time | Boolean to indicate if the mobility matrix is to be generated for the initial position of the structures. |
rho | Fluid density |
mu | Fluid viscosity. |
scale | Scale for improving the conditioning number of dense mobility matrix. The matrix is scaled as: |
managing_rank | Rank of the processor managing this dense matrix. |
Reimplemented in IBAMR::CIBMethod.
|
virtual |
Copy data from array defined on a single processor for specified stucture indices to distributed PETScVec. A default empty implementation is provided.
b | Copy to PETSc Vec. |
array | Copy from data pointer. |
struct_ids | Vector of structure indices. |
data_depth | Depth of the data stored at each Lagrangian node. |
array_rank | Rank of the processor on which the array is located. |
Reimplemented in IBAMR::CIBMethod.
|
virtual |
Copy data from array defined on a single processor for specified stucture indices to distributed PETScVec. A default implementation is provided.
b | Copy to PETSc Vec. The Vec stores only free DOFs of all the structures. |
array | Copy from data pointer. It is a linear array of maximum free DOFs of the passed structure IDs. |
struct_ids | Vector of structure indices. |
array_rank | Rank of the processor on which the array is located. |
|
virtual |
Copy data from distributed PETSc Vec for specified stucture indices to an array defined on a single processor. A default implementation is provided.
b | PETSc Vec to copy from. The Vec stores only free DOFs of all the structures. |
array | Data pointer to copy to. It is a linear array of maximum free DOFs of the passed structure IDs. |
struct_ids | Vector of structure indices. |
array_rank | Rank of the processor on which the array is located. |
|
virtual |
Copy data from distributed PETSc Vec for specified stucture indices to an array defined on a single processor. A default empty implementation is provided.
b | PETSc Vec to copy from. The Vec stores data for nodal/marker points. |
array | Data pointer to copy to. |
struct_ids | Vector of structure indices. |
data_depth | Depth of the data stored at each Lagrangian node. |
array_rank | Rank of the processor on which the array is located. |
Reimplemented in IBAMR::CIBMethod.
|
pure virtual |
Get the constraint rigid body force at the specified time within the current time interval. Generally, the implementation class maintains and stores the constraint force. This routine is called by constraint solver to update the contraint force after the (converged) solution is obtained.
time | Time (current_time or new_time) at which constraint force is required. |
Implemented in IBAMR::CIBMethod.
|
virtual |
Get the free rigid velocities (DOFs) at the specified time within the current time interval. This routine is called by constraint solver to update the free rigid DOFs after the (converged) solution is obtained.
time | Time (current_time or new_time) at which constraint force is required. |
Reimplemented in IBAMR::CIBMethod.
|
pure virtual |
Get the interpolated velocity from the Eulerian grid at the specified time.
V | Vector that should contain the interpolated velocity. |
data_time | Time at which Eulerian velocity is to be interpolated. |
scale | Scales the velocity vector after interpolating from the Eulerian grid. |
Implemented in IBAMR::CIBMethod.
|
virtual |
Get net external force and torque at the specified time within the current time interval. This routine is called by constraint solver to form the appropriate RHS.
time | Time (current_time or new_time) at which external force and torque is required. |
Reimplemented in IBAMR::CIBMethod.
const RigidDOFVector & IBAMR::CIBStrategy::getNetRigidGeneralizedForce | ( | const unsigned int | part | ) |
Get total torque and force on the structure at new_time within the current time interval.
part | The rigid part. |
|
virtual |
Rotate vector using rotation matrix to/from the reference frame of the structures (which is at the initial time of the simulation).
array | Raw data pointer containing the vector enteries. |
struct_ids | Structure ID indices. |
use_transpose | Use transpose of rotation matrix to rotate the vector. |
managing_rank | Rank of the processor managing the matrix. |
depth | Depth of the data array components. |
Reimplemented in IBAMR::CIBMethod.
|
pure virtual |
Prepare the implementation class for sprading constraint force. In particular, set the constraint Lagrangian force in the internal data structure of the class.
L | Vec containing the constraint force for all structures. |
data_time | Time at which constraint force is to be spread. |
scale | Scales the constraint force before spreading. |
Implemented in IBAMR::CIBMethod.
void IBAMR::CIBStrategy::setInitialCenterOfMass | ( | const unsigned int | part, |
const Eigen::Vector3d & | XCOM_0 | ||
) |
Set the initial center of mass location for the structures.
part | The rigid body for which we are setting the initial center of mass position. |
|
virtual |
Prepare the implementation class for getting the interpolated fluid velocity on the Lagrangian vector V
.
V | Vector that should contain the interpolated velocity. |
data_time | Time at which Eulerian velocity is to be interpolated. |
Reimplemented in IBAMR::CIBMethod.
|
pure virtual |
Set the rigid body velocity at the nodal/marker points contained in the Vec V.
part | The rigid part for which velocity needs to be set. |
U | RDV contains the rigid component of velocities. For two-dimensions the vector contains the values and for three-dimensions the vector values are . |
Implemented in IBAMR::CIBMethod.
|
virtual |
Set the rigid body velocity at the nodal/marker points contained in the Vec V.
part | The rigid part for which velocity needs to be set. |
U | Vec contains the rigid component of velocities. For two-dimensions the vector contains the values and for three-dimensions the vector values are . |
|
virtual |
Set the rigid body velocity at the nodal/marker points contained in the Vec V.
U | Vec that contains the rigid component of velocities for the required components. For two-dimensions each sub Vec contains the values and for three-dimensions the vector values are . |
only_free_dofs | Boolean indicating if the rigid body velocity is to be set only for free DOFS for all parts. The corresponding size of U_sub would be . |
only_imposed_dofs | Boolean indicating if the rigid body velocity is to be set only for prescribed kinematics dofs for all parts. The corresponding size of U_sub would be . |
all_dofs | Boolean indicating if the rigid body velocity is to be set for all parts. The corresponding size of U_sub would be . |
|
protected |
Fill the rotation matrix.
q_old | Previous applied quaternion. |
q_new | New quaternion to set. |
rot_mat | Matrix to set. |
dt | Time interval of rotation. |
void IBAMR::CIBStrategy::setSolveRigidBodyVelocity | ( | const unsigned int | part, |
const IBTK::FreeRigidDOFVector & | solve_rigid_dofs | ||
) |
Set what rigid DOFs need to be solved for this particular structure.
part | The rigid body for which we are setting the free DOFs. |
|
pure virtual |
Subtract the mean of constraint force from the background Eulerian grid. This is required for certain cases like periodic steady Stokes.
L | Vec containing the constraint force. |
f_data_idx | Patch data index of Eulerian body force. |
scale | Factor by which L is scaled. |
Implemented in IBAMR::CIBMethod.
|
protected |
Center of mass.
|
protected |
Quaternion of the body.
|
protected |
Rigid body velocity of the structures.