v0.13.2
Loading...
Searching...
No Matches
Classes | Typedefs | Enumerations | Functions
Forms Integrators

Classes and functions used to evaluate fields at integration pts, jacobians, etc.. More...

Collaboration diagram for Forms Integrators:

Classes

struct  MoFEM::EssentialBC< EleOp >
 Essential boundary conditions. More...
 
struct  MoFEM::EssentialBC< EleOp >::Assembly< A >
 Assembly methods. More...
 
struct  MoFEM::EssentialBC< EleOp >::Assembly< A >::LinearForm< I >
 
struct  MoFEM::EssentialBC< EleOp >::Assembly< A >::BiLinearForm< I >
 
struct  MoFEM::NaturalBC< EleOp >
 Natural boundary conditions. More...
 
struct  MoFEM::NaturalBC< EleOp >::Assembly< A >
 Assembly methods. More...
 
struct  MoFEM::NaturalBC< EleOp >::Assembly< A >::LinearForm< I >
 
struct  MoFEM::NaturalBC< EleOp >::Assembly< A >::BiLinearForm< I >
 
struct  MoFEM::OpSetBc
 Set indices on entities on finite element. More...
 
struct  MoFEM::FormsIntegrators< EleOp >
 Integrator forms. More...
 
struct  MoFEM::FormsIntegrators< EleOp >::Assembly< A >
 Assembly methods. More...
 
struct  MoFEM::FormsIntegrators< EleOp >::Assembly< A >::LinearForm< I >
 Linear form. More...
 
struct  MoFEM::FormsIntegrators< EleOp >::Assembly< A >::BiLinearForm< I >
 Bi linear form. More...
 
struct  MoFEM::FormsIntegrators< EleOp >::Assembly< A >
 Bilinear integrator form. More...
 

Typedefs

template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradGrad = OpGradGradImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase >
 Integrate \((v_{,i},\beta(\mathbf{x}) u_{,j}))_\Omega\). More...
 
template<int BASE_DIM, int FIELD_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpMass = OpMassImpl< BASE_DIM, FIELD_DIM, I, OpBase >
 Integrate \((v_i,\beta(\mathbf{x}) u_j)_\Omega\). More...
 
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradSymTensorGrad = OpGradSymTensorGradImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase >
 Integrate \((v_k,D_{ijkl} u_{,l})_\Omega\). More...
 
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradGradSymTensorGradGrad = OpGradGradSymTensorGradGradImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase >
 Integrate \((v_{,ij},D_{ijkl} u_{,kl})_\Omega\). More...
 
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradTensorGrad = OpGradTensorGradImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase >
 Integrate \((v_{,j},D_{ijkl} u_{,l})_\Omega\). More...
 
using MoFEM::ScalarFun = boost::function< double(const double, const double, const double)>
 Scalar function type. More...
 
template<int DIM>
using MoFEM::VectorFun = boost::function< FTensor::Tensor1< double, DIM >(const double, const double, const double)>
 Vector function type. More...
 
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradTimesTensor = OpGradTimesTensorImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase >
 Integrate \((v,f(\mathbf{x}))_\Omega\), f is a scalar. More...
 
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpConvectiveTermRhs = OpConvectiveTermRhsImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase >
 Integrate \((v_{,i},f_{ij})_\Omega\), f is symmetric tensor. More...
 

Enumerations

enum  MoFEM::AssemblyType { MoFEM::PETSC , MoFEM::SCHUR , MoFEM::USER_ASSEMBLE , MoFEM::LAST_ASSEMBLE }
 [Storage and set boundary conditions] More...
 
enum  MoFEM::IntegrationType { MoFEM::GAUSS , MoFEM::USER_INTEGRATION , MoFEM::LAST_INTEGRATION }
 Form integrator integration types. More...
 

Functions

template<>
MoFEMErrorCode MoFEM::VecSetValues< EssentialBcStorage > (Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
 Set values to vector in operator. More...
 

Detailed Description

Classes and functions used to evaluate fields at integration pts, jacobians, etc..

Typedef Documentation

◆ OpConvectiveTermRhs

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpConvectiveTermRhs = OpConvectiveTermRhsImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase>

Integrate \((v_{,i},f_{ij})_\Omega\), f is symmetric tensor.

Note
\(f_{ij}\$ is tensor at integration points @tparam BASE_DIM @tparam FIELD_DIM @tparam SPACE_DIM */ template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1> using OpGradTimesSymTensor = OpGradTimesSymTensorImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase>; /** @brief Integrate \)(\lambda_{ij,j},u_{i})_\Omega\f$
Template Parameters
BASE_DIM
FIELD_DIM
SPACE_DIM*/ template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM> using OpMixDivTimesU = OpMixDivTimesUImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase>;

/**

Integrate \((\lambda_{ij},u_{i,j})_\Omega\)

Template Parameters
SPACE_DIM*/ template <int SPACE_DIM> using OpMixTensorTimesGradU = OpMixTensorTimesGradUImpl<SPACE_DIM, I, OpBase>;

/**

Integrate \((u_{i},\lambda_{ij,j})_\Omega\)

Template Parameters
SPACE_DIM*/ template <int SPACE_DIM> using OpMixVecTimesDivLambda = OpMixVecTimesDivLambdaImpl<SPACE_DIM, I, OpBase>;

/**

Multiply vector times normal on the face times scalar function

This operator typically will be used to evaluate natural boundary conditions for mixed formulation.

Template Parameters
BASE_DIM
SPACE_DIM
OpBase*/ template <int SPACE_DIM> using OpNormalMixVecTimesScalar = OpNormalMixVecTimesScalarImpl<SPACE_DIM, I, OpBase>;

/**

Convective term

\[ (v, u_i \mathbf{y}_{,i}) \]

where

Definition at line 560 of file LinearFormsIntegrators.hpp.

◆ OpGradGrad

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradGrad = OpGradGradImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase>

Integrate \((v_{,i},\beta(\mathbf{x}) u_{,j}))_\Omega\).

Template Parameters
SPACE_DIM

Definition at line 456 of file BiLinearFormsIntegrators.hpp.

◆ OpGradGradSymTensorGradGrad

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradGradSymTensorGradGrad = OpGradGradSymTensorGradGradImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase>

Integrate \((v_{,ij},D_{ijkl} u_{,kl})_\Omega\).

Note
\(D_{ijkl}\) is * tensor with no symmetries
Template Parameters
SPACE_DIM
S

Definition at line 492 of file BiLinearFormsIntegrators.hpp.

◆ OpGradSymTensorGrad

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradSymTensorGrad = OpGradSymTensorGradImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase>

Integrate \((v_k,D_{ijkl} u_{,l})_\Omega\).

Note
\(D_{ijkl}\) is * tensor with minor & major symmetry
Template Parameters
SPACE_DIM
S

Definition at line 478 of file BiLinearFormsIntegrators.hpp.

◆ OpGradTensorGrad

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradTensorGrad = OpGradTensorGradImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase>

Integrate \((v_{,j},D_{ijkl} u_{,l})_\Omega\).

Note
\(D_{ijkl}\) is * tensor with no symmetries
Template Parameters
SPACE_DIM
S

Definition at line 507 of file BiLinearFormsIntegrators.hpp.

◆ OpGradTimesTensor

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpGradTimesTensor = OpGradTimesTensorImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I, OpBase>

Integrate \((v,f(\mathbf{x}))_\Omega\), f is a scalar.

Note
\(f(\mathbf{x})\$ is scalar function of coordinates @tparam BASE_DIM @tparam FIELD_DIM */ template <int BASE_DIM, int FIELD_DIM, typename S = SourceFunctionSpecialization> using OpSource = OpSourceImpl<BASE_DIM, FIELD_DIM, I, typename S::template S<OpBase>>; /** @brief Vector field integrator \)(v_i,f)_\Omega\f$, f is a vector
Template Parameters
BASE_DIM*/ template <int BASE_DIM, int S = 1> using OpBaseTimesScalar = OpBaseTimesScalarImpl<BASE_DIM, S, I, OpBase>;

/**

Deprecated:
use instead OpBaseTimesScalar */ template <int BASE_DIM, int S = 1> using OpBaseTimesScalarField = OpBaseTimesScalar<BASE_DIM, S>;

/**

Vector field integrator \((v,f_i)_\Omega\), f is a vector

Template Parameters
BASE_DIM
FIELD_DIM
0*/ template <int BASE_DIM, int FIELD_DIM, int S> using OpBaseTimesVector = OpBaseTimesVectorImpl<BASE_DIM, FIELD_DIM, S, I, OpBase>; [Source operator]

[Grad times tensor]

/**

Integrate \((v_{,i},f_{ij})_\Omega\), f is a vector

Note

Definition at line 488 of file LinearFormsIntegrators.hpp.

◆ OpMass

template<typename EleOp >
template<typename EleOp >
template<int BASE_DIM, int FIELD_DIM>
using MoFEM::FormsIntegrators< EleOp >::Assembly< A >::OpMass = OpMassImpl<BASE_DIM, FIELD_DIM, I, OpBase>

Integrate \((v_i,\beta(\mathbf{x}) u_j)_\Omega\).

Template Parameters

Definition at line 465 of file BiLinearFormsIntegrators.hpp.

◆ ScalarFun

using MoFEM::ScalarFun = typedef boost::function<double(const double, const double, const double)>

Scalar function type.

Definition at line 135 of file FormsIntegrators.hpp.

◆ VectorFun

template<int DIM>
using MoFEM::VectorFun = typedef boost::function<FTensor::Tensor1<double, DIM>( const double, const double, const double)>

Vector function type.

Template Parameters
DIMdimension of the return

Definition at line 167 of file FormsIntegrators.hpp.

Enumeration Type Documentation

◆ AssemblyType

[Storage and set boundary conditions]

Form integrator assembly types

Enumerator
PETSC 
SCHUR 
USER_ASSEMBLE 
LAST_ASSEMBLE 

Definition at line 104 of file FormsIntegrators.hpp.

◆ IntegrationType

Form integrator integration types.

Enumerator
GAUSS 
USER_INTEGRATION 
LAST_INTEGRATION 

Definition at line 128 of file FormsIntegrators.hpp.

Function Documentation

◆ VecSetValues< EssentialBcStorage >()

template<>
MoFEMErrorCode MoFEM::VecSetValues< EssentialBcStorage > ( Vec  V,
const EntitiesFieldData::EntData data,
const double ptr,
InsertMode  iora 
)

Set values to vector in operator.

MoFEM::FieldEntity provides MoFEM::FieldEntity::getWeakStoragePtr() storage function which allows to transfer data between FEs or operators processing the same entities.

When MoFEM::OpSetBc is pushed in weak storage indices taking in account indices which are skip to take boundary conditions are stored. Those entities are used by VecSetValues.

Parameters
V
data
ptr
iora
Returns
MoFEMErrorCode
Parameters
V
data
ptr
iora
Returns
MoFEMErrorCode

Definition at line 75 of file FormsIntegrators.cpp.

78 {
80
81 if(!V)
82 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
83 "Pointer to PETSc vector is null");
84
85 CHKERR VecSetOption(V, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
86 if (!data.getFieldEntities().empty()) {
87 if (auto e_ptr = data.getFieldEntities()[0]) {
88 if (auto stored_data_ptr =
89 e_ptr->getSharedStoragePtr<EssentialBcStorage>()) {
90 return ::VecSetValues(V, stored_data_ptr->entityIndices.size(),
91 &*stored_data_ptr->entityIndices.begin(), ptr,
92 iora);
93 }
94 }
95 }
96 return ::VecSetValues(V, data.getIndices().size(),
97 &*data.getIndices().begin(), ptr, iora);
99}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_INVALID_DATA
Definition: definitions.h:36
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535