v0.13.1
Classes | Static Public Member Functions | List of all members
MetaSpringBC Struct Reference

Set of functions declaring elements and setting operators to apply spring boundary condition. More...

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

Classes

struct  BlockOptionDataSprings
 
struct  DataAtIntegrationPtsSprings
 
struct  OpCalculateDeformation
 Operator for computing deformation gradients in side volumes. More...
 
struct  OpGetNormalSpEle
 Computes, for material configuration, normal to face that lies on a surface with springs. More...
 
struct  OpGetTangentSpEle
 Computes, for material configuration, tangent vectors to face that lie on a surface with springs. More...
 
struct  OpSpringALEMaterialLhs
 Base class for LHS-operators for pressure element (material configuration) More...
 
struct  OpSpringALEMaterialLhs_dX_dx
 LHS-operator for the pressure element (material configuration) More...
 
struct  OpSpringALEMaterialLhs_dX_dX
 LHS-operator for the pressure element (material configuration) More...
 
struct  OpSpringFs
 RHS-operator for the spring boundary condition element. More...
 
struct  OpSpringFsMaterial
 RHS-operator for the spring boundary condition element for ALE formulation. More...
 
struct  OpSpringKs
 LHS-operator for the springs element. More...
 
struct  OpSpringKs_dX
 
struct  SpringALEMaterialVolOnSideLhs
 Base class for LHS-operators (material) on side volumes. More...
 
struct  SpringALEMaterialVolOnSideLhs_dX_dx
 LHS-operator (material configuration) on the side volume for spring element. More...
 
struct  SpringALEMaterialVolOnSideLhs_dX_dX
 LHS-operator (material configuration) on the side volume. More...
 

Static Public Member Functions

static MoFEMErrorCode addSpringElements (MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 Declare spring element. More...
 
static MoFEMErrorCode addSpringElementsALE (MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions, Range &spring_triangles)
 Declare spring element in ALE. More...
 
static MoFEMErrorCode setSpringOperators (MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", double stiffness_scale=1.)
 Implementation of spring element. Set operators to calculate LHS and RHS. More...
 
static MoFEMErrorCode setSpringOperatorsMaterial (MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr_dx, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr_dX, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_integration_pts, const std::string field_name, const std::string mesh_nodals_positions, std::string side_fe_name)
 Implementation of spring element. Set operators to calculate LHS and RHS. More...
 

Detailed Description

Set of functions declaring elements and setting operators to apply spring boundary condition.

Definition at line 31 of file SpringElement.hpp.

Member Function Documentation

◆ addSpringElements()

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

Declare spring element.

Search cubit sidesets and blocksets with spring bc and declare surface element

Blockset has to have name “SPRING_BC”. The first two attributes of the blockset are spring stiffness value.

Parameters
m_fieldInterface instance
field_nameField name (e.g. SPATIAL_POSITION)
mesh_nodals_positionsName of field on which ho-geometry is defined
Returns
Error code
Examples
elasticity.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, simple_contact.cpp, and simple_contact_thermal.cpp.

Definition at line 1170 of file SpringElement.cpp.

1172 {
1174
1175 // Define boundary element that operates on rows, columns and data of a
1176 // given field
1177 CHKERR m_field.add_finite_element("SPRING", MF_ZERO);
1181 if (m_field.check_field(mesh_nodals_positions)) {
1183 mesh_nodals_positions);
1184 }
1185 // Add entities to that element, here we add all triangles with SPRING_BC
1186 // from cubit
1188 if (bit->getName().compare(0, 9, "SPRING_BC") == 0) {
1189 CHKERR m_field.add_ents_to_finite_element_by_dim(bit->getMeshset(), 2,
1190 "SPRING");
1191 }
1192 }
1193
1195}
@ MF_ZERO
Definition: definitions.h:111
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ BLOCKSET
Definition: definitions.h:161
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
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_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 MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add 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
#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.
auto bit
set bit
constexpr auto field_name

◆ addSpringElementsALE()

MoFEMErrorCode MetaSpringBC::addSpringElementsALE ( MoFEM::Interface m_field,
const std::string  field_name,
const std::string  mesh_nodals_positions,
Range &  spring_triangles 
)
static

Declare spring element in ALE.

Search cubit sidesets and blocksets with spring bc and declare surface element

Blockset has to have name “SPRING_BC”. The first two attributes of the blockset are spring stiffness value.

Parameters
m_fieldInterface instance
field_nameField name (e.g. SPATIAL_POSITION)
mesh_nodals_positionsName of field on which ho-geometry is defined
Returns
Error code

Definition at line 1197 of file SpringElement.cpp.

1199 {
1201
1202 // Define boundary element that operates on rows, columns and data of a
1203 // given field
1204 CHKERR m_field.add_finite_element("SPRING_ALE", MF_ZERO);
1208 CHKERR m_field.modify_finite_element_add_field_row("SPRING_ALE",
1209 mesh_nodals_positions);
1210 CHKERR m_field.modify_finite_element_add_field_col("SPRING_ALE",
1211 mesh_nodals_positions);
1212 CHKERR m_field.modify_finite_element_add_field_data("SPRING_ALE",
1213 mesh_nodals_positions);
1214
1215 CHKERR m_field.add_ents_to_finite_element_by_dim(spring_triangles, 2,
1216 "SPRING_ALE");
1217
1219}

◆ setSpringOperators()

MoFEMErrorCode MetaSpringBC::setSpringOperators ( MoFEM::Interface m_field,
boost::shared_ptr< FaceElementForcesAndSourcesCore fe_spring_lhs_ptr,
boost::shared_ptr< FaceElementForcesAndSourcesCore fe_spring_rhs_ptr,
const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS",
double  stiffness_scale = 1. 
)
static

Implementation of spring element. Set operators to calculate LHS and RHS.

Parameters
m_fieldInterface instance
fe_spring_lhs_ptrPointer to the FE instance for LHS
fe_spring_rhs_ptrPointer to the FE instance for RHS
field_nameField name (e.g. SPATIAL_POSITION)
mesh_nodals_positionsName of field on which ho-geometry is defined
Returns
Error code
Examples
elasticity.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, simple_contact.cpp, and simple_contact_thermal.cpp.

Definition at line 1221 of file SpringElement.cpp.

1225 {
1227
1228 // Push operators to instances for springs
1229 // loop over blocks
1230 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> commonDataPtr =
1231 boost::make_shared<MetaSpringBC::DataAtIntegrationPtsSprings>(m_field, stiffness_scale);
1232 CHKERR commonDataPtr->getParameters();
1233
1234 if (m_field.check_field(mesh_nodals_positions)) {
1235 CHKERR addHOOpsFace3D(mesh_nodals_positions, *fe_spring_lhs_ptr, false,
1236 false);
1237 CHKERR addHOOpsFace3D(mesh_nodals_positions, *fe_spring_rhs_ptr, false,
1238 false);
1239 }
1240
1241 for (auto &sitSpring : commonDataPtr->mapSpring) {
1242
1243 fe_spring_lhs_ptr->getOpPtrVector().push_back(
1244 new OpGetTangentSpEle(mesh_nodals_positions, commonDataPtr));
1245 fe_spring_lhs_ptr->getOpPtrVector().push_back(
1246 new OpGetNormalSpEle(mesh_nodals_positions, commonDataPtr));
1247 fe_spring_lhs_ptr->getOpPtrVector().push_back(
1248 new OpSpringKs(commonDataPtr, sitSpring.second, field_name));
1249
1250 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1251 new OpCalculateVectorFieldValues<3>(field_name, commonDataPtr->xAtPts));
1252 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1253 new OpCalculateVectorFieldValues<3>(mesh_nodals_positions,
1254 commonDataPtr->xInitAtPts));
1255 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1256 new OpGetTangentSpEle(mesh_nodals_positions, commonDataPtr));
1257 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1258 new OpGetNormalSpEle(mesh_nodals_positions, commonDataPtr));
1259 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1260 new OpSpringFs(commonDataPtr, sitSpring.second, field_name));
1261 }
1262
1264}
MoFEMErrorCode addHOOpsFace3D(const std::string field, E &e, bool hcurl, bool hdiv)
Get values at integration pts for tensor filed rank 1, i.e. vector field.

◆ setSpringOperatorsMaterial()

MoFEMErrorCode MetaSpringBC::setSpringOperatorsMaterial ( MoFEM::Interface m_field,
boost::shared_ptr< FaceElementForcesAndSourcesCore fe_spring_lhs_ptr_dx,
boost::shared_ptr< FaceElementForcesAndSourcesCore fe_spring_lhs_ptr_dX,
boost::shared_ptr< FaceElementForcesAndSourcesCore fe_spring_rhs_ptr,
boost::shared_ptr< DataAtIntegrationPtsSprings data_at_integration_pts,
const std::string  field_name,
const std::string  mesh_nodals_positions,
std::string  side_fe_name 
)
static

Implementation of spring element. Set operators to calculate LHS and RHS.

Parameters
m_fieldInterface instance
fe_spring_lhs_ptrPointer to the FE instance for LHS
fe_spring_rhs_ptrPointer to the FE instance for RHS
field_nameField name (e.g. SPATIAL_POSITION)
mesh_nodals_positionsName of field on which ho-geometry is defined
Returns
Error code

Definition at line 1266 of file SpringElement.cpp.

1274 {
1276
1277 // Push operators to instances for springs
1278 // loop over blocks
1279 CHKERR data_at_integration_pts->getParameters();
1280
1281 if (m_field.check_field(mesh_nodals_positions)) {
1282 CHKERR addHOOpsFace3D(mesh_nodals_positions, *fe_spring_lhs_ptr_dx, false,
1283 false);
1284 CHKERR addHOOpsFace3D(mesh_nodals_positions, *fe_spring_lhs_ptr_dX, false,
1285 false);
1286 CHKERR addHOOpsFace3D(mesh_nodals_positions, *fe_spring_rhs_ptr, false,
1287 false);
1288 }
1289
1290 for (auto &sitSpring : data_at_integration_pts->mapSpring) {
1291
1292 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> feMatSideRhs =
1293 boost::make_shared<VolumeElementForcesAndSourcesCoreOnSide>(m_field);
1294
1295 feMatSideRhs->getOpPtrVector().push_back(
1297 mesh_nodals_positions, data_at_integration_pts->HMat));
1298 feMatSideRhs->getOpPtrVector().push_back(
1300 field_name, data_at_integration_pts->hMat));
1301
1302 feMatSideRhs->getOpPtrVector().push_back(new OpCalculateDeformation(
1303 mesh_nodals_positions, data_at_integration_pts));
1304
1305 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1306 new OpGetTangentSpEle(mesh_nodals_positions, data_at_integration_pts));
1307 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1308 new OpGetNormalSpEle(mesh_nodals_positions, data_at_integration_pts));
1309
1310 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1312 data_at_integration_pts->xAtPts));
1313 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1315 mesh_nodals_positions, data_at_integration_pts->xInitAtPts));
1316
1317 fe_spring_rhs_ptr->getOpPtrVector().push_back(
1318 new OpSpringFsMaterial(mesh_nodals_positions, data_at_integration_pts,
1319 feMatSideRhs, side_fe_name, sitSpring.second));
1320
1321 fe_spring_lhs_ptr_dx->getOpPtrVector().push_back(
1323 data_at_integration_pts->xAtPts));
1324 fe_spring_lhs_ptr_dx->getOpPtrVector().push_back(
1326 mesh_nodals_positions, data_at_integration_pts->xInitAtPts));
1327 fe_spring_lhs_ptr_dx->getOpPtrVector().push_back(
1328 new OpGetTangentSpEle(mesh_nodals_positions, data_at_integration_pts));
1329 fe_spring_lhs_ptr_dx->getOpPtrVector().push_back(
1330 new OpGetNormalSpEle(mesh_nodals_positions, data_at_integration_pts));
1331
1332 fe_spring_lhs_ptr_dx->getOpPtrVector().push_back(
1333 new OpSpringKs_dX(data_at_integration_pts, sitSpring.second, field_name,
1334 mesh_nodals_positions));
1335
1336 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> feMatSideLhs_dx =
1337 boost::make_shared<VolumeElementForcesAndSourcesCoreOnSide>(m_field);
1338
1339 feMatSideLhs_dx->getOpPtrVector().push_back(
1341 mesh_nodals_positions, data_at_integration_pts->HMat));
1342 feMatSideLhs_dx->getOpPtrVector().push_back(
1344 field_name, data_at_integration_pts->hMat));
1345
1346 feMatSideLhs_dx->getOpPtrVector().push_back(new OpCalculateDeformation(
1347 mesh_nodals_positions, data_at_integration_pts));
1348
1349 feMatSideLhs_dx->getOpPtrVector().push_back(
1350 new SpringALEMaterialVolOnSideLhs_dX_dx(
1351 mesh_nodals_positions, field_name, data_at_integration_pts));
1352
1353 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1355 data_at_integration_pts->xAtPts));
1356 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1358 mesh_nodals_positions, data_at_integration_pts->xInitAtPts));
1359 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1360 new OpGetTangentSpEle(mesh_nodals_positions, data_at_integration_pts));
1361 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1362 new OpGetNormalSpEle(mesh_nodals_positions, data_at_integration_pts));
1363
1364 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1365 new OpSpringALEMaterialLhs_dX_dx(mesh_nodals_positions, field_name,
1366 data_at_integration_pts,
1367 feMatSideLhs_dx, side_fe_name));
1368
1369 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> feMatSideLhs_dX =
1370 boost::make_shared<VolumeElementForcesAndSourcesCoreOnSide>(m_field);
1371
1372 feMatSideLhs_dX->getOpPtrVector().push_back(
1374 mesh_nodals_positions, data_at_integration_pts->HMat));
1375 feMatSideLhs_dX->getOpPtrVector().push_back(
1377 field_name, data_at_integration_pts->hMat));
1378
1379 feMatSideLhs_dX->getOpPtrVector().push_back(new OpCalculateDeformation(
1380 mesh_nodals_positions, data_at_integration_pts));
1381
1382 feMatSideLhs_dX->getOpPtrVector().push_back(
1383 new SpringALEMaterialVolOnSideLhs_dX_dX(mesh_nodals_positions,
1384 mesh_nodals_positions,
1385 data_at_integration_pts));
1386
1387 fe_spring_lhs_ptr_dX->getOpPtrVector().push_back(
1388 new OpSpringALEMaterialLhs_dX_dX(
1389 mesh_nodals_positions, mesh_nodals_positions,
1390 data_at_integration_pts, feMatSideLhs_dX, side_fe_name));
1391 }
1392
1394}
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.

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