![Logo](MoFEMLogo.png) |
| v0.14.0
|
#include <src/finite_elements/ForcesAndSourcesCore.hpp>
|
enum | OpType {
OPROW = 1 << 0,
OPCOL = 1 << 1,
OPROWCOL = 1 << 2,
OPSPACE = 1 << 3,
OPLAST = 1 << 3
} |
| Controls loop over entities on element. More...
|
|
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)> |
|
|
| 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 NumeredEntFiniteElement > | getNumeredEntFiniteElementPtr () 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< 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 () 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 FEMethod * | getFEMethod () 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...
|
|
ForcesAndSourcesCore * | getPtrFE () const |
|
ForcesAndSourcesCore * | getSidePtrFE () const |
|
ForcesAndSourcesCore * | getRefinePtrFE () const |
|
|
const PetscData::Switches & | getDataCtx () 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 |
|
|
Vec | getSNESf () const |
|
Vec | getSNESx () const |
|
Mat | getSNESA () const |
|
Mat | getSNESB () const |
|
|
Vec | getTSu () const |
|
Vec | getTSu_t () const |
|
Vec | getTSu_tt () const |
|
Vec | getTSf () const |
|
Mat | getTSA () const |
|
Mat | getTSB () const |
|
int | getTSstep () const |
|
double | getTStime () const |
|
double | getTStimeStep () const |
|
double | getTSa () const |
|
double | getTSaa () const |
|
|
MatrixDouble & | getGaussPts () |
| matrix of integration (Gauss) points for Volume Element More...
|
|
auto | getFTensor0IntegrationWeight () |
| Get integration weights. More...
|
|
|
MatrixDouble & | getCoordsAtGaussPts () |
| Gauss points and weight, matrix (nb. of points x 3) More...
|
|
auto | getFTensor1CoordsAtGaussPts () |
| Get coordinates at integration points assuming linear geometry. More...
|
|
|
double | getMeasure () const |
| get measure of element More...
|
|
double & | getMeasure () |
| get measure of element More...
|
|
| 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...
|
|
|
}
|
using | AdjCache = std::map< EntityHandle, std::vector< boost::weak_ptr< NumeredEntFiniteElement > >> |
|
MoFEMErrorCode | loopSide (const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy, AdjCache *adj_cache=nullptr) |
| User calls this function to loop over elements on the side of face. This function calls finite element with its operator to do calculations. More...
|
|
MoFEMErrorCode | loopThis (const string &fe_name, ForcesAndSourcesCore *this_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy) |
| User calls this function to loop over the same element using a different set of integration points. This function calls finite element with its operator to do calculations. More...
|
|
MoFEMErrorCode | loopParent (const string &fe_name, ForcesAndSourcesCore *parent_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy) |
| User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations. More...
|
|
MoFEMErrorCode | loopChildren (const string &fe_name, ForcesAndSourcesCore *child_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy) |
| User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations. More...
|
|
- Examples
- analytical_poisson_field_split.cpp, continuity_check_on_contact_prism_side_ele.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, dynamic_first_order_con_law.cpp, ElasticityMixedFormulation.hpp, EshelbianPlasticity.cpp, forces_and_sources_testing_users_base.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, hello_world.cpp, HookeElement.hpp, level_set.cpp, MagneticElement.hpp, mesh_smoothing.cpp, mortar_contact_thermal.cpp, NavierStokesElement.hpp, prism_elements_from_surface.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, schur_test_diag_mat.cpp, simple_contact_thermal.cpp, simple_elasticity.cpp, and simple_interface.cpp.
Definition at line 549 of file ForcesAndSourcesCore.hpp.
◆ AdjCache
◆ OpType
Controls loop over entities on element.
- OPROW is used if row vector is assembled
- OPCOL is usually used if column vector is assembled
- OPROWCOL is usually used for assemble matrices.
- OPSPACE no field is defined for such operator. Is usually used to modify base
For typical problem like Bubnov-Galerkin OPROW and OPCOL are the same. In more general case for example for non-square matrices columns and rows could have different numeration and/or different set of fields.
Enumerator |
---|
OPROW | operator doWork function is executed on FE rows
|
OPCOL | operator doWork function is executed on FE columns
|
OPROWCOL | operator doWork is executed on FE rows &columns
|
OPSPACE | operator do Work is execute on space data
|
OPLAST | - Deprecated:
- would be removed
|
Definition at line 566 of file ForcesAndSourcesCore.hpp.
◆ UserDataOperator() [1/3]
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
- EshelbianPlasticity.cpp.
Definition at line 1974 of file ForcesAndSourcesCore.cpp.
◆ UserDataOperator() [2/3]
UserDataOperator::UserDataOperator |
( |
const std::string |
field_name, |
|
|
const char |
type, |
|
|
const bool |
symm = true |
|
) |
| |
◆ UserDataOperator() [3/3]
UserDataOperator::UserDataOperator |
( |
const std::string |
row_field_name, |
|
|
const std::string |
col_field_name, |
|
|
const char |
type, |
|
|
const bool |
symm = true |
|
) |
| |
◆ addOpType()
void UserDataOperator::addOpType |
( |
const OpType |
type | ) |
|
|
inline |
◆ getCoordsAtGaussPts()
◆ getDataCtx()
◆ getFEDim()
int UserDataOperator::getFEDim |
( |
| ) |
const |
|
inline |
◆ getFEEntityHandle()
◆ getFEMethod()
const FEMethod * UserDataOperator::getFEMethod |
( |
| ) |
const |
|
inline |
◆ getFEName()
std::string UserDataOperator::getFEName |
( |
| ) |
const |
|
inline |
◆ getFEType()
EntityType UserDataOperator::getFEType |
( |
| ) |
const |
|
inline |
◆ getFTensor0IntegrationWeight()
auto UserDataOperator::getFTensor0IntegrationWeight |
( |
| ) |
|
|
inline |
◆ getFTensor1CoordsAtGaussPts()
auto UserDataOperator::getFTensor1CoordsAtGaussPts |
( |
| ) |
|
|
inline |
◆ 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
- analytical_poisson_field_split.cpp, MagneticElement.hpp, and UnsaturatedFlow.hpp.
Definition at line 1236 of file ForcesAndSourcesCore.hpp.
◆ getKSPA()
Mat UserDataOperator::getKSPA |
( |
| ) |
const |
|
inline |
◆ getKSPB()
Mat UserDataOperator::getKSPB |
( |
| ) |
const |
|
inline |
◆ getKSPCtx()
◆ getKSPf()
Vec UserDataOperator::getKSPf |
( |
| ) |
const |
|
inline |
◆ getLoopSize()
int UserDataOperator::getLoopSize |
( |
| ) |
const |
|
inline |
◆ getMeasure() [1/2]
double& MoFEM::ForcesAndSourcesCore::UserDataOperator::getMeasure |
( |
| ) |
|
|
inline |
get measure of element
- Returns
- volume
◆ getMeasure() [2/2]
double & UserDataOperator::getMeasure |
( |
| ) |
const |
|
inline |
◆ getNinTheLoop()
int UserDataOperator::getNinTheLoop |
( |
| ) |
const |
|
inline |
◆ getNumberOfNodesOnElement()
int UserDataOperator::getNumberOfNodesOnElement |
( |
| ) |
const |
|
inline |
◆ getNumeredEntFiniteElementPtr()
◆ getOpType()
int UserDataOperator::getOpType |
( |
| ) |
const |
|
inline |
◆ getProblemColIndices()
MoFEMErrorCode 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 | |
type | entity type |
side | side 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 1669 of file ForcesAndSourcesCore.cpp.
◆ getProblemRowIndices()
MoFEMErrorCode 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 | |
type | entity type |
side | side 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 1652 of file ForcesAndSourcesCore.cpp.
◆ getPtrFE()
◆ getRefinePtrFE()
◆ getSideEntity()
EntityHandle UserDataOperator::getSideEntity |
( |
const int |
side_number, |
|
|
const EntityType |
type |
|
) |
| |
|
inline |
Get the side entity.
- Note
- For vertex is expectation. Side basses in argument of function doWork is zero. For other entity types side can be used as argument of this function.
for (
int n = 0;
n != number_of_nodes; ++
n)
} else {
}
}
- Parameters
-
Definition at line 1030 of file ForcesAndSourcesCore.hpp.
1033 return side_ptr->ent;
◆ getSideNumberPtr()
boost::weak_ptr< SideNumber > UserDataOperator::getSideNumberPtr |
( |
const int |
side_number, |
|
|
const EntityType |
type |
|
) |
| |
|
inline |
Get the side number pointer.
- Note
- For vertex is expectation. Side basses in argument of function doWork is zero. For other entity types side can be used as argument of this function.
- Parameters
-
- Returns
- boost::weak_ptr<SideNumber>
Definition at line 1017 of file ForcesAndSourcesCore.hpp.
1019 auto &side_table_by_side_and_type =
1022 side_table_by_side_and_type.find(boost::make_tuple(
type, side_number));
1023 if (side_it != side_table_by_side_and_type.end())
1026 return boost::weak_ptr<SideNumber>();
◆ getSidePtrFE()
◆ getSNESA()
Mat UserDataOperator::getSNESA |
( |
| ) |
const |
|
inline |
◆ getSNESB()
Mat UserDataOperator::getSNESB |
( |
| ) |
const |
|
inline |
◆ getSNESCtx()
◆ getSNESf()
Vec UserDataOperator::getSNESf |
( |
| ) |
const |
|
inline |
◆ getSNESx()
Vec UserDataOperator::getSNESx |
( |
| ) |
const |
|
inline |
◆ getTSA()
Mat UserDataOperator::getTSA |
( |
| ) |
const |
|
inline |
◆ getTSa()
double UserDataOperator::getTSa |
( |
| ) |
const |
|
inline |
◆ getTSaa()
double UserDataOperator::getTSaa |
( |
| ) |
const |
|
inline |
◆ getTSB()
Mat UserDataOperator::getTSB |
( |
| ) |
const |
|
inline |
◆ getTSCtx()
◆ getTSf()
Vec UserDataOperator::getTSf |
( |
| ) |
const |
|
inline |
◆ getTSstep()
int UserDataOperator::getTSstep |
( |
| ) |
const |
|
inline |
◆ getTStime()
double UserDataOperator::getTStime |
( |
| ) |
const |
|
inline |
◆ getTStimeStep()
double UserDataOperator::getTStimeStep |
( |
| ) |
const |
|
inline |
◆ getTSu()
Vec UserDataOperator::getTSu |
( |
| ) |
const |
|
inline |
◆ getTSu_t()
Vec UserDataOperator::getTSu_t |
( |
| ) |
const |
|
inline |
◆ getTSu_tt()
Vec UserDataOperator::getTSu_tt |
( |
| ) |
const |
|
inline |
◆ loopChildren()
User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations.
- Parameters
-
fe_name | |
child_fe | |
verb | |
sev | |
- Returns
- MoFEMErrorCode
Definition at line 1871 of file ForcesAndSourcesCore.cpp.
1877 ->get<FiniteElement_name_mi_tag>()
1880 ->get<FiniteElement_name_mi_tag>()
1886 problem_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>();
1891 ref_ents.get<Composite_ParentEnt_And_EntType_mi_tag>().equal_range(
1892 boost::make_tuple(parent_type, parent_ent));
1894 if (
auto size = std::distance(range.first, range.second)) {
1896 std::vector<EntityHandle> childs_vec;
1897 childs_vec.reserve(size);
1898 for (; range.first != range.second; ++range.first)
1899 childs_vec.emplace_back((*range.first)->getEnt());
1903 if ((childs_vec.back() - childs_vec.front() + 1) == size)
1904 childs =
Range(childs_vec.front(), childs_vec.back());
1906 childs.insert_list(childs_vec.begin(), childs_vec.end());
1908 child_fe->feName = fe_name;
1919 CHKERR child_fe->preProcess();
1922 child_fe->loopSize = size;
1924 for (
auto p = childs.pair_begin(); p != childs.pair_end(); ++p) {
1928 p->first, (*fe_miit)->getFEUId()));
1931 p->second, (*fe_miit)->getFEUId()));
1933 for (; miit != hi_miit; ++miit) {
1937 <<
"Child finite element " <<
"(" << nn <<
"): " << **miit;
1939 child_fe->nInTheLoop = nn++;
1940 child_fe->numeredEntFiniteElementPtr = *miit;
1946 CHKERR child_fe->postProcess();
◆ loopParent()
User calls this function to loop over parent elements. This function calls finite element with its operator to do calculations.
- Parameters
-
fe_name | |
parent_fe | |
verb | |
sev | |
- Returns
- MoFEMErrorCode
Definition at line 1822 of file ForcesAndSourcesCore.cpp.
1829 auto fe_miit = fes.find(fe_name);
1830 if (fe_miit != fes.end()) {
1836 parent_fe->feName = fe_name;
1848 parent_ent, (*fe_miit)->getFEUId()));
1849 if (miit != numered_fe.end()) {
1851 MOFEM_LOG(
"SELF", sev) <<
"Parent finite element: " << **miit;
1852 parent_fe->loopSize = 1;
1853 parent_fe->nInTheLoop = 0;
1854 parent_fe->numeredEntFiniteElementPtr = *miit;
1855 CHKERR parent_fe->preProcess();
1857 CHKERR parent_fe->postProcess();
1860 MOFEM_LOG(
"SELF", sev) <<
"Parent finite element: no parent";
1861 parent_fe->loopSize = 0;
1862 parent_fe->nInTheLoop = 0;
1863 CHKERR parent_fe->preProcess();
1864 CHKERR parent_fe->postProcess();
◆ loopSide()
User calls this function to loop over elements on the side of face. This function calls finite element with its operator to do calculations.
- Parameters
-
fe_name | name of the side element |
side_fe | pointer to the side element instance |
dim | dimension the of side element |
ent_for_side | entity handle for which adjacent volume or face will be accessed |
verb | |
sev | |
adj_cache | |
- Returns
- MoFEMErrorCode
Definition at line 1700 of file ForcesAndSourcesCore.cpp.
1711 ->get<FiniteElement_name_mi_tag>()
1714 ->get<FiniteElement_name_mi_tag>()
1719 side_fe->feName = fe_name;
1730 CHKERR side_fe->preProcess();
1732 std::vector<boost::weak_ptr<NumeredEntFiniteElement>> fe_vec;
1733 auto get_numered_fe_ptr = [&](
auto &fe_uid,
Range &&adjacent_ents)
1734 -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1735 fe_vec.reserve(adjacent_ents.size());
1736 for (
auto fe_ent : adjacent_ents) {
1737 auto miit = numered_fe.find(
1739 if (miit != numered_fe.end()) {
1740 fe_vec.emplace_back(*miit);
1746 auto get_bit_entity_adjacency = [&]() {
1747 Range adjacent_ents;
1749 ent, side_dim, adjacent_ents);
1750 return adjacent_ents;
1753 auto get_adj = [&](
auto &fe_uid)
1754 -> std::vector<boost::weak_ptr<NumeredEntFiniteElement>> & {
1757 return (*adj_cache).at(ent);
1758 }
catch (
const std::out_of_range &) {
1759 return (*adj_cache)[ent] =
1760 get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1763 return get_numered_fe_ptr(fe_uid, get_bit_entity_adjacency());
1767 auto adj = get_adj((*fe_miit)->getFEUId());
1770 side_fe->loopSize = adj.size();
1771 for (
auto fe_weak_ptr : adj) {
1772 if (
auto fe_ptr = fe_weak_ptr.lock()) {
1775 <<
"Side finite element " <<
"(" << nn <<
"): " << *fe_ptr;
1776 side_fe->nInTheLoop = nn;
1777 side_fe->numeredEntFiniteElementPtr = fe_ptr;
1783 CHKERR side_fe->postProcess();
◆ loopThis()
User calls this function to loop over the same element using a different set of integration points. This function calls finite element with its operator to do calculations.
- Parameters
-
- Returns
- MoFEMErrorCode
Definition at line 1789 of file ForcesAndSourcesCore.cpp.
1795 MOFEM_LOG(
"SELF", sev) <<
"This finite element: "
1798 this_fe->feName = fe_name;
1809 CHKERR this_fe->preProcess();
1817 CHKERR this_fe->postProcess();
◆ setOpType()
void UserDataOperator::setOpType |
( |
const OpType |
type | ) |
|
|
inline |
◆ setPtrFE()
◆ ContactPrismElementForcesAndSourcesCore
◆ EdgeElementForcesAndSourcesCore
◆ FaceElementForcesAndSourcesCore
◆ ForcesAndSourcesCore
◆ colFieldName
std::string MoFEM::ForcesAndSourcesCore::UserDataOperator::colFieldName |
◆ opType
char MoFEM::ForcesAndSourcesCore::UserDataOperator::opType |
◆ OpTypeNames
const char *const UserDataOperator::OpTypeNames |
|
static |
◆ ptrFE
◆ rowFieldName
std::string MoFEM::ForcesAndSourcesCore::UserDataOperator::rowFieldName |
◆ sPace
FieldSpace MoFEM::ForcesAndSourcesCore::UserDataOperator::sPace |
The documentation for this struct was generated from the following files:
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
auto dimension_from_handle(const EntityHandle h)
get entity dimension form handle
EntityHandle getSideEntity(const int side_number, const EntityType type)
Get the side entity.
int getNinTheLoop() const
get number of evaluated element in the loop
boost::weak_ptr< SideNumber > getSideNumberPtr(const int side_number, const EntityType type)
Get the side number pointer.
PetscReal ts_a
shift for U_t (see PETSc Time Solver)
@ OPSPACE
operator do Work is execute on space data
static constexpr Switches CtxSetB
Vec & ts_F
residual vector
static constexpr Switches CtxSetA
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Mat & snes_B
preconditioner of jacobian matrix
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
Mat & snes_A
jacobian matrix
@ OPROWCOL
operator doWork is executed on FE rows &columns
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.
MoFEMErrorCode getProblemTypeRowIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
MoFEMErrorCode getProblemTypeColIndices(const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
DataOperator(const bool symm=true)
auto getFTensor0IntegrationWeight()
Get integration weights.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
friend class ForcesAndSourcesCore
#define CHKERR
Inline error check.
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.
int getLoopSize() const
get loop size
MatrixDouble coordsAtGaussPts
coordinated at gauss points
auto getNumberOfNodes() const
@ OPCOL
operator doWork function is executed on FE columns
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.
Vec & ts_u_t
time derivative of state vector
int getNinTheLoop() const
get number of finite element in the loop
int getLoopSize() const
get size of elements in the loop
static constexpr Switches CtxSetX_T
ForcesAndSourcesCore * getPtrFE() const
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
auto getFEName() const
get finite element name
MatrixDouble & getCoordsAtGaussPts()
Gauss points and weight, matrix (nb. of points x 3)
auto type_from_handle(const EntityHandle h)
get type from entity handle
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElementsPtr
store finite elements
EntitiesFieldData::EntData EntData
constexpr auto field_name
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
static constexpr Switches CtxSetX_TT
PetscInt ts_step
time step number
#define MOFEM_LOG(channel, severity)
Log.
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
PetscReal ts_aa
shift for U_tt shift for U_tt
MoFEMErrorCode getProblemNodesColIndices(const std::string &field_name, VectorInt &nodes_indices) const
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
boost::weak_ptr< CacheTuple > cacheWeakPtr
@ MOFEM_DATA_INCONSISTENCY
const Problem * problemPtr
raw pointer to problem
MatrixDouble gaussPts
Matrix of integration points.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscReal ts_dt
time step size
Mat & ts_B
Preconditioner for ts_A.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode getProblemNodesRowIndices(const std::string &field_name, VectorInt &nodes_indices) const
ForcesAndSourcesCore * ptrFE
Vec & ts_u_tt
second time derivative of state vector
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
@ OPROW
operator doWork function is executed on FE rows
KSPContext ksp_ctx
Context.