v0.13.1
NavierStokesElement.hpp
Go to the documentation of this file.
1/**
2 * \file NavierStokesElement.hpp
3 * \example NavierStokesElement.hpp
4 *
5 * \brief Implementation of operators for fluid flow
6 *
7 * Implementation of operators for computations of the fluid flow, governed by
8 * for Stokes and Navier-Stokes equations, and computation of the drag force
9 */
10
11/* This file is part of MoFEM.
12 * MoFEM is free software: you can redistribute it and/or modify it under
13 * the terms of the GNU Lesser General Public License as published by the
14 * Free Software Foundation, either version 3 of the License, or (at your
15 * option) any later version.
16 *
17 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
18 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
19 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
20 * License for more details.
21 *
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
24
25#ifndef __NAVIERSTOKESELEMENT_HPP__
26#define __NAVIERSTOKESELEMENT_HPP__
27
28#ifndef __BASICFINITEELEMENTS_HPP__
30#endif // __BASICFINITEELEMENTS_HPP__
31
32using namespace boost::numeric;
33
34/**
35 * @brief Element for simulating viscous fluid flow
36*/
38
44
45 struct BlockData {
46 int iD;
51 Range eNts;
53 : iD(-1), fluidViscosity(-1), fluidDensity(-1), inertiaCoef(-1),
54 viscousCoef(-1) {}
55 };
56
57 struct CommonData {
58
59 boost::shared_ptr<MatrixDouble> gradVelPtr;
60 boost::shared_ptr<MatrixDouble> velPtr;
61 boost::shared_ptr<VectorDouble> pressPtr;
62
63 boost::shared_ptr<MatrixDouble> pressureDragTract;
64 boost::shared_ptr<MatrixDouble> shearDragTract;
65 boost::shared_ptr<MatrixDouble> totalDragTract;
66
70
72
73 std::map<int, BlockData> setOfBlocksData;
74 std::map<int, BlockData> setOfFacesData;
75
77
78 gradVelPtr = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
79 velPtr = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
80 pressPtr = boost::shared_ptr<VectorDouble>(new VectorDouble());
81
82 pressureDragTract = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
83 shearDragTract = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
84 totalDragTract = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
85
86 int vec_size;
87 if (m_field.get_comm_rank() == 0)
88 vec_size = 3;
89 else
90 vec_size = 0;
91
92 pressureDragForceVec = createSmartVectorMPI(m_field.get_comm(), vec_size, 3);
93 shearDragForceVec = createSmartVectorMPI(m_field.get_comm(), vec_size, 3);
94 totalDragForceVec = createSmartVectorMPI(m_field.get_comm(), vec_size, 3);
95
96 volumeFluxVec = createSmartVectorMPI(m_field.get_comm(), vec_size, 3);
97 }
98
101 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Problem", "none");
102
103 ierr = PetscOptionsEnd();
104 CHKERRQ(ierr);
106 }
107 };
108
110
111 static double lambda;
112
115 nf *= lambda;
117 }
118 };
119
120 /**
121 * @brief Setting up elements
122 *
123 * This functions adds element to the database, adds provided fields to rows
124 * and columns of the element, provides access of the element to the fields
125 * data and adds entities of particular dimension (or a given range of
126 * entities to the element)
127 *
128 * @param m_field MoFEM interface
129 * @param element_name Name of the element
130 * @param velocity_field_name Name of the velocity field
131 * @param pressure_field_name Name of the pressure field
132 * @param mesh_field_name Name for mesh node positions field
133 * @param ents Range of entities to be added to element
134 * @return Error code
135 */
137 const string element_name,
138 const string velocity_field_name,
139 const string pressure_field_name,
140 const string mesh_field_name,
141 const int dim = 3,
142 Range *ents = nullptr) {
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 }
173
174 /**
175 * @brief Setting up operators for solving Navier-Stokes equations
176 *
177 * Pushes operators for solving Navier-Stokes equations to pipelines of RHS
178 * and LHS element instances
179 *
180 * @param feRhs pointer to RHS element instance
181 * @param feLhs pointer to LHS element instance
182 * @param velocity_field name of the velocity field
183 * @param pressure_field name of the pressure field
184 * @param common_data pointer to common data object
185 * @param type type of entities in the domain
186 * @return Error code
187 */
189 boost::shared_ptr<VolumeElementForcesAndSourcesCore> feRhs,
190 boost::shared_ptr<VolumeElementForcesAndSourcesCore> feLhs,
191 const std::string velocity_field, const std::string pressure_field,
192 boost::shared_ptr<CommonData> common_data, const EntityType type = MBTET);
193
194 /**
195 * @brief Setting up operators for solving Stokes equations
196 *
197 * Pushes operators for solving Stokes equations to pipelines of RHS
198 * and LHS element instances
199 *
200 * @param feRhs pointer to RHS element instance
201 * @param feLhs pointer to LHS element instance
202 * @param velocity_field name of the velocity field
203 * @param pressure_field name of the pressure field
204 * @param common_data pointer to common data object
205 * @param type type of entities in the domain
206 * @return Error code
207 */
209 boost::shared_ptr<VolumeElementForcesAndSourcesCore> feRhs,
210 boost::shared_ptr<VolumeElementForcesAndSourcesCore> feLhs,
211 const std::string velocity_field, const std::string pressure_field,
212 boost::shared_ptr<CommonData> common_data, const EntityType type = MBTET);
213
214 /**
215 * @brief Setting up operators for calculating drag force on the solid surface
216 *
217 * Pushes operators for caluclating drag force components on the fluid-solid
218 * interface
219 *
220 * @param dragFe pointer to face element instance
221 * @param sideDragFe pointer to volume on side element instance
222 * @param velocity_field name of the velocity field
223 * @param pressure_field name of the pressure field
224 * @param common_data pointer to common data object
225 * @return Error code
226 */
228 boost::shared_ptr<FaceElementForcesAndSourcesCore> dragFe,
229 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideDragFe,
230 std::string side_fe_name, const std::string velocity_field,
231 const std::string pressure_field,
232 boost::shared_ptr<CommonData> common_data);
233
234 /**
235 * @brief Setting up operators for post processing output of drag traction
236 *
237 * Pushes operators for post processing ouput of drag traction components on
238 * the fluid-solid interface
239 *
240 * @param dragFe pointer to face element instance
241 * @param sideDragFe pointer to volume on side element instance
242 * @param velocity_field name of the velocity field
243 * @param pressure_field name of the pressure field
244 * @param common_data pointer to common data object
245 * @return Error code
246 */
248 boost::shared_ptr<PostProcFaceOnRefinedMesh> postProcDragPtr,
249 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideDragFe,
250 std::string side_fe_name, const std::string velocity_field,
251 const std::string pressure_field,
252 boost::shared_ptr<CommonData> common_data);
253
254 /**
255 * @brief Setting up operators for calculation of volume flux
256 */
258 boost::shared_ptr<VolumeElementForcesAndSourcesCore> fe_flux_ptr,
259 const std::string velocity_field,
260 boost::shared_ptr<CommonData> common_data, const EntityType type = MBTET);
261
262 /**
263 * \brief Set integration rule to volume elements
264 *
265 * Integration rule is order of polynomial which is calculated exactly.
266 * Finite element selects integration method based on return of this
267 * function.
268 *
269 */
270 struct VolRule {
271 int operator()(int order_row, int order_col, int order_data) const {
272 return 2 * order_data;
273 }
274 };
275
276 struct FaceRule {
277 int operator()(int order_row, int order_col, int order_data) const {
278 return order_data + 2;
279 }
280 };
281
282 /**
283 * \brief Base class for operators assembling LHS
284 */
286
287 boost::shared_ptr<CommonData> commonData;
290
292
296
298
299 OpAssembleLhs(const string field_name_row, const string field_name_col,
300 boost::shared_ptr<CommonData> common_data,
301 BlockData &block_data)
302 : UserDataOperator(field_name_row, field_name_col,
304 commonData(common_data), blockData(block_data) {
305 sYmm = false;
306 diagonalBlock = false;
307 };
308
309 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
310 EntityType col_type, EntData &row_data,
311 EntData &col_data);
312
313 virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
316 };
317
318 MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
319 };
320
321 /**
322 * @brief Assemble off-diagonal block of the LHS
323 * Operator for assembling off-diagonal block of the LHS
324 */
326
327 OpAssembleLhsOffDiag(const string field_name_row,
328 const string field_name_col,
329 boost::shared_ptr<CommonData> common_data,
330 BlockData &block_data)
331 : OpAssembleLhs(field_name_row, field_name_col, common_data,
332 block_data) {
333 sYmm = false;
334 diagonalBlock = false;
335 };
336
337 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
338 };
339
340 /**
341 * @brief Assemble linear (symmetric) part of the diagonal block of the LHS
342 * Operator for assembling linear (symmetric) part of the diagonal block of
343 * the LHS
344 */
346
348
349 OpAssembleLhsDiagLin(const string field_name_row,
350 const string field_name_col,
351 boost::shared_ptr<CommonData> common_data,
352 BlockData &block_data)
353 : OpAssembleLhs(field_name_row, field_name_col, common_data,
354 block_data) {
355 sYmm = true;
356 diagonalBlock = true;
357 };
358
359 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
360 };
361
362 /**
363 * @brief Assemble non-linear (non-symmetric) part of the diagonal block of
364 * the LHS Operator for assembling non-linear (non-symmetric) part of the
365 * diagonal block of the LHS
366 */
368
370
371 OpAssembleLhsDiagNonLin(const string field_name_row,
372 const string field_name_col,
373 boost::shared_ptr<CommonData> common_data,
374 BlockData &block_data)
375 : OpAssembleLhs(field_name_row, field_name_col, common_data,
376 block_data) {
377 sYmm = false;
378 diagonalBlock = true;
379 };
380
381 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
382 };
383
384 /**
385 * \brief Base class for operators assembling RHS
386 */
388
389 boost::shared_ptr<CommonData> commonData;
391 int nbRows; ///< number of dofs on row
393
395
397 boost::shared_ptr<CommonData> common_data,
398 BlockData &block_data)
400 commonData(common_data), blockData(block_data){};
401
402 MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
403
407 };
408
410 };
411
412 /**
413 * @brief Assemble linear part of the velocity component of the RHS vector
414 *
415 * Operator for assembling linear part of the velocity component of the RHS
416 * vector:
417 * \f[ \mathbf{R}^{\textrm{S}}_{\mathbf{u}} =
418 * C_\text{V}\int\limits_{\Omega}\nabla\mathbf{u}\mathbin{:}\nabla\mathbf{v}\,
419 * d\Omega
420 * \f]
421 * where \f$C_\text{V}\f$ is the viscosity
422 * coefficient: \f$C_\text{V}=\mu\f$ in the dimensional case and
423 * \f$C_\text{V}=1\f$ in the non-dimensional case.
424 */
426
428 boost::shared_ptr<CommonData> common_data,
429 BlockData &block_data)
430 : OpAssembleRhs(field_name, common_data, block_data){};
431
432 /**
433 * \brief Integrate local entity vector
434 * @param data entity data on element row
435 * @return error code
436 */
438 };
439
440 /**
441 * @brief Assemble non-linear part of the velocity component of the RHS vector
442 *
443 * Operator for assembling non-linear part of the velocity component of the
444 * RHS vector:
445 * \f[
446 * \mathbf{R}^{\textrm{NS}}_{\mathbf{u}} =
447 * C_\text{I}\int\limits_\Omega \left(\mathbf{u}\cdot\nabla\right)\mathbf{u}
448 * \cdot\mathbf{v} \,d\Omega,
449 * \f]
450 * where \f$C_\text{I}\f$ is the inertia
451 * coefficient: \f$C_\text{V}=\rho\f$ in the dimensional case and
452 * \f$C_\text{V}=\mathcal{R}\f$ in the non-dimensional case.
453 */
455
457 boost::shared_ptr<CommonData> common_data,
458 BlockData &block_data)
459 : OpAssembleRhs(field_name, common_data, block_data){};
460
462 };
463
464 /**
465 * @brief Assemble the pressure component of the RHS vector
466 *
467 * Operator for assembling pressure component of the RHS vector:
468 * \f[
469 * \mathbf{R}^{\textrm{S}}_{p} = -\int\limits_{\Omega}p\,
470 * \nabla\cdot\mathbf{v} \, d\Omega
471 * \f]
472 */
474
476 boost::shared_ptr<CommonData> common_data,
477 BlockData &block_data)
478 : OpAssembleRhs(field_name, common_data, block_data){};
479
480 /**
481 * \brief Integrate local constrains vector
482 */
484 };
485
486 /**
487 * @brief Calculate drag force on the fluid-solid interface
488 *
489 * Operator fo calculating drag force components on the fluid-solid interface.
490 * Integrates components of the drag traction:
491 * \f[
492 * \mathbf{F}_{\textrm{D}} =
493 * -\int\limits_{\Gamma_{\textrm{S}}}\left(-p\mathbf{I} +
494 * \mu\left(\nabla\mathbf{u}+\mathbf{u}^{\intercal}\right)\right) \, d\Gamma
495 * \f]
496 */
498
499 boost::shared_ptr<CommonData> commonData;
501
503 boost::shared_ptr<CommonData> &common_data,
504 BlockData &block_data)
507 commonData(common_data), blockData(block_data) {
508 doVertices = true;
509 doEdges = false;
510 doQuads = false;
511 doTris = false;
512 doTets = false;
513 doPrisms = false;
514 };
515
516 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
517 };
518
519 /**
520 * @brief Calculate drag traction on the fluid-solid interface
521 *
522 * Operator fo calculating drag traction on the fluid-solid interface
523 * \f[
524 * \mathbf{t} = -p\mathbf{I} +
525 * \mu\left(\nabla\mathbf{u}+\mathbf{u}^{\intercal}\right)
526 * \f]
527 */
529
530 boost::shared_ptr<CommonData> commonData;
532 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
533 std::string sideFeName;
534
536 const string field_name,
537 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &side_fe,
538 std::string side_fe_name, boost::shared_ptr<CommonData> &common_data,
539 BlockData &block_data)
542 sideFe(side_fe), sideFeName(side_fe_name), commonData(common_data),
543 blockData(block_data) {
544 doVertices = true;
545 doEdges = false;
546 doQuads = false;
547 doTris = false;
548 doTets = false;
549 doPrisms = false;
550 };
551
552 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
553 };
554
555 /**
556 * @brief Post processing output of drag traction on the fluid-solid interface
557 *
558 */
560
561 boost::shared_ptr<CommonData> commonData;
563 std::vector<EntityHandle> &mapGaussPts;
565
566 OpPostProcDrag(const string field_name, moab::Interface &post_proc_mesh,
567 std::vector<EntityHandle> &map_gauss_pts,
568 boost::shared_ptr<CommonData> &common_data,
569 BlockData &block_data)
572 commonData(common_data), postProcMesh(post_proc_mesh),
573 mapGaussPts(map_gauss_pts), blockData(block_data) {
574 doVertices = true;
575 doEdges = false;
576 doQuads = false;
577 doTris = false;
578 doTets = false;
579 doPrisms = false;
580 };
581
582 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
583 };
584
585 /**
586 * @brief Post processing output of the vorticity criterion levels
587 *
588 */
590
591 boost::shared_ptr<CommonData> commonData;
593 std::vector<EntityHandle> &mapGaussPts;
595
597 std::vector<EntityHandle> &map_gauss_pts,
598 boost::shared_ptr<CommonData> &common_data,
599 BlockData &block_data)
601 commonData(common_data), postProcMesh(post_proc_mesh),
602 mapGaussPts(map_gauss_pts), blockData(block_data) {
603 doVertices = true;
604 doEdges = false;
605 doQuads = false;
606 doTris = false;
607 doTets = false;
608 doPrisms = false;
609 };
610
611 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
612 };
613
614 /**
615 * @brief calculating volumetric flux
616 *
617 */
619
620 boost::shared_ptr<CommonData> commonData;
622 int nbRows; ///< number of dofs on row
624
626 boost::shared_ptr<CommonData> common_data,
627 BlockData &block_data)
629 commonData(common_data), blockData(block_data){};
630
631 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
632 };
633};
634
635#endif //__NAVIERSTOKESELEMENT_HPP__
EntitiesFieldData::EntData EntData
#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
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
auto createSmartVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
constexpr auto field_name
Class used to scale loads, f.e. in arc-length control.
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
bool & doTris
\deprectaed
bool & doPrisms
\deprectaed
bool sYmm
If true assume that matrix is symmetric structure.
bool & doVertices
\deprectaed If false skip vertices
bool & doTets
\deprectaed
bool & doEdges
\deprectaed If false skip edges
bool & doQuads
\deprectaed
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
structure for User Loop Methods on finite elements
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
boost::shared_ptr< MatrixDouble > pressureDragTract
boost::shared_ptr< MatrixDouble > gradVelPtr
std::map< int, BlockData > setOfBlocksData
boost::shared_ptr< MatrixDouble > shearDragTract
std::map< int, BlockData > setOfFacesData
boost::shared_ptr< VectorDouble > pressPtr
boost::shared_ptr< MatrixDouble > totalDragTract
CommonData(MoFEM::Interface &m_field)
SmartPetscObj< Vec > pressureDragForceVec
boost::shared_ptr< MatrixDouble > velPtr
int operator()(int order_row, int order_col, int order_data) const
MoFEMErrorCode scaleNf(const FEMethod *fe, VectorDouble &nf)
Assemble linear (symmetric) part of the diagonal block of the LHS Operator for assembling linear (sym...
OpAssembleLhsDiagLin(const string field_name_row, const string field_name_col, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
FTensor::Tensor2< double, 3, 3 > diffDiff
Assemble non-linear (non-symmetric) part of the diagonal block of the LHS Operator for assembling non...
OpAssembleLhsDiagNonLin(const string field_name_row, const string field_name_col, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Base class for operators assembling LHS.
MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
boost::shared_ptr< CommonData > commonData
OpAssembleLhs(const string field_name_row, const string field_name_col, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
Assemble off-diagonal block of the LHS Operator for assembling off-diagonal block of the LHS.
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
OpAssembleLhsOffDiag(const string field_name_row, const string field_name_col, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
Base class for operators assembling RHS.
OpAssembleRhs(const string field_name, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
boost::shared_ptr< CommonData > commonData
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode aSsemble(EntData &data)
virtual MoFEMErrorCode iNtegrate(EntData &data)
Assemble the pressure component of the RHS vector.
MoFEMErrorCode iNtegrate(EntData &data)
Integrate local constrains vector.
OpAssembleRhsPressure(const string field_name, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
Assemble linear part of the velocity component of the RHS vector.
MoFEMErrorCode iNtegrate(EntData &data)
Integrate local entity vector.
OpAssembleRhsVelocityLin(const string field_name, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
Assemble non-linear part of the velocity component of the RHS vector.
OpAssembleRhsVelocityNonLin(const string field_name, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
Calculate drag force on the fluid-solid interface.
boost::shared_ptr< CommonData > commonData
OpCalcDragForce(const string field_name, boost::shared_ptr< CommonData > &common_data, BlockData &block_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Calculate drag traction on the fluid-solid interface.
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< CommonData > commonData
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
OpCalcDragTraction(const string field_name, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > &side_fe, std::string side_fe_name, boost::shared_ptr< CommonData > &common_data, BlockData &block_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalcVolumeFlux(const string field_name, boost::shared_ptr< CommonData > common_data, BlockData &block_data)
boost::shared_ptr< CommonData > commonData
Post processing output of drag traction on the fluid-solid interface.
boost::shared_ptr< CommonData > commonData
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpPostProcDrag(const string field_name, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > &common_data, BlockData &block_data)
moab::Interface & postProcMesh
BlockData & blockData
std::vector< EntityHandle > & mapGaussPts
Post processing output of the vorticity criterion levels.
moab::Interface & postProcMesh
BlockData & blockData
std::vector< EntityHandle > & mapGaussPts
OpPostProcVorticity(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > &common_data, BlockData &block_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< CommonData > commonData
Set integration rule to volume elements.
int operator()(int order_row, int order_col, int order_data) const
Element for simulating viscous fluid flow.
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< 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.
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 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.