v0.14.0
Loading...
Searching...
No Matches
Classes | Public Types | Public Member Functions | Public Attributes | List of all members
FluidPressure Struct Reference

Fluid pressure forces. More...

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

Collaboration diagram for FluidPressure:
[legend]

Classes

struct  FluidData
 
struct  MyTriangleFE
 
struct  OpCalculatePressure
 

Public Types

typedef int MeshSetId
 

Public Member Functions

MyTriangleFEgetLoopFe ()
 
 FluidPressure (MoFEM::Interface &m_field)
 
MoFEMErrorCode addNeumannFluidPressureBCElements (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 
MoFEMErrorCode setNeumannFluidPressureFiniteElementOperators (string field_name, Vec F, bool allow_negative_pressure=true, bool ho_geometry=false)
 

Public Attributes

MoFEM::InterfacemField
 
MyTriangleFE fe
 
std::map< MeshSetId, FluidDatasetOfFluids
 
boost::ptr_vector< MethodForForceScalingmethodsOp
 

Detailed Description

Fluid pressure forces.

Todo:
Implementation for large displacements
Examples
elasticity.cpp, and nonlinear_dynamics.cpp.

Definition at line 20 of file FluidPressure.hpp.

Member Typedef Documentation

◆ MeshSetId

Definition at line 39 of file FluidPressure.hpp.

Constructor & Destructor Documentation

◆ FluidPressure()

FluidPressure::FluidPressure ( MoFEM::Interface m_field)
inline

Definition at line 37 of file FluidPressure.hpp.

37: mField(m_field), fe(mField) {}
MyTriangleFE fe
MoFEM::Interface & mField

Member Function Documentation

◆ addNeumannFluidPressureBCElements()

MoFEMErrorCode FluidPressure::addNeumannFluidPressureBCElements ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

Definition at line 67 of file FluidPressure.cpp.

68 {
70
71 CHKERR mField.add_finite_element("FLUID_PRESSURE_FE", MF_ZERO);
78 if (mField.check_field(mesh_nodals_positions)) {
80 mesh_nodals_positions);
81 }
82
83 // takes skin of block of entities
84 Skinner skin(&mField.get_moab());
85 // loop over all blocksets and get data which name is FluidPressure
87
88 if (bit->getName().compare(0, 14, "FLUID_PRESSURE") == 0) {
89
90 // get block attributes
91 std::vector<double> attributes;
92 CHKERR bit->getAttributes(attributes);
93 if (attributes.size() < 7) {
94 SETERRQ1(PETSC_COMM_SELF, 1,
95 "not enough block attributes to deffine fluid pressure "
96 "element, attributes.size() = %d ",
97 attributes.size());
98 }
99 setOfFluids[bit->getMeshsetId()].dEnsity = attributes[0];
100 setOfFluids[bit->getMeshsetId()].aCCeleration.resize(3);
101 setOfFluids[bit->getMeshsetId()].aCCeleration[0] = attributes[1];
102 setOfFluids[bit->getMeshsetId()].aCCeleration[1] = attributes[2];
103 setOfFluids[bit->getMeshsetId()].aCCeleration[2] = attributes[3];
104 setOfFluids[bit->getMeshsetId()].zEroPressure.resize(3);
105 setOfFluids[bit->getMeshsetId()].zEroPressure[0] = attributes[4];
106 setOfFluids[bit->getMeshsetId()].zEroPressure[1] = attributes[5];
107 setOfFluids[bit->getMeshsetId()].zEroPressure[2] = attributes[6];
108 // get blok tetrahedron and triangles
109 Range tets;
110 CHKERR mField.get_moab().get_entities_by_type(bit->meshset, MBTET, tets,
111 true);
112 Range tris;
113 CHKERR mField.get_moab().get_entities_by_type(
114 bit->meshset, MBTRI, setOfFluids[bit->getMeshsetId()].tRis, true);
115 // this get triangles only on block surfaces
116 Range tets_skin_tris;
117 CHKERR skin.find_skin(0, tets, false, tets_skin_tris);
118 setOfFluids[bit->getMeshsetId()].tRis.merge(tets_skin_tris);
119 std::ostringstream ss;
120 ss << setOfFluids[bit->getMeshsetId()] << std::endl;
121 PetscPrintf(mField.get_comm(), ss.str().c_str());
122
124 setOfFluids[bit->getMeshsetId()].tRis, MBTRI, "FLUID_PRESSURE_FE");
125 }
126 }
127
129}
@ MF_ZERO
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ BLOCKSET
Definition: definitions.h:148
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
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_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_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_row(const std::string &fe_name, const std::string name_row)=0
set field row 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
std::map< MeshSetId, FluidData > setOfFluids
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0

◆ getLoopFe()

MyTriangleFE & FluidPressure::getLoopFe ( )
inline

Definition at line 35 of file FluidPressure.hpp.

35{ return fe; }

◆ setNeumannFluidPressureFiniteElementOperators()

MoFEMErrorCode FluidPressure::setNeumannFluidPressureFiniteElementOperators ( string  field_name,
Vec  F,
bool  allow_negative_pressure = true,
bool  ho_geometry = false 
)

Definition at line 131 of file FluidPressure.cpp.

132 {
134 std::map<MeshSetId, FluidData>::iterator sit = setOfFluids.begin();
135 for (; sit != setOfFluids.end(); sit++) {
136 // add finite element
137 fe.getOpPtrVector().push_back(
138 new OpCalculatePressure(field_name, F, sit->second, methodsOp,
139 allow_negative_pressure, ho_geometry));
140 }
142}
@ F
boost::ptr_vector< MethodForForceScaling > methodsOp
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.

Member Data Documentation

◆ fe

MyTriangleFE FluidPressure::fe

Definition at line 34 of file FluidPressure.hpp.

◆ methodsOp

boost::ptr_vector<MethodForForceScaling> FluidPressure::methodsOp

Definition at line 51 of file FluidPressure.hpp.

◆ mField

MoFEM::Interface& FluidPressure::mField

Definition at line 22 of file FluidPressure.hpp.

◆ setOfFluids

std::map<MeshSetId, FluidData> FluidPressure::setOfFluids

Definition at line 49 of file FluidPressure.hpp.


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