Class CCLaplaceOperator is a concrete LaplaceOperator which implements a globally second-order accurate cell-centered finite difference discretization of a scalar elliptic operator of the form . More...
Class CellCellNoCornersFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the cell-centered ghost region surrounding a patch, excluding all corners. In 3D, it is also possible to configure this fill pattern object also to exclude all edges. More...
Class CopyToRootTransaction is a concrete implementation of the abstract base class SAMRAI::tbox::Transaction. It is used to communicate distributed patch data to a unified patch data object on a root MPI process. More...
Class EdgeCellSynchCopyFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the edge-centered ghost region surrounding a patch appropriate for "synchronizing" edge-centered values in an axis-by-axis manner at patch boundaries. More...
Class FaceCellSynchCopyFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the face-centered ghost region surrounding a patch appropriate for "synchronizing" face-centered values at patch boundaries. More...
Class FACPreconditionerStrategy provides an interface for specifying the problem-specific operations needed to implement a specific FAC preconditioner. More...
Class GeneralSolver provides an abstract interface for the implementation of linear or nonlinear solvers for systems of equations defined on an AMR patch hierarchy. More...
Class HierarchyGhostCellInterpolation encapsulates the operations required to set ghost cell values at physical and coarse-fine boundaries across a range of levels of a locally refined patch hierarchy. More...
Class HierarchyIntegrator provides an abstract interface for a time integrator for a system of equations defined on an AMR grid hierarchy, along with basic data management for variables defined on that hierarchy. More...
Class IndexUtilities is a utility class that defines simple functions such as conversion routines between physical coordinates and Cartesian index space. More...
Class LEInteractor is a utility class that defines several functions to interpolate data from Eulerian grid patches onto Lagrangian meshes and to spread values (not densities) from Lagrangian meshes to Eulerian grid patches. More...
Class LibMeshSystemIBVectors is a convenience class that manages access to libMesh vectors for the same system defined on multiple parts. It extends the base class LibMeshSystemVectors to provide access to vectors ghosted with both the Lagrangian partitioning (i.e., libMesh's computed partitioning, as in the base class) as well as the IB partitioning (i.e., the partitioning based on the distribution of SAMRAI data). More...
Class LIndexSetData is a specialization of the templated class LSetData that is intended to be used with Lagrangian data objects that provide Lagrangian and PETSc indexing information. More...
Class LSetData is a specialization of the templated class SAMRAI::pdat::IndexData that provides access to Lagrangian objects that are embedded in the a Cartesian grid patch. More...
Class LSetDataIterator is an iterator class which may be used to iterate through LSet objects associated with a specified box in cell-centered index space. More...
Class LTransaction represents a communication transaction between two processors or a local data copy for communicating or copying Lagrangian objects. More...
Class MergingLoadBalancer merges the boxes generated by a load balancer in a final step to decrease the total number of boxes. In essence, it postprocesses the list of boxes generated by its parent class to try and coalesce the set of boxes on each process. More...
Class NewtonKrylovSolver provides an abstract interface for the implementation of inexact Newton-Krylov solvers for nonlinear problems of the form . More...
Class NodeCellSynchCopyFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the node-centered ghost region surrounding a patch appropriate for "synchronizing" node-centered values in an axis-by-axis manner at patch boundaries. More...
Class PartitioningBox implements an NDIM-dimensional bounding box defined by two points. Unlike a standard bounding box, a PartitioningBox is an NDIM-dimensional tensor product of half-open intervals: i.e., it is a half-open box. This property allows one to partition a domain into a set of boxes. More...
Class PartitioningBoxes stores a set of bounding boxes and can check if a point is in the set of partitioning boxes or not in a more optimized way than just looping over a std::vector<PartitioningBox>. More...
Class PETScMFFDJacobianOperator provides a method for computing Jacobian-vector products, i.e., , via a matrix-free finite-difference approach. More...
Class RobinPhysBdryPatchStrategy is an abstract strategy class that extends the SAMRAI::xfer::RefinePatchStrategy interface to provide routines specific for setting Robin boundary conditions at physical boundaries. This class also provides default implementations of some methods defined in SAMRAI::xfer::RefinePatchStrategy that are generally not needed for filling ghost cell values at physical boundaries. More...
Class SAMRAIDataCache is a utility class for caching cloned SAMRAI patch data. Patch data are allocated as needed and should not be deallocated by the caller. More...
Class SCLaplaceOperator is a concrete LaplaceOperator which implements a globally second-order accurate side-centered finite difference discretization of a scalar elliptic operator of the form . More...
Class SideCellNoCornersFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the cell-centered ghost region surrounding a patch, excluding all corners. In 3D, it is also possible to configure this fill pattern object also to exclude all edges. More...
Class SideCellSynchCopyFillPattern is a concrete implementation of the abstract base class SAMRAI::xfer::VariableFillPattern. It is used to calculate overlaps according to a pattern which limits overlaps to the side-centered ghost region surrounding a patch appropriate for "synchronizing" side-centered values at patch boundaries. More...
Class VCSCViscousOperator is a subclass of SCLaplaceOperator which implements a globally second-order accurate side-centered finite difference discretization of a vector elliptic operator of the form . More...
Compute the product of the finite element representation of the force and (optionally) a set of weights (e.g., JxW values) at all quadrature points on an element. See integrate_elem_rhs for details on the first two template arguments.
Defines "using" declarations for all IBTK namespaces. This header file may be included in application codes, but it MUST NOT be included in any other header (.h) or inline (-inl.h) file in the library.
Eigen types have special alignment requirements and require a specific memory allocator. This is a convenience type alias for a std::vector with the correct allocator.
Utility function which asserts that the SAMRAI pointer is not null.
This is useful for writing generic code in which we might want to assert that a pointer is not null in the initialization list where we cannot use the normal assertion macro.
This function interpolates in time between two snapshots. We find the two closest time points stored in the snapshot cache and linear interpolate in time between them. Each snapshot will be transferred onto the provided hierarchy. A scratch patch index with the same patch layout as u_idx must be provided, and data for both must be already allocated. The refine_type must be a valid refinement operation for the data layout.
This assumes that the time is between two snapshots stored in the cache. An error will occur if two time points can not be found. If there is a single snapshot, this function returns that snapshot.
The hier_data_ops object, if provided, must match the same variable type used in for the snapshots.
If the optional period parameter is provided, this function will treat the first snapshot time point t_1 as also the last snapshot time point t_1 + period.
Transfer the snapshot at the specified time value within a given tolerance to the patch index u_idx on the supplied patch hierarchy. The patch index u_idx must contain sufficient ghost cell width to perform the operations used by the refinement operator.
Note this can require a significant amount of communication if the supplied patch hierarchy has a different configuration of patches than the snapshot patch hierarchy.
This function does not synchronize the data on the current hierarchy (i.e. no coarsening is performed). It transfers the snapshot using the refinement operator.
Assemble the element contribution to a load vector.
The three functions below allow compile-time selection of the number of basis functions and number of variables. Fixing the loop bounds via templates and switch statements makes most IBFE codes a few percent faster.
This function (the last of the three) takes the number of variables and number of basis functions as template argument which allows the compiler to unroll the two inner loops when appropriate. Setting either value to -1 results in run-time selection of loop bounds (which is considerably slower).
Parameters
[in]
qp_offset
Offset used to calculate an index into F_qp. The relevant part of the array is assumed to start at n_vars * qp_offset.
[in]
phi_F
Values of test functions evaluated at quadrature points, indexed by test function number and then quadrature point number.
[in]
JxW_F
Products of Jacobian and quadrature weight at each quadrature point.
[in]
F_qp
Globally indexed array containing function values at quadrature points. The array is indexed by quadrature point and then by variable: i.e., if n_vars is greater than one then components of the vector-valued function being projected at a certain quadrature point are contiguous.
[out]
F_rhs_concatenated
Vector containing element integrals of products of test functions and values of the interpolated function. Unlike F_qp, the values in this vector are first indexed by variable and then by quadrature point.
Given a set of boxes (which describe a region in index space), return another set of boxes whose union covers the same index space. The returned set of boxes is formed by merging boxes in boxes along their longest edges.
Compute the product of the finite element representation of the force and (optionally) a set of weights (e.g., JxW values) at all quadrature points on an element. See integrate_elem_rhs for details on the first two template arguments.
Template Parameters
weights_are_unity
Assume that all values in weights are 1. This value should be true when we are only interpolating values at quadrature points and otherwise false.
Parameters
[in]
qp_offset
Offset used to calculate an index into F_w_qp. The relevant part of the array is assumed to start at n_vars * qp_offset.
[in]
phi_F
Values of test functions evaluated at quadrature points, indexed by test function number and then quadrature point number.
[in]
weights
Vector containing weights at quadrature points: in the case of force spreading this is the standard JxW array. If weights_are_unity is true then its values are assumed to be 1.
[in]
F_node
Values of the finite element solution (i.e., the multipliers on each trial function) on the current element.
[out]
F_w_qp
Globally indexed array containing the product of the finite element representation of the force and weights at the quadrature points on the current element. The array is indexed by quadrature point and then by variable: i.e., if n_vars is greater than one then variables of the vector-valued force being evaluated at a certain quadrature point are contiguous.
Update a snapshot at the given time. This function copies the data stored in u_idx into the snapshot index in the cache object. The snapshotted hierarchy is also updated to the current hierarchy.
This outputs an error if a snapshot with that time point within the specified tolerance can not be found.
Calls to this function may use a different patch index than the one used in setSnapshot(), but the underlying data layout must be consistent. This data layout is specified by the variable used to construct the cache.