![]() |
MFEM v4.8.0
Finite element discretization library
|
#include <bilininteg.hpp>
Public Member Functions | |
ElasticityIntegrator (Coefficient &l, Coefficient &m) | |
ElasticityIntegrator (Coefficient &m, real_t q_l, real_t q_m) | |
void | AssembleElementMatrix (const FiniteElement &el, ElementTransformation &Tr, DenseMatrix &elmat) override |
Given a particular Finite Element computes the element matrix elmat. | |
void | AssemblePA (const FiniteElementSpace &fes) override |
Method defining partial assembly. | |
void | AssembleDiagonalPA (Vector &diag) override |
Assemble diagonal and add it to Vector diag. | |
void | AddMultPA (const Vector &x, Vector &y) const override |
Method for partially assembled action. | |
void | AddMultTransposePA (const Vector &x, Vector &y) const override |
Method for partially assembled transposed action. | |
void | ComputeElementFlux (const FiniteElement &el, ElementTransformation &Trans, Vector &u, const FiniteElement &fluxelem, Vector &flux, bool with_coef=true, const IntegrationRule *ir=NULL) override |
real_t | ComputeFluxEnergy (const FiniteElement &fluxelem, ElementTransformation &Trans, Vector &flux, Vector *d_energy=NULL) override |
void | AssemblePA (const FiniteElementSpace &trial_fes, const FiniteElementSpace &test_fes) override |
![]() | |
virtual void | AssembleNURBSPA (const FiniteElementSpace &fes) |
Method defining partial assembly on NURBS patches. | |
virtual void | AssemblePABoundary (const FiniteElementSpace &fes) |
virtual void | AssemblePAInteriorFaces (const FiniteElementSpace &fes) |
virtual void | AssemblePABoundaryFaces (const FiniteElementSpace &fes) |
virtual void | AssembleDiagonalPA_ADAt (const Vector &D, Vector &diag) |
Assemble diagonal of \(A D A^T\) ( \(A\) is this integrator) and add it to diag. | |
virtual void | AddMultNURBSPA (const Vector &x, Vector &y) const |
Method for partially assembled action on NURBS patches. | |
virtual void | AssembleEA (const FiniteElementSpace &fes, Vector &emat, const bool add=true) |
Method defining element assembly. | |
void | AssembleMF (const FiniteElementSpace &fes) override |
Method defining matrix-free assembly. | |
void | AddMultMF (const Vector &x, Vector &y) const override |
virtual void | AddMultTransposeMF (const Vector &x, Vector &y) const |
virtual void | AssembleDiagonalMF (Vector &diag) |
Assemble diagonal and add it to Vector diag. | |
virtual void | AssembleEABoundary (const FiniteElementSpace &fes, Vector &ea_data_bdr, const bool add=true) |
virtual void | AssembleEAInteriorFaces (const FiniteElementSpace &fes, Vector &ea_data_int, Vector &ea_data_ext, const bool add=true) |
virtual void | AssembleEAInteriorFaces (const FiniteElementSpace &trial_fes, const FiniteElementSpace &test_fes, Vector &emat, const bool add=true) |
Method defining element assembly for mixed trace integrators. | |
virtual void | AssembleEABoundaryFaces (const FiniteElementSpace &fes, Vector &ea_data_bdr, const bool add=true) |
virtual void | AssembleElementMatrix2 (const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat) |
virtual void | AssemblePatchMatrix (const int patch, const FiniteElementSpace &fes, SparseMatrix *&smat) |
virtual void | AssembleFaceMatrix (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat) |
virtual void | AssembleFaceMatrix (const FiniteElement &trial_fe1, const FiniteElement &test_fe1, const FiniteElement &trial_fe2, const FiniteElement &test_fe2, FaceElementTransformations &Trans, DenseMatrix &elmat) |
virtual void | AssembleFaceMatrix (const FiniteElement &trial_face_fe, const FiniteElement &test_fe1, const FiniteElement &test_fe2, FaceElementTransformations &Trans, DenseMatrix &elmat) |
virtual void | AssembleTraceFaceMatrix (int elem, const FiniteElement &trial_face_fe, const FiniteElement &test_fe, FaceElementTransformations &Trans, DenseMatrix &elmat) |
void | AssembleElementVector (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect) override |
Perform the local action of the BilinearFormIntegrator. Note that the default implementation in the base class is general but not efficient. | |
void | AssembleFaceVector (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr, const Vector &elfun, Vector &elvect) override |
Perform the local action of the BilinearFormIntegrator resulting from a face integral term. Note that the default implementation in the base class is general but not efficient. | |
void | AssembleElementGrad (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, DenseMatrix &elmat) override |
Assemble the local gradient matrix. | |
void | AssembleFaceGrad (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr, const Vector &elfun, DenseMatrix &elmat) override |
Assemble the local action of the gradient of the NonlinearFormIntegrator resulting from a face integral term. | |
virtual bool | RequiresFaceNormalDerivatives () const |
For bilinear forms on element faces, specifies if the normal derivatives are needed on the faces or just the face restriction. | |
virtual void | AddMultPAFaceNormalDerivatives (const Vector &x, const Vector &dxdn, Vector &y, Vector &dydn) const |
Method for partially assembled action. | |
virtual | ~BilinearFormIntegrator () |
![]() | |
void | SetIntegrationMode (Mode m) |
bool | Patchwise () const |
void | SetPAMemoryType (MemoryType mt) |
virtual real_t | GetElementEnergy (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun) |
Compute the local energy. | |
virtual void | AssembleGradPA (const Vector &x, const FiniteElementSpace &fes) |
Prepare the integrator for partial assembly (PA) gradient evaluations on the given FE space fes at the state x. | |
virtual real_t | GetLocalStateEnergyPA (const Vector &x) const |
Compute the local (to the MPI rank) energy with partial assembly. | |
virtual void | AddMultGradPA (const Vector &x, Vector &y) const |
Method for partially assembled gradient action. | |
virtual void | AssembleGradDiagonalPA (Vector &diag) const |
Method for computing the diagonal of the gradient with partial assembly. | |
virtual bool | SupportsCeed () const |
Indicates whether this integrator can use a Ceed backend. | |
ceed::Operator & | GetCeedOp () |
virtual | ~NonlinearFormIntegrator () |
![]() | |
Integrator (const IntegrationRule *ir=NULL) | |
Create a new Integrator, optionally providing a prescribed quadrature rule to use in assembly. | |
virtual void | SetIntRule (const IntegrationRule *ir) |
Prescribe a fixed IntegrationRule to use, or set to null to let the integrator choose an appropriate rule. | |
void | SetIntegrationRule (const IntegrationRule &ir) |
Prescribe a fixed IntegrationRule to use. Sets the NURBS patch integration rule to null. | |
void | SetNURBSPatchIntRule (NURBSMeshRules *pr) |
Sets an integration rule for use on NURBS patches. | |
bool | HasNURBSPatchIntRule () const |
Check if a NURBS patch integration rule has been set. | |
const IntegrationRule * | GetIntRule () const |
Directly return the IntRule pointer (possibly null) without checking for NURBS patch rules or falling back on a default. | |
const IntegrationRule * | GetIntegrationRule () const |
Equivalent to GetIntRule, but retained for backward compatibility with applications. | |
Protected Attributes | |
real_t | q_lambda |
real_t | q_mu |
Coefficient * | lambda |
Coefficient * | mu |
![]() | |
Mode | integrationMode = Mode::ELEMENTWISE |
ceed::Operator * | ceedOp |
MemoryType | pa_mt = MemoryType::DEFAULT |
![]() | |
const IntegrationRule * | IntRule |
NURBSMeshRules * | patchRules = nullptr |
Friends | |
class | ElasticityComponentIntegrator |
Additional Inherited Members | |
![]() | |
enum | Mode { ELEMENTWISE = 0 , PATCHWISE = 1 , PATCHWISE_REDUCED = 2 } |
![]() | |
BilinearFormIntegrator (const IntegrationRule *ir=NULL) | |
![]() | |
NonlinearFormIntegrator (const IntegrationRule *ir=NULL) | |
![]() | |
const IntegrationRule * | GetIntegrationRule (const FiniteElement &trial_fe, const FiniteElement &test_fe, const ElementTransformation &trans) const |
Returns an integration rule based on the arguments and internal state of the Integrator object. | |
const IntegrationRule * | GetIntegrationRule (const FiniteElement &el, const ElementTransformation &trans) const |
Returns an integration rule based on the arguments and internal state. (Version for identical trial_fe and test_fe) | |
virtual const IntegrationRule * | GetDefaultIntegrationRule (const FiniteElement &trial_fe, const FiniteElement &test_fe, const ElementTransformation &trans) const |
Subclasses should override to choose a default integration rule. | |
Integrator for the linear elasticity form:
\[ a(u,v) = (\lambda \mathrm{div}(u), \mathrm{div}(v)) + (2 \mu \varepsilon(u), \varepsilon(v)), \]
where \(\varepsilon(v) = \frac{1}{2} (\mathrm{grad}(v) + \mathrm{grad}(v)^{\mathrm{T}})\). This is a 'Vector' integrator, i.e. defined for FE spaces using multiple copies of a scalar FE space.
Definition at line 3148 of file bilininteg.hpp.
|
inline |
Definition at line 3180 of file bilininteg.hpp.
|
inline |
With this constructor \(\lambda = q_l m\) and \(\mu = q_m m\) if \(dim q_l + 2 q_m = 0\) then \(tr(\sigma) = 0\).
Definition at line 3184 of file bilininteg.hpp.
Method for partially assembled action.
Perform the action of integrator on the input x and add the result to the output y. Both x and y are E-vectors, i.e. they represent the element-wise discontinuous version of the FE space.
This method can be called only after the method AssemblePA() has been called.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 67 of file bilininteg_elasticity_pa.cpp.
|
overridevirtual |
Method for partially assembled transposed action.
Perform the transpose action of integrator on the input x and add the result to the output y. Both x and y are E-vectors, i.e. they represent the element-wise discontinuous version of the FE space.
This method can be called only after the method AssemblePA() has been called.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 73 of file bilininteg_elasticity_pa.cpp.
|
overridevirtual |
Assemble diagonal and add it to Vector diag.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 60 of file bilininteg_elasticity_pa.cpp.
|
overridevirtual |
Given a particular Finite Element computes the element matrix elmat.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 3180 of file bilininteg.cpp.
|
overridevirtual |
Method defining partial assembly.
The result of the partial assembly is stored internally so that it can be used later in the methods AddMultPA() and AddMultTransposePA().
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 40 of file bilininteg_elasticity_pa.cpp.
|
overridevirtual |
Used with BilinearFormIntegrators that have different spaces.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 52 of file bilininteg.cpp.
|
overridevirtual |
Compute the stress corresponding to the local displacement \(u\) and interpolate it at the nodes of the given fluxelem. Only the symmetric part of the stress is stored, so that the size of flux is equal to the number of DOFs in fluxelem times dim*(dim+1)/2. In 2D, the order of the stress components is: \(s_xx, s_yy, s_xy\). In 3D, it is: \(s_xx, s_yy, s_zz, s_xy, s_xz, s_yz\). In other words, flux is the local vector for a FE space with dim*(dim+1)/2 vector components, based on the finite element fluxelem. The integration rule is taken from fluxelem. ir exists to specific an alternative integration rule.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 3262 of file bilininteg.cpp.
|
overridevirtual |
Compute the element energy (integral of the strain energy density) corresponding to the stress represented by flux which is a vector of coefficients multiplying the basis functions defined by fluxelem. In other words, flux is the local vector for a FE space with dim*(dim+1)/2 vector components, based on the finite element fluxelem. The number of components, dim*(dim+1)/2 is such that it represents the symmetric part of the (symmetric) stress tensor. The order of the components is: \(s_xx, s_yy, s_xy\) in 2D, and \(s_xx, s_yy, s_zz, s_xy, s_xz, s_yz\) in 3D.
Reimplemented from mfem::BilinearFormIntegrator.
Definition at line 3341 of file bilininteg.cpp.
|
friend |
Definition at line 3150 of file bilininteg.hpp.
|
protected |
Definition at line 3154 of file bilininteg.hpp.
|
protected |
Definition at line 3154 of file bilininteg.hpp.
|
protected |
Definition at line 3153 of file bilininteg.hpp.
|
protected |
Definition at line 3153 of file bilininteg.hpp.