v0.15.0
Loading...
Searching...
No Matches
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 PostProcVol
 
using PostProcFace
 
using UserDataOperator
 
using FaceUserDataOperator
 
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.
 
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.
 
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.
 
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.
 
static MoFEMErrorCode setPostProcDragOperators (boost::shared_ptr< PostProcFace > 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.
 
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.
 

Detailed Description

Element for simulating viscous fluid flow.

Examples
NavierStokesElement.hpp.

Definition at line 25 of file NavierStokesElement.hpp.

Member Typedef Documentation

◆ EntData

◆ FaceUserDataOperator

◆ PostProcFace

Initial value:
PostProcBrokenMeshInMoab<FaceElementForcesAndSourcesCore>

Definition at line 29 of file NavierStokesElement.hpp.

◆ PostProcVol

Initial value:
PostProcBrokenMeshInMoab<VolumeElementForcesAndSourcesCore>

Definition at line 27 of file NavierStokesElement.hpp.

◆ UserDataOperator

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 )
inlinestatic

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 127 of file NavierStokesElement.hpp.

133 {
135
136 CHKERR m_field.add_finite_element(element_name);
137
138 CHKERR m_field.modify_finite_element_add_field_row(element_name,
139 velocity_field_name);
140 CHKERR m_field.modify_finite_element_add_field_col(element_name,
141 velocity_field_name);
142 CHKERR m_field.modify_finite_element_add_field_data(element_name,
143 velocity_field_name);
144
145 CHKERR m_field.modify_finite_element_add_field_row(element_name,
146 pressure_field_name);
147 CHKERR m_field.modify_finite_element_add_field_col(element_name,
148 pressure_field_name);
149 CHKERR m_field.modify_finite_element_add_field_data(element_name,
150 pressure_field_name);
151
152 CHKERR m_field.modify_finite_element_add_field_data(element_name,
153 mesh_field_name);
154
155 if (ents != nullptr) {
156 CHKERR m_field.add_ents_to_finite_element_by_dim(*ents, dim,
157 element_name);
158 } else {
159 CHKERR m_field.add_ents_to_finite_element_by_dim(0, dim, element_name);
160 }
161
163 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
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_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
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 modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data

◆ 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 150 of file NavierStokesElement.cpp.

155 {
157
158 auto det_ptr = boost::make_shared<VectorDouble>();
159 auto jac_ptr = boost::make_shared<MatrixDouble>();
160 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
161
162 for (auto &sit : common_data->setOfFacesData) {
163 sideDragFe->getOpPtrVector().push_back(
164 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
165 common_data->gradVelPtr));
166 dragFe->getOpPtrVector().push_back(new OpCalculateHOJac<2>(jac_ptr));
167 dragFe->getOpPtrVector().push_back(
168 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
169 dragFe->getOpPtrVector().push_back(
170 new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
171
172 dragFe->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
173 pressure_field, common_data->pressPtr));
174
175 dragFe->getOpPtrVector().push_back(
177 velocity_field, sideDragFe, side_fe_name, common_data, sit.second));
178 dragFe->getOpPtrVector().push_back(new NavierStokesElement::OpCalcDragForce(
179 velocity_field, common_data, sit.second));
180 }
182};
@ H1
continuous field
Definition definitions.h:85
Get value at integration points for scalar field.
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Set inverse jacobian to base functions.
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 124 of file NavierStokesElement.cpp.

127 {
129
130 for (auto &sit : common_data->setOfBlocksData) {
131
132 if (type == MBPRISM) {
133 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
134 fe_flux_ptr->getOpPtrVector().push_back(
136 fe_flux_ptr->getOpPtrVector().push_back(
137 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
138 fe_flux_ptr->getOpPtrVector().push_back(
139 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
140 }
141 fe_flux_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
142 velocity_field, common_data->velPtr));
143 fe_flux_ptr->getOpPtrVector().push_back(
144 new OpCalcVolumeFlux(velocity_field, common_data, sit.second));
145 }
146
148};
Calculate inverse of jacobian for face element.
Get values at integration pts for tensor field 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 22 of file NavierStokesElement.cpp.

26 {
28
29 for (auto &sit : common_data->setOfBlocksData) {
30
31 if (type == MBPRISM) {
32 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
33 fe_lhs_ptr->getOpPtrVector().push_back(
35 fe_lhs_ptr->getOpPtrVector().push_back(
36 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
37 fe_lhs_ptr->getOpPtrVector().push_back(
38 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
39 fe_rhs_ptr->getOpPtrVector().push_back(
41 fe_rhs_ptr->getOpPtrVector().push_back(
42 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
43 fe_rhs_ptr->getOpPtrVector().push_back(
44 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
45 }
46
47 fe_lhs_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
48 velocity_field, common_data->velPtr));
49 fe_lhs_ptr->getOpPtrVector().push_back(
50 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
51 common_data->gradVelPtr));
52
53 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagLin(
54 velocity_field, velocity_field, common_data, sit.second));
55 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagNonLin(
56 velocity_field, velocity_field, common_data, sit.second));
57 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsOffDiag(
58 velocity_field, pressure_field, common_data, sit.second));
59
60 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
61 velocity_field, common_data->velPtr));
62 fe_rhs_ptr->getOpPtrVector().push_back(
63 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
64 common_data->gradVelPtr));
65 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
66 pressure_field, common_data->pressPtr));
67
68 fe_rhs_ptr->getOpPtrVector().push_back(
69 new OpAssembleRhsVelocityLin(velocity_field, common_data, sit.second));
70 fe_rhs_ptr->getOpPtrVector().push_back(new OpAssembleRhsVelocityNonLin(
71 velocity_field, common_data, sit.second));
72 fe_rhs_ptr->getOpPtrVector().push_back(
73 new OpAssembleRhsPressure(pressure_field, common_data, sit.second));
74 }
75
77};

◆ setPostProcDragOperators()

MoFEMErrorCode NavierStokesElement::setPostProcDragOperators ( boost::shared_ptr< PostProcFace > 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 184 of file NavierStokesElement.cpp.

189 {
191
192 auto det_ptr = boost::make_shared<VectorDouble>();
193 auto jac_ptr = boost::make_shared<MatrixDouble>();
194 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
195
196 for (auto &sit : common_data->setOfFacesData) {
197 sideDragFe->getOpPtrVector().push_back(
198 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
199 common_data->gradVelPtr));
200 postProcDragPtr->getOpPtrVector().push_back(
201 new OpCalculateHOJac<2>(jac_ptr));
202 postProcDragPtr->getOpPtrVector().push_back(
203 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
204 postProcDragPtr->getOpPtrVector().push_back(
205 new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
206
207 postProcDragPtr->getOpPtrVector().push_back(
208 new OpCalculateVectorFieldValues<3>(velocity_field,
209 common_data->velPtr));
210 postProcDragPtr->getOpPtrVector().push_back(
211 new OpCalculateScalarFieldValues(pressure_field,
212 common_data->pressPtr));
213
214 postProcDragPtr->getOpPtrVector().push_back(
216 velocity_field, sideDragFe, side_fe_name, common_data, sit.second));
217 postProcDragPtr->getOpPtrVector().push_back(
219 velocity_field, postProcDragPtr->getPostProcMesh(),
220 postProcDragPtr->getMapGaussPts(), common_data, sit.second));
221 }
223};
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 79 of file NavierStokesElement.cpp.

83 {
85
86 for (auto &sit : common_data->setOfBlocksData) {
87
88 if (type == MBPRISM) {
89 boost::shared_ptr<MatrixDouble> inv_jac_ptr(new MatrixDouble);
90 fe_lhs_ptr->getOpPtrVector().push_back(
92 fe_lhs_ptr->getOpPtrVector().push_back(
93 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
94 fe_lhs_ptr->getOpPtrVector().push_back(
95 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
96 fe_rhs_ptr->getOpPtrVector().push_back(
98 fe_rhs_ptr->getOpPtrVector().push_back(
99 new OpCalculateInvJacForFatPrism(inv_jac_ptr));
100 fe_rhs_ptr->getOpPtrVector().push_back(
101 new OpSetInvJacH1ForFatPrism(inv_jac_ptr));
102 }
103
104 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsDiagLin(
105 velocity_field, velocity_field, common_data, sit.second));
106 fe_lhs_ptr->getOpPtrVector().push_back(new OpAssembleLhsOffDiag(
107 velocity_field, pressure_field, common_data, sit.second));
108
109 fe_rhs_ptr->getOpPtrVector().push_back(
110 new OpCalculateVectorFieldGradient<3, 3>(velocity_field,
111 common_data->gradVelPtr));
112
113 fe_rhs_ptr->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
114 pressure_field, common_data->pressPtr));
115 fe_rhs_ptr->getOpPtrVector().push_back(
116 new OpAssembleRhsVelocityLin(velocity_field, common_data, sit.second));
117 fe_rhs_ptr->getOpPtrVector().push_back(
118 new OpAssembleRhsPressure(pressure_field, common_data, sit.second));
119 }
120
122};

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