v0.13.0
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 }
 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...
 
const std::string & getFEName () const
 Get name of the element. More...
 
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...
 
double & getMeasure ()
 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)
 
ForcesAndSourcesCoregetPtrFE () const
 
ForcesAndSourcesCoregetSidePtrFE () const
 
ForcesAndSourcesCoregetRefinePtrFE () const
 

Protected Attributes

ForcesAndSourcesCoreptrFE
 

Friends

class ForcesAndSourcesCore
 
class EdgeElementForcesAndSourcesCoreBase
 
class FaceElementForcesAndSourcesCoreBase
 
class ContactPrismElementForcesAndSourcesCore
 

Detailed Description

Examples
mesh_smoothing.cpp.

Definition at line 553 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 
OPCOL 
OPROWCOL 
OPSPACE 

Definition at line 570 of file ForcesAndSourcesCore.hpp.

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 1835 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 1840 of file ForcesAndSourcesCore.cpp.

1842  : DataOperator(symm), opType(type), rowFieldName(field_name),
1843  colFieldName(field_name), sPace(LASTSPACE), ptrFE(nullptr) {}
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:102

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

1848  : DataOperator(symm), opType(type), rowFieldName(row_field_name),
1849  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 1066 of file ForcesAndSourcesCore.hpp.

1066  {
1067  opType |= type;
1068 }

◆ getCoordsAtGaussPts()

MatrixDouble & ForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts ( )

Gauss points and weight, matrix (nb. of points x 3)

Column 0-2 integration points coordinate x, y and z, respectively. At rows are integration points.

Examples
UnsaturatedFlow.hpp.

Definition at line 1271 of file ForcesAndSourcesCore.hpp.

1271  {
1272  return static_cast<ForcesAndSourcesCore *>(ptrFE)->coordsAtGaussPts;
1273 }
MatrixDouble coordsAtGaussPts
coordinated at gauss points

◆ getDataCtx()

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

Definition at line 1083 of file ForcesAndSourcesCore.hpp.

1083  {
1084  return getFEMethod()->data_ctx;
1085 }
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 1022 of file ForcesAndSourcesCore.hpp.

1022  {
1024 };
auto dimension_from_handle
get entity dimension form handle
Definition: Templates.hpp:1449
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.

Definition at line 1018 of file ForcesAndSourcesCore.hpp.

1018  {
1019  return getNumeredEntFiniteElementPtr()->getEnt();
1020 }
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.

◆ getFEMethod()

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

Return raw pointer to Finite Element Method object.

Examples
UnsaturatedFlow.hpp.

Definition at line 1056 of file ForcesAndSourcesCore.hpp.

1056  {
1057  return ptrFE;
1058 }

◆ getFEName()

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

Get name of the element.

Definition at line 1078 of file ForcesAndSourcesCore.hpp.

1078  {
1079  return getFEMethod()->feName;
1080 }
std::string feName
Name of finite element.

◆ getFEType()

EntityType ForcesAndSourcesCore::UserDataOperator::getFEType ( ) const

Get dimension of finite element.

Returns
int

Definition at line 1026 of file ForcesAndSourcesCore.hpp.

1026  {
1028 };
auto type_from_handle
get type from entity handle
Definition: Templates.hpp:1441

◆ getFTensor0IntegrationWeight()

auto ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight ( )

Get integration weights.

for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
// integrate something
++t_w;
}
auto getFTensor0IntegrationWeight()
Get integration weights.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Returns
FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
Examples
PoissonDiscontinousGalerkin.hpp, UnsaturatedFlow.hpp, heat_method.cpp, plate.cpp, prism_polynomial_approximation.cpp, and quad_polynomial_approximation.cpp.

Definition at line 1246 of file ForcesAndSourcesCore.hpp.

1246  {
1248  &(getGaussPts()(getGaussPts().size1() - 1, 0)));
1249 }

◆ 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, prism_polynomial_approximation.cpp, and quad_polynomial_approximation.cpp.

Definition at line 1275 of file ForcesAndSourcesCore.hpp.

1275  {
1277  &getCoordsAtGaussPts()(0, 0), &getCoordsAtGaussPts()(0, 1),
1278  &getCoordsAtGaussPts()(0, 2));
1279 }
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, UnsaturatedFlow.hpp, heat_method.cpp, and plate.cpp.

Definition at line 1242 of file ForcesAndSourcesCore.hpp.

1242  {
1243  return static_cast<ForcesAndSourcesCore *>(ptrFE)->gaussPts;
1244 }
MatrixDouble gaussPts
Matrix of integration points.

◆ getKSPA()

Mat ForcesAndSourcesCore::UserDataOperator::getKSPA ( ) const

Definition at line 1110 of file ForcesAndSourcesCore.hpp.

1110  {
1111 #ifndef NDEBUG
1112  if (getFEMethod()->ksp_A == PETSC_NULL)
1113  THROW_MESSAGE("KSP not set A vector");
1114 #endif
1115  return getFEMethod()->ksp_A;
1116 }
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:574

◆ getKSPB()

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

Definition at line 1118 of file ForcesAndSourcesCore.hpp.

1118  {
1119 #ifndef NDEBUG
1120  if (getFEMethod()->ksp_B == PETSC_NULL)
1121  THROW_MESSAGE("KSP not set B vector");
1122 #endif
1123  return getFEMethod()->ksp_B;
1124 }

◆ getKSPCtx()

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

Definition at line 1088 of file ForcesAndSourcesCore.hpp.

1088  {
1089  return getFEMethod()->ksp_ctx;
1090 }
KSPContext ksp_ctx
Context.
Definition: LoopMethods.hpp:96

◆ getKSPf()

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

Definition at line 1102 of file ForcesAndSourcesCore.hpp.

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

◆ getLoopSize()

int ForcesAndSourcesCore::UserDataOperator::getLoopSize ( ) const

get size of elements in the loop

Returns
loop size

Definition at line 1074 of file ForcesAndSourcesCore.hpp.

1074  {
1075  return getFEMethod()->getLoopSize();
1076 }
int getLoopSize() const
get loop size

◆ getMeasure() [1/2]

double & ForcesAndSourcesCore::UserDataOperator::getMeasure ( )

get measure of element

Returns
volume

Definition at line 1285 of file ForcesAndSourcesCore.hpp.

1285  {
1286  return static_cast<ForcesAndSourcesCore *>(ptrFE)->elementMeasure;
1287 }

◆ getMeasure() [2/2]

double ForcesAndSourcesCore::UserDataOperator::getMeasure ( ) const

get measure of element

Returns
volume
Examples
heat_method.cpp.

Definition at line 1281 of file ForcesAndSourcesCore.hpp.

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

◆ getNinTheLoop()

int ForcesAndSourcesCore::UserDataOperator::getNinTheLoop ( ) const

get number of finite element in the loop

Returns
number of finite element

Definition at line 1070 of file ForcesAndSourcesCore.hpp.

1070  {
1071  return getFEMethod()->getNinTheLoop();
1072 }
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

Definition at line 1052 of file ForcesAndSourcesCore.hpp.

1052  {
1053  return ptrFE->getNumberOfNodes();
1054 }
auto getNumberOfNodes() const

◆ getNumeredEntFiniteElementPtr()

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

Return raw pointer to NumeredEntFiniteElement.

Examples
ElasticityMixedFormulation.hpp, and UnsaturatedFlow.hpp.

Definition at line 1014 of file ForcesAndSourcesCore.hpp.

1014  {
1016 };
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getOpType()

int ForcesAndSourcesCore::UserDataOperator::getOpType ( ) const

Get operator types.

Returns
Return operator type

Definition at line 1060 of file ForcesAndSourcesCore.hpp.

1060 { 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 1597 of file ForcesAndSourcesCore.cpp.

1599  {
1601  if (ptrFE == NULL)
1602  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1603 
1604  switch (type) {
1605  case MBVERTEX:
1606  CHKERR ptrFE->getProblemNodesColIndices(field_name, indices);
1607  break;
1608  default:
1609  CHKERR ptrFE->getProblemTypeColIndices(field_name, type, side, indices);
1610  }
1612 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
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 1580 of file ForcesAndSourcesCore.cpp.

1582  {
1584  if (ptrFE == NULL)
1585  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1586 
1587  switch (type) {
1588  case MBVERTEX:
1589  CHKERR ptrFE->getProblemNodesRowIndices(field_name, indices);
1590  break;
1591  default:
1592  CHKERR ptrFE->getProblemTypeRowIndices(field_name, type, side, indices);
1593  }
1595 }
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
protected

Definition at line 1257 of file ForcesAndSourcesCore.hpp.

1257  {
1258  return ptrFE;
1259 }

◆ getRefinePtrFE()

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

Definition at line 1267 of file ForcesAndSourcesCore.hpp.

1267  {
1268  return ptrFE->refinePtrFE;
1269 }
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.
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
for (int n = 0; n != number_of_nodes; ++n)
// Do somthing
} else {
// Do somthing
}
}
static Index< 'n', 3 > n
EntitiesFieldData::EntData EntData
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
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.
EntityHandle getSideEntity(const int side_number, const EntityType type)
Get the side entity.
Parameters
side_number
type

Definition at line 1044 of file ForcesAndSourcesCore.hpp.

1045  {
1046  if (auto side_ptr = getSideNumberPtr(side_number, type).lock())
1047  return side_ptr->ent;
1048  else
1049  return 0;
1050 }
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 1031 of file ForcesAndSourcesCore.hpp.

1032  {
1033  auto &side_table_by_side_and_type =
1034  ptrFE->numeredEntFiniteElementPtr->getSideNumberTable().get<1>();
1035  auto side_it =
1036  side_table_by_side_and_type.find(boost::make_tuple(type, side_number));
1037  if (side_it != side_table_by_side_and_type.end())
1038  return *side_it;
1039  else
1040  return boost::weak_ptr<SideNumber>();
1041 }

◆ getSidePtrFE()

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

Definition at line 1262 of file ForcesAndSourcesCore.hpp.

1262  {
1263  return ptrFE->sidePtrFE;
1264 }
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.

◆ getSNESA()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESA ( ) const

Definition at line 1142 of file ForcesAndSourcesCore.hpp.

1142  {
1143 #ifndef NDEBUG
1144  if (getFEMethod()->snes_A == PETSC_NULL)
1145  THROW_MESSAGE("SNES not set A vector");
1146 #endif
1147  return getFEMethod()->snes_A;
1148 }
Mat & snes_A
jacobian matrix

◆ getSNESB()

Mat ForcesAndSourcesCore::UserDataOperator::getSNESB ( ) const

Definition at line 1150 of file ForcesAndSourcesCore.hpp.

1150  {
1151 #ifndef NDEBUG
1152  if (getFEMethod()->snes_B == PETSC_NULL)
1153  THROW_MESSAGE("SNES not set A matrix");
1154 #endif
1155  return getFEMethod()->snes_B;
1156 }
Mat & snes_B
preconditioner of jacobian matrix

◆ getSNESCtx()

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

Definition at line 1093 of file ForcesAndSourcesCore.hpp.

1093  {
1094  return getFEMethod()->snes_ctx;
1095 }
SNESContext snes_ctx

◆ getSNESf()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESf ( ) const

Definition at line 1126 of file ForcesAndSourcesCore.hpp.

1126  {
1127 #ifndef NDEBUG
1128  if (getFEMethod()->snes_f == PETSC_NULL)
1129  THROW_MESSAGE("SNES not set F vector");
1130 #endif
1131  return getFEMethod()->snes_f;
1132 }
Vec & snes_f
residual

◆ getSNESx()

Vec ForcesAndSourcesCore::UserDataOperator::getSNESx ( ) const

Definition at line 1134 of file ForcesAndSourcesCore.hpp.

1134  {
1135 #ifndef NDEBUG
1136  if (getFEMethod()->snes_x == PETSC_NULL)
1137  THROW_MESSAGE("SNESnot set X vector");
1138 #endif
1139  return getFEMethod()->snes_x;
1140 }
Vec & snes_x
state vector

◆ getTSA()

Mat ForcesAndSourcesCore::UserDataOperator::getTSA ( ) const

Definition at line 1190 of file ForcesAndSourcesCore.hpp.

1190  {
1191 #ifndef NDEBUG
1192  if (getFEMethod()->ts_A == PETSC_NULL)
1193  THROW_MESSAGE("TS not set A matrix");
1194 #endif
1195  return getFEMethod()->ts_A;
1196 }

◆ getTSa()

double ForcesAndSourcesCore::UserDataOperator::getTSa ( ) const

Definition at line 1222 of file ForcesAndSourcesCore.hpp.

1222  {
1223 #ifndef NDEBUG
1225  .none() ||
1226  (getFEMethod()->data_ctx & (PetscData::CtxSetX_T)).none())
1227  THROW_MESSAGE("TS not set B matrix");
1228 #endif
1229  return getFEMethod()->ts_a;
1230 }
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:47
static constexpr Switches CtxSetX_T
Definition: LoopMethods.hpp:50
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:48
PetscReal ts_a
shift for U_t (see PETSc Time Solver)

◆ getTSaa()

double ForcesAndSourcesCore::UserDataOperator::getTSaa ( ) const

Definition at line 1232 of file ForcesAndSourcesCore.hpp.

1232  {
1233 #ifndef NDEBUG
1235  .none() ||
1236  (getFEMethod()->data_ctx & (PetscData::CtxSetX_TT)).none())
1237  THROW_MESSAGE("TS not set B matrix");
1238 #endif
1239  return getFEMethod()->ts_aa;
1240 }
static constexpr Switches CtxSetX_TT
Definition: LoopMethods.hpp:51
PetscReal ts_aa
shift for U_tt shift for U_tt

◆ getTSB()

Mat ForcesAndSourcesCore::UserDataOperator::getTSB ( ) const

Definition at line 1198 of file ForcesAndSourcesCore.hpp.

1198  {
1199 #ifndef NDEBUG
1200  if (getFEMethod()->ts_B == PETSC_NULL)
1201  THROW_MESSAGE("TS not set B matrix");
1202 #endif
1203  return getFEMethod()->ts_B;
1204 }
Mat & ts_B
Preconditioner for ts_A.

◆ getTSCtx()

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

Definition at line 1098 of file ForcesAndSourcesCore.hpp.

1098  {
1099  return getFEMethod()->ts_ctx;
1100 }
TSContext ts_ctx

◆ getTSf()

Vec ForcesAndSourcesCore::UserDataOperator::getTSf ( ) const

Definition at line 1182 of file ForcesAndSourcesCore.hpp.

1182  {
1183 #ifndef NDEBUG
1184  if (getFEMethod()->ts_F == PETSC_NULL)
1185  THROW_MESSAGE("TS not set F vector");
1186 #endif
1187  return getFEMethod()->ts_F;
1188 }
Vec & ts_F
residual vector

◆ getTSstep()

int ForcesAndSourcesCore::UserDataOperator::getTSstep ( ) const

Definition at line 1206 of file ForcesAndSourcesCore.hpp.

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

◆ getTStime()

double ForcesAndSourcesCore::UserDataOperator::getTStime ( ) const

Definition at line 1214 of file ForcesAndSourcesCore.hpp.

1214  {
1215 #ifndef NDEBUG
1216  if ((getFEMethod()->data_ctx & PetscData::PetscData::CtxSetTime).none())
1217  THROW_MESSAGE("TS not set time");
1218 #endif
1219  return getFEMethod()->ts_t;
1220 }
PetscReal ts_t
time

◆ getTSu()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu ( ) const

Definition at line 1158 of file ForcesAndSourcesCore.hpp.

1158  {
1159 #ifndef NDEBUG
1160  if (getFEMethod()->ts_u == PETSC_NULL)
1161  THROW_MESSAGE("TS not set U vector");
1162 #endif
1163  return getFEMethod()->ts_u;
1164 }
Vec & ts_u
state vector

◆ getTSu_t()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_t ( ) const

Definition at line 1166 of file ForcesAndSourcesCore.hpp.

1166  {
1167 #ifndef NDEBUG
1168  if (getFEMethod()->ts_u_t == PETSC_NULL)
1169  THROW_MESSAGE("TS not set U_t vector");
1170 #endif
1171  return getFEMethod()->ts_u_t;
1172 }
Vec & ts_u_t
time derivative of state vector

◆ getTSu_tt()

Vec ForcesAndSourcesCore::UserDataOperator::getTSu_tt ( ) const

Definition at line 1174 of file ForcesAndSourcesCore.hpp.

1174  {
1175 #ifndef NDEBUG
1176  if (getFEMethod()->ts_u_tt == PETSC_NULL)
1177  THROW_MESSAGE("TS not set U_tt vector");
1178 #endif
1179  return getFEMethod()->ts_u_tt;
1180 }
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 1745 of file ForcesAndSourcesCore.cpp.

1747  {
1749 
1750  const auto *problem_ptr = getFEMethod()->problemPtr;
1751  auto &ref_ents = *getPtrFE()->mField.get_ref_ents();
1752  auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1753  ->get<Composite_Name_And_Ent_mi_tag>();
1754 
1755  const auto parent_ent = getNumeredEntFiniteElementPtr()->getEnt();
1756  const auto parent_type = getNumeredEntFiniteElementPtr()->getEntType();
1757  auto range =
1758  ref_ents.get<Composite_ParentEnt_And_EntType_mi_tag>().equal_range(
1759  boost::make_tuple(parent_type, parent_ent));
1760 
1761  if (auto size = std::distance(range.first, range.second)) {
1762 
1763  std::vector<EntityHandle> childs_vec;
1764  childs_vec.reserve(size);
1765  for (; range.first != range.second; ++range.first)
1766  childs_vec.emplace_back((*range.first)->getEnt());
1767 
1768  Range childs;
1769 
1770  if ((childs_vec.back() - childs_vec.front() + 1) == size)
1771  childs = Range(childs_vec.front(), childs_vec.back());
1772  else
1773  childs.insert_list(childs_vec.begin(), childs_vec.end());
1774 
1775  child_fe->feName = fe_name;
1776 
1777  CHKERR child_fe->setRefineFEPtr(ptrFE);
1778  CHKERR child_fe->copyBasicMethod(*getFEMethod());
1779  CHKERR child_fe->copyPetscData(*getFEMethod());
1780  CHKERR child_fe->copyKsp(*getFEMethod());
1781  CHKERR child_fe->copySnes(*getFEMethod());
1782  CHKERR child_fe->copyTs(*getFEMethod());
1783 
1784  CHKERR child_fe->preProcess();
1785 
1786  int nn = 0;
1787  child_fe->loopSize = size;
1788 
1789  for (auto p = childs.pair_begin(); p != childs.pair_end(); ++p) {
1790 
1791  auto miit = numered_fe.lower_bound(boost::make_tuple(fe_name, p->first));
1792  auto hi_miit =
1793  numered_fe.upper_bound(boost::make_tuple(fe_name, p->second));
1794 
1795  for (; miit != hi_miit; ++miit) {
1796 
1797  if (verb >= VERBOSE)
1798  MOFEM_LOG("SELF", sev) << "Child finite element "
1799  << "(" << nn << "): " << **miit;
1800 
1801  child_fe->nInTheLoop = nn++;
1802  child_fe->numeredEntFiniteElementPtr = *miit;
1803  CHKERR (*child_fe)();
1804  }
1805  }
1806 
1807  CHKERR child_fe->postProcess();
1808  };
1809 
1811 }
static Index< 'p', 3 > p
@ VERBOSE
Definition: definitions.h:222
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
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 1705 of file ForcesAndSourcesCore.cpp.

1707  {
1709 
1710  const auto *problem_ptr = getFEMethod()->problemPtr;
1711  auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1712  ->get<Composite_Name_And_Ent_mi_tag>();
1713 
1714  const auto parent_ent = getNumeredEntFiniteElementPtr()->getParentEnt();
1715  auto miit = numered_fe.find(boost::make_tuple(fe_name, parent_ent));
1716  if (miit != numered_fe.end()) {
1717 
1718  if (verb >= VERBOSE)
1719  MOFEM_LOG("SELF", sev) << "Parent finite element: " << **miit;
1720 
1721  parent_fe->feName = fe_name;
1722 
1723  CHKERR parent_fe->setRefineFEPtr(ptrFE);
1724  CHKERR parent_fe->copyBasicMethod(*getFEMethod());
1725  CHKERR parent_fe->copyPetscData(*getFEMethod());
1726  CHKERR parent_fe->copyKsp(*getFEMethod());
1727  CHKERR parent_fe->copySnes(*getFEMethod());
1728  CHKERR parent_fe->copyTs(*getFEMethod());
1729 
1730  CHKERR parent_fe->preProcess();
1731 
1732  int nn = 0;
1733  parent_fe->loopSize = 1;
1734 
1735  parent_fe->nInTheLoop = nn++;
1736  parent_fe->numeredEntFiniteElementPtr = *miit;
1737  CHKERR (*parent_fe)();
1738 
1739  CHKERR parent_fe->postProcess();
1740  }
1741 
1743 }
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 1628 of file ForcesAndSourcesCore.cpp.

1631  {
1633  const auto ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1634  const auto *problem_ptr = getFEMethod()->problemPtr;
1635 
1636  side_fe->feName = fe_name;
1637 
1638  CHKERR side_fe->setSideFEPtr(ptrFE);
1639  CHKERR side_fe->copyBasicMethod(*getFEMethod());
1640  CHKERR side_fe->copyPetscData(*getFEMethod());
1641  CHKERR side_fe->copyKsp(*getFEMethod());
1642  CHKERR side_fe->copySnes(*getFEMethod());
1643  CHKERR side_fe->copyTs(*getFEMethod());
1644 
1645  CHKERR side_fe->preProcess();
1646 
1647  Range adjacent_ents;
1648  CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1649  ent, side_dim, adjacent_ents);
1650  auto &numered_fe = problem_ptr->numeredFiniteElementsPtr
1651  ->get<Composite_Name_And_Ent_mi_tag>();
1652 
1653  int nn = 0;
1654  side_fe->loopSize = adjacent_ents.size();
1655  for (auto fe_ent : adjacent_ents) {
1656  auto miit = numered_fe.find(boost::make_tuple(fe_name, fe_ent));
1657  if (miit != numered_fe.end()) {
1658 
1659  if (verb >= VERBOSE)
1660  MOFEM_LOG("SELF", sev) << "Side finite element "
1661  << "(" << nn << "): " << **miit;
1662 
1663  side_fe->nInTheLoop = nn++;
1664  side_fe->numeredEntFiniteElementPtr = *miit;
1665  CHKERR (*side_fe)();
1666  }
1667  }
1668 
1669  CHKERR side_fe->postProcess();
1670 
1672 }
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 1674 of file ForcesAndSourcesCore.cpp.

1676  {
1678 
1679  if (verb >= VERBOSE)
1680  MOFEM_LOG("SELF", sev) << "This finite element: "
1682 
1683  const auto *problem_ptr = getFEMethod()->problemPtr;
1684  this_fe->feName = fe_name;
1685 
1686  CHKERR this_fe->setRefineFEPtr(ptrFE);
1687  CHKERR this_fe->copyBasicMethod(*getFEMethod());
1688  CHKERR this_fe->copyPetscData(*getFEMethod());
1689  CHKERR this_fe->copyKsp(*getFEMethod());
1690  CHKERR this_fe->copySnes(*getFEMethod());
1691  CHKERR this_fe->copyTs(*getFEMethod());
1692 
1693  CHKERR this_fe->preProcess();
1694 
1695  this_fe->nInTheLoop = getNinTheLoop();
1696  this_fe->numeredEntFiniteElementPtr = getNumeredEntFiniteElementPtr();
1697 
1698  CHKERR (*this_fe)();
1699 
1700  CHKERR this_fe->postProcess();
1701 
1703 }
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 1062 of file ForcesAndSourcesCore.hpp.

1062  {
1063  opType = type;
1064 }

◆ setPtrFE()

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

Reimplemented in MoFEM::VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator, MoFEM::VolumeElementForcesAndSourcesCoreBase::UserDataOperator, MoFEM::VertexElementForcesAndSourcesCore::UserDataOperator, MoFEM::FlatPrismElementForcesAndSourcesCore::UserDataOperator, MoFEM::FatPrismElementForcesAndSourcesCore::UserDataOperator, MoFEM::FaceElementForcesAndSourcesCoreOnSideBase::UserDataOperator, MoFEM::FaceElementForcesAndSourcesCoreBase::UserDataOperator, and MoFEM::EdgeElementForcesAndSourcesCoreBase::UserDataOperator.

Definition at line 1884 of file ForcesAndSourcesCore.cpp.

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

Friends And Related Function Documentation

◆ ContactPrismElementForcesAndSourcesCore

Definition at line 1007 of file ForcesAndSourcesCore.hpp.

◆ EdgeElementForcesAndSourcesCoreBase

Definition at line 1005 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCoreBase

Definition at line 1006 of file ForcesAndSourcesCore.hpp.

◆ ForcesAndSourcesCore

friend class ForcesAndSourcesCore
friend

Definition at line 1004 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ colFieldName

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

Definition at line 581 of file ForcesAndSourcesCore.hpp.

◆ opType

char MoFEM::ForcesAndSourcesCore::UserDataOperator::opType

Definition at line 579 of file ForcesAndSourcesCore.hpp.

◆ OpTypeNames

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

Definition at line 577 of file ForcesAndSourcesCore.hpp.

◆ ptrFE

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

Definition at line 993 of file ForcesAndSourcesCore.hpp.

◆ rowFieldName

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

Definition at line 580 of file ForcesAndSourcesCore.hpp.

◆ sPace

FieldSpace MoFEM::ForcesAndSourcesCore::UserDataOperator::sPace

Definition at line 582 of file ForcesAndSourcesCore.hpp.


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