v0.10.0
Static Public Member Functions | List of all members
MetaNeumannForces Struct Reference

Set of high-level function declaring elements and setting operators to apply forces/fluxes. More...

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

Static Public Member Functions

static MoFEMErrorCode addNeumannBCElements (MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", Range *intersect_ptr=NULL)
 Declare finite element. More...
 
static MoFEMErrorCode setMomentumFluxOperators (MoFEM::Interface &m_field, boost::ptr_map< std::string, NeumannForcesSurface > &neumann_forces, Vec F, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 Set operators to finite elements calculating right hand side vector. More...
 
static MoFEMErrorCode addNeumannFluxBCElements (MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 
static MoFEMErrorCode setMassFluxOperators (MoFEM::Interface &m_field, boost::ptr_map< std::string, NeumannForcesSurface > &neumann_forces, Vec F, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 

Detailed Description

Set of high-level function declaring elements and setting operators to apply forces/fluxes.

Definition at line 794 of file SurfacePressure.hpp.

Member Function Documentation

◆ addNeumannBCElements()

MoFEMErrorCode MetaNeumannForces::addNeumannBCElements ( MoFEM::Interface m_field,
const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS",
Range *  intersect_ptr = NULL 
)
static

Declare finite element.

Search cubit sidesets and blocksets with pressure bc and declare surface elemen

Block set has to have name “PRESSURE”. Can have name “PRESSURE_01” or any other name with prefix. The first attribute of block set is pressure value.

Parameters
m_fieldInterface insurance
field_nameField name (f.e. DISPLACEMENT)
mesh_nodals_positionsName of field on which ho-geometry is defined
intersect_ptrPointer to range to interect meshset entities
Returns
Error code
Examples
elasticity.cpp, elasticity_mixed_formulation.cpp, and simple_contact.cpp.

Definition at line 1272 of file SurfacePressure.cpp.

1274  {
1276 
1277  // Define boundary element that operates on rows, columns and data of a
1278  // given field
1279  CHKERR m_field.add_finite_element("FORCE_FE", MF_ZERO);
1280  CHKERR m_field.modify_finite_element_add_field_row("FORCE_FE", field_name);
1281  CHKERR m_field.modify_finite_element_add_field_col("FORCE_FE", field_name);
1282  CHKERR m_field.modify_finite_element_add_field_data("FORCE_FE", field_name);
1283  if (m_field.check_field(mesh_nodals_positions)) {
1284  CHKERR m_field.modify_finite_element_add_field_data("FORCE_FE",
1285  mesh_nodals_positions);
1286  }
1287  // Add entities to that element, here we add all triangles with FORCESET
1288  // from cubit
1290  it)) {
1291  Range tris;
1292  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1293  true);
1294  if (intersect_ptr)
1295  tris = intersect(tris, *intersect_ptr);
1296  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "FORCE_FE");
1297  }
1298 
1299  CHKERR m_field.add_finite_element("PRESSURE_FE", MF_ZERO);
1300  CHKERR m_field.modify_finite_element_add_field_row("PRESSURE_FE", field_name);
1301  CHKERR m_field.modify_finite_element_add_field_col("PRESSURE_FE", field_name);
1302  CHKERR m_field.modify_finite_element_add_field_data("PRESSURE_FE",
1303  field_name);
1304  if (m_field.check_field(mesh_nodals_positions)) {
1305  CHKERR m_field.modify_finite_element_add_field_data("PRESSURE_FE",
1306  mesh_nodals_positions);
1307  }
1308 
1310  SIDESET | PRESSURESET, it)) {
1311  Range tris;
1312  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1313  true);
1314  if (intersect_ptr)
1315  tris = intersect(tris, *intersect_ptr);
1316  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "PRESSURE_FE");
1317  }
1318 
1319  // Reading forces from BLOCKSET
1320 
1321  const string block_set_force_name("FORCE");
1322  // search for block named FORCE and add its attributes to FORCE_FE element
1324  if (it->getName().compare(0, block_set_force_name.length(),
1325  block_set_force_name) == 0) {
1326  Range tris;
1327  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1328  true);
1329  if (intersect_ptr)
1330  tris = intersect(tris, *intersect_ptr);
1331  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "FORCE_FE");
1332  }
1333  }
1334  // search for block named PRESSURE and add its attributes to PRESSURE_FE
1335  // element
1336  const string block_set_pressure_name("PRESSURE");
1338  if (it->getName().compare(0, block_set_pressure_name.length(),
1339  block_set_pressure_name) == 0) {
1340  Range tris;
1341  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1342  true);
1343  if (intersect_ptr)
1344  tris = intersect(tris, *intersect_ptr);
1345  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "PRESSURE_FE");
1346  }
1347  }
1348 
1349  // search for block named LINEAR_PRESSURE and add its attributes to
1350  // PRESSURE_FE element
1351  const string block_set_linear_pressure_name("LINEAR_PRESSURE");
1353  if (it->getName().compare(0, block_set_linear_pressure_name.length(),
1354  block_set_linear_pressure_name) == 0) {
1355  Range tris;
1356  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1357  true);
1358  if (intersect_ptr)
1359  tris = intersect(tris, *intersect_ptr);
1360  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "PRESSURE_FE");
1361  }
1362  }
1363 
1365 }
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 MoFEMErrorCode add_ents_to_finite_element_by_dim(const EntityHandle entities, const int dim, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:604
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 bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ addNeumannFluxBCElements()

MoFEMErrorCode MetaNeumannForces::addNeumannFluxBCElements ( MoFEM::Interface m_field,
const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)
static

Definition at line 1423 of file SurfacePressure.cpp.

1425  {
1427 
1428  CHKERR m_field.add_finite_element("FLUX_FE", MF_ZERO);
1429  CHKERR m_field.modify_finite_element_add_field_row("FLUX_FE", field_name);
1430  CHKERR m_field.modify_finite_element_add_field_col("FLUX_FE", field_name);
1431  CHKERR m_field.modify_finite_element_add_field_data("FLUX_FE", field_name);
1432  if (m_field.check_field(mesh_nodals_positions)) {
1433  CHKERR m_field.modify_finite_element_add_field_data("FLUX_FE",
1434  mesh_nodals_positions);
1435  }
1436 
1438  SIDESET | PRESSURESET, it)) {
1439  Range tris;
1440  CHKERR m_field.get_moab().get_entities_by_dimension(it->meshset, 2, tris,
1441  true);
1442  CHKERR m_field.add_ents_to_finite_element_by_dim(tris, 2, "FLUX_FE");
1443  }
1444 
1446 }
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 MoFEMErrorCode add_ents_to_finite_element_by_dim(const EntityHandle entities, const int dim, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define CHKERR
Inline error check.
Definition: definitions.h:604
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 bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ setMassFluxOperators()

MoFEMErrorCode MetaNeumannForces::setMassFluxOperators ( MoFEM::Interface m_field,
boost::ptr_map< std::string, NeumannForcesSurface > &  neumann_forces,
Vec  F,
const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)
static

Definition at line 1448 of file SurfacePressure.cpp.

1451  {
1453 
1454  string fe_name;
1455  fe_name = "FLUX_FE";
1456  neumann_forces.insert(fe_name, new NeumannForcesSurface(m_field));
1458  SIDESET | PRESSURESET, it)) {
1459  bool ho_geometry = m_field.check_field(mesh_nodals_positions);
1460  CHKERR neumann_forces.at(fe_name).addFlux(field_name, F, it->getMeshsetId(),
1461  ho_geometry);
1462  }
1464 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:604
virtual bool check_field(const std::string &name) const =0
check if field is in database
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
Finite element and operators to apply force/pressures applied to surfaces.

◆ setMomentumFluxOperators()

MoFEMErrorCode MetaNeumannForces::setMomentumFluxOperators ( MoFEM::Interface m_field,
boost::ptr_map< std::string, NeumannForcesSurface > &  neumann_forces,
Vec  F,
const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)
static

Set operators to finite elements calculating right hand side vector.

Parameters
m_fieldInterface
neumann_forcesMap of pointers to force/pressure elements
FRight hand side vector
field_nameField name (f.e. DISPLACEMENT)
mesh_nodals_positionsName of field on which ho-geometry is defined
Returns
Error code
Examples
elasticity.cpp, elasticity_mixed_formulation.cpp, and simple_contact.cpp.

Definition at line 1367 of file SurfacePressure.cpp.

1370  {
1372 
1373  string fe_name;
1374  fe_name = "FORCE_FE";
1375  neumann_forces.insert(fe_name, new NeumannForcesSurface(m_field));
1376  bool ho_geometry = m_field.check_field(mesh_nodals_positions);
1378  it)) {
1379  CHKERR neumann_forces.at(fe_name).addForce(
1380  field_name, F, it->getMeshsetId(), ho_geometry, false);
1381  }
1382  // Reading forces from BLOCKSET
1383  const string block_set_force_name("FORCE");
1385  if (it->getName().compare(0, block_set_force_name.length(),
1386  block_set_force_name) == 0) {
1387  CHKERR neumann_forces.at(fe_name).addForce(
1388  field_name, F, it->getMeshsetId(), ho_geometry, true);
1389  }
1390  }
1391 
1392  fe_name = "PRESSURE_FE";
1393  neumann_forces.insert(fe_name, new NeumannForcesSurface(m_field));
1395  SIDESET | PRESSURESET, it)) {
1396  CHKERR neumann_forces.at(fe_name).addPressure(
1397  field_name, F, it->getMeshsetId(), ho_geometry, false);
1398  }
1399 
1400  // Reading pressures from BLOCKSET
1401  const string block_set_pressure_name("PRESSURE");
1403  if (it->getName().compare(0, block_set_pressure_name.length(),
1404  block_set_pressure_name) == 0) {
1405  CHKERR neumann_forces.at(fe_name).addPressure(
1406  field_name, F, it->getMeshsetId(), ho_geometry, true);
1407  }
1408  }
1409 
1410  // Reading pressures from BLOCKSET
1411  const string block_set_linear_pressure_name("LINEAR_PRESSURE");
1413  if (it->getName().compare(0, block_set_linear_pressure_name.length(),
1414  block_set_linear_pressure_name) == 0) {
1415  CHKERR neumann_forces.at(fe_name).addLinearPressure(
1416  field_name, F, it->getMeshsetId(), ho_geometry);
1417  }
1418  }
1419 
1421 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:604
virtual bool check_field(const std::string &name) const =0
check if field is in database
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
Finite element and operators to apply force/pressures applied to surfaces.

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