v0.10.0
FluidPressure.hpp
Go to the documentation of this file.
1 /* \file FluidPressure.hpp
2  *
3  * \brief Implementation of fluid pressure element
4  *
5  * \todo Implement nonlinear case (consrvative force, i.e. normal follows
6  * surface normal)
7  *
8  */
9 
10 /* This file is part of MoFEM.
11  * MoFEM is free software: you can redistribute it and/or modify it under
12  * the terms of the GNU Lesser General Public License as published by the
13  * Free Software Foundation, either version 3 of the License, or (at your
14  * option) any later version.
15  *
16  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19  * License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23 
24 #ifndef __FLUID_PRESSURE_HPP
25 #define __FLUID_PRESSURE_HPP
26 
27 /** \brief Fluid pressure forces
28 
29 \todo Implementation for large displacements
30 
31 */
32 struct FluidPressure {
33 
36 
39  int getRule(int order) { return order + 1; };
40 
44  }
45  };
47  MyTriangleFE &getLoopFe() { return fe; }
48 
49  FluidPressure(MoFEM::Interface &m_field) : mField(m_field), fe(mField) {}
50 
51  typedef int MeshSetId;
52  struct FluidData {
53  double dEnsity; ///< fluid density [kg/m^2] or any consistent unit
54  VectorDouble aCCeleration; ///< acceleration [m/s^2]
55  VectorDouble zEroPressure; ///< fluid level of reference zero pressure.
56  Range
57  tRis; ///< range of surface element to which fluid pressure is applied
58  friend std::ostream &operator<<(std::ostream &os,
59  const FluidPressure::FluidData &e);
60  };
61  std::map<MeshSetId, FluidData> setOfFluids;
62 
63  boost::ptr_vector<MethodForForceScaling> methodsOp;
64 
67 
68  Vec F;
70  boost::ptr_vector<MethodForForceScaling> &methodsOp;
71  bool allowNegativePressure; ///< allows for negative pressures
72  bool hoGeometry;
73 
74  OpCalculatePressure(const std::string field_name, Vec _F, FluidData &data,
75  boost::ptr_vector<MethodForForceScaling> &methods_op,
76  bool allow_negative_pressure, bool ho_geometry)
78  field_name, UserDataOperator::OPROW),
79  F(_F), dAta(data), methodsOp(methods_op),
80  allowNegativePressure(allow_negative_pressure),
81  hoGeometry(ho_geometry) {}
82 
84 
85  MoFEMErrorCode doWork(int side, EntityType type,
87 
88  };
89 
91  const std::string field_name,
92  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
93 
95  string field_name, Vec F, bool allow_negative_pressure = true,
96  bool ho_geometry = false);
97 };
98 
99 std::ostream &operator<<(std::ostream &os, const FluidPressure::FluidData &e);
100 
101 #endif //__FLUID_PRESSSURE_HPP
Fluid pressure forces.
OpCalculatePressure(const std::string field_name, Vec _F, FluidData &data, boost::ptr_vector< MethodForForceScaling > &methods_op, bool allow_negative_pressure, bool ho_geometry)
Deprecated interface functions.
FluidPressure(MoFEM::Interface &m_field)
MoFEMErrorCode setNeumannFluidPressureFiniteElementOperators(string field_name, Vec F, bool allow_negative_pressure=true, bool ho_geometry=false)
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
double dEnsity
fluid density [kg/m^2] or any consistent unit
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MyTriangleFE & getLoopFe()
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
bool allowNegativePressure
allows for negative pressures
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
MoFEM::Interface & mField
boost::ptr_vector< MethodForForceScaling > methodsOp
MoFEMErrorCode addNeumannFluidPressureBCElements(const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
DataForcesAndSourcesCore::EntData EntData
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
Range tRis
range of surface element to which fluid pressure is applied
friend std::ostream & operator<<(std::ostream &os, const FluidPressure::FluidData &e)
constexpr int order
boost::ptr_vector< MethodForForceScaling > & methodsOp
std::ostream & operator<<(std::ostream &os, const FluidPressure::FluidData &e)
VectorDouble zEroPressure
fluid level of reference zero pressure.
std::map< MeshSetId, FluidData > setOfFluids
MyTriangleFE fe
MyTriangleFE(MoFEM::Interface &m_field)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
VectorDouble aCCeleration
acceleration [m/s^2]
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...
#define _F(n)
Definition: quad.c:25