|
| FiniteElement ()=delete |
| Default constructor. Disabled. More...
|
|
| FiniteElement (numerics::Matrix< double > &M, CellType cellType, bool useExternal=false) |
| Constructs a FiniteElement instance from a matrix consisting of the element coordinates and a cell type. More...
|
|
| ~FiniteElement () |
| Destructor. More...
|
|
bool | usesExternalBuffer () const |
| Checks if this instance is pointing to an external, user-supplied, buffer for the element coordinate information. More...
|
|
void | setMaxSolverIterations (int numIters) |
| Overrides the max number of iterations for the Newton-Raphson, used for the isoparametric inverse mapping from physical coordinates to reference coordinates. More...
|
|
int | getMaxSolverIterations () const |
| Returns the max number of iterations used for the Newton-Raphson. More...
|
|
CellType | getCellType () const |
| Returns the cell type associated with this FiniteElement instance. More...
|
|
int | getBasisType () const |
| Returns the Basis associated with this FiniteElement instance. More...
|
|
int | getPhysicalDimension () const |
| Returns the physical dimension for this FiniteElement instance. More...
|
|
int | getReferenceDimension () const |
| Returns the dimension of the element in reference space. More...
|
|
int | getNumNodes () const |
| Returns the number of nodes, i.e., degrees-of-freedom (dofs) of the FiniteElement. More...
|
|
int | getNumDofs () const |
| Returns the number of degrees of freedom associated with the Finite Element basis bound to this FiniteElement instance. More...
|
|
int | computeReferenceCoords (const double *xp, double *xr, double TOL=1.e-12) |
| Given a point in physical space, \( \bar{x_p} \), this method computes the corresponding reference coordinates, \( \bar{xi} \) in the reference space of the element. More...
|
|
void | computePhysicalCoords (const double *xr, double *xp) |
| Given reference coordinates, \( \xi \in \bar{\Omega}^e \), this method computes the corresponding physical coordinates, given by \( \mathbf{x}(\xi)=\sum\limits_{i=0}^n{N}_i^e({\xi}){x_i} \). More...
|
|
void | jacobian (const double *xr, numerics::Matrix< double > &J) |
| Given reference coordinates, \( \xi \in \bar{\Omega}^e \), this method computes the jacobian, \( \mathcal{J}(\xi) \). More...
|
|
void | evaluateShapeFunctions (const double *xr, double *phi) |
| Given reference coordinates, \( \xi \in \bar{\Omega}^e \), this method evaluates the shape functions at each degree of freedom, \( \phi_i(\xi)=N_i^e(\xi) \). More...
|
|
void | evaluateDerivatives (const double *xr, double *phidot) |
| Given reference coordinates, \( \xi \in \bar{\Omega}^e \), this method evaluates the first derivatives of the shape function at each degree of freedom, \( \partial\phi_i(\xi)=N_i^e(\xi) \). More...
|
|
|
double * | getPhysicalNodes () |
| Returns a pointer to the nodes of the element in physical space, \( \forall\hat{x_i} \in \Omega^e \). More...
|
|
const double * | getPhysicalNodes () const |
|
|
double * | getReferenceNodes () |
| Returns a pointer to the nodes of the element in reference space, \( \forall\hat{\xi_i} \in \bar{\Omega}^e \). More...
|
|
const double * | getReferenceNodes () const |
|
|
double * | getReferenceCenter () |
| Returns a pointer to the centroid of the reference element \( \xi_c \in \bar{\Omega}^e \). More...
|
|
const double * | getReferenceCenter () const |
|
The FiniteElement object is used to represent a mesh element \( \Omega^e \) corresponding to a mesh \( \mathcal{M} \) .
The FiniteElement is associated with a Finite Element Basis (FEBasis), which defines a set of shape functions \( \mathcal{N}_i( \hat{\xi} ) \), on a corresponding element, \( \bar{\Omega}^e \) in reference space. This class provides the following key operations:
-
Evaluate the shape function and its derivatives
-
Forward/Inverse mapping between physical/reference coordinates
-
Compute the jacobian \( \mathcal{J}(\hat{\xi} ) \)
The FiniteElement class is the primary object that application codes will use to perform these operations. Once a FiniteElement object is instantiated it can then be bound to a Finite Element Basis by invoking the bind_basis( ) method on the target FiniteElement instance, which is templated on two enum values: (a) the BasisType, defined in (FEBasisTypes.hpp) and (b) the CellType defined in(CellType.hpp). The rationale behind this is to insulate the user from the low-level classes and provide a simple and unified interface to Finite Element operations.
A simple example illustrating how to use the FiniteElement class is given below:
double xp[2];
double xr[2];
double wgts[4];
...
mint::Mesh*
m = get_mesh();
const int N =
m->getMeshNumberOfCells();
const bool zero_copy = true;
for ( int idx=0; idx < N; ++idx ) {
...
double coords[ ] = {
x1[ idx ], y1[ idx ],
x2[ idx ], y2[ idx ],
x3[ idx ], y3[ idx ],
x4[ idx ], y4[ idx ],
};
numeric::Matrix< double >
m( 2, 4, coords, zero_copy );
mint::FiniteElement fe( m,
mint::QUAD, zero_copy );
mint::bind_basis< MINT_LAGRANGE_BASIS, mint::QUAD >( fe );
int status = fe.computeReferenceCoords( xp, xr );
fe.evaluateShapeFunction( xr, wgts );
....
}
}
...
@ INSIDE_ELEMENT
Definition: FiniteElement.hpp:32
constexpr CellType QUAD
Definition: CellTypes.hpp:47
- See also
- FEBasis
-
CellType.hpp
-
FEBasisTypes.hpp
-
ShapeFunction