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

Element for simulating viscous fluid flow. More...

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

Classes

struct  BlockData
 
struct  CommonData
 
struct  FaceRule
 
struct  LoadScale
 
struct  OpAssembleLhs
 Base class for operators assembling LHS. More...
 
struct  OpAssembleLhsDiagLin
 Assemble linear (symmetric) part of the diagonal block of the LHS Operator for assembling linear (symmetric) part of the diagonal block of the LHS. More...
 
struct  OpAssembleLhsDiagNonLin
 Assemble non-linear (non-symmetric) part of the diagonal block of the LHS Operator for assembling non-linear (non-symmetric) part of the diagonal block of the LHS. More...
 
struct  OpAssembleLhsOffDiag
 Assemble off-diagonal block of the LHS Operator for assembling off-diagonal block of the LHS. More...
 
struct  OpAssembleRhs
 Base class for operators assembling RHS. More...
 
struct  OpAssembleRhsPressure
 Assemble the pressure component of the RHS vector. More...
 
struct  OpAssembleRhsVelocityLin
 Assemble linear part of the velocity component of the RHS vector. More...
 
struct  OpAssembleRhsVelocityNonLin
 Assemble non-linear part of the velocity component of the RHS vector. More...
 
struct  OpCalcDragForce
 Calculate drag force on the fluid-solid interface. More...
 
struct  OpCalcDragTraction
 Calculate drag traction on the fluid-solid interface. More...
 
struct  OpCalcVolumeFlux
 calculating volumetric flux More...
 
struct  OpPostProcDrag
 Post processing output of drag traction on the fluid-solid interface. More...
 
struct  OpPostProcVorticity
 Post processing output of the vorticity criterion levels. More...
 
struct  VolRule
 Set integration rule to volume elements. More...
 

Public Types

using UserDataOperator = MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator
 
using FaceUserDataOperator = MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator
 
using EntData = EntitiesFieldData::EntData
 

Static Public Member Functions

static MoFEMErrorCode addElement (MoFEM::Interface &m_field, const string element_name, const string velocity_field_name, const string pressure_field_name, const string mesh_field_name, const int dim=3, Range *ents=nullptr)
 Setting up elements. More...
 
static MoFEMErrorCode setNavierStokesOperators (boost::shared_ptr< VolumeElementForcesAndSourcesCore > feRhs, boost::shared_ptr< VolumeElementForcesAndSourcesCore > feLhs, const std::string velocity_field, const std::string pressure_field, boost::shared_ptr< CommonData > common_data, const EntityType type=MBTET)
 Setting up operators for solving Navier-Stokes equations. More...
 
static MoFEMErrorCode setStokesOperators (boost::shared_ptr< VolumeElementForcesAndSourcesCore > feRhs, boost::shared_ptr< VolumeElementForcesAndSourcesCore > feLhs, const std::string velocity_field, const std::string pressure_field, boost::shared_ptr< CommonData > common_data, const EntityType type=MBTET)
 Setting up operators for solving Stokes equations. More...
 
static MoFEMErrorCode setCalcDragOperators (boost::shared_ptr< FaceElementForcesAndSourcesCore > dragFe, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideDragFe, std::string side_fe_name, const std::string velocity_field, const std::string pressure_field, boost::shared_ptr< CommonData > common_data)
 Setting up operators for calculating drag force on the solid surface. More...
 
static MoFEMErrorCode setPostProcDragOperators (boost::shared_ptr< PostProcFaceOnRefinedMesh > postProcDragPtr, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideDragFe, std::string side_fe_name, const std::string velocity_field, const std::string pressure_field, boost::shared_ptr< CommonData > common_data)
 Setting up operators for post processing output of drag traction. More...
 
static MoFEMErrorCode setCalcVolumeFluxOperators (boost::shared_ptr< VolumeElementForcesAndSourcesCore > fe_flux_ptr, const std::string velocity_field, boost::shared_ptr< CommonData > common_data, const EntityType type=MBTET)
 Setting up operators for calculation of volume flux. More...
 

Detailed Description

Element for simulating viscous fluid flow.

Examples
NavierStokesElement.hpp.

Definition at line 37 of file NavierStokesElement.hpp.

Member Typedef Documentation

◆ EntData

Definition at line 43 of file NavierStokesElement.hpp.

◆ FaceUserDataOperator

Examples
NavierStokesElement.hpp.

Definition at line 41 of file NavierStokesElement.hpp.

◆ UserDataOperator

Definition at line 39 of file NavierStokesElement.hpp.

Member Function Documentation

◆ addElement()

static MoFEMErrorCode NavierStokesElement::addElement ( MoFEM::Interface m_field,
const string  element_name,
const string  velocity_field_name,
const string  pressure_field_name,
const string  mesh_field_name,
const int  dim = 3,
Range *  ents = nullptr 
)
static

Setting up elements.

This functions adds element to the database, adds provided fields to rows and columns of the element, provides access of the element to the fields data and adds entities of particular dimension (or a given range of entities to the element)

Parameters
m_fieldMoFEM interface
element_nameName of the element
velocity_field_nameName of the velocity field
pressure_field_nameName of the pressure field
mesh_field_nameName for mesh node positions field
entsRange of entities to be added to element
Returns
Error code
Examples
NavierStokesElement.hpp, and navier_stokes.cpp.

Definition at line 136 of file NavierStokesElement.hpp.

142 {
144
145 CHKERR m_field.add_finite_element(element_name);
146
147 CHKERR m_field.modify_finite_element_add_field_row(element_name,
148 velocity_field_name);
149 CHKERR m_field.modify_finite_element_add_field_col(element_name,
150 velocity_field_name);
151 CHKERR m_field.modify_finite_element_add_field_data(element_name,
152 velocity_field_name);
153
154 CHKERR m_field.modify_finite_element_add_field_row(element_name,
155 pressure_field_name);
156 CHKERR m_field.modify_finite_element_add_field_col(element_name,
157 pressure_field_name);
158 CHKERR m_field.modify_finite_element_add_field_data(element_name,
159 pressure_field_name);
160
161 CHKERR m_field.modify_finite_element_add_field_data(element_name,
162 mesh_field_name);
163
164 if (ents != nullptr) {
166 element_name);
167 } else {
168 CHKERR m_field.add_ents_to_finite_element_by_dim(0, dim, element_name);
169 }
170
172 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#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
const int dim
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

◆ setCalcDragOperators()

MoFEMErrorCode NavierStokesElement::setCalcDragOperators ( boost::shared_ptr< FaceElementForcesAndSourcesCore dragFe,
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide sideDragFe,
std::string  side_fe_name,
const std::string  velocity_field,
const std::string  pressure_field,
boost::shared_ptr< CommonData common_data 
)
static

Setting up operators for calculating drag force on the solid surface.

Pushes operators for caluclating drag force components on the fluid-solid interface

Parameters
dragFepointer to face element instance
sideDragFepointer to volume on side element instance
velocity_fieldname of the velocity field
pressure_fieldname of the pressure field
common_datapointer to common data object
Returns
Error code
Examples
NavierStokesElement.cpp, NavierStokesElement.hpp, and navier_stokes.cpp.

Definition at line 162 of file NavierStokesElement.cpp.

167 {
169
170 auto det_ptr = boost::make_shared<VectorDouble>();
171 auto jac_ptr = boost::make_shared<MatrixDouble>();
172 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
173
174 for (auto &sit : common_data->setOfFacesData) {
175 sideDragFe->getOpPtrVector().push_back(
176 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
177 common_data->gradVelPtr));
178 dragFe->getOpPtrVector().push_back(new OpCalculateHOJac<2>(jac_ptr));
179 dragFe->getOpPtrVector().push_back(
180 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
181 dragFe->getOpPtrVector().push_back(
182 new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
183
184 dragFe->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
185 pressure_field, common_data->pressPtr));
186
187 dragFe->getOpPtrVector().push_back(
189 velocity_field, sideDragFe, side_fe_name, common_data, sit.second));
190 dragFe->getOpPtrVector().push_back(new NavierStokesElement::OpCalcDragForce(
191 velocity_field, common_data, sit.second));
192 }
194};
@ H1
continuous field
Definition: definitions.h:98
Get value at integration points for scalar field.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Calculate drag force on the fluid-solid interface.
Calculate drag traction on the fluid-solid interface.

◆ setCalcVolumeFluxOperators()

MoFEMErrorCode NavierStokesElement::setCalcVolumeFluxOperators ( boost::shared_ptr< VolumeElementForcesAndSourcesCore fe_flux_ptr,
const std::string  velocity_field,
boost::shared_ptr< CommonData common_data,
const EntityType  type = MBTET 
)
static

Setting up operators for calculation of volume flux.

Examples
NavierStokesElement.cpp, and NavierStokesElement.hpp.

Definition at line 136 of file NavierStokesElement.cpp.

139 {
141
142 for (auto &sit : common_data->setOfBlocksData) {
143
144 if (type == MBPRISM) {
145 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
146 fe_flux_ptr->getOpPtrVector().push_back(
148 fe_flux_ptr->getOpPtrVector().push_back(
149 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
150 fe_flux_ptr->getOpPtrVector().push_back(
151 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
152 }
153 fe_flux_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
154 velocity_field, common_data->velPtr));
155 fe_flux_ptr->getOpPtrVector().push_back(
156 new OpCalcVolumeFlux(velocity_field, common_data, sit.second));
157 }
158
160};
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
Calculate inverse of jacobian for face element.
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Operator for fat prism element updating integration weights in the volume.
Transform local reference derivatives of shape functions to global derivatives.

◆ setNavierStokesOperators()

MoFEMErrorCode NavierStokesElement::setNavierStokesOperators ( boost::shared_ptr< VolumeElementForcesAndSourcesCore feRhs,
boost::shared_ptr< VolumeElementForcesAndSourcesCore feLhs,
const std::string  velocity_field,
const std::string  pressure_field,
boost::shared_ptr< CommonData common_data,
const EntityType  type = MBTET 
)
static

Setting up operators for solving Navier-Stokes equations.

Pushes operators for solving Navier-Stokes equations to pipelines of RHS and LHS element instances

Parameters
feRhspointer to RHS element instance
feLhspointer to LHS element instance
velocity_fieldname of the velocity field
pressure_fieldname of the pressure field
common_datapointer to common data object
typetype of entities in the domain
Returns
Error code
Examples
NavierStokesElement.cpp, NavierStokesElement.hpp, and navier_stokes.cpp.

Definition at line 34 of file NavierStokesElement.cpp.

38 {
40
41 for (auto &sit : common_data->setOfBlocksData) {
42
43 if (type == MBPRISM) {
44 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
45 fe_lhs_ptr->getOpPtrVector().push_back(
47 fe_lhs_ptr->getOpPtrVector().push_back(
48 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
49 fe_lhs_ptr->getOpPtrVector().push_back(
50 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
51 fe_rhs_ptr->getOpPtrVector().push_back(
53 fe_rhs_ptr->getOpPtrVector().push_back(
54 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
55 fe_rhs_ptr->getOpPtrVector().push_back(
56 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
57 }
58
59 fe_lhs_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
60 velocity_field, common_data->velPtr));
61 fe_lhs_ptr->getOpPtrVector().push_back(
62 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
63 common_data->gradVelPtr));
64
65 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagLin(
66 velocity_field, velocity_field, common_data, sit.second));
67 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagNonLin(
68 velocity_field, velocity_field, common_data, sit.second));
69 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsOffDiag(
70 velocity_field, pressure_field, common_data, sit.second));
71
72 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
73 velocity_field, common_data->velPtr));
74 fe_rhs_ptr->getOpPtrVector().push_back(
75 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
76 common_data->gradVelPtr));
77 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
78 pressure_field, common_data->pressPtr));
79
80 fe_rhs_ptr->getOpPtrVector().push_back(
81 new OpAssembleRhsVelocityLin(velocity_field, common_data, sit.second));
82 fe_rhs_ptr->getOpPtrVector().push_back(new OpAssembleRhsVelocityNonLin(
83 velocity_field, common_data, sit.second));
84 fe_rhs_ptr->getOpPtrVector().push_back(
85 new OpAssembleRhsPressure(pressure_field, common_data, sit.second));
86 }
87
89};

◆ setPostProcDragOperators()

MoFEMErrorCode NavierStokesElement::setPostProcDragOperators ( boost::shared_ptr< PostProcFaceOnRefinedMesh postProcDragPtr,
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide sideDragFe,
std::string  side_fe_name,
const std::string  velocity_field,
const std::string  pressure_field,
boost::shared_ptr< CommonData common_data 
)
static

Setting up operators for post processing output of drag traction.

Pushes operators for post processing ouput of drag traction components on the fluid-solid interface

Parameters
dragFepointer to face element instance
sideDragFepointer to volume on side element instance
velocity_fieldname of the velocity field
pressure_fieldname of the pressure field
common_datapointer to common data object
Returns
Error code
Examples
NavierStokesElement.cpp, NavierStokesElement.hpp, and navier_stokes.cpp.

Definition at line 196 of file NavierStokesElement.cpp.

201 {
203
204 auto det_ptr = boost::make_shared<VectorDouble>();
205 auto jac_ptr = boost::make_shared<MatrixDouble>();
206 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
207
208 for (auto &sit : common_data->setOfFacesData) {
209 sideDragFe->getOpPtrVector().push_back(
210 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
211 common_data->gradVelPtr));
212 postProcDragPtr->getOpPtrVector().push_back(
213 new OpCalculateHOJac<2>(jac_ptr));
214 postProcDragPtr->getOpPtrVector().push_back(
215 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
216 postProcDragPtr->getOpPtrVector().push_back(
217 new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
218
219 postProcDragPtr->getOpPtrVector().push_back(
220 new OpCalculateVectorFieldValues<3>(velocity_field,
221 common_data->velPtr));
222 postProcDragPtr->getOpPtrVector().push_back(
223 new OpCalculateScalarFieldValues(pressure_field,
224 common_data->pressPtr));
225
226 postProcDragPtr->getOpPtrVector().push_back(
228 velocity_field, sideDragFe, side_fe_name, common_data, sit.second));
229 postProcDragPtr->getOpPtrVector().push_back(
231 velocity_field, postProcDragPtr->postProcMesh,
232 postProcDragPtr->mapGaussPts, common_data, sit.second));
233 }
235};
Post processing output of drag traction on the fluid-solid interface.

◆ setStokesOperators()

MoFEMErrorCode NavierStokesElement::setStokesOperators ( boost::shared_ptr< VolumeElementForcesAndSourcesCore feRhs,
boost::shared_ptr< VolumeElementForcesAndSourcesCore feLhs,
const std::string  velocity_field,
const std::string  pressure_field,
boost::shared_ptr< CommonData common_data,
const EntityType  type = MBTET 
)
static

Setting up operators for solving Stokes equations.

Pushes operators for solving Stokes equations to pipelines of RHS and LHS element instances

Parameters
feRhspointer to RHS element instance
feLhspointer to LHS element instance
velocity_fieldname of the velocity field
pressure_fieldname of the pressure field
common_datapointer to common data object
typetype of entities in the domain
Returns
Error code
Examples
NavierStokesElement.cpp, NavierStokesElement.hpp, and navier_stokes.cpp.

Definition at line 91 of file NavierStokesElement.cpp.

95 {
97
98 for (auto &sit : common_data->setOfBlocksData) {
99
100 if (type == MBPRISM) {
101 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
102 fe_lhs_ptr->getOpPtrVector().push_back(
104 fe_lhs_ptr->getOpPtrVector().push_back(
105 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
106 fe_lhs_ptr->getOpPtrVector().push_back(
107 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
108 fe_rhs_ptr->getOpPtrVector().push_back(
110 fe_rhs_ptr->getOpPtrVector().push_back(
111 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
112 fe_rhs_ptr->getOpPtrVector().push_back(
113 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
114 }
115
116 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagLin(
117 velocity_field, velocity_field, common_data, sit.second));
118 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsOffDiag(
119 velocity_field, pressure_field, common_data, sit.second));
120
121 fe_rhs_ptr->getOpPtrVector().push_back(
122 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
123 common_data->gradVelPtr));
124
125 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
126 pressure_field, common_data->pressPtr));
127 fe_rhs_ptr->getOpPtrVector().push_back(
128 new OpAssembleRhsVelocityLin(velocity_field, common_data, sit.second));
129 fe_rhs_ptr->getOpPtrVector().push_back(
130 new OpAssembleRhsPressure(pressure_field, common_data, sit.second));
131 }
132
134};

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