v0.13.1
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 42 of file ConvectiveMassElement.hpp.

Constructor & Destructor Documentation

◆ ConvectiveMassElement()

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

Definition at line 91 of file ConvectiveMassElement.cpp.

93 : feMassRhs(m_field), feMassLhs(m_field), feMassAuxLhs(m_field),
94 feVelRhs(m_field), feVelLhs(m_field), feTRhs(m_field), feTLhs(m_field),
95 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 1854 of file ConvectiveMassElement.cpp.

1857 {
1859
1860 //
1861
1862 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1864 velocity_field_name);
1866 velocity_field_name);
1868 velocity_field_name);
1870 element_name, spatial_position_field_name);
1872 element_name, spatial_position_field_name);
1874 element_name, spatial_position_field_name);
1875 if (mField.check_field(material_position_field_name)) {
1876 if (ale) {
1878 element_name, material_position_field_name);
1880 element_name, material_position_field_name);
1882 element_name, "DOT_" + material_position_field_name);
1883 }
1885 element_name, material_position_field_name);
1886 }
1888 element_name, "DOT_" + velocity_field_name);
1890 element_name, "DOT_" + spatial_position_field_name);
1891
1892 Range tets;
1893 if (bit.any()) {
1894 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1895 bit, BitRefLevel().set(), MBTET, tets);
1896 }
1897
1898 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1899 for (; sit != setOfBlocks.end(); sit++) {
1900 Range add_tets = sit->second.tEts;
1901 if (!tets.empty()) {
1902 add_tets = intersect(add_tets, tets);
1903 }
1905 element_name);
1906 }
1907
1909}
@ MF_ZERO
Definition: definitions.h:111
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
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:51
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 1964 of file ConvectiveMassElement.cpp.

1967 {
1969
1970 //
1971
1972 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1974 velocity_field_name);
1976 velocity_field_name);
1978 element_name, spatial_position_field_name);
1980 element_name, spatial_position_field_name);
1981 if (mField.check_field(material_position_field_name)) {
1982 if (ale) {
1984 element_name, material_position_field_name);
1986 element_name, material_position_field_name);
1988 element_name, "DOT_" + material_position_field_name);
1989 }
1991 element_name, material_position_field_name);
1992 }
1994 element_name, "DOT_" + velocity_field_name);
1996 element_name, "DOT_" + spatial_position_field_name);
1997
1998 Range tets;
1999 if (bit.any()) {
2000 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
2001 bit, BitRefLevel().set(), MBTET, tets);
2002 }
2003 if (intersected != NULL) {
2004 if (tets.empty()) {
2005 tets = *intersected;
2006 } else {
2007 tets = intersect(*intersected, tets);
2008 }
2009 }
2010
2011 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2012 for (; sit != setOfBlocks.end(); sit++) {
2013 Range add_tets = sit->second.tEts;
2014 if (!tets.empty()) {
2015 add_tets = intersect(add_tets, tets);
2016 }
2018 element_name);
2019 }
2020
2022}

◆ addHOOpsVol()

MoFEMErrorCode ConvectiveMassElement::addHOOpsVol ( )

Definition at line 106 of file ConvectiveMassElement.hpp.

106 {
108 auto add_ops = [&](auto &fe) {
109 return MoFEM::addHOOpsVol("MESH_NODE_POSITIONS", fe, true, false, false,
110 false);
111 };
112 CHKERR add_ops(feMassRhs);
113 CHKERR add_ops(feMassLhs);
114 CHKERR add_ops(feMassAuxLhs);
115 CHKERR add_ops(feVelRhs);
116 CHKERR add_ops(feTRhs);
117 CHKERR add_ops(feTLhs);
118 CHKERR add_ops(feEnergy);
120 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
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 1911 of file ConvectiveMassElement.cpp.

1914 {
1916
1917 //
1918
1919 CHKERR mField.add_finite_element(element_name, MF_ZERO);
1921 velocity_field_name);
1923 velocity_field_name);
1925 velocity_field_name);
1927 element_name, spatial_position_field_name);
1929 element_name, spatial_position_field_name);
1930 if (mField.check_field(material_position_field_name)) {
1931 if (ale) {
1933 element_name, material_position_field_name);
1935 element_name, "DOT_" + material_position_field_name);
1936 }
1938 element_name, material_position_field_name);
1939 }
1941 element_name, "DOT_" + velocity_field_name);
1943 element_name, "DOT_" + spatial_position_field_name);
1944
1945 Range tets;
1946 if (bit.any()) {
1947 CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1948 bit, BitRefLevel().set(), MBTET, tets);
1949 }
1950
1951 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1952 for (; sit != setOfBlocks.end(); sit++) {
1953 Range add_tets = sit->second.tEts;
1954 if (!tets.empty()) {
1955 add_tets = intersect(add_tets, tets);
1956 }
1958 element_name);
1959 }
1960
1962}

◆ getLoopFeEnergy()

MyVolumeFE & ConvectiveMassElement::getLoopFeEnergy ( )

get kinetic energy element

Definition at line 102 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassAuxLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassAuxLhs ( )

get lhs volume element for Kuu shell matrix

Definition at line 87 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassLhs ( )

get lhs volume element

Definition at line 82 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeMassRhs ( )

get rhs volume element

Definition at line 77 of file ConvectiveMassElement.hpp.

◆ getLoopFeTLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeTLhs ( )

get lhs volume element

Definition at line 99 of file ConvectiveMassElement.hpp.

◆ getLoopFeTRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeTRhs ( )

get rhs volume element

Definition at line 97 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelLhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeVelLhs ( )

get lhs volume element

Definition at line 94 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelRhs()

MyVolumeFE & ConvectiveMassElement::getLoopFeVelRhs ( )

get rhs volume element

Definition at line 92 of file ConvectiveMassElement.hpp.

◆ MultOpA()

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

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 560 of file ConvectiveMassElement.hpp.

560 {
562
563 void *void_ctx;
564 CHKERR MatShellGetContext(A, &void_ctx);
565 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
566 if (!ctx->iNitialized) {
567 CHKERR ctx->iNit();
568 }
569 CHKERR VecZeroEntries(f);
570 // Mult Ku
571 CHKERR VecScatterBegin(ctx->scatterU, x, ctx->u, INSERT_VALUES,
572 SCATTER_FORWARD);
573 CHKERR VecScatterEnd(ctx->scatterU, x, ctx->u, INSERT_VALUES,
574 SCATTER_FORWARD);
575 CHKERR MatMult(ctx->K, ctx->u, ctx->Ku);
576 CHKERR VecScatterBegin(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
577 SCATTER_REVERSE);
578 CHKERR VecScatterEnd(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
579 SCATTER_REVERSE);
580 // Mult Mv
581 CHKERR VecScatterBegin(ctx->scatterV, x, ctx->v, INSERT_VALUES,
582 SCATTER_FORWARD);
583 CHKERR VecScatterEnd(ctx->scatterV, x, ctx->v, INSERT_VALUES,
584 SCATTER_FORWARD);
585 CHKERR MatMult(ctx->M, ctx->v, ctx->Mv);
586 CHKERR VecScatterBegin(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
587 SCATTER_REVERSE);
588 CHKERR VecScatterEnd(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
589 SCATTER_REVERSE);
590 // Velocities
591 CHKERR VecAXPY(ctx->v, -ctx->ts_a, ctx->u);
592 // CHKERR VecScale(ctx->v,ctx->scale);
593 CHKERR VecScatterBegin(ctx->scatterV, ctx->v, f, INSERT_VALUES,
594 SCATTER_REVERSE);
595 CHKERR VecScatterEnd(ctx->scatterV, ctx->v, f, INSERT_VALUES,
596 SCATTER_REVERSE);
597 // Assemble
598 CHKERR VecAssemblyBegin(f);
599 CHKERR VecAssemblyEnd(f);
601 }
double A

◆ PCShellApplyOp()

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

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 685 of file ConvectiveMassElement.hpp.

685 {
687
688 void *void_ctx;
689 CHKERR PCShellGetContext(pc, &void_ctx);
690 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
691 MatShellCtx *shell_mat_ctx;
692 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
693 // forward
694 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
695 INSERT_VALUES, SCATTER_FORWARD);
696 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
697 INSERT_VALUES, SCATTER_FORWARD);
698 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
699 INSERT_VALUES, SCATTER_FORWARD);
700 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
701 INSERT_VALUES, SCATTER_FORWARD);
702 // CHKERR VecScale(shell_mat_ctx->v,1/shell_mat_ctx->scale);
703 // apply pre-conditioner and calculate u
704 CHKERR MatMult(shell_mat_ctx->M, shell_mat_ctx->v,
705 shell_mat_ctx->Mv); // Mrv
706 CHKERR VecAXPY(shell_mat_ctx->Ku, -1, shell_mat_ctx->Mv); // f-Mrv
707 CHKERR PCApply(ctx->pC, shell_mat_ctx->Ku,
708 shell_mat_ctx->u); // u = (aM+K)^(-1)(ru-Mrv)
709 // VecView(shell_mat_ctx->u,PETSC_VIEWER_STDOUT_WORLD);
710 // calculate velocities
711 CHKERR VecAXPY(shell_mat_ctx->v, shell_mat_ctx->ts_a,
712 shell_mat_ctx->u); // v = v + a*u
713 // VecView(shell_mat_ctx->v,PETSC_VIEWER_STDOUT_WORLD);
714 // reverse
715 CHKERR VecZeroEntries(x);
716 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
717 INSERT_VALUES, SCATTER_REVERSE);
718 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
719 INSERT_VALUES, SCATTER_REVERSE);
720 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
721 INSERT_VALUES, SCATTER_REVERSE);
722 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
723 INSERT_VALUES, SCATTER_REVERSE);
724 CHKERR VecAssemblyBegin(x);
725 CHKERR VecAssemblyEnd(x);
727 }

◆ PCShellDestroy()

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

Definition at line 647 of file ConvectiveMassElement.hpp.

647 {
649
650 void *void_ctx;
651 CHKERR PCShellGetContext(pc, &void_ctx);
652 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
653 CHKERR ctx->dEstroy();
655 }

◆ PCShellSetUpOp()

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

Definition at line 632 of file ConvectiveMassElement.hpp.

632 {
634
635 void *void_ctx;
636 CHKERR PCShellGetContext(pc, &void_ctx);
637 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
638 CHKERR ctx->iNit();
639 MatShellCtx *shell_mat_ctx;
640 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
641 CHKERR PCSetFromOptions(ctx->pC);
642 CHKERR PCSetOperators(ctx->pC, shell_mat_ctx->barK, shell_mat_ctx->barK);
643 CHKERR PCSetUp(ctx->pC);
645 }

◆ setBlocks() [1/2]

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

Definition at line 1780 of file ConvectiveMassElement.cpp.

1780 {
1782
1783 Range added_tets;
1785 mField, BLOCKSET | BODYFORCESSET, it)) {
1786 int id = it->getMeshsetId();
1787 EntityHandle meshset = it->getMeshset();
1788 CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET,
1789 setOfBlocks[id].tEts, true);
1790 added_tets.merge(setOfBlocks[id].tEts);
1791 Block_BodyForces mydata;
1792 CHKERR it->getAttributeDataStructure(mydata);
1793 setOfBlocks[id].rho0 = mydata.data.density;
1794 setOfBlocks[id].a0.resize(3);
1795 setOfBlocks[id].a0[0] = mydata.data.acceleration_x;
1796 setOfBlocks[id].a0[1] = mydata.data.acceleration_y;
1797 setOfBlocks[id].a0[2] = mydata.data.acceleration_z;
1798 // std::cerr << setOfBlocks[id].tEts << std::endl;
1799 }
1800
1802 mField, BLOCKSET | MAT_ELASTICSET, it)) {
1803 Mat_Elastic mydata;
1804 CHKERR it->getAttributeDataStructure(mydata);
1805 if (mydata.data.User1 == 0)
1806 continue;
1807 Range tets;
1808 EntityHandle meshset = it->getMeshset();
1809 CHKERR mField.get_moab().get_entities_by_type(meshset, MBTET, tets, true);
1810 tets = subtract(tets, added_tets);
1811 if (tets.empty())
1812 continue;
1813 int id = it->getMeshsetId();
1814 setOfBlocks[-id].tEts = tets;
1815 setOfBlocks[-id].rho0 = mydata.data.User1;
1816 setOfBlocks[-id].a0.resize(3);
1817 setOfBlocks[-id].a0[0] = mydata.data.User2;
1818 setOfBlocks[-id].a0[1] = mydata.data.User3;
1819 setOfBlocks[-id].a0[2] = mydata.data.User4;
1820 // std::cerr << setOfBlocks[id].tEts << std::endl;
1821 }
1822
1824}
@ BODYFORCESSET
block name is "BODY_FORCES"
Definition: definitions.h:175
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
Definition: definitions.h:172
@ BLOCKSET
Definition: definitions.h:161
#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 1826 of file ConvectiveMassElement.cpp.

1828 {
1830
1831 if (!block_sets_ptr)
1832 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1833 "Pointer to block of sets is null");
1834
1836 m_field, BLOCKSET | BODYFORCESSET, it)) {
1837 Block_BodyForces mydata;
1838 CHKERR it->getAttributeDataStructure(mydata);
1839 int id = it->getMeshsetId();
1840 auto &block_data = (*block_sets_ptr)[id];
1841 EntityHandle meshset = it->getMeshset();
1842 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, 3,
1843 block_data.tEts, true);
1844 block_data.rho0 = mydata.data.density;
1845 block_data.a0.resize(3);
1846 block_data.a0[0] = mydata.data.acceleration_x;
1847 block_data.a0[1] = mydata.data.acceleration_y;
1848 block_data.a0[2] = mydata.data.acceleration_z;
1849 }
1850
1852}
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44

◆ 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 2024 of file ConvectiveMassElement.cpp.

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

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

◆ 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 2262 of file ConvectiveMassElement.cpp.

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

◆ 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 2122 of file ConvectiveMassElement.cpp.

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

◆ ZeroEntriesOp()

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

Definition at line 603 of file ConvectiveMassElement.hpp.

603 {
605
606 void *void_ctx;
607 CHKERR MatShellGetContext(A, &void_ctx);
608 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
609 CHKERR MatZeroEntries(ctx->K);
610 CHKERR MatZeroEntries(ctx->M);
612 }

Member Data Documentation

◆ commonData

CommonData ConvectiveMassElement::commonData

Definition at line 162 of file ConvectiveMassElement.hpp.

◆ feEnergy

MyVolumeFE ConvectiveMassElement::feEnergy

calculate kinetic energy

Definition at line 101 of file ConvectiveMassElement.hpp.

◆ feMassAuxLhs

MyVolumeFE ConvectiveMassElement::feMassAuxLhs

calculate left hand side for tetrahedral elements for Kuu shell matrix

Definition at line 85 of file ConvectiveMassElement.hpp.

◆ feMassLhs

MyVolumeFE ConvectiveMassElement::feMassLhs

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

Definition at line 80 of file ConvectiveMassElement.hpp.

◆ feMassRhs

MyVolumeFE ConvectiveMassElement::feMassRhs

calculate right hand side for tetrahedral elements

Definition at line 76 of file ConvectiveMassElement.hpp.

◆ feTLhs

MyVolumeFE ConvectiveMassElement::feTLhs

calculate left hand side for tetrahedral elements

Definition at line 98 of file ConvectiveMassElement.hpp.

◆ feTRhs

MyVolumeFE ConvectiveMassElement::feTRhs

calculate right hand side for tetrahedral elements

Definition at line 96 of file ConvectiveMassElement.hpp.

◆ feVelLhs

MyVolumeFE ConvectiveMassElement::feVelLhs

calculate left hand side for tetrahedral elements

Definition at line 93 of file ConvectiveMassElement.hpp.

◆ feVelRhs

MyVolumeFE ConvectiveMassElement::feVelRhs

calculate right hand side for tetrahedral elements

Definition at line 91 of file ConvectiveMassElement.hpp.

◆ methodsOp

boost::ptr_vector<MethodForForceScaling> ConvectiveMassElement::methodsOp

Definition at line 164 of file ConvectiveMassElement.hpp.

◆ mField

MoFEM::Interface& ConvectiveMassElement::mField

Definition at line 122 of file ConvectiveMassElement.hpp.

◆ setOfBlocks

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

maps block set id with appropriate BlockData

Definition at line 136 of file ConvectiveMassElement.hpp.

◆ tAg

short int ConvectiveMassElement::tAg

Definition at line 123 of file ConvectiveMassElement.hpp.


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