v0.10.0
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Friends | List of all members
MoFEM::ForcesAndSourcesCore::UserDataOperator Struct Reference

Data operator to do calculations at integration points.Is inherited and implemented by user to do calculations. It can be used in many different ways but typically is used to integrate matrices (f.e. stiffness matrix) and the right hand vector (f.e. force vector). More...

#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, OPLAST = 1 << 3 }
 Controls loop over entities on element. More...
 

Public Member Functions

 UserDataOperator (const FieldSpace space, const char type=OPLAST, 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...
 
boost::weak_ptr< SideNumber > getSideNumberPtr (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 ()
 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
 
DEPRECATED Vec getSnesF () const
 
DEPRECATED Vec getSnesX () const
 
DEPRECATED Mat getSnesA () const
 
DEPRECATED 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
 
Base funtions and integration points
MatrixDoublegetGaussPts ()
 matrix of integration (Gauss) points for Volume Element More...
 
auto getFTensor0IntegrationWeight ()
 Get integration weights. More...
 
Deprecated (do not use)
DEPRECATED MoFEMErrorCode getPorblemRowIndices (const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
 
- 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, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 Operator for bi-linear form, usually to calculate values on left hand side. More...
 
virtual MoFEMErrorCode opLhs (DataForcesAndSourcesCore &row_data, DataForcesAndSourcesCore &col_data)
 
virtual MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 Operator for linear form, usually to calculate values on right hand side. More...
 
virtual MoFEMErrorCode opRhs (DataForcesAndSourcesCore &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
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...
 

Protected Member Functions

virtual MoFEMErrorCode setPtrFE (ForcesAndSourcesCore *ptr)
 
ForcesAndSourcesCoregetPtrFE () const
 
ForcesAndSourcesCoregetSidePtrFE () const
 

Protected Attributes

ForcesAndSourcesCoreptrFE
 

Private Member Functions

MoFEMErrorCode loopSide (const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0)
 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...
 

Friends

class ForcesAndSourcesCore
 
class EdgeElementForcesAndSourcesCoreBase
 
class FaceElementForcesAndSourcesCoreBase
 
class ContactPrismElementForcesAndSourcesCore
 

Detailed Description

Data operator to do calculations at integration points.

Is inherited and implemented by user to do calculations. It can be used in many different ways but typically is used to integrate matrices (f.e. stiffness matrix) and the right hand vector (f.e. force vector).

Note: It is assumed that operator is executed for symmetric problem. That means that is executed for not repeating entities on finite element. For example on triangle we have nodes, 3 edges and one face. Because of symmetry calculations are for: nodes-nodes, nodes-edge0, nodes-edge_1, nodes-edge_2, nodes-face, edges_1-edges_1, edges_1-edge_1, edge_1-edge_2, edge_1-edge_1, edge_1-edge_2, edge_2-edge_2, edge_1-face, edge_1-face, edges_3-face, face - face

In case of non symmetric problem in addition calculations of lower off diagonal terms. F.e. edge_1-edge_0, esges_3-edge_0, edge_3-edge_1,

In that case class variable UserDataOperator::sYmm = false;

NoteL: By default sYmm is set for symmetric problems

Examples
continuity_check_on_contact_prism_side_ele.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, forces_and_sources_testing_users_base.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, hello_world.cpp, mesh_smoothing.cpp, prism_elements_from_surface.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, simple_elasticity.cpp, and simple_interface.cpp.

Definition at line 84 of file ForcesAndSourcesCore.hpp.

Member Enumeration Documentation

◆ 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.

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 
OPLAST 

Definition at line 98 of file ForcesAndSourcesCore.hpp.

Constructor & Destructor Documentation

◆ UserDataOperator() [1/3]

MoFEM::ForcesAndSourcesCore::UserDataOperator::UserDataOperator ( const FieldSpace  space,
const char  type = OPLAST,
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 1602 of file ForcesAndSourcesCore.cpp.

◆ UserDataOperator() [2/3]

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

◆ UserDataOperator() [3/3]

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

Member Function Documentation

◆ addOpType()

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

Add operator type.

Definition at line 910 of file ForcesAndSourcesCore.hpp.

◆ getDataCtx()

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

Definition at line 927 of file ForcesAndSourcesCore.hpp.

927  {
928  return getFEMethod()->data_ctx;
929 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getFEEntityHandle()

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

Return finite element entity handle.

Returns
Finite element entity handle
Examples
UnsaturatedFlow.hpp.

Definition at line 868 of file ForcesAndSourcesCore.hpp.

868  {
869  return getNumeredEntFiniteElementPtr()->getEnt();
870 }
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.

◆ getFEMethod()

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

Return raw pointer to Finite Element Method object.

Examples
heat_equation.hpp, UnsaturatedFlow.hpp, and wave_equation.hpp.

Definition at line 900 of file ForcesAndSourcesCore.hpp.

900  {
901  return ptrFE;
902 }

◆ getFEName()

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

Get name of the element.

Definition at line 922 of file ForcesAndSourcesCore.hpp.

922  {
923  return getFEMethod()->feName;
924 }
std::string feName
Name of finite element.
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getFTensor0IntegrationWeight()

auto MoFEM::ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight ( )

Get integration weights.

for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
// integrate something
++t_w;
}
Returns
FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
Examples
heat_equation.hpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, UnsaturatedFlow.hpp, and wave_equation.hpp.

Definition at line 1014 of file ForcesAndSourcesCore.hpp.

1014  {
1016  &(getGaussPts()(getGaussPts().size1() - 1, 0)));
1017 }
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element

◆ getGaussPts()

MatrixDouble & MoFEM::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
heat_equation.hpp, MagneticElement.hpp, UnsaturatedFlow.hpp, and wave_equation.hpp.

Definition at line 1010 of file ForcesAndSourcesCore.hpp.

1010  {
1011  return static_cast<ForcesAndSourcesCore *>(ptrFE)->gaussPts;
1012 }

◆ getKSPA()

Mat MoFEM::ForcesAndSourcesCore::UserDataOperator::getKSPA ( ) const

Definition at line 950 of file ForcesAndSourcesCore.hpp.

950  {
951  return getFEMethod()->ksp_A;
952 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getKSPB()

Mat MoFEM::ForcesAndSourcesCore::UserDataOperator::getKSPB ( ) const

Definition at line 954 of file ForcesAndSourcesCore.hpp.

954  {
955  return getFEMethod()->ksp_B;
956 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getKSPCtx()

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

Definition at line 932 of file ForcesAndSourcesCore.hpp.

932  {
933  return getFEMethod()->ksp_ctx;
934 }
KSPContext ksp_ctx
Context.
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getKSPf()

Vec MoFEM::ForcesAndSourcesCore::UserDataOperator::getKSPf ( ) const

Definition at line 946 of file ForcesAndSourcesCore.hpp.

946  {
947  return getFEMethod()->ksp_f;
948 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getLoopSize()

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

get size of elements in the loop

Returns
loop size

Definition at line 918 of file ForcesAndSourcesCore.hpp.

918  {
919  return getFEMethod()->getLoopSize();
920 }
int getLoopSize() const
get loop size
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getNinTheLoop()

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

get number of finite element in the loop

Returns
number of finite element

Definition at line 914 of file ForcesAndSourcesCore.hpp.

914  {
915  return getFEMethod()->getNinTheLoop();
916 }
int getNinTheLoop() const
get number of evaluated element in the loop
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getNumberOfNodesOnElement()

int MoFEM::ForcesAndSourcesCore::UserDataOperator::getNumberOfNodesOnElement ( )

Get the number of nodes on finite element.

Returns
int

Definition at line 894 of file ForcesAndSourcesCore.hpp.

894  {
895  int num_nodes;
896  CHKERR ptrFE->getNumberOfNodes(num_nodes);
897  return num_nodes;
898 }
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.

◆ getNumeredEntFiniteElementPtr()

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

Return raw pointer to NumeredEntFiniteElement.

Examples
ElasticityMixedFormulation.hpp, and UnsaturatedFlow.hpp.

Definition at line 864 of file ForcesAndSourcesCore.hpp.

864  {
866 };
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getOpType()

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

Get operator types.

Returns
Return operator type

Definition at line 904 of file ForcesAndSourcesCore.hpp.

◆ getPorblemRowIndices()

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

Definition at line 1019 of file ForcesAndSourcesCore.hpp.

1021  {
1022  return getProblemRowIndices(filed_name, type, side, indices);
1023 }
MoFEMErrorCode getProblemRowIndices(const std::string filed_name, const EntityType type, const int side, VectorInt &indices) const
Get row indices.

◆ getProblemColIndices()

MoFEMErrorCode MoFEM::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 1513 of file ForcesAndSourcesCore.cpp.

1515  {
1517  if (ptrFE == NULL)
1518  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1519 
1520  switch (type) {
1521  case MBVERTEX:
1522  CHKERR ptrFE->getProblemNodesColIndices(field_name, indices);
1523  break;
1524  default:
1525  CHKERR ptrFE->getProblemTypeColIndices(field_name, type, side, indices);
1526  }
1528 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getProblemRowIndices()

MoFEMErrorCode MoFEM::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 1496 of file ForcesAndSourcesCore.cpp.

1498  {
1500  if (ptrFE == NULL)
1501  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
1502 
1503  switch (type) {
1504  case MBVERTEX:
1505  CHKERR ptrFE->getProblemNodesRowIndices(field_name, indices);
1506  break;
1507  default:
1508  CHKERR ptrFE->getProblemTypeRowIndices(field_name, type, side, indices);
1509  }
1511 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
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
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getPtrFE()

ForcesAndSourcesCore * MoFEM::ForcesAndSourcesCore::UserDataOperator::getPtrFE ( ) const
protected

Definition at line 1025 of file ForcesAndSourcesCore.hpp.

1025  {
1026  return ptrFE;
1027 }

◆ getSideEntity()

EntityHandle MoFEM::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
}
}
Parameters
side_number
type

Definition at line 886 of file ForcesAndSourcesCore.hpp.

887  {
888  if (auto side_ptr = getSideNumberPtr(side_number, type).lock())
889  return side_ptr->ent;
890  else
891  return 0;
892 }
boost::weak_ptr< SideNumber > getSideNumberPtr(const int side_number, const EntityType type)
Get the side number pointer.

◆ getSideNumberPtr()

boost::weak_ptr< SideNumber > MoFEM::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 873 of file ForcesAndSourcesCore.hpp.

874  {
875  auto &side_table_by_side_and_type =
876  ptrFE->numeredEntFiniteElementPtr->getSideNumberTable().get<1>();
877  auto side_it =
878  side_table_by_side_and_type.find(boost::make_tuple(type, side_number));
879  if (side_it != side_table_by_side_and_type.end())
880  return *side_it;
881  else
882  return boost::weak_ptr<SideNumber>();
883 }
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getSidePtrFE()

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

Definition at line 1030 of file ForcesAndSourcesCore.hpp.

1030  {
1031  return ptrFE->sidePtrFE;
1032 }
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.

◆ getSNESA()

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

Definition at line 966 of file ForcesAndSourcesCore.hpp.

966  {
967  return getFEMethod()->snes_A;
968 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Mat & snes_A
jacobian matrix

◆ getSnesA()

DEPRECATED Mat MoFEM::ForcesAndSourcesCore::UserDataOperator::getSnesA ( ) const
Deprecated:
Use getSNESA intead

Definition at line 310 of file ForcesAndSourcesCore.hpp.

◆ getSNESB()

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

Definition at line 970 of file ForcesAndSourcesCore.hpp.

970  {
971  return getFEMethod()->snes_B;
972 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Mat & snes_B
preconditioner of jacobian matrix

◆ getSnesB()

DEPRECATED Mat MoFEM::ForcesAndSourcesCore::UserDataOperator::getSnesB ( ) const
Deprecated:
Use getSNESB intead

Definition at line 313 of file ForcesAndSourcesCore.hpp.

◆ getSNESCtx()

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

Definition at line 937 of file ForcesAndSourcesCore.hpp.

937  {
938  return getFEMethod()->snes_ctx;
939 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
SNESContext snes_ctx

◆ getSNESf()

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

Definition at line 958 of file ForcesAndSourcesCore.hpp.

958  {
959  return getFEMethod()->snes_f;
960 }
Vec & snes_f
residual
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getSnesF()

DEPRECATED Vec MoFEM::ForcesAndSourcesCore::UserDataOperator::getSnesF ( ) const
Deprecated:
Use getSNESF intead

Definition at line 304 of file ForcesAndSourcesCore.hpp.

◆ getSNESx()

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

Definition at line 962 of file ForcesAndSourcesCore.hpp.

962  {
963  return getFEMethod()->snes_x;
964 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Vec & snes_x
state vector

◆ getSnesX()

DEPRECATED Vec MoFEM::ForcesAndSourcesCore::UserDataOperator::getSnesX ( ) const
Deprecated:
Use getSNESX intead

Definition at line 307 of file ForcesAndSourcesCore.hpp.

◆ getTSA()

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

Definition at line 990 of file ForcesAndSourcesCore.hpp.

990  {
991  return getFEMethod()->ts_A;
992 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getTSa()

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

Definition at line 1006 of file ForcesAndSourcesCore.hpp.

1006  {
1007  return getFEMethod()->ts_a;
1008 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
PetscReal ts_a
shift for U_tt (see PETSc Time Solver)

◆ getTSB()

Mat MoFEM::ForcesAndSourcesCore::UserDataOperator::getTSB ( ) const
Examples
heat_equation.hpp, and wave_equation.hpp.

Definition at line 994 of file ForcesAndSourcesCore.hpp.

994  {
995  return getFEMethod()->ts_B;
996 }
Mat & ts_B
Preconditioner for ts_A.
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getTSCtx()

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

Definition at line 942 of file ForcesAndSourcesCore.hpp.

942  {
943  return getFEMethod()->ts_ctx;
944 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
TSContext ts_ctx

◆ getTSf()

Vec MoFEM::ForcesAndSourcesCore::UserDataOperator::getTSf ( ) const
Examples
heat_equation.hpp, and wave_equation.hpp.

Definition at line 986 of file ForcesAndSourcesCore.hpp.

986  {
987  return getFEMethod()->ts_F;
988 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Vec & ts_F
residual vector

◆ getTSstep()

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

Definition at line 998 of file ForcesAndSourcesCore.hpp.

998  {
999  return getFEMethod()->ts_step;
1000 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
PetscInt ts_step
time step

◆ getTStime()

double MoFEM::ForcesAndSourcesCore::UserDataOperator::getTStime ( ) const

Definition at line 1002 of file ForcesAndSourcesCore.hpp.

1002  {
1003  return getFEMethod()->ts_t;
1004 }
PetscReal ts_t
time
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.

◆ getTSu()

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

Definition at line 974 of file ForcesAndSourcesCore.hpp.

974  {
975  return getFEMethod()->ts_u;
976 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Vec & ts_u
state vector

◆ getTSu_t()

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

Definition at line 978 of file ForcesAndSourcesCore.hpp.

978  {
979  return getFEMethod()->ts_u_t;
980 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Vec & ts_u_t
time derivative of state vector

◆ getTSu_tt()

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

Definition at line 982 of file ForcesAndSourcesCore.hpp.

982  {
983  return getFEMethod()->ts_u_tt;
984 }
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Vec & ts_u_tt
second time derivative of state vector

◆ loopSide()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::UserDataOperator::loopSide ( const string &  fe_name,
ForcesAndSourcesCore side_fe,
const size_t  dim,
const EntityHandle  ent_for_side = 0 
)
private

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
Returns
MoFEMErrorCode

Definition at line 1537 of file ForcesAndSourcesCore.cpp.

1539  {
1541  const EntityHandle ent = ent_for_side ? ent_for_side : getFEEntityHandle();
1542 
1543  const Problem *problem_ptr = getFEMethod()->problemPtr;
1544  Range adjacent_ents;
1545  CHKERR ptrFE->mField.getInterface<BitRefManager>()->getAdjacenciesAny(
1546  ent, side_dim, adjacent_ents);
1547  typedef NumeredEntFiniteElement_multiIndex::index<
1548  Composite_Name_And_Ent_mi_tag>::type FEByComposite;
1549  FEByComposite &numered_fe =
1550  problem_ptr->numeredFiniteElements->get<Composite_Name_And_Ent_mi_tag>();
1551 
1552  side_fe->feName = fe_name;
1553 
1554  CHKERR side_fe->setSideFEPtr(ptrFE);
1555  CHKERR side_fe->copyBasicMethod(*getFEMethod());
1556  CHKERR side_fe->copyKsp(*getFEMethod());
1557  CHKERR side_fe->copySnes(*getFEMethod());
1558  CHKERR side_fe->copyTs(*getFEMethod());
1559 
1560  CHKERR side_fe->preProcess();
1561 
1562  int nn = 0;
1563  side_fe->loopSize = adjacent_ents.size();
1564  for (Range::iterator vit = adjacent_ents.begin(); vit != adjacent_ents.end();
1565  vit++) {
1566  FEByComposite::iterator miit =
1567  numered_fe.find(boost::make_tuple(fe_name, *vit));
1568  if (miit != numered_fe.end()) {
1569  side_fe->nInTheLoop = nn++;
1570  side_fe->numeredEntFiniteElementPtr = *miit;
1571  CHKERR (*side_fe)();
1572  }
1573  }
1574 
1575  CHKERR side_fe->postProcess();
1576 
1578 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const Problem * problemPtr
raw pointer to problem
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ setOpType()

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

Set operator type.

Parameters
Operatortype

Definition at line 906 of file ForcesAndSourcesCore.hpp.

◆ setPtrFE()

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

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

Definition at line 1644 of file ForcesAndSourcesCore.cpp.

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

Friends And Related Function Documentation

◆ ContactPrismElementForcesAndSourcesCore

Definition at line 415 of file ForcesAndSourcesCore.hpp.

◆ EdgeElementForcesAndSourcesCoreBase

Definition at line 413 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCoreBase

Definition at line 414 of file ForcesAndSourcesCore.hpp.

◆ ForcesAndSourcesCore

friend class ForcesAndSourcesCore
friend

Definition at line 412 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ colFieldName

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

Definition at line 107 of file ForcesAndSourcesCore.hpp.

◆ opType

char MoFEM::ForcesAndSourcesCore::UserDataOperator::opType

Definition at line 105 of file ForcesAndSourcesCore.hpp.

◆ ptrFE

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

Definition at line 387 of file ForcesAndSourcesCore.hpp.

◆ rowFieldName

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

Definition at line 106 of file ForcesAndSourcesCore.hpp.

◆ sPace

FieldSpace MoFEM::ForcesAndSourcesCore::UserDataOperator::sPace

Definition at line 108 of file ForcesAndSourcesCore.hpp.


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