MFEM  v3.3
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
mfem::DGElasticityIntegrator Class Reference

#include <bilininteg.hpp>

Inheritance diagram for mfem::DGElasticityIntegrator:
[legend]
Collaboration diagram for mfem::DGElasticityIntegrator:
[legend]

Public Member Functions

 DGElasticityIntegrator (double alpha_, double kappa_)
 
 DGElasticityIntegrator (Coefficient &lambda_, Coefficient &mu_, double alpha_, double kappa_)
 
virtual void AssembleFaceMatrix (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
 
- Public Member Functions inherited from mfem::BilinearFormIntegrator
virtual void AssembleElementMatrix (const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
 Given a particular Finite Element computes the element matrix elmat. More...
 
virtual void AssembleElementMatrix2 (const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &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 AssembleElementVector (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect)
 Perform the local action of the BilinearFormIntegrator. More...
 
virtual void AssembleElementGrad (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, DenseMatrix &elmat)
 Assemble the local gradient matrix. More...
 
virtual void ComputeElementFlux (const FiniteElement &el, ElementTransformation &Trans, Vector &u, const FiniteElement &fluxelem, Vector &flux, int with_coef=1)
 
virtual double ComputeFluxEnergy (const FiniteElement &fluxelem, ElementTransformation &Trans, Vector &flux, Vector *d_energy=NULL)
 
void SetIntRule (const IntegrationRule *ir)
 
virtual ~BilinearFormIntegrator ()
 
- Public Member Functions inherited from mfem::NonlinearFormIntegrator
virtual double GetElementEnergy (const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun)
 Compute the local energy. More...
 
virtual ~NonlinearFormIntegrator ()
 

Static Protected Member Functions

static void AssembleBlock (const int dim, const int row_ndofs, const int col_ndofs, const int row_offset, const int col_offset, const double jmatcoef, const Vector &col_nL, const Vector &col_nM, const Vector &row_shape, const Vector &col_shape, const Vector &col_dshape_dnM, const DenseMatrix &col_dshape, DenseMatrix &elmat, DenseMatrix &jmat)
 

Protected Attributes

Coefficientlambda
 
Coefficientmu
 
double alpha
 
double kappa
 
Vector shape1
 
Vector shape2
 
DenseMatrix dshape1
 
DenseMatrix dshape2
 
DenseMatrix adjJ
 
DenseMatrix dshape1_ps
 
DenseMatrix dshape2_ps
 
Vector nor
 
Vector nL1
 
Vector nL2
 
Vector nM1
 
Vector nM2
 
Vector dshape1_dnM
 
Vector dshape2_dnM
 
DenseMatrix jmat
 
- Protected Attributes inherited from mfem::BilinearFormIntegrator
const IntegrationRuleIntRule
 

Additional Inherited Members

- Protected Member Functions inherited from mfem::BilinearFormIntegrator
 BilinearFormIntegrator (const IntegrationRule *ir=NULL)
 

Detailed Description

Integrator for the DG elasticity form, for the formulations see:

where $ \left<u, v\right> = \int_{F} u \cdot v $, and $ F $ is a face which is either a boundary face $ F_b $ of an element $ K $ or an interior face $ F_i $ separating elements $ K_1 $ and $ K_2 $.

In the bilinear form above $ \tau(u) $ is traction, and it's also $ \tau(u) = \sigma(u) \cdot \vec{n} $, where $ \sigma(u) $ is stress, and $ \vec{n} $ is the unit normal vector w.r.t. to $ F $.

In other words, we have

\[ - \left< \{ \sigma(u) \cdot \vec{n} \}, [v] \right> + \alpha \left< \{ \sigma(v) \cdot \vec{n} \}, [u] \right> + \kappa \left< h^{-1} \{ \lambda + 2 \mu \} [u], [v] \right> \]

For isotropic media

\[ \begin{split} \sigma(u) &= \lambda \nabla \cdot u I + 2 \mu \varepsilon(u) \\ &= \lambda \nabla \cdot u I + 2 \mu \frac{1}{2} (\nabla u + \nabla u^T) \\ &= \lambda \nabla \cdot u I + \mu (\nabla u + \nabla u^T) \end{split} \]

where $ I $ is identity matrix, $ \lambda $ and $ \mu $ are Lame coefficients (see ElasticityIntegrator), $ u, v $ are the trial and test functions, respectively.

The parameters $ \alpha $ and $ \kappa $ determine the DG method to use (when this integrator is added to the "broken" ElasticityIntegrator):

This is a 'Vector' integrator, i.e. defined for FE spaces using multiple copies of a scalar FE space.

Definition at line 2160 of file bilininteg.hpp.

Constructor & Destructor Documentation

mfem::DGElasticityIntegrator::DGElasticityIntegrator ( double  alpha_,
double  kappa_ 
)
inline

Definition at line 2163 of file bilininteg.hpp.

mfem::DGElasticityIntegrator::DGElasticityIntegrator ( Coefficient lambda_,
Coefficient mu_,
double  alpha_,
double  kappa_ 
)
inline

Definition at line 2166 of file bilininteg.hpp.

Member Function Documentation

void mfem::DGElasticityIntegrator::AssembleBlock ( const int  dim,
const int  row_ndofs,
const int  col_ndofs,
const int  row_offset,
const int  col_offset,
const double  jmatcoef,
const Vector col_nL,
const Vector col_nM,
const Vector row_shape,
const Vector col_shape,
const Vector col_dshape_dnM,
const DenseMatrix col_dshape,
DenseMatrix elmat,
DenseMatrix jmat 
)
staticprotected

Definition at line 2601 of file bilininteg.cpp.

void mfem::DGElasticityIntegrator::AssembleFaceMatrix ( const FiniteElement el1,
const FiniteElement el2,
FaceElementTransformations Trans,
DenseMatrix elmat 
)
virtual

Reimplemented from mfem::BilinearFormIntegrator.

Definition at line 2644 of file bilininteg.cpp.

Member Data Documentation

DenseMatrix mfem::DGElasticityIntegrator::adjJ
protected

Definition at line 2188 of file bilininteg.hpp.

double mfem::DGElasticityIntegrator::alpha
protected

Definition at line 2178 of file bilininteg.hpp.

DenseMatrix mfem::DGElasticityIntegrator::dshape1
protected

Definition at line 2186 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::dshape1_dnM
protected

Definition at line 2196 of file bilininteg.hpp.

DenseMatrix mfem::DGElasticityIntegrator::dshape1_ps
protected

Definition at line 2192 of file bilininteg.hpp.

DenseMatrix mfem::DGElasticityIntegrator::dshape2
protected

Definition at line 2186 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::dshape2_dnM
protected

Definition at line 2196 of file bilininteg.hpp.

DenseMatrix mfem::DGElasticityIntegrator::dshape2_ps
protected

Definition at line 2192 of file bilininteg.hpp.

DenseMatrix mfem::DGElasticityIntegrator::jmat
protected

Definition at line 2198 of file bilininteg.hpp.

double mfem::DGElasticityIntegrator::kappa
protected

Definition at line 2178 of file bilininteg.hpp.

Coefficient* mfem::DGElasticityIntegrator::lambda
protected

Definition at line 2177 of file bilininteg.hpp.

Coefficient * mfem::DGElasticityIntegrator::mu
protected

Definition at line 2177 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::nL1
protected

Definition at line 2194 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::nL2
protected

Definition at line 2194 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::nM1
protected

Definition at line 2195 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::nM2
protected

Definition at line 2195 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::nor
protected

Definition at line 2193 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::shape1
protected

Definition at line 2183 of file bilininteg.hpp.

Vector mfem::DGElasticityIntegrator::shape2
protected

Definition at line 2183 of file bilininteg.hpp.


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