v0.14.0
NodeForce.cpp
Go to the documentation of this file.
1 /** \file NodeForce.cpp
2  \ingroup mofem_static_boundary_conditions
3 */
4 
5 
6 
7 #include <MoFEM.hpp>
8 using namespace MoFEM;
9 #include <MethodForForceScaling.hpp>
10 #include <SurfacePressure.hpp>
11 #include <NodalForce.hpp>
12 
13 using namespace boost::numeric;
14 
17 
19  const std::string field_name, Vec _F, bCForce &data,
20  boost::ptr_vector<MethodForForceScaling> &methods_op, bool use_snes_f)
23  F(_F), useSnesF(use_snes_f), dAta(data), methodsOp(methods_op) {}
24 
29 
30  if (data.getIndices().size() == 0)
32  EntityHandle ent = getNumeredEntFiniteElementPtr()->getEnt();
33  if (dAta.nOdes.find(ent) == dAta.nOdes.end())
35 
36  const auto &dof_ptr = data.getFieldDofs()[0];
37  int rank = dof_ptr->getNbOfCoeffs();
38 
39  if (data.getIndices().size() != (unsigned int)rank) {
40  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
41  }
42 
43  Nf.resize(3);
44  for (int rr = 0; rr != rank; ++rr) {
45  if (rr == 0) {
46  Nf[0] = dAta.data.data.value3 * dAta.data.data.value1;
47  } else if (rr == 1) {
48  Nf[1] = dAta.data.data.value4 * dAta.data.data.value1;
49  } else if (rr == 2) {
50  Nf[2] = dAta.data.data.value5 * dAta.data.data.value1;
51  } else {
52  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
53  }
54  }
55 
57 
58  Vec myF = F;
59  if (useSnesF || F == PETSC_NULL)
60  myF = getKSPf();
61 
62  CHKERR VecSetValues(myF, data.getIndices().size(), &data.getIndices()[0],
63  &Nf[0], ADD_VALUES);
64 
66 }
67 
69  int ms_id, bool use_snes_f) {
70 
71  const CubitMeshSets *cubit_meshset_ptr;
72  MeshsetsManager *mmanager_ptr;
74  CHKERR mField.getInterface(mmanager_ptr);
75  CHKERR mmanager_ptr->getCubitMeshsetPtr(ms_id, NODESET, &cubit_meshset_ptr);
76  CHKERR cubit_meshset_ptr->getBcDataStructure(mapForce[ms_id].data);
77  CHKERR mField.get_moab().get_entities_by_type(
78  cubit_meshset_ptr->meshset, MBVERTEX, mapForce[ms_id].nOdes, true);
79  fe.getOpPtrVector().push_back(
80  new OpNodalForce(field_name, F, mapForce[ms_id], methodsOp, use_snes_f));
82 }
83 
85  : mField(m_field) {
86 
87  double def_scale = 1.;
88  const EntityHandle root_meshset = mField.get_moab().get_root_set();
89  rval = mField.get_moab().tag_get_handle(
90  "_LoadFactor_Scale_", 1, MB_TYPE_DOUBLE, thScale,
91  MB_TAG_CREAT | MB_TAG_EXCL | MB_TAG_MESH, &def_scale);
92  if (rval == MB_ALREADY_ALLOCATED) {
93  rval = mField.get_moab().tag_get_by_ptr(thScale, &root_meshset, 1,
94  (const void **)&sCale);
96  } else {
98  rval =
99  mField.get_moab().tag_set_data(thScale, &root_meshset, 1, &def_scale);
100  MOAB_THROW(rval);
101  rval = mField.get_moab().tag_get_by_ptr(thScale, &root_meshset, 1,
102  (const void **)&sCale);
103  MOAB_THROW(rval);
104  }
105 }
106 
108  VectorDouble &Nf) {
110  Nf *= *sCale;
112 }
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MetaNodalForces::TagForceScale::sCale
double * sCale
Definition: NodalForce.hpp:68
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
NodalForce::OpNodalForce::OpNodalForce
OpNodalForce(const std::string field_name, Vec _F, bCForce &data, boost::ptr_vector< MethodForForceScaling > &methods_op, bool use_snes_f=false)
Definition: NodeForce.cpp:18
NodalForce::mField
MoFEM::Interface & mField
Definition: NodalForce.hpp:15
EntityHandle
MoFEM::FEMethod
structure for User Loop Methods on finite elements
Definition: LoopMethods.hpp:369
_F
#define _F(n)
Definition: quad.c:25
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::CubitMeshSets
this struct keeps basic methods for moab meshset about material and boundary conditions
Definition: BCMultiIndices.hpp:19
NodalForce.hpp
NodalForce::methodsOp
boost::ptr_vector< MethodForForceScaling > methodsOp
Definition: NodalForce.hpp:31
MetaNodalForces::TagForceScale::thScale
Tag thScale
Definition: NodalForce.hpp:69
MoFEM.hpp
MOAB_THROW
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:554
MoFEM::CubitMeshSets::getBcDataStructure
MoFEMErrorCode getBcDataStructure(CUBIT_BC_DATA_TYPE &data) const
Definition: BCMultiIndices.hpp:296
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
NODESET
@ NODESET
Definition: definitions.h:159
MoFEM::Exceptions::rval
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
NodalForce::bCForce
Definition: NodalForce.hpp:25
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
NodalForce::fe
MyFE fe
Definition: NodalForce.hpp:22
convert.type
type
Definition: convert.py:64
MetaNodalForces::TagForceScale::scaleNf
MoFEMErrorCode scaleNf(const FEMethod *fe, VectorDouble &Nf)
Definition: NodeForce.cpp:107
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: EntitiesFieldData.hpp:1214
MoFEM::CubitMeshSets::meshset
EntityHandle meshset
Definition: BCMultiIndices.hpp:21
MoFEM::EntitiesFieldData::EntData::getFieldDofs
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity
Definition: EntitiesFieldData.hpp:1269
MetaNodalForces::TagForceScale::TagForceScale
TagForceScale(MoFEM::Interface &m_field)
Definition: NodeForce.cpp:84
NodalForce::OpNodalForce
Operator to assemble nodal force into right hand side vector.
Definition: NodalForce.hpp:34
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
NodalForce::mapForce
std::map< int, bCForce > mapForce
Definition: NodalForce.hpp:29
NodalForce::OpNodalForce::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: NodeForce.cpp:26
MoFEM::ForcesAndSourcesCore::getOpPtrVector
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
Definition: ForcesAndSourcesCore.hpp:83
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MetaNodalForces::TagForceScale::mField
MoFEM::Interface & mField
Definition: NodalForce.hpp:67
MoFEM::MeshsetsManager
Interface for managing meshsets containing materials and boundary conditions.
Definition: MeshsetsManager.hpp:104
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
NodalForce::MyFE::MyFE
MyFE(MoFEM::Interface &m_field)
Definition: NodeForce.cpp:15
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MethodForForceScaling::applyScale
static MoFEMErrorCode applyScale(const FEMethod *fe, boost::ptr_vector< MethodForForceScaling > &methods_op, VectorDouble &nf)
Definition: MethodForForceScaling.hpp:21
MoFEM::MeshsetsManager::getCubitMeshsetPtr
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
Definition: MeshsetsManager.cpp:578
NodalForce::addForce
MoFEMErrorCode addForce(const std::string field_name, Vec F, int ms_id, bool use_snes_f=false)
Definition: NodeForce.cpp:68
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::VertexElementForcesAndSourcesCore
Vertex finite element.
Definition: VertexElementForcesAndSourcesCore.hpp:28
SurfacePressure.hpp
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
F
@ F
Definition: free_surface.cpp:394