v0.14.0
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Protected Attributes | Friends | List of all members
MoFEM::ForcesAndSourcesCore::UserDataOperator Struct Reference

#include <src/finite_elements/ForcesAndSourcesCore.hpp>

Inheritance diagram for MoFEM::ForcesAndSourcesCore::UserDataOperator:
[legend]
Collaboration diagram for MoFEM::ForcesAndSourcesCore::UserDataOperator:
[legend]

Public Types

enum  OpType {
  OPROW = 1 << 0 , OPCOL = 1 << 1 , OPROWCOL = 1 << 2 , OPSPACE = 1 << 3 ,
  OPLAST = 1 << 3
}
 Controls loop over entities on element. More...
 
- Public Types inherited from MoFEM::DataOperator
using DoWorkLhsHookFunType = boost::function< MoFEMErrorCode(DataOperator *op_ptr, int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)>
 
using DoWorkRhsHookFunType = boost::function< MoFEMErrorCode(DataOperator *op_ptr, int side, EntityType type, EntitiesFieldData::EntData &data)>
 

Public Member Functions

 UserDataOperator (const FieldSpace space, const char type=OPSPACE, const bool symm=true)
 
 UserDataOperator (const std::string field_name, const char type, const bool symm=true)
 
 UserDataOperator (const std::string row_field_name, const std::string col_field_name, const char type, const bool symm=true)
 
boost::shared_ptr< const NumeredEntFiniteElementgetNumeredEntFiniteElementPtr () const
 Return raw pointer to NumeredEntFiniteElement. More...
 
EntityHandle getFEEntityHandle () const
 Return finite element entity handle. More...
 
int getFEDim () const
 Get dimension of finite element. More...
 
EntityType getFEType () const
 Get dimension of finite element. More...
 
boost::weak_ptr< SideNumbergetSideNumberPtr (const int side_number, const EntityType type)
 Get the side number pointer. More...
 
EntityHandle getSideEntity (const int side_number, const EntityType type)
 Get the side entity. More...
 
int getNumberOfNodesOnElement () const
 Get the number of nodes on finite element. More...
 
MoFEMErrorCode getProblemRowIndices (const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
 Get row indices. More...
 
MoFEMErrorCode getProblemColIndices (const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
 Get col indices. More...
 
const FEMethodgetFEMethod () const
 Return raw pointer to Finite Element Method object. More...
 
int getOpType () const
 Get operator types. More...
 
void setOpType (const OpType type)
 Set operator type. More...
 
void addOpType (const OpType type)
 Add operator type. More...
 
int getNinTheLoop () const
 get number of finite element in the loop More...
 
int getLoopSize () const
 get size of elements in the loop More...
 
std::string getFEName () const
 Get name of the element. More...
 
ForcesAndSourcesCoregetPtrFE () const
 
ForcesAndSourcesCoregetSidePtrFE () const
 
ForcesAndSourcesCoregetRefinePtrFE () const
 
Accessing KSP
const PetscData::SwitchesgetDataCtx () const
 
const KspMethod::KSPContext getKSPCtx () const
 
const SnesMethod::SNESContext getSNESCtx () const
 
const TSMethod::TSContext getTSCtx () const
 
Vec getKSPf () const
 
Mat getKSPA () const
 
Mat getKSPB () const
 
Accessing SNES
Vec getSNESf () const
 
Vec getSNESx () const
 
Mat getSNESA () const
 
Mat getSNESB () const
 
Accessing TS
Vec getTSu () const
 
Vec getTSu_t () const
 
Vec getTSu_tt () const
 
Vec getTSf () const
 
Mat getTSA () const
 
Mat getTSB () const
 
int getTSstep () const
 
double getTStime () const
 
double getTStimeStep () const
 
double getTSa () const
 
double getTSaa () const
 
Base funtions and integration points
MatrixDoublegetGaussPts ()
 matrix of integration (Gauss) points for Volume Element More...
 
auto getFTensor0IntegrationWeight ()
 Get integration weights. More...
 
Coordinates and access to internal data
MatrixDoublegetCoordsAtGaussPts ()
 Gauss points and weight, matrix (nb. of points x 3) More...
 
auto getFTensor1CoordsAtGaussPts ()
 Get coordinates at integration points assuming linear geometry. More...
 
Measures (area, volume, length, etc.)
double getMeasure () const
 get measure of element More...
 
doublegetMeasure ()
 get measure of element More...
 
- Public Member Functions inherited from MoFEM::DataOperator
 DataOperator (const bool symm=true)
 
virtual ~DataOperator ()=default
 
virtual MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 Operator for bi-linear form, usually to calculate values on left hand side. More...
 
virtual MoFEMErrorCode opLhs (EntitiesFieldData &row_data, EntitiesFieldData &col_data)
 
virtual MoFEMErrorCode doWork (int side, EntityType type, EntitiesFieldData::EntData &data)
 Operator for linear form, usually to calculate values on right hand side. More...
 
virtual MoFEMErrorCode opRhs (EntitiesFieldData &data, const bool error_if_no_base=false)
 
bool getSymm () const
 Get if operator uses symmetry of DOFs or not. More...
 
void setSymm ()
 set if operator is executed taking in account symmetry More...
 
void unSetSymm ()
 unset if operator is executed for non symmetric problem More...
 

Public Attributes

char opType
 
std::string rowFieldName
 
std::string colFieldName
 
FieldSpace sPace
 
- Public Attributes inherited from MoFEM::DataOperator
DoWorkLhsHookFunType doWorkLhsHook
 
DoWorkRhsHookFunType doWorkRhsHook
 
bool sYmm
 If true assume that matrix is symmetric structure. More...
 
std::array< bool, MBMAXTYPE > doEntities
 If true operator is executed for entity. More...
 
booldoVertices
 \deprectaed If false skip vertices More...
 
booldoEdges
 \deprectaed If false skip edges More...
 
booldoQuads
 \deprectaed More...
 
booldoTris
 \deprectaed More...
 
booldoTets
 \deprectaed More...
 
booldoPrisms
 \deprectaed More...
 

Static Public Attributes

static const char *const OpTypeNames []
 

Protected Member Functions

virtual MoFEMErrorCode setPtrFE (ForcesAndSourcesCore *ptr)
 

Protected Attributes

ForcesAndSourcesCoreptrFE
 

Friends

class ForcesAndSourcesCore
 
class EdgeElementForcesAndSourcesCore
 
class FaceElementForcesAndSourcesCore
 
class ContactPrismElementForcesAndSourcesCore
 

Loops

}

using AdjCache = std::map< EntityHandle, std::vector< boost::weak_ptr< NumeredEntFiniteElement > > >
 
MoFEMErrorCode loopSide (const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy, AdjCache *adj_cache=nullptr)
 User calls this function to loop over elements on the side of face. This function calls finite element with its operator to do calculations. More...
 
MoFEMErrorCode loopThis (const string &fe_name, ForcesAndSourcesCore *this_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
 User calls this function to loop over the same element using a different set of integration points. This function calls finite element with its operator to do calculations. More...
 
MoFEMErrorCode loopParent (const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
 User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations. More...
 
MoFEMErrorCode loopChildren (const string &fe_name, ForcesAndSourcesCore *child_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
 User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations. More...
 

Detailed Description

Examples
EshelbianPlasticity.cpp, dynamic_first_order_con_law.cpp, free_surface.cpp, hanging_node_approx.cpp, hello_world.cpp, level_set.cpp, mesh_smoothing.cpp, and plastic.cpp.

Definition at line 549 of file ForcesAndSourcesCore.hpp.

Member Typedef Documentation

◆ AdjCache

using MoFEM::ForcesAndSourcesCore::UserDataOperator::AdjCache = std::map<EntityHandle, std::vector<boost::weak_ptr<NumeredEntFiniteElement> >>

Definition at line 902 of file ForcesAndSourcesCore.hpp.

Member Enumeration Documentation

◆ OpType

enum MoFEM::ForcesAndSourcesCore::UserDataOperator::OpType

Controls loop over entities on element.

  • OPROW is used if row vector is assembled
  • OPCOL is usually used if column vector is assembled
  • OPROWCOL is usually used for assemble matrices.
  • OPSPACE no field is defined for such operator. Is usually used to modify base

For typical problem like Bubnov-Galerkin OPROW and OPCOL are the same. In more general case for example for non-square matrices columns and rows could have different numeration and/or different set of fields.

Enumerator
OPROW 

operator doWork function is executed on FE rows

OPCOL 

operator doWork function is executed on FE columns

OPROWCOL 

operator doWork is executed on FE rows &columns

OPSPACE 

operator do Work is execute on space data

OPLAST 
Deprecated:
would be removed

Definition at line 566 of file ForcesAndSourcesCore.hpp.

566 {
567 OPROW = 1 << 0, ///< operator doWork function is executed on FE rows
568 OPCOL = 1 << 1, ///< operator doWork function is executed on FE columns
569 OPROWCOL = 1 << 2, ///< operator doWork is executed on FE rows &columns
570 OPSPACE = 1 << 3, ///< operator do Work is execute on space data
571 OPLAST = 1 << 3 ///< @deprecated would be removed
572 };
@ OPCOL
operator doWork function is executed on FE columns
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
@ OPSPACE
operator do Work is execute on space data

Constructor & Destructor Documentation

◆ UserDataOperator() [1/3]

ForcesAndSourcesCore::UserDataOperator::UserDataOperator ( const FieldSpace  space,
const char  type = OPSPACE,
const bool  symm = true 
)

This Constructor is used typically when some modification base shape functions on some approx. space is applied. Operator is run for all data on space.

User has no access to field data from this operator.

Definition at line 1942 of file ForcesAndSourcesCore.cpp.

◆ UserDataOperator() [2/3]

ForcesAndSourcesCore::UserDataOperator::UserDataOperator ( const std::string  field_name,
const char  type,
const bool  symm = true 
)

Definition at line 1947 of file ForcesAndSourcesCore.cpp.

1949 : DataOperator(symm), opType(type), rowFieldName(field_name),
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:89
constexpr auto field_name

◆ UserDataOperator() [3/3]

ForcesAndSourcesCore::UserDataOperator::UserDataOperator ( const std::string  row_field_name,
const std::string  col_field_name,
const char  type,
const bool  symm = true 
)

Definition at line 1952 of file ForcesAndSourcesCore.cpp.

1955 : DataOperator(symm), opType(type), rowFieldName(row_field_name),
1956 colFieldName(col_field_name), sPace(LASTSPACE), ptrFE(nullptr) {}

Member Function Documentation

◆ addOpType()

void ForcesAndSourcesCore::UserDataOperator::addOpType ( const OpType  type)
inline

Add operator type.

Definition at line 1051 of file ForcesAndSourcesCore.hpp.

1051 {
1052 opType |= type;
1053}

◆ getCoordsAtGaussPts()

MatrixDouble & ForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts ( )
inline

◆ getDataCtx()

const PetscData::Switches & ForcesAndSourcesCore::UserDataOperator::getDataCtx ( ) const
inline

Definition at line 1068 of file ForcesAndSourcesCore.hpp.

1068 {
1069 return getFEMethod()->data_ctx;
1070}
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getFEDim()

int ForcesAndSourcesCore::UserDataOperator::getFEDim ( ) const
inline

Get dimension of finite element.

Returns
int

Definition at line 1007 of file ForcesAndSourcesCore.hpp.

1007 {
1009};
auto dimension_from_handle(const EntityHandle h)
get entity dimension form handle
Definition: Templates.hpp:1885
EntityHandle getFEEntityHandle() const
Return finite element entity handle.

◆ getFEEntityHandle()

EntityHandle ForcesAndSourcesCore::UserDataOperator::getFEEntityHandle ( ) const
inline

Return finite element entity handle.

Returns
Finite element entity handle
Examples
HookeInternalStressElement.hpp, UnsaturatedFlow.hpp, and mixed_poisson.cpp.

Definition at line 1003 of file ForcesAndSourcesCore.hpp.

1003 {
1004 return getNumeredEntFiniteElementPtr()->getEnt();
1005}
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.

◆ getFEMethod()

const FEMethod * ForcesAndSourcesCore::UserDataOperator::getFEMethod ( ) const
inline

◆ getFEName()

std::string ForcesAndSourcesCore::UserDataOperator::getFEName ( ) const
inline

Get name of the element.

Definition at line 1063 of file ForcesAndSourcesCore.hpp.

1063 {
1064 return getFEMethod()->getFEName();
1065}
auto getFEName() const
get finite element name

◆ getFEType()

EntityType ForcesAndSourcesCore::UserDataOperator::getFEType ( ) const
inline

Get dimension of finite element.

Returns
int

Definition at line 1011 of file ForcesAndSourcesCore.hpp.

1011 {
1013};
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1869

◆ getFTensor0IntegrationWeight()

auto ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight ( )
inline

◆ getFTensor1CoordsAtGaussPts()

auto ForcesAndSourcesCore::UserDataOperator::getFTensor1CoordsAtGaussPts ( )
inline

Get coordinates at integration points assuming linear geometry.

auto t_coords = getFTensor1CoordsAtGaussPts();
for(int gg = 0;gg!=nb_int_ptrs;gg++) {
// do something
++t_coords;
}
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
Examples
ContactOps.hpp, PoissonDiscontinousGalerkin.hpp, UnsaturatedFlow.hpp, child_and_parent.cpp, dg_projection.cpp, hanging_node_approx.cpp, higher_derivatives.cpp, mixed_poisson.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, and simple_interface.cpp.

Definition at line 1268 of file ForcesAndSourcesCore.hpp.

1268 {
1270 &getCoordsAtGaussPts()(0, 0), &getCoordsAtGaussPts()(0, 1),
1271 &getCoordsAtGaussPts()(0, 2));
1272}
MatrixDouble & getCoordsAtGaussPts()
Gauss points and weight, matrix (nb. of points x 3)

◆ getGaussPts()

MatrixDouble & ForcesAndSourcesCore::UserDataOperator::getGaussPts ( )
inline

matrix of integration (Gauss) points for Volume Element

For triangle: columns 0,1 are x,y coordinates respectively and column 2 is a weight value for example getGaussPts()(1,13) returns y coordinate of 13th Gauss point on particular volume element

For tetrahedron: columns 0,1,2 are x,y,z coordinates respectively and column 3 is a weight value for example getGaussPts()(1,13) returns y coordinate of 13th Gauss point on particular volume element

Examples
ContactOps.hpp, MagneticElement.hpp, PlasticOpsGeneric.hpp, PoissonDiscontinousGalerkin.hpp, ThermoElasticOps.hpp, UnsaturatedFlow.hpp, approx_sphere.cpp, child_and_parent.cpp, continuity_check_on_skeleton_3d.cpp, dg_projection.cpp, free_surface.cpp, hanging_node_approx.cpp, higher_derivatives.cpp, level_set.cpp, mixed_poisson.cpp, phase.cpp, photon_diffusion.cpp, plate.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, and simple_interface.cpp.

Definition at line 1235 of file ForcesAndSourcesCore.hpp.

1235 {
1236 return static_cast<ForcesAndSourcesCore *>(ptrFE)->gaussPts;
1237}
MatrixDouble gaussPts
Matrix of integration points.

◆ getKSPA()

Mat ForcesAndSourcesCore::UserDataOperator::getKSPA ( ) const
inline
Examples
plastic.cpp.

Definition at line 1095 of file ForcesAndSourcesCore.hpp.

1095 {
1096#ifndef NDEBUG
1097 if (getFEMethod()->ksp_A == PETSC_NULL)
1098 THROW_MESSAGE("KSP not set A vector");
1099#endif
1100 return getFEMethod()->ksp_A;
1101}
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561

◆ getKSPB()

Mat ForcesAndSourcesCore::UserDataOperator::getKSPB ( ) const
inline
Examples
EshelbianPlasticity.cpp, PoissonDiscontinousGalerkin.hpp, plastic.cpp, and plate.cpp.

Definition at line 1103 of file ForcesAndSourcesCore.hpp.

1103 {
1104#ifndef NDEBUG
1105 if (getFEMethod()->ksp_B == PETSC_NULL)
1106 THROW_MESSAGE("KSP not set B vector");
1107#endif
1108 return getFEMethod()->ksp_B;
1109}

◆ getKSPCtx()

const KspMethod::KSPContext ForcesAndSourcesCore::UserDataOperator::getKSPCtx ( ) const
inline

Definition at line 1073 of file ForcesAndSourcesCore.hpp.

1073 {
1074 return getFEMethod()->ksp_ctx;
1075}
KSPContext ksp_ctx
Context.
Definition: LoopMethods.hpp:86

◆ getKSPf()

Vec ForcesAndSourcesCore::UserDataOperator::getKSPf ( ) const
inline
Examples
PoissonDiscontinousGalerkin.hpp.

Definition at line 1087 of file ForcesAndSourcesCore.hpp.

1087 {
1088#ifndef NDEBUG
1089 if (getFEMethod()->ksp_f == PETSC_NULL)
1090 THROW_MESSAGE("KSP not set F vector");
1091#endif
1092 return getFEMethod()->ksp_f;
1093}

◆ getLoopSize()

int ForcesAndSourcesCore::UserDataOperator::getLoopSize ( ) const
inline

get size of elements in the loop

Returns
loop size

Definition at line 1059 of file ForcesAndSourcesCore.hpp.

1059 {
1060 return getFEMethod()->getLoopSize();
1061}
int getLoopSize() const
get loop size

◆ getMeasure() [1/2]

double & ForcesAndSourcesCore::UserDataOperator::getMeasure ( )
inline

get measure of element

Returns
volume

Definition at line 1278 of file ForcesAndSourcesCore.hpp.

1278 {
1279 return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1280}

◆ getMeasure() [2/2]

double ForcesAndSourcesCore::UserDataOperator::getMeasure ( ) const
inline

◆ getNinTheLoop()

int ForcesAndSourcesCore::UserDataOperator::getNinTheLoop ( ) const
inline

get number of finite element in the loop

Returns
number of finite element

Definition at line 1055 of file ForcesAndSourcesCore.hpp.

1055 {
1056 return getFEMethod()->getNinTheLoop();
1057}
int getNinTheLoop() const
get number of evaluated element in the loop

◆ getNumberOfNodesOnElement()

int ForcesAndSourcesCore::UserDataOperator::getNumberOfNodesOnElement ( ) const
inline

Get the number of nodes on finite element.

Returns
int
Examples
boundary_marker.cpp.

Definition at line 1037 of file ForcesAndSourcesCore.hpp.

1037 {
1038 return ptrFE->getNumberOfNodes();
1039}
auto getNumberOfNodes() const

◆ getNumeredEntFiniteElementPtr()

boost::shared_ptr< const NumeredEntFiniteElement > ForcesAndSourcesCore::UserDataOperator::getNumeredEntFiniteElementPtr ( ) const
inline

Return raw pointer to NumeredEntFiniteElement.

Examples
ElasticityMixedFormulation.hpp, NavierStokesElement.cpp, and UnsaturatedFlow.hpp.

Definition at line 999 of file ForcesAndSourcesCore.hpp.

999 {
1001};
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getOpType()

int ForcesAndSourcesCore::UserDataOperator::getOpType ( ) const
inline

Get operator types.

Returns
Return operator type

Definition at line 1045 of file ForcesAndSourcesCore.hpp.

1045{ return opType; }

◆ getProblemColIndices()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::getProblemColIndices ( const std::string  filed_name,
const EntityType  type,
const int  side,
VectorInt indices 
) const

Get col indices.

Field could be or not declared for this element but is declared for problem

Parameters
field_name
typeentity type
sideside number, any number if type is MBVERTEX
Returns
indices

NOTE: Using those indices to assemble matrix will result in error if new non-zero values need to be created.

Definition at line 1637 of file ForcesAndSourcesCore.cpp.

1639 {
1641 if (ptrFE == NULL)
1642 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1643
1644 switch (type) {
1645 case MBVERTEX:
1647 break;
1648 default:
1649 CHKERR ptrFE->getProblemTypeColIndices(field_name, type, side, indices);
1650 }
1652}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#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
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const

◆ getProblemRowIndices()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::getProblemRowIndices ( const std::string  filed_name,
const EntityType  type,
const int  side,
VectorInt indices 
) const

Get row indices.

Field could be or not declared for this element but is declared for problem

Parameters
field_name
typeentity type
sideside number, any number if type is MBVERTEX
Returns
indices

NOTE: Using those indices to assemble matrix will result in error if new non-zero values need to be created.

Definition at line 1620 of file ForcesAndSourcesCore.cpp.

1622 {
1624 if (ptrFE == NULL)
1625 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1626
1627 switch (type) {
1628 case MBVERTEX:
1630 break;
1631 default:
1632 CHKERR ptrFE->getProblemTypeRowIndices(field_name, type, side, indices);
1633 }
1635}
MoFEMErrorCode getProblemNodesRowIndices(const std::string &field_name, VectorInt &nodes_indices) const
MoFEMErrorCode getProblemTypeRowIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const

◆ getPtrFE()

ForcesAndSourcesCore * ForcesAndSourcesCore::UserDataOperator::getPtrFE ( ) const
inline

Definition at line 1250 of file ForcesAndSourcesCore.hpp.

1250 {
1251 return ptrFE;
1252}

◆ getRefinePtrFE()

ForcesAndSourcesCore * ForcesAndSourcesCore::UserDataOperator::getRefinePtrFE ( ) const
inline

Definition at line 1260 of file ForcesAndSourcesCore.hpp.

1260 {
1261 return ptrFE->refinePtrFE;
1262}
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.

◆ getSideEntity()

EntityHandle ForcesAndSourcesCore::UserDataOperator::getSideEntity ( const int  side_number,
const EntityType  type 
)
inline

Get the side entity.

Note
For vertex is expectation. Side basses in argument of function doWork is zero. For other entity types side can be used as argument of this function.
if (type == MBVERTEX) {
for (int n = 0; n != number_of_nodes; ++n)
// Do somthing
} else {
EntityHandle ent = getSideEntity(side, type);
// Do somthing
}
}
FTensor::Index< 'n', SPACE_DIM > n
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
Data on single entity (This is passed as argument to DataOperator::doWork)
EntityHandle getSideEntity(const int side_number, const EntityType type)
Get the side entity.
Parameters
side_number
type
Examples
boundary_marker.cpp.

Definition at line 1029 of file ForcesAndSourcesCore.hpp.

1030 {
1031 if (auto side_ptr = getSideNumberPtr(side_number, type).lock())
1032 return side_ptr->ent;
1033 else
1034 return 0;
1035}
boost::weak_ptr< SideNumber > getSideNumberPtr(const int side_number, const EntityType type)
Get the side number pointer.

◆ getSideNumberPtr()

boost::weak_ptr< SideNumber > ForcesAndSourcesCore::UserDataOperator::getSideNumberPtr ( const int  side_number,
const EntityType  type 
)
inline

Get the side number pointer.

Note
For vertex is expectation. Side basses in argument of function doWork is zero. For other entity types side can be used as argument of this function.
Parameters
side_number
type
Returns
boost::weak_ptr<SideNumber>

Definition at line 1016 of file ForcesAndSourcesCore.hpp.

1017 {
1018 auto &side_table_by_side_and_type =
1019 ptrFE->numeredEntFiniteElementPtr->getSideNumberTable().get<1>();
1020 auto side_it =
1021 side_table_by_side_and_type.find(boost::make_tuple(type, side_number));
1022 if (side_it != side_table_by_side_and_type.end())
1023 return *side_it;
1024 else
1025 return boost::weak_ptr<SideNumber>();
1026}

◆ getSidePtrFE()

ForcesAndSourcesCore * ForcesAndSourcesCore::UserDataOperator::getSidePtrFE ( ) const
inline

Definition at line 1255 of file ForcesAndSourcesCore.hpp.

1255 {
1256 return ptrFE->sidePtrFE;
1257}
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.

◆ getSNESA()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESA ( ) const
inline

Definition at line 1127 of file ForcesAndSourcesCore.hpp.

1127 {
1128#ifndef NDEBUG
1129 if (getFEMethod()->snes_A == PETSC_NULL)
1130 THROW_MESSAGE("SNES not set A vector");
1131#endif
1132 return getFEMethod()->snes_A;
1133}
Mat & snes_A
jacobian matrix

◆ getSNESB()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESB ( ) const
inline

Definition at line 1135 of file ForcesAndSourcesCore.hpp.

1135 {
1136#ifndef NDEBUG
1137 if (getFEMethod()->snes_B == PETSC_NULL)
1138 THROW_MESSAGE("SNES not set A matrix");
1139#endif
1140 return getFEMethod()->snes_B;
1141}
Mat & snes_B
preconditioner of jacobian matrix

◆ getSNESCtx()

const SnesMethod::SNESContext ForcesAndSourcesCore::UserDataOperator::getSNESCtx ( ) const
inline

Definition at line 1078 of file ForcesAndSourcesCore.hpp.

1078 {
1079 return getFEMethod()->snes_ctx;
1080}
SNESContext snes_ctx

◆ getSNESf()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESf ( ) const
inline

Definition at line 1111 of file ForcesAndSourcesCore.hpp.

1111 {
1112#ifndef NDEBUG
1113 if (getFEMethod()->snes_f == PETSC_NULL)
1114 THROW_MESSAGE("SNES not set F vector");
1115#endif
1116 return getFEMethod()->snes_f;
1117}
Vec & snes_f
residual

◆ getSNESx()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESx ( ) const
inline

Definition at line 1119 of file ForcesAndSourcesCore.hpp.

1119 {
1120#ifndef NDEBUG
1121 if (getFEMethod()->snes_x == PETSC_NULL)
1122 THROW_MESSAGE("SNESnot set X vector");
1123#endif
1124 return getFEMethod()->snes_x;
1125}
Vec & snes_x
state vector

◆ getTSA()

Mat ForcesAndSourcesCore::UserDataOperator::getTSA ( ) const
inline

Definition at line 1175 of file ForcesAndSourcesCore.hpp.

1175 {
1176#ifndef NDEBUG
1177 if (getFEMethod()->ts_A == PETSC_NULL)
1178 THROW_MESSAGE("TS not set A matrix");
1179#endif
1180 return getFEMethod()->ts_A;
1181}

◆ getTSa()

double ForcesAndSourcesCore::UserDataOperator::getTSa ( ) const
inline

Definition at line 1215 of file ForcesAndSourcesCore.hpp.

1215 {
1216#ifndef NDEBUG
1218 .none() ||
1220 THROW_MESSAGE("TS not set B matrix");
1221#endif
1222 return getFEMethod()->ts_a;
1223}
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:37
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:40
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:38
PetscReal ts_a
shift for U_t (see PETSc Time Solver)

◆ getTSaa()

double ForcesAndSourcesCore::UserDataOperator::getTSaa ( ) const
inline

Definition at line 1225 of file ForcesAndSourcesCore.hpp.

1225 {
1226#ifndef NDEBUG
1228 .none() ||
1230 THROW_MESSAGE("TS not set B matrix");
1231#endif
1232 return getFEMethod()->ts_aa;
1233}
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:41
PetscReal ts_aa
shift for U_tt shift for U_tt

◆ getTSB()

Mat ForcesAndSourcesCore::UserDataOperator::getTSB ( ) const
inline

Definition at line 1183 of file ForcesAndSourcesCore.hpp.

1183 {
1184#ifndef NDEBUG
1185 if (getFEMethod()->ts_B == PETSC_NULL)
1186 THROW_MESSAGE("TS not set B matrix");
1187#endif
1188 return getFEMethod()->ts_B;
1189}
Mat & ts_B
Preconditioner for ts_A.

◆ getTSCtx()

const TSMethod::TSContext ForcesAndSourcesCore::UserDataOperator::getTSCtx ( ) const
inline
Examples
PlasticOpsGeneric.hpp.

Definition at line 1083 of file ForcesAndSourcesCore.hpp.

1083 {
1084 return getFEMethod()->ts_ctx;
1085}
TSContext ts_ctx

◆ getTSf()

Vec ForcesAndSourcesCore::UserDataOperator::getTSf ( ) const
inline

Definition at line 1167 of file ForcesAndSourcesCore.hpp.

1167 {
1168#ifndef NDEBUG
1169 if (getFEMethod()->ts_F == PETSC_NULL)
1170 THROW_MESSAGE("TS not set F vector");
1171#endif
1172 return getFEMethod()->ts_F;
1173}
Vec & ts_F
residual vector

◆ getTSstep()

int ForcesAndSourcesCore::UserDataOperator::getTSstep ( ) const
inline

Definition at line 1191 of file ForcesAndSourcesCore.hpp.

1191 {
1192#ifndef NDEBUG
1193 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1194 THROW_MESSAGE("TS not set time");
1195#endif
1196 return getFEMethod()->ts_step;
1197}
PetscInt ts_step
time step number

◆ getTStime()

double ForcesAndSourcesCore::UserDataOperator::getTStime ( ) const
inline
Examples
ContactOps.hpp.

Definition at line 1199 of file ForcesAndSourcesCore.hpp.

1199 {
1200#ifndef NDEBUG
1201 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1202 THROW_MESSAGE("TS not set time");
1203#endif
1204 return getFEMethod()->ts_t;
1205}
PetscReal ts_t
time

◆ getTStimeStep()

double ForcesAndSourcesCore::UserDataOperator::getTStimeStep ( ) const
inline

Definition at line 1207 of file ForcesAndSourcesCore.hpp.

1207 {
1208#ifndef NDEBUG
1209 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1210 THROW_MESSAGE("TS not set time");
1211#endif
1212 return getFEMethod()->ts_dt;
1213}
PetscReal ts_dt
time step size

◆ getTSu()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu ( ) const
inline

Definition at line 1143 of file ForcesAndSourcesCore.hpp.

1143 {
1144#ifndef NDEBUG
1145 if (getFEMethod()->ts_u == PETSC_NULL)
1146 THROW_MESSAGE("TS not set U vector");
1147#endif
1148 return getFEMethod()->ts_u;
1149}
Vec & ts_u
state vector

◆ getTSu_t()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_t ( ) const
inline

Definition at line 1151 of file ForcesAndSourcesCore.hpp.

1151 {
1152#ifndef NDEBUG
1153 if (getFEMethod()->ts_u_t == PETSC_NULL)
1154 THROW_MESSAGE("TS not set U_t vector");
1155#endif
1156 return getFEMethod()->ts_u_t;
1157}
Vec & ts_u_t
time derivative of state vector

◆ getTSu_tt()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_tt ( ) const
inline

Definition at line 1159 of file ForcesAndSourcesCore.hpp.

1159 {
1160#ifndef NDEBUG
1161 if (getFEMethod()->ts_u_tt == PETSC_NULL)
1162 THROW_MESSAGE("TS not set U_tt vector");
1163#endif
1164 return getFEMethod()->ts_u_tt;
1165}
Vec & ts_u_tt
second time derivative of state vector

◆ loopChildren()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::loopChildren ( const string &  fe_name,
ForcesAndSourcesCore child_fe,
const int  verb = QUIET,
const LogManager::SeverityLevel  sev = Sev::noisy 
)

User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations.

Parameters
fe_name
child_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1839 of file ForcesAndSourcesCore.cpp.

1841 {
1843
1844 auto fe_miit = ptrFE->mField.get_finite_elements()
1845 ->get<FiniteElement_name_mi_tag>()
1846 .find(fe_name);
1847 if (fe_miit != ptrFE->mField.get_finite_elements()
1848 ->get<FiniteElement_name_mi_tag>()
1849 .end()) {
1850
1851 const auto *problem_ptr = getFEMethod()->problemPtr;
1852 auto &ref_ents = *getPtrFE()->mField.get_ref_ents();
1853 auto &numered_fe =
1854 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1855
1856 const auto parent_ent = getNumeredEntFiniteElementPtr()->getEnt();
1857 const auto parent_type = getNumeredEntFiniteElementPtr()->getEntType();
1858 auto range =
1859 ref_ents.get<Composite_ParentEnt_And_EntType_mi_tag>().equal_range(
1860 boost::make_tuple(parent_type, parent_ent));
1861
1862 if (auto size = std::distance(range.first, range.second)) {
1863
1864 std::vector<EntityHandle> childs_vec;
1865 childs_vec.reserve(size);
1866 for (; range.first != range.second; ++range.first)
1867 childs_vec.emplace_back((*range.first)->getEnt());
1868
1869 Range childs;
1870
1871 if ((childs_vec.back() - childs_vec.front() + 1) == size)
1872 childs = Range(childs_vec.front(), childs_vec.back());
1873 else
1874 childs.insert_list(childs_vec.begin(), childs_vec.end());
1875
1876 child_fe->feName = fe_name;
1877
1878 CHKERR child_fe->setRefineFEPtr(ptrFE);
1879 CHKERR child_fe->copyBasicMethod(*getFEMethod());
1880 CHKERR child_fe->copyPetscData(*getFEMethod());
1881 CHKERR child_fe->copyKsp(*getFEMethod());
1882 CHKERR child_fe->copySnes(*getFEMethod());
1883 CHKERR child_fe->copyTs(*getFEMethod());
1884
1885 child_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1886
1887 CHKERR child_fe->preProcess();
1888
1889 int nn = 0;
1890 child_fe->loopSize = size;
1891
1892 for (auto p = childs.pair_begin(); p != childs.pair_end(); ++p) {
1893
1894 auto miit =
1895 numered_fe.lower_bound(EntFiniteElement::getLocalUniqueIdCalculate(
1896 p->first, (*fe_miit)->getFEUId()));
1897 auto hi_miit =
1898 numered_fe.upper_bound(EntFiniteElement::getLocalUniqueIdCalculate(
1899 p->second, (*fe_miit)->getFEUId()));
1900
1901 for (; miit != hi_miit; ++miit) {
1902
1903 if (verb >= VERBOSE)
1904 MOFEM_LOG("SELF", sev) << "Child finite element "
1905 << "(" << nn << "): " << **miit;
1906
1907 child_fe->nInTheLoop = nn++;
1908 child_fe->numeredEntFiniteElementPtr = *miit;
1909 CHKERR (*child_fe)();
1910 }
1911 }
1912 }
1913
1914 CHKERR child_fe->postProcess();
1915 }
1916
1918}
static Index< 'p', 3 > p
@ VERBOSE
Definition: definitions.h:209
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
const Problem * problemPtr
raw pointer to problem
boost::weak_ptr< CacheTuple > cacheWeakPtr
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.

◆ loopParent()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::loopParent ( const string &  fe_name,
ForcesAndSourcesCore parent_fe,
const int  verb = QUIET,
const LogManager::SeverityLevel  sev = Sev::noisy 
)

User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations.

Parameters
fe_name
parent_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1789 of file ForcesAndSourcesCore.cpp.

1791 {
1793
1794 auto &fes =
1795 ptrFE->mField.get_finite_elements()->get<FiniteElement_name_mi_tag>();
1796 auto fe_miit = fes.find(fe_name);
1797 if (fe_miit != fes.end()) {
1798
1799 const auto *problem_ptr = getFEMethod()->problemPtr;
1800 auto &numered_fe =
1801 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1802
1803
1804 parent_fe->feName = fe_name;
1805 CHKERR parent_fe->setRefineFEPtr(ptrFE);
1806 CHKERR parent_fe->copyBasicMethod(*getFEMethod());
1807 CHKERR parent_fe->copyPetscData(*getFEMethod());
1808 CHKERR parent_fe->copyKsp(*getFEMethod());
1809 CHKERR parent_fe->copySnes(*getFEMethod());
1810 CHKERR parent_fe->copyTs(*getFEMethod());
1811
1812 parent_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1813
1814 const auto parent_ent = getNumeredEntFiniteElementPtr()->getParentEnt();
1815 auto miit = numered_fe.find(EntFiniteElement::getLocalUniqueIdCalculate(
1816 parent_ent, (*fe_miit)->getFEUId()));
1817 if (miit != numered_fe.end()) {
1818 if (verb >= VERBOSE)
1819 MOFEM_LOG("SELF", sev) << "Parent finite element: " << **miit;
1820 parent_fe->loopSize = 1;
1821 parent_fe->nInTheLoop = 0;
1822 parent_fe->numeredEntFiniteElementPtr = *miit;
1823 CHKERR parent_fe->preProcess();
1824 CHKERR (*parent_fe)();
1825 CHKERR parent_fe->postProcess();
1826 } else {
1827 if (verb >= VERBOSE)
1828 MOFEM_LOG("SELF", sev) << "Parent finite element: no parent";
1829 parent_fe->loopSize = 0;
1830 parent_fe->nInTheLoop = 0;
1831 CHKERR parent_fe->preProcess();
1832 CHKERR parent_fe->postProcess();
1833 }
1834 }
1835
1837}
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElementsPtr
store finite elements

◆ loopSide()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::loopSide ( const string &  fe_name,
ForcesAndSourcesCore side_fe,
const size_t  dim,
const EntityHandle  ent_for_side = 0,
const int  verb = QUIET,
const LogManager::SeverityLevel  sev = Sev::noisy,
AdjCache adj_cache = nullptr 
)

User calls this function to loop over elements on the side of face. This function calls finite element with its operator to do calculations.

Parameters
fe_namename of the side element
side_fepointer to the side element instance
dimdimension the of side element
ent_for_sideentity handle for which adjacent volume or face will be accessed
verb
sev
Returns
MoFEMErrorCode

Definition at line 1668 of file ForcesAndSourcesCore.cpp.

1671 {
1673
1674 const auto *problem_ptr = getFEMethod()->problemPtr;
1675 auto &numered_fe =
1676 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1677
1678 auto fe_miit = ptrFE->mField.get_finite_elements()
1679 ->get<FiniteElement_name_mi_tag>()
1680 .find(fe_name);
1681 if (fe_miit != ptrFE->mField.get_finite_elements()
1682 ->get<FiniteElement_name_mi_tag>()
1683 .end()) {
1684
1685 const auto ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1686
1687 side_fe->feName = fe_name;
1688
1689 CHKERR side_fe->setSideFEPtr(ptrFE);
1690 CHKERR side_fe->copyBasicMethod(*getFEMethod());
1691 CHKERR side_fe->copyPetscData(*getFEMethod());
1692 CHKERR side_fe->copyKsp(*getFEMethod());
1693 CHKERR side_fe->copySnes(*getFEMethod());
1694 CHKERR side_fe->copyTs(*getFEMethod());
1695
1696 side_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1697
1698 CHKERR side_fe->preProcess();
1699
1700 std::vector<boost::weak_ptr<NumeredEntFiniteElement>> fe_vec;
1701 auto get_numered_fe_ptr = [&](auto &fe_uid, Range &&adjacent_ents)
1702 -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1703 fe_vec.reserve(adjacent_ents.size());
1704 for (auto fe_ent : adjacent_ents) {
1705 auto miit = numered_fe.find(
1707 if (miit != numered_fe.end()) {
1708 fe_vec.emplace_back(*miit);
1709 }
1710 }
1711 return fe_vec;
1712 };
1713
1714 auto get_bit_entity_adjacency = [&]() {
1715 Range adjacent_ents;
1716 CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1717 ent, side_dim, adjacent_ents);
1718 return adjacent_ents;
1719 };
1720
1721 auto get_adj = [&](auto &fe_uid)
1722 -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1723 if (adj_cache) {
1724 try {
1725 return (*adj_cache).at(ent);
1726 } catch (const std::out_of_range &) {
1727 return (*adj_cache)[ent] =
1728 get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1729 }
1730 } else {
1731 return get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1732 }
1733 };
1734
1735 auto adj = get_adj((*fe_miit)->getFEUId());
1736
1737 int nn = 0;
1738 side_fe->loopSize = adj.size();
1739 for (auto fe_weak_ptr : adj) {
1740 if (auto fe_ptr = fe_weak_ptr.lock()) {
1741 if (verb >= VERBOSE)
1742 MOFEM_LOG("SELF", sev) << "Side finite element "
1743 << "(" << nn << "): " << *fe_ptr;
1744 side_fe->nInTheLoop = nn;
1745 side_fe->numeredEntFiniteElementPtr = fe_ptr;
1746 CHKERR (*side_fe)();
1747 }
1748 ++nn;
1749 }
1750
1751 CHKERR side_fe->postProcess();
1752 }
1753
1755}
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ loopThis()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::loopThis ( const string &  fe_name,
ForcesAndSourcesCore this_fe,
const int  verb = QUIET,
const LogManager::SeverityLevel  sev = Sev::noisy 
)

User calls this function to loop over the same element using a different set of integration points. This function calls finite element with its operator to do calculations.

Parameters
fe_name
this_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1757 of file ForcesAndSourcesCore.cpp.

1759 {
1761
1762 if (verb >= VERBOSE)
1763 MOFEM_LOG("SELF", sev) << "This finite element: "
1765
1766 this_fe->feName = fe_name;
1767
1768 CHKERR this_fe->setRefineFEPtr(ptrFE);
1769 CHKERR this_fe->copyBasicMethod(*getFEMethod());
1770 CHKERR this_fe->copyPetscData(*getFEMethod());
1771 CHKERR this_fe->copyKsp(*getFEMethod());
1772 CHKERR this_fe->copySnes(*getFEMethod());
1773 CHKERR this_fe->copyTs(*getFEMethod());
1774
1775 this_fe->cacheWeakPtr = getFEMethod()->cacheWeakPtr;
1776
1777 CHKERR this_fe->preProcess();
1778
1779 this_fe->nInTheLoop = getNinTheLoop();
1780 this_fe->numeredEntFiniteElementPtr = getNumeredEntFiniteElementPtr();
1781
1782 CHKERR (*this_fe)();
1783
1784 CHKERR this_fe->postProcess();
1785
1787}
int getNinTheLoop() const
get number of finite element in the loop

◆ setOpType()

void ForcesAndSourcesCore::UserDataOperator::setOpType ( const OpType  type)
inline

Set operator type.

Parameters
Operatortype

Definition at line 1047 of file ForcesAndSourcesCore.hpp.

1047 {
1048 opType = type;
1049}

◆ setPtrFE()

MoFEMErrorCode ForcesAndSourcesCore::UserDataOperator::setPtrFE ( ForcesAndSourcesCore ptr)
protectedvirtual

Reimplemented in MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator, MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator, MoFEM::FaceElementForcesAndSourcesCoreOnSide::UserDataOperator, MoFEM::FatPrismElementForcesAndSourcesCore::UserDataOperator, MoFEM::FlatPrismElementForcesAndSourcesCore::UserDataOperator, MoFEM::VertexElementForcesAndSourcesCore::UserDataOperator, MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator, and MoFEM::VolumeElementForcesAndSourcesCoreOnSide::UserDataOperator.

Definition at line 1991 of file ForcesAndSourcesCore.cpp.

1991 {
1993 if (!(ptrFE = dynamic_cast<ForcesAndSourcesCore *>(ptr)))
1994 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1995 "User operator and finite element do not work together");
1997}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440

Friends And Related Function Documentation

◆ ContactPrismElementForcesAndSourcesCore

Definition at line 992 of file ForcesAndSourcesCore.hpp.

◆ EdgeElementForcesAndSourcesCore

friend class EdgeElementForcesAndSourcesCore
friend

Definition at line 990 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCore

friend class FaceElementForcesAndSourcesCore
friend

Definition at line 991 of file ForcesAndSourcesCore.hpp.

◆ ForcesAndSourcesCore

friend class ForcesAndSourcesCore
friend

Definition at line 989 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ colFieldName

std::string MoFEM::ForcesAndSourcesCore::UserDataOperator::colFieldName

Definition at line 578 of file ForcesAndSourcesCore.hpp.

◆ opType

char MoFEM::ForcesAndSourcesCore::UserDataOperator::opType

Definition at line 576 of file ForcesAndSourcesCore.hpp.

◆ OpTypeNames

const char *const ForcesAndSourcesCore::UserDataOperator::OpTypeNames
static
Initial value:
= {
"OPROW", " OPCOL", "OPROWCOL", "OPSPACE", "OPLAST"}

Definition at line 574 of file ForcesAndSourcesCore.hpp.

◆ ptrFE

ForcesAndSourcesCore* MoFEM::ForcesAndSourcesCore::UserDataOperator::ptrFE
protected

Definition at line 984 of file ForcesAndSourcesCore.hpp.

◆ rowFieldName

std::string MoFEM::ForcesAndSourcesCore::UserDataOperator::rowFieldName

Definition at line 577 of file ForcesAndSourcesCore.hpp.

◆ sPace

FieldSpace MoFEM::ForcesAndSourcesCore::UserDataOperator::sPace

Definition at line 579 of file ForcesAndSourcesCore.hpp.


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