v0.5.86
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)
 
PetscErrorCode addNeumannFluidPressureBCElements (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 
PetscErrorCode 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.

Definition at line 32 of file FluidPressure.hpp.

Member Typedef Documentation

◆ MeshSetId

Definition at line 53 of file FluidPressure.hpp.

Constructor & Destructor Documentation

◆ FluidPressure()

FluidPressure::FluidPressure ( MoFEM::Interface m_field)

Definition at line 51 of file FluidPressure.hpp.

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

Member Function Documentation

◆ addNeumannFluidPressureBCElements()

PetscErrorCode FluidPressure::addNeumannFluidPressureBCElements ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)
Examples:
elasticity.cpp.

Definition at line 23 of file FluidPressure.cpp.

25  {
26  PetscFunctionBegin;
27 
28  ierr = mField.add_finite_element("FLUID_PRESSURE_FE",MF_ZERO); CHKERRQ(ierr);
29  ierr = mField.modify_finite_element_add_field_row("FLUID_PRESSURE_FE",field_name); CHKERRQ(ierr);
30  ierr = mField.modify_finite_element_add_field_col("FLUID_PRESSURE_FE",field_name); CHKERRQ(ierr);
31  ierr = mField.modify_finite_element_add_field_data("FLUID_PRESSURE_FE",field_name); CHKERRQ(ierr);
32  if(mField.check_field(mesh_nodals_positions)) {
33  ierr = mField.modify_finite_element_add_field_data("FLUID_PRESSURE_FE",mesh_nodals_positions); CHKERRQ(ierr);
34  }
35 
36  //takes skin of block of entities
37  Skinner skin(&mField.get_moab());
38  // loop over all blocksets and get data which name is FluidPressure
40 
41  if(bit->getName().compare(0,14,"FLUID_PRESSURE") == 0) {
42 
43  //get block attributes
44  std::vector<double> attributes;
45  ierr = bit->getAttributes(attributes); CHKERRQ(ierr);
46  if(attributes.size()<7) {
47  SETERRQ1(PETSC_COMM_SELF,1,"not enough block attributes to deffine fluid pressure element, attributes.size() = %d ",attributes.size());
48  }
49  setOfFluids[bit->getMeshsetId()].dEnsity = attributes[0];
50  setOfFluids[bit->getMeshsetId()].aCCeleration.resize(3);
51  setOfFluids[bit->getMeshsetId()].aCCeleration[0] = attributes[1];
52  setOfFluids[bit->getMeshsetId()].aCCeleration[1] = attributes[2];
53  setOfFluids[bit->getMeshsetId()].aCCeleration[2] = attributes[3];
54  setOfFluids[bit->getMeshsetId()].zEroPressure.resize(3);
55  setOfFluids[bit->getMeshsetId()].zEroPressure[0] = attributes[4];
56  setOfFluids[bit->getMeshsetId()].zEroPressure[1] = attributes[5];
57  setOfFluids[bit->getMeshsetId()].zEroPressure[2] = attributes[6];
58  //get blok tetrahedron and triangles
59  Range tets;
60  rval = mField.get_moab().get_entities_by_type(bit->meshset,MBTET,tets,true); CHKERRQ_MOAB(rval);
61  Range tris;
62  rval = mField.get_moab().get_entities_by_type(bit->meshset,MBTRI,setOfFluids[bit->getMeshsetId()].tRis,true); CHKERRQ_MOAB(rval);
63  //this get triangles only on block surfaces
64  Range tets_skin_tris;
65  rval = skin.find_skin(0,tets,false,tets_skin_tris); CHKERRQ_MOAB(rval);
66  setOfFluids[bit->getMeshsetId()].tRis.merge(tets_skin_tris);
67  std::ostringstream ss;
68  ss << setOfFluids[bit->getMeshsetId()] << std::endl;
69  PetscPrintf(mField.get_comm(),ss.str().c_str());
70 
71  ierr = mField.add_ents_to_finite_element_by_type(setOfFluids[bit->getMeshsetId()].tRis,MBTRI,"FLUID_PRESSURE_FE"); CHKERRQ(ierr);
72 
73  }
74 
75  }
76 
77  PetscFunctionReturn(0);
78 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:362
#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...
static MoABErrorCode rval
Definition: Common.hpp:25
std::map< MeshSetId, FluidData > setOfFluids
virtual moab::Interface & get_moab()=0
virtual PetscErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
MoFEM::Interface & mField
CHKERRQ(ierr)
virtual PetscErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL)=0
add finite element
virtual PetscErrorCode 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 bool check_field(const std::string &name) const =0
check if field is in database
virtual PetscErrorCode 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 PetscErrorCode 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 MPI_Comm & get_comm() const =0

◆ getLoopFe()

MyTriangleFE& FluidPressure::getLoopFe ( )
Examples:
elasticity.cpp.

Definition at line 49 of file FluidPressure.hpp.

49 { return fe; }
MyTriangleFE fe

◆ setNeumannFluidPressureFiniteElementOperators()

PetscErrorCode FluidPressure::setNeumannFluidPressureFiniteElementOperators ( string  field_name,
Vec  F,
bool  allow_negative_pressure = true,
bool  ho_geometry = false 
)
Examples:
elasticity.cpp.

Definition at line 81 of file FluidPressure.cpp.

83  {
84  PetscFunctionBegin;
85  std::map<MeshSetId,FluidData>::iterator sit = setOfFluids.begin();
86  for(;sit!=setOfFluids.end();sit++) {
87  //add finite element
88  fe.getOpPtrVector().push_back(new OpCalculatePressure(
89  field_name,F,sit->second,methodsOp,allow_negative_pressure,ho_geometry
90  ));
91  }
92  PetscFunctionReturn(0);
93 }
std::map< MeshSetId, FluidData > setOfFluids
boost::ptr_vector< MethodForForceScaling > methodsOp
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
MyTriangleFE fe

Member Data Documentation

◆ fe

MyTriangleFE FluidPressure::fe

Definition at line 48 of file FluidPressure.hpp.

◆ methodsOp

boost::ptr_vector<MethodForForceScaling> FluidPressure::methodsOp

Definition at line 63 of file FluidPressure.hpp.

◆ mField

MoFEM::Interface& FluidPressure::mField

Definition at line 34 of file FluidPressure.hpp.

◆ setOfFluids

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

Definition at line 61 of file FluidPressure.hpp.


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