v0.13.1
Loading...
Searching...
No Matches
Classes | Public Member Functions | Static Public Member Functions | Public Attributes | List of all members
ConvectiveMassElement Struct Reference

structure grouping operators and data used for calculation of mass (convective) element \ nonlinear_elastic_elem More...

#include <users_modules/basic_finite_elements/src/ConvectiveMassElement.hpp>

Collaboration diagram for ConvectiveMassElement:
[legend]

Classes

struct  BlockData
 data for calculation inertia forces More...
 
struct  CommonData
 common data used by volume elements More...
 
struct  CommonFunctions
 
struct  MatShellCtx
 
struct  MyVolumeFE
 definition of volume element More...
 
struct  OpEnergy
 
struct  OpEshelbyDynamicMaterialMomentumJacobian
 
struct  OpEshelbyDynamicMaterialMomentumLhs_dv
 
struct  OpEshelbyDynamicMaterialMomentumLhs_dx
 
struct  OpEshelbyDynamicMaterialMomentumLhs_dX
 
struct  OpEshelbyDynamicMaterialMomentumRhs
 
struct  OpGetCommonDataAtGaussPts
 
struct  OpGetDataAtGaussPts
 
struct  OpMassJacobian
 
struct  OpMassLhs_dM_dv
 
struct  OpMassLhs_dM_dx
 
struct  OpMassLhs_dM_dX
 
struct  OpMassRhs
 
struct  OpVelocityJacobian
 
struct  OpVelocityLhs_dV_dv
 
struct  OpVelocityLhs_dV_dx
 
struct  OpVelocityLhs_dV_dX
 
struct  OpVelocityRhs
 
struct  PCShellCtx
 
struct  ShellResidualElement
 
struct  UpdateAndControl
 Set fields DOT_. More...
 

Public Member Functions

MyVolumeFEgetLoopFeMassRhs ()
 get rhs volume element More...
 
MyVolumeFEgetLoopFeMassLhs ()
 get lhs volume element More...
 
MyVolumeFEgetLoopFeMassAuxLhs ()
 get lhs volume element for Kuu shell matrix More...
 
MyVolumeFEgetLoopFeVelRhs ()
 get rhs volume element More...
 
MyVolumeFEgetLoopFeVelLhs ()
 get lhs volume element More...
 
MyVolumeFEgetLoopFeTRhs ()
 get rhs volume element More...
 
MyVolumeFEgetLoopFeTLhs ()
 get lhs volume element More...
 
MyVolumeFEgetLoopFeEnergy ()
 get kinetic energy element More...
 
MoFEMErrorCode addHOOpsVol ()
 
 ConvectiveMassElement (MoFEM::Interface &m_field, short int tag)
 
MoFEMErrorCode setBlocks ()
 
MoFEMErrorCode addConvectiveMassElement (string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel())
 
MoFEMErrorCode addVelocityElement (string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel())
 
MoFEMErrorCode addEshelbyDynamicMaterialMomentum (string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel(), Range *intersected=NULL)
 
MoFEMErrorCode setConvectiveMassOperators (string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, bool linear=false)
 
MoFEMErrorCode setVelocityOperators (string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false)
 
MoFEMErrorCode setKinematicEshelbyOperators (string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", Range *forces_on_entities_ptr=NULL)
 
MoFEMErrorCode setShellMatrixMassOperators (string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool linear=false)
 

Static Public Member Functions

static MoFEMErrorCode setBlocks (MoFEM::Interface &m_field, boost::shared_ptr< map< int, BlockData > > &block_sets_ptr)
 
static MoFEMErrorCode MultOpA (Mat A, Vec x, Vec f)
 Mult operator for shell matrix. More...
 
static MoFEMErrorCode ZeroEntriesOp (Mat A)
 
static MoFEMErrorCode PCShellSetUpOp (PC pc)
 
static MoFEMErrorCode PCShellDestroy (PC pc)
 
static MoFEMErrorCode PCShellApplyOp (PC pc, Vec f, Vec x)
 apply pre-conditioner for shell matrix More...
 

Public Attributes

MyVolumeFE feMassRhs
 calculate right hand side for tetrahedral elements More...
 
MyVolumeFE feMassLhs
 
MyVolumeFE feMassAuxLhs
 
MyVolumeFE feVelRhs
 calculate right hand side for tetrahedral elements More...
 
MyVolumeFE feVelLhs
 calculate left hand side for tetrahedral elements More...
 
MyVolumeFE feTRhs
 calculate right hand side for tetrahedral elements More...
 
MyVolumeFE feTLhs
 calculate left hand side for tetrahedral elements More...
 
MyVolumeFE feEnergy
 calculate kinetic energy More...
 
MoFEM::InterfacemField
 
short int tAg
 
std::map< int, BlockDatasetOfBlocks
 maps block set id with appropriate BlockData More...
 
CommonData commonData
 
boost::ptr_vector< MethodForForceScalingmethodsOp
 

Detailed Description

structure grouping operators and data used for calculation of mass (convective) element \ nonlinear_elastic_elem

structure grouping operators and data used for calculation of nonlinear elastic element

In order to assemble matrices and right hand vectors, the loops over elements, entities over that elements and finally loop over integration points are executed.

Following implementation separate those three celeries of loops and to each loop attach operator.

In order to assemble matrices and right hand vectors, the loops over elements, entities over that elements and finally loop over integration points are executed.

Following implementation separate those three categories of loops and to each loop attach operator.

Examples
HookeElement.hpp, and nonlinear_dynamics.cpp.

Definition at line 28 of file ConvectiveMassElement.hpp.

Constructor & Destructor Documentation

◆ ConvectiveMassElement()

ConvectiveMassElement::ConvectiveMassElement ( MoFEM::Interface m_field,
short int  tag 
)

Definition at line 77 of file ConvectiveMassElement.cpp.

79 : feMassRhs(m_field), feMassLhs(m_field), feMassAuxLhs(m_field),
80 feVelRhs(m_field), feVelLhs(m_field), feTRhs(m_field), feTLhs(m_field),
81 feEnergy(m_field), mField(m_field), tAg(tag) {}
MyVolumeFE feEnergy
calculate kinetic energy
MyVolumeFE feVelRhs
calculate right hand side for tetrahedral elements
MyVolumeFE feTRhs
calculate right hand side for tetrahedral elements
MyVolumeFE feMassRhs
calculate right hand side for tetrahedral elements
MyVolumeFE feTLhs
calculate left hand side for tetrahedral elements
MyVolumeFE feVelLhs
calculate left hand side for tetrahedral elements

Member Function Documentation

◆ addConvectiveMassElement()

MoFEMErrorCode ConvectiveMassElement::addConvectiveMassElement ( string  element_name,
string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  ale = false,
BitRefLevel  bit = BitRefLevel() 
)

Definition at line 1843 of file ConvectiveMassElement.cpp.

1846 {
1848
1849 //
1850
1851 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1853 velocity_field_name);
1855 velocity_field_name);
1857 velocity_field_name);
1859 element_name, spatial_position_field_name);
1861 element_name, spatial_position_field_name);
1863 element_name, spatial_position_field_name);
1864 if (mField.check_field(material_position_field_name)) {
1865 if (ale) {
1867 element_name, material_position_field_name);
1869 element_name, material_position_field_name);
1871 element_name, "DOT_" + material_position_field_name);
1872 }
1874 element_name, material_position_field_name);
1875 }
1877 element_name, "DOT_" + velocity_field_name);
1879 element_name, "DOT_" + spatial_position_field_name);
1880
1881 Range tets;
1882 if (bit.any()) {
1883 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1884 bit, BitRefLevel().set(), MBTET, tets);
1885 }
1886
1887 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1888 for (; sit != setOfBlocks.end(); sit++) {
1889 Range add_tets = sit->second.tEts;
1890 if (!tets.empty()) {
1891 add_tets = intersect(add_tets, tets);
1892 }
1894 element_name);
1895 }
1896
1898}
@ MF_ZERO
Definition: definitions.h:98
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual bool check_field(const std::string &name) const =0
check if field is in database
auto bit
set bit
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
Managing BitRefLevels.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ addEshelbyDynamicMaterialMomentum()

MoFEMErrorCode ConvectiveMassElement::addEshelbyDynamicMaterialMomentum ( string  element_name,
string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  ale = false,
BitRefLevel  bit = BitRefLevel(),
Range intersected = NULL 
)

Definition at line 1953 of file ConvectiveMassElement.cpp.

1956 {
1958
1959 //
1960
1961 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1963 velocity_field_name);
1965 velocity_field_name);
1967 element_name, spatial_position_field_name);
1969 element_name, spatial_position_field_name);
1970 if (mField.check_field(material_position_field_name)) {
1971 if (ale) {
1973 element_name, material_position_field_name);
1975 element_name, material_position_field_name);
1977 element_name, "DOT_" + material_position_field_name);
1978 }
1980 element_name, material_position_field_name);
1981 }
1983 element_name, "DOT_" + velocity_field_name);
1985 element_name, "DOT_" + spatial_position_field_name);
1986
1987 Range tets;
1988 if (bit.any()) {
1989 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1990 bit, BitRefLevel().set(), MBTET, tets);
1991 }
1992 if (intersected != NULL) {
1993 if (tets.empty()) {
1994 tets = *intersected;
1995 } else {
1996 tets = intersect(*intersected, tets);
1997 }
1998 }
1999
2000 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2001 for (; sit != setOfBlocks.end(); sit++) {
2002 Range add_tets = sit->second.tEts;
2003 if (!tets.empty()) {
2004 add_tets = intersect(add_tets, tets);
2005 }
2007 element_name);
2008 }
2009
2011}

◆ addHOOpsVol()

MoFEMErrorCode ConvectiveMassElement::addHOOpsVol ( )
inline

Definition at line 92 of file ConvectiveMassElement.hpp.

92 {
94 auto add_ops = [&](auto &fe) {
95 return MoFEM::addHOOpsVol("MESH_NODE_POSITIONS", fe, true, false, false,
96 false);
97 };
98 CHKERR add_ops(feMassRhs);
99 CHKERR add_ops(feMassLhs);
100 CHKERR add_ops(feMassAuxLhs);
101 CHKERR add_ops(feVelRhs);
102 CHKERR add_ops(feTRhs);
103 CHKERR add_ops(feTLhs);
104 CHKERR add_ops(feEnergy);
106 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)

◆ addVelocityElement()

MoFEMErrorCode ConvectiveMassElement::addVelocityElement ( string  element_name,
string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  ale = false,
BitRefLevel  bit = BitRefLevel() 
)

Definition at line 1900 of file ConvectiveMassElement.cpp.

1903 {
1905
1906 //
1907
1908 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1910 velocity_field_name);
1912 velocity_field_name);
1914 velocity_field_name);
1916 element_name, spatial_position_field_name);
1918 element_name, spatial_position_field_name);
1919 if (mField.check_field(material_position_field_name)) {
1920 if (ale) {
1922 element_name, material_position_field_name);
1924 element_name, "DOT_" + material_position_field_name);
1925 }
1927 element_name, material_position_field_name);
1928 }
1930 element_name, "DOT_" + velocity_field_name);
1932 element_name, "DOT_" + spatial_position_field_name);
1933
1934 Range tets;
1935 if (bit.any()) {
1936 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1937 bit, BitRefLevel().set(), MBTET, tets);
1938 }
1939
1940 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1941 for (; sit != setOfBlocks.end(); sit++) {
1942 Range add_tets = sit->second.tEts;
1943 if (!tets.empty()) {
1944 add_tets = intersect(add_tets, tets);
1945 }
1947 element_name);
1948 }
1949
1951}

◆ getLoopFeEnergy()

MyVolumeFE & ConvectiveMassElement::getLoopFeEnergy ( )
inline

get kinetic energy element

Definition at line 88 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassAuxLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassAuxLhs ( )
inline

get lhs volume element for Kuu shell matrix

Definition at line 73 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassLhs ( )
inline

get lhs volume element

Definition at line 68 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassRhs ( )
inline

get rhs volume element

Definition at line 63 of file ConvectiveMassElement.hpp.

◆ getLoopFeTLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeTLhs ( )
inline

get lhs volume element

Definition at line 85 of file ConvectiveMassElement.hpp.

◆ getLoopFeTRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeTRhs ( )
inline

get rhs volume element

Definition at line 83 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeVelLhs ( )
inline

get lhs volume element

Definition at line 80 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeVelRhs ( )
inline

get rhs volume element

Definition at line 78 of file ConvectiveMassElement.hpp.

◆ MultOpA()

static MoFEMErrorCode ConvectiveMassElement::MultOpA ( Mat  A,
Vec  x,
Vec  f 
)
inlinestatic

Mult operator for shell matrix.

\[ \left[ \begin{array}{cc} \mathbf{M} & \mathbf{K} \\ \mathbf{I} & -\mathbf{I}a \end{array} \right] \left[ \begin{array}{c} \mathbf{v} \\ \mathbf{u} \end{array} \right] = \left[ \begin{array}{c} \mathbf{r}_u \\ \mathbf{r}_v \end{array} \right] \]

Examples
nonlinear_dynamics.cpp.

Definition at line 546 of file ConvectiveMassElement.hpp.

546 {
548
549 void *void_ctx;
550 CHKERR MatShellGetContext(A, &void_ctx);
551 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
552 if (!ctx->iNitialized) {
553 CHKERR ctx->iNit();
554 }
555 CHKERR VecZeroEntries(f);
556 // Mult Ku
557 CHKERR VecScatterBegin(ctx->scatterU, x, ctx->u, INSERT_VALUES,
558 SCATTER_FORWARD);
559 CHKERR VecScatterEnd(ctx->scatterU, x, ctx->u, INSERT_VALUES,
560 SCATTER_FORWARD);
561 CHKERR MatMult(ctx->K, ctx->u, ctx->Ku);
562 CHKERR VecScatterBegin(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
563 SCATTER_REVERSE);
564 CHKERR VecScatterEnd(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
565 SCATTER_REVERSE);
566 // Mult Mv
567 CHKERR VecScatterBegin(ctx->scatterV, x, ctx->v, INSERT_VALUES,
568 SCATTER_FORWARD);
569 CHKERR VecScatterEnd(ctx->scatterV, x, ctx->v, INSERT_VALUES,
570 SCATTER_FORWARD);
571 CHKERR MatMult(ctx->M, ctx->v, ctx->Mv);
572 CHKERR VecScatterBegin(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
573 SCATTER_REVERSE);
574 CHKERR VecScatterEnd(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
575 SCATTER_REVERSE);
576 // Velocities
577 CHKERR VecAXPY(ctx->v, -ctx->ts_a, ctx->u);
578 // CHKERR VecScale(ctx->v,ctx->scale);
579 CHKERR VecScatterBegin(ctx->scatterV, ctx->v, f, INSERT_VALUES,
580 SCATTER_REVERSE);
581 CHKERR VecScatterEnd(ctx->scatterV, ctx->v, f, INSERT_VALUES,
582 SCATTER_REVERSE);
583 // Assemble
584 CHKERR VecAssemblyBegin(f);
585 CHKERR VecAssemblyEnd(f);
587 }
constexpr AssemblyType A
Definition: plastic.cpp:46

◆ PCShellApplyOp()

static MoFEMErrorCode ConvectiveMassElement::PCShellApplyOp ( PC  pc,
Vec  f,
Vec  x 
)
inlinestatic

apply pre-conditioner for shell matrix

\[ \left[ \begin{array}{cc} \mathbf{M} & \mathbf{K} \\ \mathbf{I} & -\mathbf{I}a \end{array} \right] \left[ \begin{array}{c} \mathbf{v} \\ \mathbf{u} \end{array} \right] = \left[ \begin{array}{c} \mathbf{r}_u \\ \mathbf{r}_v \end{array} \right] \]

where \(\mathbf{v} = \mathbf{r}_v + a\mathbf{u}\) and \(\mathbf{u}=(a\mathbf{M}+\mathbf{K})^{-1}(\mathbf{r}_u - \mathbf{M}\mathbf{r}_v\).

Examples
nonlinear_dynamics.cpp.

Definition at line 671 of file ConvectiveMassElement.hpp.

671 {
673
674 void *void_ctx;
675 CHKERR PCShellGetContext(pc, &void_ctx);
676 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
677 MatShellCtx *shell_mat_ctx;
678 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
679 // forward
680 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
681 INSERT_VALUES, SCATTER_FORWARD);
682 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
683 INSERT_VALUES, SCATTER_FORWARD);
684 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
685 INSERT_VALUES, SCATTER_FORWARD);
686 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
687 INSERT_VALUES, SCATTER_FORWARD);
688 // CHKERR VecScale(shell_mat_ctx->v,1/shell_mat_ctx->scale);
689 // apply pre-conditioner and calculate u
690 CHKERR MatMult(shell_mat_ctx->M, shell_mat_ctx->v,
691 shell_mat_ctx->Mv); // Mrv
692 CHKERR VecAXPY(shell_mat_ctx->Ku, -1, shell_mat_ctx->Mv); // f-Mrv
693 CHKERR PCApply(ctx->pC, shell_mat_ctx->Ku,
694 shell_mat_ctx->u); // u = (aM+K)^(-1)(ru-Mrv)
695 // VecView(shell_mat_ctx->u,PETSC_VIEWER_STDOUT_WORLD);
696 // calculate velocities
697 CHKERR VecAXPY(shell_mat_ctx->v, shell_mat_ctx->ts_a,
698 shell_mat_ctx->u); // v = v + a*u
699 // VecView(shell_mat_ctx->v,PETSC_VIEWER_STDOUT_WORLD);
700 // reverse
701 CHKERR VecZeroEntries(x);
702 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
703 INSERT_VALUES, SCATTER_REVERSE);
704 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
705 INSERT_VALUES, SCATTER_REVERSE);
706 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
707 INSERT_VALUES, SCATTER_REVERSE);
708 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
709 INSERT_VALUES, SCATTER_REVERSE);
710 CHKERR VecAssemblyBegin(x);
711 CHKERR VecAssemblyEnd(x);
713 }

◆ PCShellDestroy()

static MoFEMErrorCode ConvectiveMassElement::PCShellDestroy ( PC  pc)
inlinestatic
Examples
nonlinear_dynamics.cpp.

Definition at line 633 of file ConvectiveMassElement.hpp.

633 {
635
636 void *void_ctx;
637 CHKERR PCShellGetContext(pc, &void_ctx);
638 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
639 CHKERR ctx->dEstroy();
641 }

◆ PCShellSetUpOp()

static MoFEMErrorCode ConvectiveMassElement::PCShellSetUpOp ( PC  pc)
inlinestatic
Examples
nonlinear_dynamics.cpp.

Definition at line 618 of file ConvectiveMassElement.hpp.

618 {
620
621 void *void_ctx;
622 CHKERR PCShellGetContext(pc, &void_ctx);
623 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
624 CHKERR ctx->iNit();
625 MatShellCtx *shell_mat_ctx;
626 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
627 CHKERR PCSetFromOptions(ctx->pC);
628 CHKERR PCSetOperators(ctx->pC, shell_mat_ctx->barK, shell_mat_ctx->barK);
629 CHKERR PCSetUp(ctx->pC);
631 }

◆ setBlocks() [1/2]

MoFEMErrorCode ConvectiveMassElement::setBlocks ( )
Examples
elasticity.cpp.

Definition at line 1769 of file ConvectiveMassElement.cpp.

1769 {
1771
1772 Range added_tets;
1774 mField, BLOCKSET | BODYFORCESSET, it)) {
1775 int id = it->getMeshsetId();
1776 EntityHandle meshset = it->getMeshset();
1777 CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET,
1778 setOfBlocks[id].tEts, true);
1779 added_tets.merge(setOfBlocks[id].tEts);
1780 Block_BodyForces mydata;
1781 CHKERR it->getAttributeDataStructure(mydata);
1782 setOfBlocks[id].rho0 = mydata.data.density;
1783 setOfBlocks[id].a0.resize(3);
1784 setOfBlocks[id].a0[0] = mydata.data.acceleration_x;
1785 setOfBlocks[id].a0[1] = mydata.data.acceleration_y;
1786 setOfBlocks[id].a0[2] = mydata.data.acceleration_z;
1787 // std::cerr << setOfBlocks[id].tEts << std::endl;
1788 }
1789
1791 mField, BLOCKSET | MAT_ELASTICSET, it)) {
1792 Mat_Elastic mydata;
1793 CHKERR it->getAttributeDataStructure(mydata);
1794 if (mydata.data.User1 == 0)
1795 continue;
1796 Range tets;
1797 EntityHandle meshset = it->getMeshset();
1798 CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET, tets, true);
1799 tets = subtract(tets, added_tets);
1800 if (tets.empty())
1801 continue;
1802 int id = it->getMeshsetId();
1803 setOfBlocks[-id].tEts = tets;
1804 setOfBlocks[-id].rho0 = mydata.data.User1;
1805 setOfBlocks[-id].a0.resize(3);
1806 setOfBlocks[-id].a0[0] = mydata.data.User2;
1807 setOfBlocks[-id].a0[1] = mydata.data.User3;
1808 setOfBlocks[-id].a0[2] = mydata.data.User4;
1809 // std::cerr << setOfBlocks[id].tEts << std::endl;
1810 }
1811
1813}
@ BODYFORCESSET
block name is "BODY_FORCES"
Definition: definitions.h:162
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
Definition: definitions.h:159
@ BLOCKSET
Definition: definitions.h:148
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
Body force data structure.
virtual moab::Interface & get_moab()=0
Elastic material data structure.

◆ setBlocks() [2/2]

MoFEMErrorCode ConvectiveMassElement::setBlocks ( MoFEM::Interface m_field,
boost::shared_ptr< map< int, BlockData > > &  block_sets_ptr 
)
static

Definition at line 1815 of file ConvectiveMassElement.cpp.

1817 {
1819
1820 if (!block_sets_ptr)
1821 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1822 "Pointer to block of sets is null");
1823
1825 m_field, BLOCKSET | BODYFORCESSET, it)) {
1826 Block_BodyForces mydata;
1827 CHKERR it->getAttributeDataStructure(mydata);
1828 int id = it->getMeshsetId();
1829 auto &block_data = (*block_sets_ptr)[id];
1830 EntityHandle meshset = it->getMeshset();
1831 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, 3,
1832 block_data.tEts, true);
1833 block_data.rho0 = mydata.data.density;
1834 block_data.a0.resize(3);
1835 block_data.a0[0] = mydata.data.acceleration_x;
1836 block_data.a0[1] = mydata.data.acceleration_y;
1837 block_data.a0[2] = mydata.data.acceleration_z;
1838 }
1839
1841}
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31

◆ setConvectiveMassOperators()

MoFEMErrorCode ConvectiveMassElement::setConvectiveMassOperators ( string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  ale = false,
bool  linear = false 
)

Definition at line 2013 of file ConvectiveMassElement.cpp.

2015 {
2017
2018 commonData.spatialPositions = spatial_position_field_name;
2019 commonData.meshPositions = material_position_field_name;
2020 commonData.spatialVelocities = velocity_field_name;
2021 commonData.lInear = linear;
2022
2023 // Rhs
2024 feMassRhs.getOpPtrVector().push_back(
2025 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2026 feMassRhs.getOpPtrVector().push_back(
2027 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2028 feMassRhs.getOpPtrVector().push_back(
2029 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2030 feMassRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2031 "DOT_" + spatial_position_field_name, commonData));
2032 if (mField.check_field(material_position_field_name)) {
2033 feMassRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2034 material_position_field_name, commonData));
2035 if (ale) {
2036 feMassRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2037 "DOT_" + material_position_field_name, commonData));
2038 } else {
2039 feMassRhs.meshPositionsFieldName = material_position_field_name;
2040 }
2041 }
2042 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2043 for (; sit != setOfBlocks.end(); sit++) {
2044 feMassRhs.getOpPtrVector().push_back(
2045 new OpMassJacobian(spatial_position_field_name, sit->second, commonData,
2046 methodsOp, tAg, false));
2047 feMassRhs.getOpPtrVector().push_back(
2048 new OpMassRhs(spatial_position_field_name, sit->second, commonData));
2049 }
2050
2051 // Lhs
2052 feMassLhs.getOpPtrVector().push_back(
2053 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2054 feMassLhs.getOpPtrVector().push_back(
2055 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2056 feMassLhs.getOpPtrVector().push_back(
2057 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2058 feMassLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2059 "DOT_" + spatial_position_field_name, commonData));
2060 if (mField.check_field(material_position_field_name)) {
2061 feMassLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2062 material_position_field_name, commonData));
2063 if (ale) {
2064 feMassLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2065 "DOT_" + material_position_field_name, commonData));
2066 } else {
2067 feMassLhs.meshPositionsFieldName = material_position_field_name;
2068 }
2069 }
2070 sit = setOfBlocks.begin();
2071 for (; sit != setOfBlocks.end(); sit++) {
2072 feMassLhs.getOpPtrVector().push_back(
2073 new OpMassJacobian(spatial_position_field_name, sit->second, commonData,
2074 methodsOp, tAg, true));
2075 feMassLhs.getOpPtrVector().push_back(
2076 new OpMassLhs_dM_dv(spatial_position_field_name, velocity_field_name,
2077 sit->second, commonData));
2078 feMassLhs.getOpPtrVector().push_back(new OpMassLhs_dM_dx(
2079 spatial_position_field_name, spatial_position_field_name, sit->second,
2080 commonData));
2081 if (mField.check_field(material_position_field_name)) {
2082 if (ale) {
2083 feMassLhs.getOpPtrVector().push_back(new OpMassLhs_dM_dX(
2084 spatial_position_field_name, material_position_field_name,
2085 sit->second, commonData));
2086 } else {
2087 feMassLhs.meshPositionsFieldName = material_position_field_name;
2088 }
2089 }
2090 }
2091
2092 // Energy
2093 feEnergy.getOpPtrVector().push_back(
2094 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2095 feEnergy.getOpPtrVector().push_back(
2096 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2097 if (mField.check_field(material_position_field_name)) {
2098 feEnergy.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2099 material_position_field_name, commonData));
2100 feEnergy.meshPositionsFieldName = material_position_field_name;
2101 }
2102 sit = setOfBlocks.begin();
2103 for (; sit != setOfBlocks.end(); sit++) {
2104 feEnergy.getOpPtrVector().push_back(new OpEnergy(
2105 spatial_position_field_name, sit->second, commonData, feEnergy.V));
2106 }
2107
2109}
boost::ptr_vector< MethodForForceScaling > methodsOp

◆ setKinematicEshelbyOperators()

MoFEMErrorCode ConvectiveMassElement::setKinematicEshelbyOperators ( string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
Range forces_on_entities_ptr = NULL 
)

Definition at line 2189 of file ConvectiveMassElement.cpp.

2191 {
2193
2194 commonData.spatialPositions = spatial_position_field_name;
2195 commonData.meshPositions = material_position_field_name;
2196 commonData.spatialVelocities = velocity_field_name;
2197
2198 // Rhs
2199 feTRhs.getOpPtrVector().push_back(
2200 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2201 feTRhs.getOpPtrVector().push_back(
2202 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2203 feTRhs.getOpPtrVector().push_back(
2204 new OpGetCommonDataAtGaussPts(material_position_field_name, commonData));
2205 feTRhs.getOpPtrVector().push_back(
2206 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2207
2208 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2209 for (; sit != setOfBlocks.end(); sit++) {
2210 feTRhs.getOpPtrVector().push_back(
2211 new OpEshelbyDynamicMaterialMomentumJacobian(
2212 material_position_field_name, sit->second, commonData, tAg, false));
2213 feTRhs.getOpPtrVector().push_back(new OpEshelbyDynamicMaterialMomentumRhs(
2214 material_position_field_name, sit->second, commonData,
2215 forces_on_entities_ptr));
2216 }
2217
2218 // Lhs
2219 feTLhs.getOpPtrVector().push_back(
2220 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2221 feTLhs.getOpPtrVector().push_back(
2222 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2223 feTLhs.getOpPtrVector().push_back(
2224 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2225 if (mField.check_field(material_position_field_name)) {
2226 feTLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2227 material_position_field_name, commonData));
2228 }
2229 sit = setOfBlocks.begin();
2230 for (; sit != setOfBlocks.end(); sit++) {
2231 feTLhs.getOpPtrVector().push_back(
2232 new OpEshelbyDynamicMaterialMomentumJacobian(
2233 material_position_field_name, sit->second, commonData, tAg));
2234 feTLhs.getOpPtrVector().push_back(
2235 new OpEshelbyDynamicMaterialMomentumLhs_dv(
2236 material_position_field_name, velocity_field_name, sit->second,
2237 commonData, forces_on_entities_ptr));
2238 feTLhs.getOpPtrVector().push_back(
2239 new OpEshelbyDynamicMaterialMomentumLhs_dx(
2240 material_position_field_name, spatial_position_field_name,
2241 sit->second, commonData, forces_on_entities_ptr));
2242 feTLhs.getOpPtrVector().push_back(
2243 new OpEshelbyDynamicMaterialMomentumLhs_dX(
2244 material_position_field_name, material_position_field_name,
2245 sit->second, commonData, forces_on_entities_ptr));
2246 }
2247
2249}

◆ setShellMatrixMassOperators()

MoFEMErrorCode ConvectiveMassElement::setShellMatrixMassOperators ( string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  linear = false 
)

Definition at line 2251 of file ConvectiveMassElement.cpp.

2253 {
2255
2256 commonData.spatialPositions = spatial_position_field_name;
2257 commonData.meshPositions = material_position_field_name;
2258 commonData.spatialVelocities = velocity_field_name;
2259 commonData.lInear = linear;
2260
2261 // Rhs
2262 feMassRhs.getOpPtrVector().push_back(
2263 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2264 feMassRhs.getOpPtrVector().push_back(
2265 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2266 feMassRhs.getOpPtrVector().push_back(
2267 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2268 if (mField.check_field(material_position_field_name)) {
2269 feMassRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2270 material_position_field_name, commonData));
2271 feMassRhs.meshPositionsFieldName = material_position_field_name;
2272 }
2273 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2274 for (; sit != setOfBlocks.end(); sit++) {
2275 feMassRhs.getOpPtrVector().push_back(
2276 new OpMassJacobian(spatial_position_field_name, sit->second, commonData,
2277 methodsOp, tAg, false));
2278 feMassRhs.getOpPtrVector().push_back(
2279 new OpMassRhs(spatial_position_field_name, sit->second, commonData));
2280 }
2281
2282 // Lhs
2283 feMassLhs.getOpPtrVector().push_back(
2284 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2285 feMassLhs.getOpPtrVector().push_back(
2286 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2287 feMassLhs.getOpPtrVector().push_back(
2288 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2289 if (mField.check_field(material_position_field_name)) {
2290 feMassLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2291 material_position_field_name, commonData));
2292 feMassLhs.meshPositionsFieldName = material_position_field_name;
2293 }
2294 sit = setOfBlocks.begin();
2295 for (; sit != setOfBlocks.end(); sit++) {
2296 feMassLhs.getOpPtrVector().push_back(
2297 new OpMassJacobian(spatial_position_field_name, sit->second, commonData,
2298 methodsOp, tAg, true));
2299 feMassLhs.getOpPtrVector().push_back(new OpMassLhs_dM_dv(
2300 spatial_position_field_name, spatial_position_field_name, sit->second,
2301 commonData));
2302 if (mField.check_field(material_position_field_name)) {
2303 feMassLhs.meshPositionsFieldName = material_position_field_name;
2304 }
2305 }
2306
2307 // Aux Lhs
2308 feMassAuxLhs.getOpPtrVector().push_back(
2309 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2310 feMassAuxLhs.getOpPtrVector().push_back(
2311 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2312 feMassAuxLhs.getOpPtrVector().push_back(
2313 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2314 if (mField.check_field(material_position_field_name)) {
2315 feMassAuxLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2316 material_position_field_name, commonData));
2317 feMassAuxLhs.meshPositionsFieldName = material_position_field_name;
2318 }
2319 sit = setOfBlocks.begin();
2320 for (; sit != setOfBlocks.end(); sit++) {
2321 feMassAuxLhs.getOpPtrVector().push_back(
2322 new OpMassJacobian(spatial_position_field_name, sit->second, commonData,
2323 methodsOp, tAg, true));
2324 feMassAuxLhs.getOpPtrVector().push_back(new OpMassLhs_dM_dx(
2325 spatial_position_field_name, spatial_position_field_name, sit->second,
2326 commonData));
2327 if (mField.check_field(material_position_field_name)) {
2328 feMassAuxLhs.meshPositionsFieldName = material_position_field_name;
2329 }
2330 }
2331
2332 // Energy E=0.5*rho*v*v
2333 feEnergy.getOpPtrVector().push_back(
2334 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2335 feEnergy.getOpPtrVector().push_back(
2336 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2337 if (mField.check_field(material_position_field_name)) {
2338 feEnergy.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2339 material_position_field_name, commonData));
2340 feEnergy.meshPositionsFieldName = material_position_field_name;
2341 }
2342 sit = setOfBlocks.begin();
2343 for (; sit != setOfBlocks.end(); sit++) {
2344 feEnergy.getOpPtrVector().push_back(new OpEnergy(
2345 spatial_position_field_name, sit->second, commonData, feEnergy.V));
2346 }
2347
2349}

◆ setVelocityOperators()

MoFEMErrorCode ConvectiveMassElement::setVelocityOperators ( string  velocity_field_name,
string  spatial_position_field_name,
string  material_position_field_name = "MESH_NODE_POSITIONS",
bool  ale = false 
)

Definition at line 2111 of file ConvectiveMassElement.cpp.

2113 {
2115
2116 commonData.spatialPositions = spatial_position_field_name;
2117 commonData.meshPositions = material_position_field_name;
2118 commonData.spatialVelocities = velocity_field_name;
2119
2120 // Rhs
2121 feVelRhs.getOpPtrVector().push_back(
2122 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2123 feVelRhs.getOpPtrVector().push_back(
2124 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2125 feVelRhs.getOpPtrVector().push_back(
2126 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2127 if (mField.check_field(material_position_field_name)) {
2128 feVelRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2129 "DOT_" + spatial_position_field_name, commonData));
2130 feVelRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2131 material_position_field_name, commonData));
2132 if (ale) {
2133 feVelRhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2134 "DOT_" + material_position_field_name, commonData));
2135 } else {
2136 feVelRhs.meshPositionsFieldName = material_position_field_name;
2137 }
2138 }
2139 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2140 for (; sit != setOfBlocks.end(); sit++) {
2141 feVelRhs.getOpPtrVector().push_back(new OpVelocityJacobian(
2142 velocity_field_name, sit->second, commonData, tAg, false));
2143 feVelRhs.getOpPtrVector().push_back(
2144 new OpVelocityRhs(velocity_field_name, sit->second, commonData));
2145 }
2146
2147 // Lhs
2148 feVelLhs.getOpPtrVector().push_back(
2149 new OpGetCommonDataAtGaussPts(velocity_field_name, commonData));
2150 feVelLhs.getOpPtrVector().push_back(
2151 new OpGetCommonDataAtGaussPts(spatial_position_field_name, commonData));
2152 feVelLhs.getOpPtrVector().push_back(
2153 new OpGetCommonDataAtGaussPts("DOT_" + velocity_field_name, commonData));
2154 if (mField.check_field(material_position_field_name)) {
2155 feVelLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2156 "DOT_" + spatial_position_field_name, commonData));
2157 feVelLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2158 material_position_field_name, commonData));
2159 if (ale) {
2160 feVelLhs.getOpPtrVector().push_back(new OpGetCommonDataAtGaussPts(
2161 "DOT_" + material_position_field_name, commonData));
2162 } else {
2163 feVelLhs.meshPositionsFieldName = material_position_field_name;
2164 }
2165 }
2166 sit = setOfBlocks.begin();
2167 for (; sit != setOfBlocks.end(); sit++) {
2168 feVelLhs.getOpPtrVector().push_back(new OpVelocityJacobian(
2169 velocity_field_name, sit->second, commonData, tAg));
2170 feVelLhs.getOpPtrVector().push_back(new OpVelocityLhs_dV_dv(
2171 velocity_field_name, velocity_field_name, sit->second, commonData));
2172 feVelLhs.getOpPtrVector().push_back(new OpVelocityLhs_dV_dx(
2173 velocity_field_name, spatial_position_field_name, sit->second,
2174 commonData));
2175 if (mField.check_field(material_position_field_name)) {
2176 if (ale) {
2177 feVelLhs.getOpPtrVector().push_back(new OpVelocityLhs_dV_dX(
2178 velocity_field_name, material_position_field_name, sit->second,
2179 commonData));
2180 } else {
2181 feVelLhs.meshPositionsFieldName = material_position_field_name;
2182 }
2183 }
2184 }
2185
2187}

◆ ZeroEntriesOp()

static MoFEMErrorCode ConvectiveMassElement::ZeroEntriesOp ( Mat  A)
inlinestatic
Examples
nonlinear_dynamics.cpp.

Definition at line 589 of file ConvectiveMassElement.hpp.

589 {
591
592 void *void_ctx;
593 CHKERR MatShellGetContext(A, &void_ctx);
594 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
595 CHKERR MatZeroEntries(ctx->K);
596 CHKERR MatZeroEntries(ctx->M);
598 }

Member Data Documentation

◆ commonData

CommonData ConvectiveMassElement::commonData

Definition at line 148 of file ConvectiveMassElement.hpp.

◆ feEnergy

MyVolumeFE ConvectiveMassElement::feEnergy

calculate kinetic energy

Definition at line 87 of file ConvectiveMassElement.hpp.

◆ feMassAuxLhs

MyVolumeFE ConvectiveMassElement::feMassAuxLhs

calculate left hand side for tetrahedral elements for Kuu shell matrix

Definition at line 71 of file ConvectiveMassElement.hpp.

◆ feMassLhs

MyVolumeFE ConvectiveMassElement::feMassLhs

calculate left hand side for tetrahedral elements,i.e. mass element

Definition at line 66 of file ConvectiveMassElement.hpp.

◆ feMassRhs

MyVolumeFE ConvectiveMassElement::feMassRhs

calculate right hand side for tetrahedral elements

Definition at line 62 of file ConvectiveMassElement.hpp.

◆ feTLhs

MyVolumeFE ConvectiveMassElement::feTLhs

calculate left hand side for tetrahedral elements

Definition at line 84 of file ConvectiveMassElement.hpp.

◆ feTRhs

MyVolumeFE ConvectiveMassElement::feTRhs

calculate right hand side for tetrahedral elements

Definition at line 82 of file ConvectiveMassElement.hpp.

◆ feVelLhs

MyVolumeFE ConvectiveMassElement::feVelLhs

calculate left hand side for tetrahedral elements

Definition at line 79 of file ConvectiveMassElement.hpp.

◆ feVelRhs

MyVolumeFE ConvectiveMassElement::feVelRhs

calculate right hand side for tetrahedral elements

Definition at line 77 of file ConvectiveMassElement.hpp.

◆ methodsOp

boost::ptr_vector<MethodForForceScaling> ConvectiveMassElement::methodsOp

Definition at line 150 of file ConvectiveMassElement.hpp.

◆ mField

MoFEM::Interface& ConvectiveMassElement::mField

Definition at line 108 of file ConvectiveMassElement.hpp.

◆ setOfBlocks

std::map<int, BlockData> ConvectiveMassElement::setOfBlocks

maps block set id with appropriate BlockData

Definition at line 122 of file ConvectiveMassElement.hpp.

◆ tAg

short int ConvectiveMassElement::tAg

Definition at line 109 of file ConvectiveMassElement.hpp.


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