v0.13.1
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 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...
 
Loops
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)
 User call this function to loop over elements on the side of face. This function calls finite element with is operator to do calculations. More...
 
MoFEMErrorCode loopThis (const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
 User call this function to loop over parent elements. This function calls finite element with is 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 call this function to loop over parent elements. This function calls finite element with is 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 call this function to loop over parent elements. This function calls finite element with is operator to do calculations. 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
 

Detailed Description

Examples
hello_world.cpp, and mesh_smoothing.cpp.

Definition at line 546 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 563 of file ForcesAndSourcesCore.hpp.

563 {
564 OPROW = 1 << 0, ///< operator doWork function is executed on FE rows
565 OPCOL = 1 << 1, ///< operator doWork function is executed on FE columns
566 OPROWCOL = 1 << 2, ///< operator doWork is executed on FE rows &columns
567 OPSPACE = 1 << 3, ///< operator do Work is execute on space data
568 OPLAST = 1 << 3 ///< @deprecated would be removed
569 };
@ 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.

Examples
NavierStokesElement.hpp.

Definition at line 1866 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 1871 of file ForcesAndSourcesCore.cpp.

@ 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 1876 of file ForcesAndSourcesCore.cpp.

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

Member Function Documentation

◆ addOpType()

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

Add operator type.

Definition at line 1061 of file ForcesAndSourcesCore.hpp.

1061 {
1062 opType |= type;
1063}

◆ getCoordsAtGaussPts()

MatrixDouble & ForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts ( )

◆ getDataCtx()

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

Definition at line 1078 of file ForcesAndSourcesCore.hpp.

1078 {
1079 return getFEMethod()->data_ctx;
1080}
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getFEDim()

int ForcesAndSourcesCore::UserDataOperator::getFEDim ( ) const

Get dimension of finite element.

Returns
int

Definition at line 1017 of file ForcesAndSourcesCore.hpp.

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

◆ getFEEntityHandle()

EntityHandle ForcesAndSourcesCore::UserDataOperator::getFEEntityHandle ( ) const

Return finite element entity handle.

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

Definition at line 1013 of file ForcesAndSourcesCore.hpp.

1013 {
1014 return getNumeredEntFiniteElementPtr()->getEnt();
1015}
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.

◆ getFEMethod()

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

◆ getFEName()

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

Get name of the element.

Definition at line 1073 of file ForcesAndSourcesCore.hpp.

1073 {
1074 return getFEMethod()->getFEName();
1075}
auto getFEName() const
get finite element name

◆ getFEType()

EntityType ForcesAndSourcesCore::UserDataOperator::getFEType ( ) const

Get dimension of finite element.

Returns
int

Definition at line 1021 of file ForcesAndSourcesCore.hpp.

1021 {
1023};
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1479

◆ getFTensor0IntegrationWeight()

auto ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight ( )

◆ getFTensor1CoordsAtGaussPts()

auto ForcesAndSourcesCore::UserDataOperator::getFTensor1CoordsAtGaussPts ( )

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
PoissonDiscontinousGalerkin.hpp, UnsaturatedFlow.hpp, child_and_parent.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 1270 of file ForcesAndSourcesCore.hpp.

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

◆ getGaussPts()

MatrixDouble & ForcesAndSourcesCore::UserDataOperator::getGaussPts ( )

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
MagneticElement.hpp, PoissonDiscontinousGalerkin.hpp, ThermoElasticOps.hpp, UnsaturatedFlow.hpp, approx_sphere.cpp, child_and_parent.cpp, continuity_check_on_skeleton_3d.cpp, hanging_node_approx.cpp, heat_method.cpp, higher_derivatives.cpp, mixed_poisson.cpp, phase.cpp, photon_diffusion.cpp, plate.cpp, poisson_2d_homogeneous.hpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, and simple_interface.cpp.

Definition at line 1237 of file ForcesAndSourcesCore.hpp.

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

◆ getKSPA()

Mat ForcesAndSourcesCore::UserDataOperator::getKSPA ( ) const

Definition at line 1105 of file ForcesAndSourcesCore.hpp.

1105 {
1106#ifndef NDEBUG
1107 if (getFEMethod()->ksp_A == PETSC_NULL)
1108 THROW_MESSAGE("KSP not set A vector");
1109#endif
1110 return getFEMethod()->ksp_A;
1111}
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561

◆ getKSPB()

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

Definition at line 1113 of file ForcesAndSourcesCore.hpp.

1113 {
1114#ifndef NDEBUG
1115 if (getFEMethod()->ksp_B == PETSC_NULL)
1116 THROW_MESSAGE("KSP not set B vector");
1117#endif
1118 return getFEMethod()->ksp_B;
1119}

◆ getKSPCtx()

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

Definition at line 1083 of file ForcesAndSourcesCore.hpp.

1083 {
1084 return getFEMethod()->ksp_ctx;
1085}
KSPContext ksp_ctx
Context.
Definition: LoopMethods.hpp:86

◆ getKSPf()

Vec ForcesAndSourcesCore::UserDataOperator::getKSPf ( ) const
Examples
PoissonDiscontinousGalerkin.hpp, heat_method.cpp, and poisson_2d_homogeneous.hpp.

Definition at line 1097 of file ForcesAndSourcesCore.hpp.

1097 {
1098#ifndef NDEBUG
1099 if (getFEMethod()->ksp_f == PETSC_NULL)
1100 THROW_MESSAGE("KSP not set F vector");
1101#endif
1102 return getFEMethod()->ksp_f;
1103}

◆ getLoopSize()

int ForcesAndSourcesCore::UserDataOperator::getLoopSize ( ) const

get size of elements in the loop

Returns
loop size

Definition at line 1069 of file ForcesAndSourcesCore.hpp.

1069 {
1070 return getFEMethod()->getLoopSize();
1071}
int getLoopSize() const
get loop size

◆ getMeasure() [1/2]

double & ForcesAndSourcesCore::UserDataOperator::getMeasure ( )

get measure of element

Returns
volume

Definition at line 1280 of file ForcesAndSourcesCore.hpp.

1280 {
1281 return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1282}

◆ getMeasure() [2/2]

double ForcesAndSourcesCore::UserDataOperator::getMeasure ( ) const

get measure of element

Returns
volume
Examples
child_and_parent.cpp, hanging_node_approx.cpp, heat_method.cpp, higher_derivatives.cpp, and mixed_poisson.cpp.

Definition at line 1276 of file ForcesAndSourcesCore.hpp.

1276 {
1277 return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1278}

◆ getNinTheLoop()

int ForcesAndSourcesCore::UserDataOperator::getNinTheLoop ( ) const

get number of finite element in the loop

Returns
number of finite element

Definition at line 1065 of file ForcesAndSourcesCore.hpp.

1065 {
1066 return getFEMethod()->getNinTheLoop();
1067}
int getNinTheLoop() const
get number of evaluated element in the loop

◆ getNumberOfNodesOnElement()

int ForcesAndSourcesCore::UserDataOperator::getNumberOfNodesOnElement ( ) const

Get the number of nodes on finite element.

Returns
int
Examples
boundary_marker.cpp.

Definition at line 1047 of file ForcesAndSourcesCore.hpp.

1047 {
1048 return ptrFE->getNumberOfNodes();
1049}
auto getNumberOfNodes() const

◆ getNumeredEntFiniteElementPtr()

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

Return raw pointer to NumeredEntFiniteElement.

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

Definition at line 1009 of file ForcesAndSourcesCore.hpp.

1009 {
1011};
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getOpType()

int ForcesAndSourcesCore::UserDataOperator::getOpType ( ) const

Get operator types.

Returns
Return operator type

Definition at line 1055 of file ForcesAndSourcesCore.hpp.

1055{ 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 1630 of file ForcesAndSourcesCore.cpp.

1632 {
1634 if (ptrFE == NULL)
1635 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1636
1637 switch (type) {
1638 case MBVERTEX:
1640 break;
1641 default:
1643 }
1645}
#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 1613 of file ForcesAndSourcesCore.cpp.

1615 {
1617 if (ptrFE == NULL)
1618 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1619
1620 switch (type) {
1621 case MBVERTEX:
1623 break;
1624 default:
1626 }
1628}
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

Definition at line 1252 of file ForcesAndSourcesCore.hpp.

1252 {
1253 return ptrFE;
1254}

◆ getRefinePtrFE()

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

Definition at line 1262 of file ForcesAndSourcesCore.hpp.

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

◆ getSideEntity()

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

Get the side entity.

Note
For vertex is expection. 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)
EntityHandle ent = getSideEntity(n, type);
// 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 1039 of file ForcesAndSourcesCore.hpp.

1040 {
1041 if (auto side_ptr = getSideNumberPtr(side_number, type).lock())
1042 return side_ptr->ent;
1043 else
1044 return 0;
1045}
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 
)

Get the side number pointer.

Note
For vertex is expection. 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 1026 of file ForcesAndSourcesCore.hpp.

1027 {
1028 auto &side_table_by_side_and_type =
1029 ptrFE->numeredEntFiniteElementPtr->getSideNumberTable().get<1>();
1030 auto side_it =
1031 side_table_by_side_and_type.find(boost::make_tuple(type, side_number));
1032 if (side_it != side_table_by_side_and_type.end())
1033 return *side_it;
1034 else
1035 return boost::weak_ptr<SideNumber>();
1036}

◆ getSidePtrFE()

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

Definition at line 1257 of file ForcesAndSourcesCore.hpp.

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

◆ getSNESA()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESA ( ) const

Definition at line 1137 of file ForcesAndSourcesCore.hpp.

1137 {
1138#ifndef NDEBUG
1139 if (getFEMethod()->snes_A == PETSC_NULL)
1140 THROW_MESSAGE("SNES not set A vector");
1141#endif
1142 return getFEMethod()->snes_A;
1143}
Mat & snes_A
jacobian matrix

◆ getSNESB()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESB ( ) const

Definition at line 1145 of file ForcesAndSourcesCore.hpp.

1145 {
1146#ifndef NDEBUG
1147 if (getFEMethod()->snes_B == PETSC_NULL)
1148 THROW_MESSAGE("SNES not set A matrix");
1149#endif
1150 return getFEMethod()->snes_B;
1151}
Mat & snes_B
preconditioner of jacobian matrix

◆ getSNESCtx()

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

Definition at line 1088 of file ForcesAndSourcesCore.hpp.

1088 {
1089 return getFEMethod()->snes_ctx;
1090}
SNESContext snes_ctx

◆ getSNESf()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESf ( ) const

Definition at line 1121 of file ForcesAndSourcesCore.hpp.

1121 {
1122#ifndef NDEBUG
1123 if (getFEMethod()->snes_f == PETSC_NULL)
1124 THROW_MESSAGE("SNES not set F vector");
1125#endif
1126 return getFEMethod()->snes_f;
1127}
Vec & snes_f
residual

◆ getSNESx()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESx ( ) const

Definition at line 1129 of file ForcesAndSourcesCore.hpp.

1129 {
1130#ifndef NDEBUG
1131 if (getFEMethod()->snes_x == PETSC_NULL)
1132 THROW_MESSAGE("SNESnot set X vector");
1133#endif
1134 return getFEMethod()->snes_x;
1135}
Vec & snes_x
state vector

◆ getTSA()

Mat ForcesAndSourcesCore::UserDataOperator::getTSA ( ) const

Definition at line 1185 of file ForcesAndSourcesCore.hpp.

1185 {
1186#ifndef NDEBUG
1187 if (getFEMethod()->ts_A == PETSC_NULL)
1188 THROW_MESSAGE("TS not set A matrix");
1189#endif
1190 return getFEMethod()->ts_A;
1191}

◆ getTSa()

double ForcesAndSourcesCore::UserDataOperator::getTSa ( ) const

Definition at line 1217 of file ForcesAndSourcesCore.hpp.

1217 {
1218#ifndef NDEBUG
1220 .none() ||
1222 THROW_MESSAGE("TS not set B matrix");
1223#endif
1224 return getFEMethod()->ts_a;
1225}
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

Definition at line 1227 of file ForcesAndSourcesCore.hpp.

1227 {
1228#ifndef NDEBUG
1230 .none() ||
1232 THROW_MESSAGE("TS not set B matrix");
1233#endif
1234 return getFEMethod()->ts_aa;
1235}
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

Definition at line 1193 of file ForcesAndSourcesCore.hpp.

1193 {
1194#ifndef NDEBUG
1195 if (getFEMethod()->ts_B == PETSC_NULL)
1196 THROW_MESSAGE("TS not set B matrix");
1197#endif
1198 return getFEMethod()->ts_B;
1199}
Mat & ts_B
Preconditioner for ts_A.

◆ getTSCtx()

const TSMethod::TSContext ForcesAndSourcesCore::UserDataOperator::getTSCtx ( ) const

Definition at line 1093 of file ForcesAndSourcesCore.hpp.

1093 {
1094 return getFEMethod()->ts_ctx;
1095}
TSContext ts_ctx

◆ getTSf()

Vec ForcesAndSourcesCore::UserDataOperator::getTSf ( ) const

Definition at line 1177 of file ForcesAndSourcesCore.hpp.

1177 {
1178#ifndef NDEBUG
1179 if (getFEMethod()->ts_F == PETSC_NULL)
1180 THROW_MESSAGE("TS not set F vector");
1181#endif
1182 return getFEMethod()->ts_F;
1183}
Vec & ts_F
residual vector

◆ getTSstep()

int ForcesAndSourcesCore::UserDataOperator::getTSstep ( ) const

Definition at line 1201 of file ForcesAndSourcesCore.hpp.

1201 {
1202#ifndef NDEBUG
1203 if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1204 THROW_MESSAGE("TS not set time");
1205#endif
1206 return getFEMethod()->ts_step;
1207}
PetscInt ts_step
time step

◆ getTStime()

double ForcesAndSourcesCore::UserDataOperator::getTStime ( ) const

Definition at line 1209 of file ForcesAndSourcesCore.hpp.

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

◆ getTSu()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu ( ) const

Definition at line 1153 of file ForcesAndSourcesCore.hpp.

1153 {
1154#ifndef NDEBUG
1155 if (getFEMethod()->ts_u == PETSC_NULL)
1156 THROW_MESSAGE("TS not set U vector");
1157#endif
1158 return getFEMethod()->ts_u;
1159}
Vec & ts_u
state vector

◆ getTSu_t()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_t ( ) const

Definition at line 1161 of file ForcesAndSourcesCore.hpp.

1161 {
1162#ifndef NDEBUG
1163 if (getFEMethod()->ts_u_t == PETSC_NULL)
1164 THROW_MESSAGE("TS not set U_t vector");
1165#endif
1166 return getFEMethod()->ts_u_t;
1167}
Vec & ts_u_t
time derivative of state vector

◆ getTSu_tt()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_tt ( ) const

Definition at line 1169 of file ForcesAndSourcesCore.hpp.

1169 {
1170#ifndef NDEBUG
1171 if (getFEMethod()->ts_u_tt == PETSC_NULL)
1172 THROW_MESSAGE("TS not set U_tt vector");
1173#endif
1174 return getFEMethod()->ts_u_tt;
1175}
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 call this function to loop over parent elements. This function calls finite element with is operator to do calculations.

Parameters
fe_name
child_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1776 of file ForcesAndSourcesCore.cpp.

1778 {
1780
1781 const auto *problem_ptr = getFEMethod()->problemPtr;
1782 auto &ref_ents = *getPtrFE()->mField.get_ref_ents();
1783 auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1784 ->get<Composite_Name_And_Ent_mi_tag>();
1785
1786 const auto parent_ent = getNumeredEntFiniteElementPtr()->getEnt();
1787 const auto parent_type = getNumeredEntFiniteElementPtr()->getEntType();
1788 auto range =
1789 ref_ents.get<Composite_ParentEnt_And_EntType_mi_tag>().equal_range(
1790 boost::make_tuple(parent_type, parent_ent));
1791
1792 if (auto size = std::distance(range.first, range.second)) {
1793
1794 std::vector<EntityHandle> childs_vec;
1795 childs_vec.reserve(size);
1796 for (; range.first != range.second; ++range.first)
1797 childs_vec.emplace_back((*range.first)->getEnt());
1798
1799 Range childs;
1800
1801 if ((childs_vec.back() - childs_vec.front() + 1) == size)
1802 childs = Range(childs_vec.front(), childs_vec.back());
1803 else
1804 childs.insert_list(childs_vec.begin(), childs_vec.end());
1805
1806 child_fe->feName = fe_name;
1807
1808 CHKERR child_fe->setRefineFEPtr(ptrFE);
1809 CHKERR child_fe->copyBasicMethod(*getFEMethod());
1810 CHKERR child_fe->copyPetscData(*getFEMethod());
1811 CHKERR child_fe->copyKsp(*getFEMethod());
1812 CHKERR child_fe->copySnes(*getFEMethod());
1813 CHKERR child_fe->copyTs(*getFEMethod());
1814
1815 CHKERR child_fe->preProcess();
1816
1817 int nn = 0;
1818 child_fe->loopSize = size;
1819
1820 for (auto p = childs.pair_begin(); p != childs.pair_end(); ++p) {
1821
1822 auto miit = numered_fe.lower_bound(boost::make_tuple(fe_name, p->first));
1823 auto hi_miit =
1824 numered_fe.upper_bound(boost::make_tuple(fe_name, p->second));
1825
1826 for (; miit != hi_miit; ++miit) {
1827
1828 if (verb >= VERBOSE)
1829 MOFEM_LOG("SELF", sev) << "Child finite element "
1830 << "(" << nn << "): " << **miit;
1831
1832 child_fe->nInTheLoop = nn++;
1833 child_fe->numeredEntFiniteElementPtr = *miit;
1834 CHKERR (*child_fe)();
1835 }
1836 }
1837
1838 CHKERR child_fe->postProcess();
1839 };
1840
1842}
static Index< 'p', 3 > p
@ VERBOSE
Definition: definitions.h:209
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
const Problem * problemPtr
raw pointer to problem

◆ loopParent()

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

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

Parameters
fe_name
parent_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1738 of file ForcesAndSourcesCore.cpp.

1740 {
1742
1743 const auto *problem_ptr = getFEMethod()->problemPtr;
1744 auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1745 ->get<Composite_Name_And_Ent_mi_tag>();
1746
1747 const auto parent_ent = getNumeredEntFiniteElementPtr()->getParentEnt();
1748 auto miit = numered_fe.find(boost::make_tuple(fe_name, parent_ent));
1749 if (miit != numered_fe.end()) {
1750
1751 if (verb >= VERBOSE)
1752 MOFEM_LOG("SELF", sev) << "Parent finite element: " << **miit;
1753
1754 parent_fe->feName = fe_name;
1755
1756 CHKERR parent_fe->setRefineFEPtr(ptrFE);
1757 CHKERR parent_fe->copyBasicMethod(*getFEMethod());
1758 CHKERR parent_fe->copyPetscData(*getFEMethod());
1759 CHKERR parent_fe->copyKsp(*getFEMethod());
1760 CHKERR parent_fe->copySnes(*getFEMethod());
1761 CHKERR parent_fe->copyTs(*getFEMethod());
1762
1763 CHKERR parent_fe->preProcess();
1764
1765 parent_fe->loopSize = 1;
1766 parent_fe->nInTheLoop = 0;
1767 parent_fe->numeredEntFiniteElementPtr = *miit;
1768 CHKERR (*parent_fe)();
1769
1770 CHKERR parent_fe->postProcess();
1771 }
1772
1774}
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 
)

User call this function to loop over elements on the side of face. This function calls finite element with is 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 1661 of file ForcesAndSourcesCore.cpp.

1664 {
1666 const auto ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1667 const auto *problem_ptr = getFEMethod()->problemPtr;
1668
1669 side_fe->feName = fe_name;
1670
1671 CHKERR side_fe->setSideFEPtr(ptrFE);
1672 CHKERR side_fe->copyBasicMethod(*getFEMethod());
1673 CHKERR side_fe->copyPetscData(*getFEMethod());
1674 CHKERR side_fe->copyKsp(*getFEMethod());
1675 CHKERR side_fe->copySnes(*getFEMethod());
1676 CHKERR side_fe->copyTs(*getFEMethod());
1677
1678 CHKERR side_fe->preProcess();
1679
1680 Range adjacent_ents;
1681 CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1682 ent, side_dim, adjacent_ents);
1683 auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1684 ->get<Composite_Name_And_Ent_mi_tag>();
1685
1686 int nn = 0;
1687 side_fe->loopSize = adjacent_ents.size();
1688 for (auto fe_ent : adjacent_ents) {
1689 auto miit = numered_fe.find(boost::make_tuple(fe_name, fe_ent));
1690 if (miit != numered_fe.end()) {
1691
1692 if (verb >= VERBOSE)
1693 MOFEM_LOG("SELF", sev) << "Side finite element "
1694 << "(" << nn << "): " << **miit;
1695
1696 side_fe->nInTheLoop = nn++;
1697 side_fe->numeredEntFiniteElementPtr = *miit;
1698 CHKERR (*side_fe)();
1699 }
1700 }
1701
1702 CHKERR side_fe->postProcess();
1703
1705}
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ loopThis()

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

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

Parameters
fe_name
parent_fe
verb
sev
Returns
MoFEMErrorCode

Definition at line 1707 of file ForcesAndSourcesCore.cpp.

1709 {
1711
1712 if (verb >= VERBOSE)
1713 MOFEM_LOG("SELF", sev) << "This finite element: "
1715
1716 const auto *problem_ptr = getFEMethod()->problemPtr;
1717 this_fe->feName = fe_name;
1718
1719 CHKERR this_fe->setRefineFEPtr(ptrFE);
1720 CHKERR this_fe->copyBasicMethod(*getFEMethod());
1721 CHKERR this_fe->copyPetscData(*getFEMethod());
1722 CHKERR this_fe->copyKsp(*getFEMethod());
1723 CHKERR this_fe->copySnes(*getFEMethod());
1724 CHKERR this_fe->copyTs(*getFEMethod());
1725
1726 CHKERR this_fe->preProcess();
1727
1728 this_fe->nInTheLoop = getNinTheLoop();
1729 this_fe->numeredEntFiniteElementPtr = getNumeredEntFiniteElementPtr();
1730
1731 CHKERR (*this_fe)();
1732
1733 CHKERR this_fe->postProcess();
1734
1736}
int getNinTheLoop() const
get number of finite element in the loop

◆ setOpType()

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

Set operator type.

Parameters
Operatortype

Definition at line 1057 of file ForcesAndSourcesCore.hpp.

1057 {
1058 opType = type;
1059}

◆ 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 1915 of file ForcesAndSourcesCore.cpp.

1915 {
1917 if (!(ptrFE = dynamic_cast<ForcesAndSourcesCore *>(ptr)))
1918 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1919 "User operator and finite element do not work together");
1921}
#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 1002 of file ForcesAndSourcesCore.hpp.

◆ EdgeElementForcesAndSourcesCore

friend class EdgeElementForcesAndSourcesCore
friend

Definition at line 1000 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCore

friend class FaceElementForcesAndSourcesCore
friend

Definition at line 1001 of file ForcesAndSourcesCore.hpp.

◆ ForcesAndSourcesCore

friend class ForcesAndSourcesCore
friend

Definition at line 999 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ colFieldName

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

Definition at line 575 of file ForcesAndSourcesCore.hpp.

◆ opType

char MoFEM::ForcesAndSourcesCore::UserDataOperator::opType

Definition at line 573 of file ForcesAndSourcesCore.hpp.

◆ OpTypeNames

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

Definition at line 571 of file ForcesAndSourcesCore.hpp.

◆ ptrFE

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

Definition at line 994 of file ForcesAndSourcesCore.hpp.

◆ rowFieldName

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

Definition at line 574 of file ForcesAndSourcesCore.hpp.

◆ sPace

FieldSpace MoFEM::ForcesAndSourcesCore::UserDataOperator::sPace

Definition at line 576 of file ForcesAndSourcesCore.hpp.


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