v0.12.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
 

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

135  : feMassRhs(m_field), feMassLhs(m_field), feMassAuxLhs(m_field),
136  feVelRhs(m_field), feVelLhs(m_field), feTRhs(m_field), feTLhs(m_field),
137  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 1844 of file ConvectiveMassElement.cpp.

1847  {
1849 
1850  //
1851 
1852  CHKERR mField.add_finite_element(element_name, MF_ZERO);
1854  velocity_field_name);
1856  velocity_field_name);
1858  velocity_field_name);
1860  element_name, spatial_position_field_name);
1862  element_name, spatial_position_field_name);
1864  element_name, spatial_position_field_name);
1865  if (mField.check_field(material_position_field_name)) {
1866  if (ale) {
1868  element_name, material_position_field_name);
1870  element_name, material_position_field_name);
1872  element_name, "DOT_" + material_position_field_name);
1873  }
1875  element_name, material_position_field_name);
1876  }
1878  element_name, "DOT_" + velocity_field_name);
1880  element_name, "DOT_" + spatial_position_field_name);
1881 
1882  Range tets;
1883  if (bit.any()) {
1884  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1885  bit, BitRefLevel().set(), MBTET, tets);
1886  }
1887 
1888  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1889  for (; sit != setOfBlocks.end(); sit++) {
1890  Range add_tets = sit->second.tEts;
1891  if (!tets.empty()) {
1892  add_tets = intersect(add_tets, tets);
1893  }
1895  element_name);
1896  }
1897 
1899 }
@ MF_ZERO
Definition: definitions.h:109
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:440
#define CHKERR
Inline error check.
Definition: definitions.h:528
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:433
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
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 1954 of file ConvectiveMassElement.cpp.

1957  {
1959 
1960  //
1961 
1962  CHKERR mField.add_finite_element(element_name, MF_ZERO);
1964  velocity_field_name);
1966  velocity_field_name);
1968  element_name, spatial_position_field_name);
1970  element_name, spatial_position_field_name);
1971  if (mField.check_field(material_position_field_name)) {
1972  if (ale) {
1974  element_name, material_position_field_name);
1976  element_name, material_position_field_name);
1978  element_name, "DOT_" + material_position_field_name);
1979  }
1981  element_name, material_position_field_name);
1982  }
1984  element_name, "DOT_" + velocity_field_name);
1986  element_name, "DOT_" + spatial_position_field_name);
1987 
1988  Range tets;
1989  if (bit.any()) {
1990  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1991  bit, BitRefLevel().set(), MBTET, tets);
1992  }
1993  if (intersected != NULL) {
1994  if (tets.empty()) {
1995  tets = *intersected;
1996  } else {
1997  tets = intersect(*intersected, tets);
1998  }
1999  }
2000 
2001  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
2002  for (; sit != setOfBlocks.end(); sit++) {
2003  Range add_tets = sit->second.tEts;
2004  if (!tets.empty()) {
2005  add_tets = intersect(add_tets, tets);
2006  }
2008  element_name);
2009  }
2010 
2012 }

◆ addHOOpsVol()

MoFEMErrorCode ConvectiveMassElement::addHOOpsVol ( )

Definition at line 108 of file ConvectiveMassElement.hpp.

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

1904  {
1906 
1907  //
1908 
1909  CHKERR mField.add_finite_element(element_name, MF_ZERO);
1911  velocity_field_name);
1913  velocity_field_name);
1915  velocity_field_name);
1917  element_name, spatial_position_field_name);
1919  element_name, spatial_position_field_name);
1920  if (mField.check_field(material_position_field_name)) {
1921  if (ale) {
1923  element_name, material_position_field_name);
1925  element_name, "DOT_" + material_position_field_name);
1926  }
1928  element_name, material_position_field_name);
1929  }
1931  element_name, "DOT_" + velocity_field_name);
1933  element_name, "DOT_" + spatial_position_field_name);
1934 
1935  Range tets;
1936  if (bit.any()) {
1937  CHKERR mField.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
1938  bit, BitRefLevel().set(), MBTET, tets);
1939  }
1940 
1941  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
1942  for (; sit != setOfBlocks.end(); sit++) {
1943  Range add_tets = sit->second.tEts;
1944  if (!tets.empty()) {
1945  add_tets = intersect(add_tets, tets);
1946  }
1948  element_name);
1949  }
1950 
1952 }

◆ getLoopFeEnergy()

MyVolumeFE& ConvectiveMassElement::getLoopFeEnergy ( )

get kinetic energy element

Definition at line 104 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassAuxLhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeMassAuxLhs ( )

get lhs volume element for Kuu shell matrix

Definition at line 89 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassLhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeMassLhs ( )

get lhs volume element

Definition at line 84 of file ConvectiveMassElement.hpp.

◆ getLoopFeMassRhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeMassRhs ( )

get rhs volume element

Definition at line 79 of file ConvectiveMassElement.hpp.

◆ getLoopFeTLhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeTLhs ( )

get lhs volume element

Definition at line 101 of file ConvectiveMassElement.hpp.

◆ getLoopFeTRhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeTRhs ( )

get rhs volume element

Definition at line 99 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelLhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeVelLhs ( )

get lhs volume element

Definition at line 96 of file ConvectiveMassElement.hpp.

◆ getLoopFeVelRhs()

MyVolumeFE& ConvectiveMassElement::getLoopFeVelRhs ( )

get rhs volume element

Definition at line 94 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 594 of file ConvectiveMassElement.hpp.

594  {
596 
597  void *void_ctx;
598  CHKERR MatShellGetContext(A, &void_ctx);
599  MatShellCtx *ctx = (MatShellCtx *)void_ctx;
600  if (!ctx->iNitialized) {
601  CHKERR ctx->iNit();
602  }
603  CHKERR VecZeroEntries(f);
604  // Mult Ku
605  CHKERR VecScatterBegin(ctx->scatterU, x, ctx->u, INSERT_VALUES,
606  SCATTER_FORWARD);
607  CHKERR VecScatterEnd(ctx->scatterU, x, ctx->u, INSERT_VALUES,
608  SCATTER_FORWARD);
609  CHKERR MatMult(ctx->K, ctx->u, ctx->Ku);
610  CHKERR VecScatterBegin(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
611  SCATTER_REVERSE);
612  CHKERR VecScatterEnd(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
613  SCATTER_REVERSE);
614  // Mult Mv
615  CHKERR VecScatterBegin(ctx->scatterV, x, ctx->v, INSERT_VALUES,
616  SCATTER_FORWARD);
617  CHKERR VecScatterEnd(ctx->scatterV, x, ctx->v, INSERT_VALUES,
618  SCATTER_FORWARD);
619  CHKERR MatMult(ctx->M, ctx->v, ctx->Mv);
620  CHKERR VecScatterBegin(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
621  SCATTER_REVERSE);
622  CHKERR VecScatterEnd(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
623  SCATTER_REVERSE);
624  // Velocities
625  CHKERR VecAXPY(ctx->v, -ctx->ts_a, ctx->u);
626  // CHKERR VecScale(ctx->v,ctx->scale);
627  CHKERR VecScatterBegin(ctx->scatterV, ctx->v, f, INSERT_VALUES,
628  SCATTER_REVERSE);
629  CHKERR VecScatterEnd(ctx->scatterV, ctx->v, f, INSERT_VALUES,
630  SCATTER_REVERSE);
631  // Assemble
632  CHKERR VecAssemblyBegin(f);
633  CHKERR VecAssemblyEnd(f);
635  }
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 719 of file ConvectiveMassElement.hpp.

719  {
721 
722  void *void_ctx;
723  CHKERR PCShellGetContext(pc, &void_ctx);
724  PCShellCtx *ctx = (PCShellCtx *)void_ctx;
725  MatShellCtx *shell_mat_ctx;
726  CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
727  // forward
728  CHKERR VecScatterBegin(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
729  INSERT_VALUES, SCATTER_FORWARD);
730  CHKERR VecScatterEnd(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
731  INSERT_VALUES, SCATTER_FORWARD);
732  CHKERR VecScatterBegin(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
733  INSERT_VALUES, SCATTER_FORWARD);
734  CHKERR VecScatterEnd(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
735  INSERT_VALUES, SCATTER_FORWARD);
736  // CHKERR VecScale(shell_mat_ctx->v,1/shell_mat_ctx->scale);
737  // apply pre-conditioner and calculate u
738  CHKERR MatMult(shell_mat_ctx->M, shell_mat_ctx->v,
739  shell_mat_ctx->Mv); // Mrv
740  CHKERR VecAXPY(shell_mat_ctx->Ku, -1, shell_mat_ctx->Mv); // f-Mrv
741  CHKERR PCApply(ctx->pC, shell_mat_ctx->Ku,
742  shell_mat_ctx->u); // u = (aM+K)^(-1)(ru-Mrv)
743  // VecView(shell_mat_ctx->u,PETSC_VIEWER_STDOUT_WORLD);
744  // calculate velocities
745  CHKERR VecAXPY(shell_mat_ctx->v, shell_mat_ctx->ts_a,
746  shell_mat_ctx->u); // v = v + a*u
747  // VecView(shell_mat_ctx->v,PETSC_VIEWER_STDOUT_WORLD);
748  // reverse
749  CHKERR VecZeroEntries(x);
750  CHKERR VecScatterBegin(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
751  INSERT_VALUES, SCATTER_REVERSE);
752  CHKERR VecScatterEnd(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
753  INSERT_VALUES, SCATTER_REVERSE);
754  CHKERR VecScatterBegin(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
755  INSERT_VALUES, SCATTER_REVERSE);
756  CHKERR VecScatterEnd(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
757  INSERT_VALUES, SCATTER_REVERSE);
758  CHKERR VecAssemblyBegin(x);
759  CHKERR VecAssemblyEnd(x);
761  }

◆ PCShellDestroy()

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

Definition at line 681 of file ConvectiveMassElement.hpp.

681  {
683 
684  void *void_ctx;
685  CHKERR PCShellGetContext(pc, &void_ctx);
686  PCShellCtx *ctx = (PCShellCtx *)void_ctx;
687  CHKERR ctx->dEstroy();
689  }

◆ PCShellSetUpOp()

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

Definition at line 666 of file ConvectiveMassElement.hpp.

666  {
668 
669  void *void_ctx;
670  CHKERR PCShellGetContext(pc, &void_ctx);
671  PCShellCtx *ctx = (PCShellCtx *)void_ctx;
672  CHKERR ctx->iNit();
673  MatShellCtx *shell_mat_ctx;
674  CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
675  CHKERR PCSetFromOptions(ctx->pC);
676  CHKERR PCSetOperators(ctx->pC, shell_mat_ctx->barK, shell_mat_ctx->barK);
677  CHKERR PCSetUp(ctx->pC);
679  }

◆ setBlocks() [1/2]

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

Definition at line 1768 of file ConvectiveMassElement.cpp.

1768  {
1770 
1771  Range added_tets;
1773  mField, BLOCKSET | BODYFORCESSET, it)) {
1774  int id = it->getMeshsetId();
1775  EntityHandle meshset = it->getMeshset();
1776  rval = mField.get_moab().get_entities_by_type(meshset, MBTET,
1777  setOfBlocks[id].tEts, true);
1778  CHKERRG(rval);
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  rval = mField.get_moab().get_entities_by_type(meshset, MBTET, tets, true);
1799  CHKERRG(rval);
1800  tets = subtract(tets, added_tets);
1801  if (tets.empty())
1802  continue;
1803  int id = it->getMeshsetId();
1804  setOfBlocks[-id].tEts = tets;
1805  setOfBlocks[-id].rho0 = mydata.data.User1;
1806  setOfBlocks[-id].a0.resize(3);
1807  setOfBlocks[-id].a0[0] = mydata.data.User2;
1808  setOfBlocks[-id].a0[1] = mydata.data.User3;
1809  setOfBlocks[-id].a0[2] = mydata.data.User4;
1810  // std::cerr << setOfBlocks[id].tEts << std::endl;
1811  }
1812 
1814 }
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:476
@ BODYFORCESSET
block name is "BODY_FORCES"
Definition: definitions.h:155
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
Definition: definitions.h:152
@ BLOCKSET
Definition: definitions.h:141
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:85
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 1816 of file ConvectiveMassElement.cpp.

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

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

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

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

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

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

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

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

◆ ZeroEntriesOp()

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

Definition at line 637 of file ConvectiveMassElement.hpp.

637  {
639 
640  void *void_ctx;
641  CHKERR MatShellGetContext(A, &void_ctx);
642  MatShellCtx *ctx = (MatShellCtx *)void_ctx;
643  CHKERR MatZeroEntries(ctx->K);
644  CHKERR MatZeroEntries(ctx->M);
646  }

Member Data Documentation

◆ commonData

CommonData ConvectiveMassElement::commonData

Definition at line 164 of file ConvectiveMassElement.hpp.

◆ feEnergy

MyVolumeFE ConvectiveMassElement::feEnergy

calculate kinetic energy

Definition at line 103 of file ConvectiveMassElement.hpp.

◆ feMassAuxLhs

MyVolumeFE ConvectiveMassElement::feMassAuxLhs

calculate left hand side for tetrahedral elements for Kuu shell matrix

Definition at line 87 of file ConvectiveMassElement.hpp.

◆ feMassLhs

MyVolumeFE ConvectiveMassElement::feMassLhs

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

Definition at line 82 of file ConvectiveMassElement.hpp.

◆ feMassRhs

MyVolumeFE ConvectiveMassElement::feMassRhs

calculate right hand side for tetrahedral elements

Definition at line 78 of file ConvectiveMassElement.hpp.

◆ feTLhs

MyVolumeFE ConvectiveMassElement::feTLhs

calculate left hand side for tetrahedral elements

Definition at line 100 of file ConvectiveMassElement.hpp.

◆ feTRhs

MyVolumeFE ConvectiveMassElement::feTRhs

calculate right hand side for tetrahedral elements

Definition at line 98 of file ConvectiveMassElement.hpp.

◆ feVelLhs

MyVolumeFE ConvectiveMassElement::feVelLhs

calculate left hand side for tetrahedral elements

Definition at line 95 of file ConvectiveMassElement.hpp.

◆ feVelRhs

MyVolumeFE ConvectiveMassElement::feVelRhs

calculate right hand side for tetrahedral elements

Definition at line 93 of file ConvectiveMassElement.hpp.

◆ methodsOp

boost::ptr_vector<MethodForForceScaling> ConvectiveMassElement::methodsOp

Definition at line 166 of file ConvectiveMassElement.hpp.

◆ mField

MoFEM::Interface& ConvectiveMassElement::mField

Definition at line 124 of file ConvectiveMassElement.hpp.

◆ setOfBlocks

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

maps block set id with appropriate BlockData

Definition at line 138 of file ConvectiveMassElement.hpp.

◆ tAg

short int ConvectiveMassElement::tAg

Definition at line 125 of file ConvectiveMassElement.hpp.


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