v0.14.0
navier_stokes.cpp
Go to the documentation of this file.
1 /** \file navier_stokes.cpp
2  * \example navier_stokes.cpp
3  *
4  * Example of viscous fluid flow problem
5  *
6  */
7 
8 #include <MoFEM.hpp>
9 using namespace MoFEM;
10 
11 using PostProcVol = PostProcBrokenMeshInMoab<VolumeElementForcesAndSourcesCore>;
12 using PostProcFace = PostProcBrokenMeshInMoab<FaceElementForcesAndSourcesCore>;
13 
14 using namespace boost::numeric;
15 using namespace std;
16 
17 static char help[] = "Navier-Stokes Example\n";
18 
19 double NavierStokesElement::LoadScale::lambda = 1;
20 
21 //! [Example Navier Stokes]
22 struct NavierStokesExample {
23 
24  NavierStokesExample(MoFEM::Interface &m_field) : mField(m_field) {
25  commonData = boost::make_shared<NavierStokesElement::CommonData>(m_field);
26  }
27  ~NavierStokesExample() {
28  CHKERR SNESDestroy(&snes);
29  CHKERR DMDestroy(&dM);
30  }
31 
32  MoFEMErrorCode runProblem();
33 
34 private:
35  MoFEM::Interface &mField;
36 
37  PetscBool isPartitioned = PETSC_FALSE;
38 
39  int orderVelocity = 2; // default approximation orderVelocity
40  int orderPressure = 1; // default approximation orderPressure
41  int numHoLevels = 0;
42  PetscBool isDiscontPressure = PETSC_FALSE;
43  PetscBool isHoGeometry = PETSC_TRUE;
44 
45  double pressureScale = 1.0;
46  double lengthScale = 1.0;
47  double velocityScale = 1.0;
48  double reNumber = 1.0;
49  double density;
50  double viscosity;
51  PetscBool isStokesFlow = PETSC_FALSE;
52 
53  int numSteps = 1;
54  double lambdaStep = 1.0;
55  double lambda = 0.0;
56  int step;
57 
58  Range solidFaces;
59  BitRefLevel bitLevel;
60 
61  DM dM;
62  SNES snes;
63 
64  boost::shared_ptr<NavierStokesElement::CommonData> commonData;
65 
66  SmartPetscObj<Vec> D;
67  SmartPetscObj<Vec> F;
68  SmartPetscObj<Mat> Aij;
69 
70  boost::shared_ptr<VolumeElementForcesAndSourcesCore> feLhsPtr;
71  boost::shared_ptr<VolumeElementForcesAndSourcesCore> feRhsPtr;
72 
73  boost::shared_ptr<FaceElementForcesAndSourcesCore> feDragPtr;
74  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> feDragSidePtr;
75 
76  boost::shared_ptr<DirichletDisplacementBc> dirichletBcPtr;
77  boost::ptr_map<std::string, NeumannForcesSurface> neumannForces;
78 
79  boost::shared_ptr<PostProcVol> fePostProcPtr;
80  boost::shared_ptr<PostProcFace> fePostProcDragPtr;
81 
82  MoFEMErrorCode readInput();
83  MoFEMErrorCode findBlocks();
84  MoFEMErrorCode setupFields();
85  MoFEMErrorCode defineFiniteElements();
86  MoFEMErrorCode setupDiscreteManager();
87  MoFEMErrorCode setupAlgebraicStructures();
88  MoFEMErrorCode setupElementInstances();
89  MoFEMErrorCode setupSNES();
90  MoFEMErrorCode solveProblem();
91  MoFEMErrorCode postProcess();
92 };
93 //! [Example Navier Stokes]
94 
95 //! [Run problem]
96 MoFEMErrorCode NavierStokesExample::runProblem() {
97  MoFEMFunctionBegin;
98  CHKERR readInput();
99  CHKERR findBlocks();
100  CHKERR setupFields();
101  CHKERR defineFiniteElements();
102  CHKERR setupDiscreteManager();
103  CHKERR setupAlgebraicStructures();
104  CHKERR setupElementInstances();
105  CHKERR setupSNES();
106  CHKERR solveProblem();
107  MoFEMFunctionReturn(0);
108 }
109 //! [Run problem]
110 
111 //! [Read input]
112 MoFEMErrorCode NavierStokesExample::readInput() {
113  MoFEMFunctionBegin;
114  char mesh_file_name[255];
115  PetscBool flg_mesh_file;
116 
117  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "NAVIER_STOKES problem",
118  "none");
119 
120  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
121  mesh_file_name, 255, &flg_mesh_file);
122  CHKERR PetscOptionsBool("-my_is_partitioned", "is partitioned?", "",
123  isPartitioned, &isPartitioned, PETSC_NULL);
124 
125  CHKERR PetscOptionsInt("-my_order_u", "approximation orderVelocity", "",
126  orderVelocity, &orderVelocity, PETSC_NULL);
127  CHKERR PetscOptionsInt("-my_order_p", "approximation orderPressure", "",
128  orderPressure, &orderPressure, PETSC_NULL);
129  CHKERR PetscOptionsInt("-my_num_ho_levels",
130  "number of higher order boundary levels", "",
131  numHoLevels, &numHoLevels, PETSC_NULL);
132  CHKERR PetscOptionsBool("-my_discont_pressure", "discontinuous pressure", "",
133  isDiscontPressure, &isDiscontPressure, PETSC_NULL);
134  CHKERR PetscOptionsBool("-my_ho_geometry", "use second order geometry", "",
135  isHoGeometry, &isHoGeometry, PETSC_NULL);
136 
137  CHKERR PetscOptionsScalar("-my_length_scale", "length scale", "", lengthScale,
138  &lengthScale, PETSC_NULL);
139  CHKERR PetscOptionsScalar("-my_velocity_scale", "velocity scale", "",
140  velocityScale, &velocityScale, PETSC_NULL);
141  CHKERR PetscOptionsBool("-my_stokes_flow", "stokes flow", "", isStokesFlow,
142  &isStokesFlow, PETSC_NULL);
143 
144  CHKERR PetscOptionsInt("-my_step_num", "number of steps", "", numSteps,
145  &numSteps, PETSC_NULL);
146 
147  ierr = PetscOptionsEnd();
148 
149  if (flg_mesh_file != PETSC_TRUE) {
150  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
151  }
152 
153  // Read whole mesh or part of it if partitioned
154  if (isPartitioned == PETSC_TRUE) {
155  // This is a case of distributed mesh and algebra. In that case each
156  // processor keeps only part of the problem.
157  const char *option;
158  option = "PARALLEL=READ_PART;"
159  "PARALLEL_RESOLVE_SHARED_ENTS;"
160  "PARTITION=PARALLEL_PARTITION;";
161  CHKERR mField.get_moab().load_file(mesh_file_name, 0, option);
162  } else {
163  // In this case, we have distributed algebra, i.e. assembly of vectors and
164  // matrices is in parallel, but whole mesh is stored on all processors.
165  // snes and matrix scales well, however problem set-up of problem is
166  // not fully parallel.
167  const char *option;
168  option = "";
169  CHKERR mField.get_moab().load_file(mesh_file_name, 0, option);
170  }
171  CHKERR mField.rebuild_database();
172 
173  bitLevel.set(0);
174  CHKERR mField.getInterface<BitRefManager>()->setBitRefLevelByDim(0, 3,
175  bitLevel);
176 
177  MoFEMFunctionReturn(0);
178 }
179 //! [Read input]
180 
181 //! [Find blocks]
182 MoFEMErrorCode NavierStokesExample::findBlocks() {
183  MoFEMFunctionBegin;
184 
185  for (_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(mField, BLOCKSET, bit)) {
186  if (bit->getName().compare(0, 9, "MAT_FLUID") == 0) {
187  const int id = bit->getMeshsetId();
188  CHKERR mField.get_moab().get_entities_by_dimension(
189  bit->getMeshset(), 3, commonData->setOfBlocksData[id].eNts, true);
190 
191  std::vector<double> attributes;
192  bit->getAttributes(attributes);
193  if (attributes.size() < 2) {
194  SETERRQ1(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
195  "should be at least 2 attributes but is %d",
196  attributes.size());
197  }
198  commonData->setOfBlocksData[id].iD = id;
199  commonData->setOfBlocksData[id].fluidViscosity = attributes[0];
200  commonData->setOfBlocksData[id].fluidDensity = attributes[1];
201 
202  viscosity = commonData->setOfBlocksData[id].fluidViscosity;
203  density = commonData->setOfBlocksData[id].fluidDensity;
204  }
205  }
206 
207  for (_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(mField, BLOCKSET, bit)) {
208  if (bit->getName().compare(0, 9, "INT_SOLID") == 0) {
209  Range tets, tet;
210  const int id = bit->getMeshsetId();
211  CHKERR mField.get_moab().get_entities_by_type(
212  bit->getMeshset(), MBTRI, commonData->setOfFacesData[id].eNts, true);
213  solidFaces.merge(commonData->setOfFacesData[id].eNts);
214  CHKERR mField.get_moab().get_adjacencies(
215  commonData->setOfFacesData[id].eNts, 3, true, tets,
216  moab::Interface::UNION);
217  tet = Range(tets.front(), tets.front());
218  for (auto &bit : commonData->setOfBlocksData) {
219  if (bit.second.eNts.contains(tet)) {
220  commonData->setOfFacesData[id].fluidViscosity =
221  bit.second.fluidViscosity;
222  commonData->setOfFacesData[id].fluidDensity = bit.second.fluidDensity;
223  commonData->setOfFacesData[id].iD = id;
224  break;
225  }
226  }
227  if (commonData->setOfFacesData[id].fluidViscosity < 0) {
228  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
229  "Cannot find a fluid block adjacent to a given solid face");
230  }
231  }
232  }
233 
234  MoFEMFunctionReturn(0);
235 }
236 //! [Find blocks]
237 
238 //! [Setup fields]
239 MoFEMErrorCode NavierStokesExample::setupFields() {
240  MoFEMFunctionBegin;
241 
242  CHKERR mField.add_field("VELOCITY", H1, AINSWORTH_LEGENDRE_BASE, 3);
243  if (isDiscontPressure) {
244  CHKERR mField.add_field("PRESSURE", L2, AINSWORTH_LEGENDRE_BASE, 1);
245  } else {
246  CHKERR mField.add_field("PRESSURE", H1, AINSWORTH_LEGENDRE_BASE, 1);
247  }
248  CHKERR mField.add_field("MESH_NODE_POSITIONS", H1, AINSWORTH_LEGENDRE_BASE,
249  3);
250 
251  CHKERR mField.add_ents_to_field_by_dim(0, 3, "VELOCITY");
252  CHKERR mField.add_ents_to_field_by_dim(0, 3, "PRESSURE");
253  CHKERR mField.add_ents_to_field_by_dim(0, 3, "MESH_NODE_POSITIONS");
254 
255  CHKERR mField.set_field_order(0, MBVERTEX, "VELOCITY", 1);
256  CHKERR mField.set_field_order(0, MBEDGE, "VELOCITY", orderVelocity);
257  CHKERR mField.set_field_order(0, MBTRI, "VELOCITY", orderVelocity);
258  CHKERR mField.set_field_order(0, MBTET, "VELOCITY", orderVelocity);
259 
260  if (!isDiscontPressure) {
261  CHKERR mField.set_field_order(0, MBVERTEX, "PRESSURE", 1);
262  CHKERR mField.set_field_order(0, MBEDGE, "PRESSURE", orderPressure);
263  CHKERR mField.set_field_order(0, MBTRI, "PRESSURE", orderPressure);
264  }
265  CHKERR mField.set_field_order(0, MBTET, "PRESSURE", orderPressure);
266 
267  if (numHoLevels > 0) {
268  std::vector<Range> levels(numHoLevels);
269  Range ents;
270  ents.merge(solidFaces);
271  for (int ll = 0; ll != numHoLevels; ll++) {
272  Range verts;
273  CHKERR mField.get_moab().get_connectivity(ents, verts, true);
274  for (auto d : {1, 2, 3}) {
275  CHKERR mField.get_moab().get_adjacencies(verts, d, false, ents,
276  moab::Interface::UNION);
277  }
278  levels[ll] = subtract(ents, ents.subset_by_type(MBVERTEX));
279  }
280  for (int ll = numHoLevels - 1; ll >= 1; ll--) {
281  levels[ll] = subtract(levels[ll], levels[ll - 1]);
282  }
283 
284  int add_order = 1;
285  for (int ll = numHoLevels - 1; ll >= 0; ll--) {
286  if (isPartitioned)
287  CHKERR mField.getInterface<CommInterface>()->synchroniseEntities(
288  levels[ll]);
289 
290  CHKERR mField.set_field_order(levels[ll], "VELOCITY",
291  orderVelocity + add_order);
292  if (!isDiscontPressure)
293  CHKERR mField.set_field_order(levels[ll], "PRESSURE",
294  orderPressure + add_order);
295  else
296  CHKERR mField.set_field_order(levels[ll].subset_by_type(MBTET),
297  "PRESSURE", orderPressure + add_order);
298  ++add_order;
299  }
300  }
301 
302  CHKERR mField.set_field_order(0, MBVERTEX, "MESH_NODE_POSITIONS", 1);
303  // Set 2nd order of approximation of geometry
304  if (isHoGeometry) {
305  Range ents, edges;
306  CHKERR mField.get_moab().get_entities_by_dimension(0, 3, ents);
307  CHKERR mField.get_moab().get_adjacencies(ents, 1, false, edges,
308  moab::Interface::UNION);
309  if (isPartitioned)
310  CHKERR mField.getInterface<CommInterface>()->synchroniseEntities(edges);
311  CHKERR mField.set_field_order(edges, "MESH_NODE_POSITIONS", 2);
312  }
313 
314  if (isPartitioned) {
315  CHKERR mField.getInterface<CommInterface>()->synchroniseFieldEntities(
316  "VELOCITY");
317  CHKERR mField.getInterface<CommInterface>()->synchroniseFieldEntities(
318  "PRESSURE");
319  CHKERR mField.getInterface<CommInterface>()->synchroniseFieldEntities(
320  "MESH_NODE_POSITIONS");
321  }
322 
323  CHKERR mField.build_fields();
324 
325  Projection10NodeCoordsOnField ent_method_material(mField,
326  "MESH_NODE_POSITIONS");
327  CHKERR mField.loop_dofs("MESH_NODE_POSITIONS", ent_method_material);
328 
329  MoFEMFunctionReturn(0);
330 }
331 //! [Setup fields]
332 
333 //! [Define finite elements]
334 MoFEMErrorCode NavierStokesExample::defineFiniteElements() {
335  MoFEMFunctionBegin;
336 
337  // Add finite element (this defines element, declaration comes later)
338  CHKERR NavierStokesElement::addElement(mField, "NAVIER_STOKES", "VELOCITY",
339  "PRESSURE", "MESH_NODE_POSITIONS");
340 
341  CHKERR NavierStokesElement::addElement(mField, "DRAG", "VELOCITY", "PRESSURE",
342  "MESH_NODE_POSITIONS", 2, &solidFaces);
343 
344  // setup elements for loading
345  CHKERR MetaNeumannForces::addNeumannBCElements(mField, "VELOCITY");
346 
347  // build finite elements
348  CHKERR mField.build_finite_elements();
349  // build adjacencies between elements and degrees of freedom
350  CHKERR mField.build_adjacencies(bitLevel);
351 
352  MoFEMFunctionReturn(0);
353 }
354 //! [Define finite elements]
355 
356 //! [Setup discrete manager]
357 MoFEMErrorCode NavierStokesExample::setupDiscreteManager() {
358  MoFEMFunctionBegin;
359  DMType dm_name = "DM_NAVIER_STOKES";
360  // Register DM problem
361  CHKERR DMRegister_MoFEM(dm_name);
362  CHKERR DMCreate(PETSC_COMM_WORLD, &dM);
363  CHKERR DMSetType(dM, dm_name);
364  // Create DM instance
365  CHKERR DMMoFEMCreateMoFEM(dM, &mField, dm_name, bitLevel);
366  // Configure DM form line command options s
367  CHKERR DMSetFromOptions(dM);
368  // Add elements to dM (only one here)
369  CHKERR DMMoFEMAddElement(dM, "NAVIER_STOKES");
370  CHKERR DMMoFEMAddElement(dM, "DRAG");
371  if (mField.check_finite_element("PRESSURE_FE"))
372  CHKERR DMMoFEMAddElement(dM, "PRESSURE_FE");
373  CHKERR DMMoFEMSetIsPartitioned(dM, isPartitioned);
374  // setup the DM
375  CHKERR DMSetUp(dM);
376  MoFEMFunctionReturn(0);
377 }
378 //! [Setup discrete manager]
379 
380 //! [Setup algebraic structures]
381 MoFEMErrorCode NavierStokesExample::setupAlgebraicStructures() {
382  MoFEMFunctionBegin;
383 
384  D = createDMVector(dM);
385  F = vectorDuplicate(D);
386  Aij = createDMMatrix(dM);
387 
388  CHKERR VecZeroEntries(F);
389  CHKERR VecGhostUpdateBegin(F, INSERT_VALUES, SCATTER_FORWARD);
390  CHKERR VecGhostUpdateEnd(F, INSERT_VALUES, SCATTER_FORWARD);
391 
392  CHKERR VecZeroEntries(D);
393  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
394  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
395 
396  CHKERR MatSetOption(Aij, MAT_SPD, PETSC_TRUE);
397  CHKERR MatZeroEntries(Aij);
398 
399  MoFEMFunctionReturn(0);
400 }
401 //! [Setup algebraic structures]
402 
403 //! [Setup element instances]
404 MoFEMErrorCode NavierStokesExample::setupElementInstances() {
405  MoFEMFunctionBegin;
406 
407  feLhsPtr = boost::shared_ptr<VolumeElementForcesAndSourcesCore>(
408  new VolumeElementForcesAndSourcesCore(mField));
409  feRhsPtr = boost::shared_ptr<VolumeElementForcesAndSourcesCore>(
410  new VolumeElementForcesAndSourcesCore(mField));
411 
412  feLhsPtr->getRuleHook = NavierStokesElement::VolRule();
413  feRhsPtr->getRuleHook = NavierStokesElement::VolRule();
414  CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *feLhsPtr, true, false, false,
415  true);
416  CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *feRhsPtr, true, false, false,
417  true);
418 
419  feDragPtr = boost::shared_ptr<FaceElementForcesAndSourcesCore>(
420  new FaceElementForcesAndSourcesCore(mField));
421  feDragSidePtr = boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide>(
422  new VolumeElementForcesAndSourcesCoreOnSide(mField));
423 
424  feDragPtr->getRuleHook = NavierStokesElement::FaceRule();
425  CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *feDragSidePtr, true, false, false,
426  true);
427 
428  if (isStokesFlow) {
429  CHKERR NavierStokesElement::setStokesOperators(
430  feRhsPtr, feLhsPtr, "VELOCITY", "PRESSURE", commonData);
431  } else {
432  CHKERR NavierStokesElement::setNavierStokesOperators(
433  feRhsPtr, feLhsPtr, "VELOCITY", "PRESSURE", commonData);
434  }
435 
436  NavierStokesElement::setCalcDragOperators(feDragPtr, feDragSidePtr,
437  "NAVIER_STOKES", "VELOCITY",
438  "PRESSURE", commonData);
439 
440  dirichletBcPtr = boost::shared_ptr<DirichletDisplacementBc>(
441  new DirichletDisplacementBc(mField, "VELOCITY", Aij, D, F));
442  dirichletBcPtr->methodsOp.push_back(new NavierStokesElement::LoadScale());
443  // dirichletBcPtr->snes_ctx = FEMethod::CTX_SNESNONE;
444  dirichletBcPtr->snes_x = D;
445 
446  // Assemble pressure and traction forces
447  CHKERR MetaNeumannForces::setMomentumFluxOperators(mField, neumannForces,
448  NULL, "VELOCITY");
449 
450  // for postprocessing:
451  fePostProcPtr = boost::make_shared<PostProcVol>(mField);
452  CHKERR addHOOpsVol("MESH_NODE_POSITIONS", *fePostProcPtr, true, false, false,
453  true);
454 
455  auto v_ptr = boost::make_shared<MatrixDouble>();
456  auto grad_ptr = boost::make_shared<MatrixDouble>();
457  auto pos_ptr = boost::make_shared<MatrixDouble>();
458  auto p_ptr = boost::make_shared<VectorDouble>();
459 
460  fePostProcPtr->getOpPtrVector().push_back(
461  new OpCalculateVectorFieldValues<3>("VELOCITY", v_ptr));
462  fePostProcPtr->getOpPtrVector().push_back(
463  new OpCalculateVectorFieldGradient<3, 3>("VELOCITY", grad_ptr));
464  fePostProcPtr->getOpPtrVector().push_back(
465  new OpCalculateScalarFieldValues("PRESSURE", p_ptr));
466  fePostProcPtr->getOpPtrVector().push_back(
467  new OpCalculateVectorFieldValues<3>("MESH_NODE_POSITIONS", pos_ptr));
468 
469  using OpPPMap = OpPostProcMapInMoab<3, 3>;
470 
471  fePostProcPtr->getOpPtrVector().push_back(
472 
473  new OpPPMap(
474 
475  fePostProcPtr->getPostProcMesh(), fePostProcPtr->getMapGaussPts(),
476 
477  {{"PRESSURE", p_ptr}},
478 
479  {{"VELOCITY", v_ptr}, {"MESH_NODE_POSITIONS", pos_ptr}},
480 
481  {{"VELOCITY_GRAD", grad_ptr}},
482 
483  {}
484 
485  )
486 
487  );
488 
489  fePostProcDragPtr = boost::make_shared<PostProcFace>(mField);
490  fePostProcDragPtr->getOpPtrVector().push_back(
491  new OpCalculateVectorFieldValues<3>("MESH_NODE_POSITIONS", pos_ptr));
492  fePostProcDragPtr->getOpPtrVector().push_back(
493 
494  new OpPPMap(
495 
496  fePostProcDragPtr->getPostProcMesh(),
497  fePostProcDragPtr->getMapGaussPts(),
498 
499  {},
500 
501  {{"MESH_NODE_POSITIONS", pos_ptr}},
502 
503  {},
504 
505  {}
506 
507  )
508 
509  );
510 
511  CHKERR NavierStokesElement::setPostProcDragOperators(
512  fePostProcDragPtr, feDragSidePtr, "NAVIER_STOKES", "VELOCITY", "PRESSURE",
513  commonData);
514 
515  MoFEMFunctionReturn(0);
516 }
517 //! [Setup element instances]
518 
519 //! [Setup SNES]
520 MoFEMErrorCode NavierStokesExample::setupSNES() {
521  MoFEMFunctionBegin;
522 
523  boost::ptr_map<std::string, NeumannForcesSurface>::iterator mit =
524  neumannForces.begin();
525  for (; mit != neumannForces.end(); mit++) {
526  CHKERR DMMoFEMSNESSetFunction(dM, mit->first.c_str(),
527  &mit->second->getLoopFe(), NULL, NULL);
528  }
529 
530  boost::shared_ptr<FEMethod> null_fe;
531  CHKERR DMMoFEMSNESSetFunction(dM, "NAVIER_STOKES", feRhsPtr, null_fe,
532  null_fe);
533  CHKERR DMMoFEMSNESSetFunction(dM, DM_NO_ELEMENT, null_fe, null_fe,
534  dirichletBcPtr);
535 
536  CHKERR DMMoFEMSNESSetJacobian(dM, DM_NO_ELEMENT, null_fe, dirichletBcPtr,
537  null_fe);
538  CHKERR DMMoFEMSNESSetJacobian(dM, "NAVIER_STOKES", feLhsPtr, null_fe,
539  null_fe);
540  CHKERR DMMoFEMSNESSetJacobian(dM, DM_NO_ELEMENT, null_fe, null_fe,
541  dirichletBcPtr);
542 
543  SnesCtx *snes_ctx;
544  // create snes nonlinear solver
545  CHKERR SNESCreate(PETSC_COMM_WORLD, &snes);
546  CHKERR DMMoFEMGetSnesCtx(dM, &snes_ctx);
547  CHKERR SNESSetFunction(snes, F, SnesRhs, snes_ctx);
548  CHKERR SNESSetJacobian(snes, Aij, Aij, SnesMat, snes_ctx);
549  CHKERR SNESSetFromOptions(snes);
550 
551  MoFEMFunctionReturn(0);
552 }
553 //! [Setup SNES]
554 
555 //! [Solve problem]
556 MoFEMErrorCode NavierStokesExample::solveProblem() {
557  MoFEMFunctionBegin;
558 
559  auto scale_problem = [&](double U, double L, double P) {
560  MoFEMFunctionBegin;
561  CHKERR mField.getInterface<FieldBlas>()->fieldScale(L,
562  "MESH_NODE_POSITIONS");
563 
564  ProjectionFieldOn10NodeTet ent_method_on_10nodeTet(
565  mField, "MESH_NODE_POSITIONS", true, true);
566  CHKERR mField.loop_dofs("MESH_NODE_POSITIONS", ent_method_on_10nodeTet);
567  ent_method_on_10nodeTet.setNodes = false;
568  CHKERR mField.loop_dofs("MESH_NODE_POSITIONS", ent_method_on_10nodeTet);
569 
570  CHKERR mField.getInterface<FieldBlas>()->fieldScale(U, "VELOCITY");
571  CHKERR mField.getInterface<FieldBlas>()->fieldScale(P, "PRESSURE");
572  MoFEMFunctionReturn(0);
573  };
574 
575  pressureScale = viscosity * velocityScale / lengthScale;
576  NavierStokesElement::LoadScale::lambda = 1.0 / velocityScale;
577  CHKERR scale_problem(1.0 / velocityScale, 1.0 / lengthScale,
578  1.0 / pressureScale);
579 
580  step = 0;
581 
582  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Viscosity: %6.4e\n", viscosity);
583  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Density: %6.4e\n", density);
584  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Length scale: %6.4e\n", lengthScale);
585  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Velocity scale: %6.4e\n",
586  velocityScale);
587  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Pressure scale: %6.4e\n",
588  pressureScale);
589  if (isStokesFlow) {
590  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Re number: 0 (Stokes flow)\n");
591  } else {
592  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Re number: %6.4e\n",
593  density / viscosity * velocityScale * lengthScale);
594  }
595 
596  lambdaStep = 1.0 / numSteps;
597 
598  while (lambda < 1.0 - 1e-12) {
599 
600  lambda += lambdaStep;
601 
602  if (isStokesFlow) {
603  reNumber = 0.0;
604  for (auto &bit : commonData->setOfBlocksData) {
605  bit.second.inertiaCoef = 0.0;
606  bit.second.viscousCoef = 1.0;
607  }
608  } else {
609  reNumber = density / viscosity * velocityScale * lengthScale * lambda;
610  for (auto &bit : commonData->setOfBlocksData) {
611  bit.second.inertiaCoef = reNumber;
612  bit.second.viscousCoef = 1.0;
613  }
614  }
615 
616  CHKERR PetscPrintf(
617  PETSC_COMM_WORLD,
618  "Step: %d | Lambda: %6.4e | Inc: %6.4e | Re number: %6.4e \n", step,
619  lambda, lambdaStep, reNumber);
620 
621  CHKERR DMoFEMPreProcessFiniteElements(dM, dirichletBcPtr.get());
622 
623  CHKERR VecAssemblyBegin(D);
624  CHKERR VecAssemblyEnd(D);
625  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
626  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
627 
628  CHKERR SNESSolve(snes, PETSC_NULL, D);
629 
630  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
631  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
632  CHKERR DMoFEMMeshToGlobalVector(dM, D, INSERT_VALUES, SCATTER_REVERSE);
633 
634  CHKERR scale_problem(velocityScale, lengthScale, pressureScale);
635 
636  CHKERR postProcess();
637 
638  CHKERR scale_problem(1.0 / velocityScale, 1.0 / lengthScale,
639  1.0 / pressureScale);
640 
641  step++;
642  }
643 
644  MoFEMFunctionReturn(0);
645 }
646 //! [Solve problem]
647 
648 //! [Post process]
649 MoFEMErrorCode NavierStokesExample::postProcess() {
650  MoFEMFunctionBegin;
651 
652  string out_file_name;
653 
654  CHKERR DMoFEMLoopFiniteElements(dM, "NAVIER_STOKES", fePostProcPtr);
655  {
656  std::ostringstream stm;
657  stm << "out_" << step << ".h5m";
658  out_file_name = stm.str();
659  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Write file %s\n",
660  out_file_name.c_str());
661  CHKERR fePostProcPtr->writeFile(out_file_name.c_str());
662  }
663 
664  CHKERR VecZeroEntries(commonData->pressureDragForceVec);
665  CHKERR VecZeroEntries(commonData->shearDragForceVec);
666  CHKERR VecZeroEntries(commonData->totalDragForceVec);
667  CHKERR DMoFEMLoopFiniteElements(dM, "DRAG", feDragPtr);
668 
669  auto get_vec_data = [&](auto vec, std::array<double, 3> &data) {
670  MoFEMFunctionBegin;
671  CHKERR VecAssemblyBegin(vec);
672  CHKERR VecAssemblyEnd(vec);
673  const double *array;
674  CHKERR VecGetArrayRead(vec, &array);
675  if (mField.get_comm_rank() == 0) {
676  for (int i : {0, 1, 2})
677  data[i] = array[i];
678  }
679  CHKERR VecRestoreArrayRead(vec, &array);
680  MoFEMFunctionReturn(0);
681  };
682 
683  std::array<double, 3> pressure_drag;
684  std::array<double, 3> shear_drag;
685  std::array<double, 3> total_drag;
686 
687  CHKERR get_vec_data(commonData->pressureDragForceVec, pressure_drag);
688  CHKERR get_vec_data(commonData->shearDragForceVec, shear_drag);
689  CHKERR get_vec_data(commonData->totalDragForceVec, total_drag);
690 
691  if (mField.get_comm_rank() == 0) {
692  MOFEM_LOG_C("SELF", Sev::inform,
693  "Pressure drag force: [ %6.4e, %6.4e, %6.4e ]",
694  pressure_drag[0], pressure_drag[1], pressure_drag[2]);
695  MOFEM_LOG_C("SELF", Sev::inform,
696  "Shear drag force: [ %6.4e, %6.4e, %6.4e ]", shear_drag[0],
697  shear_drag[1], shear_drag[2]);
698  MOFEM_LOG_C("SELF", Sev::inform,
699  "Total drag force: [ %6.4e, %6.4e, %6.4e ]", total_drag[0],
700  total_drag[1], total_drag[2]);
701  }
702 
703  CHKERR DMoFEMLoopFiniteElements(dM, "DRAG", fePostProcDragPtr);
704  {
705  std::ostringstream stm;
706  stm << "out_drag_" << step << ".h5m";
707  out_file_name = stm.str();
708  CHKERR PetscPrintf(PETSC_COMM_WORLD, "out file %s\n",
709  out_file_name.c_str());
710  CHKERR fePostProcDragPtr->writeFile(out_file_name);
711  }
712  MoFEMFunctionReturn(0);
713 }
714 //! [Post process]
715 
716 //! [Main function]
717 int main(int argc, char *argv[]) {
718 
719  const char param_file[] = "param_file.petsc";
720  // Initialise MoFEM
721  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
722 
723  try {
724 
725  // Create mesh database
726  moab::Core mb_instance; // create database
727  moab::Interface &moab = mb_instance; // create interface to database
728 
729  // Create MoFEM database and link it to MoAB
730  MoFEM::Core core(moab);
731  MoFEM::Interface &m_field = core;
732 
733  NavierStokesExample ex(m_field);
734  CHKERR ex.runProblem();
735  }
736  CATCH_ERRORS;
737 
738  // finish work cleaning memory, getting statistics, etc
739  MoFEM::Core::Finalize();
740  return 0;
741 }
742 //! [Main function]
NavierStokesExample::density
double density
Definition: navier_stokes.cpp:49
NavierStokesExample::NavierStokesExample
NavierStokesExample(MoFEM::Interface &m_field)
Definition: navier_stokes.cpp:24
NavierStokesExample::runProblem
MoFEMErrorCode runProblem()
[Example Navier Stokes]
Definition: navier_stokes.cpp:96
NavierStokesExample::F
SmartPetscObj< Vec > F
Definition: navier_stokes.cpp:67
MetaNeumannForces::addNeumannBCElements
static MoFEMErrorCode addNeumannBCElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", Range *intersect_ptr=NULL)
Declare finite element.
Definition: SurfacePressure.cpp:1974
NavierStokesExample::feDragSidePtr
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > feDragSidePtr
Definition: navier_stokes.cpp:74
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
help
static char help[]
Definition: navier_stokes.cpp:17
MoFEM::addHOOpsVol
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
Definition: HODataOperators.hpp:674
NavierStokesExample::feDragPtr
boost::shared_ptr< FaceElementForcesAndSourcesCore > feDragPtr
Definition: navier_stokes.cpp:73
NavierStokesExample::solidFaces
Range solidFaces
Definition: navier_stokes.cpp:58
NavierStokesExample
[Example Navier Stokes]
Definition: navier_stokes.cpp:22
NavierStokesExample::dirichletBcPtr
boost::shared_ptr< DirichletDisplacementBc > dirichletBcPtr
Definition: navier_stokes.cpp:76
NavierStokesElement::setPostProcDragOperators
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.
Definition: NavierStokesElement.cpp:184
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
NavierStokesExample::Aij
SmartPetscObj< Mat > Aij
Definition: navier_stokes.cpp:68
NavierStokesExample::mField
MoFEM::Interface & mField
Definition: navier_stokes.cpp:35
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::OpCalculateVectorFieldValues
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Definition: UserDataOperators.hpp:466
NavierStokesElement::FaceRule
Definition: NavierStokesElement.hpp:269
MoFEM.hpp
MoFEM::Projection10NodeCoordsOnField
Projection of edge entities with one mid-node on hierarchical basis.
Definition: Projection10NodeCoordsOnField.hpp:24
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
out_file_name
char out_file_name[255]
Definition: initial_diffusion.cpp:53
NavierStokesExample::fePostProcDragPtr
boost::shared_ptr< PostProcFace > fePostProcDragPtr
Definition: navier_stokes.cpp:80
NavierStokesExample::findBlocks
MoFEMErrorCode findBlocks()
[Read input]
Definition: navier_stokes.cpp:182
NavierStokesExample::dM
DM dM
Definition: navier_stokes.cpp:61
DirichletDisplacementBc
Set Dirichlet boundary conditions on displacements.
Definition: DirichletBC.hpp:57
MoFEM::createDMMatrix
auto createDMMatrix(DM dm)
Get smart matrix from DM.
Definition: DMMoFEM.hpp:1056
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:497
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
NavierStokesElement::VolRule
Set integration rule to volume elements.
Definition: NavierStokesElement.hpp:263
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
main
int main(int argc, char *argv[])
[Post process]
Definition: navier_stokes.cpp:717
NavierStokesExample::~NavierStokesExample
~NavierStokesExample()
Definition: navier_stokes.cpp:27
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::ProjectionFieldOn10NodeTet
Definition: Projection10NodeCoordsOnField.hpp:49
NavierStokesExample::defineFiniteElements
MoFEMErrorCode defineFiniteElements()
[Setup fields]
Definition: navier_stokes.cpp:334
MoFEM::createDMVector
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1099
NavierStokesExample::step
int step
Definition: navier_stokes.cpp:56
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MetaNeumannForces::setMomentumFluxOperators
static MoFEMErrorCode setMomentumFluxOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, NeumannForcesSurface > &neumann_forces, Vec F, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Set operators to finite elements calculating right hand side vector.
Definition: SurfacePressure.cpp:2069
NavierStokesExample::viscosity
double viscosity
Definition: navier_stokes.cpp:50
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
bit
auto bit
set bit
Definition: hanging_node_approx.cpp:75
OpPPMap
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Definition: photon_diffusion.cpp:29
DM_NO_ELEMENT
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:10
NavierStokesExample::setupElementInstances
MoFEMErrorCode setupElementInstances()
[Setup algebraic structures]
Definition: navier_stokes.cpp:404
EshelbianPlasticity::P
@ P
Definition: EshelbianContact.cpp:197
NavierStokesExample::feLhsPtr
boost::shared_ptr< VolumeElementForcesAndSourcesCore > feLhsPtr
Definition: navier_stokes.cpp:70
MoFEM::DMoFEMMeshToGlobalVector
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMoFEM.cpp:535
NavierStokesExample::setupFields
MoFEMErrorCode setupFields()
[Find blocks]
Definition: navier_stokes.cpp:239
MoFEM::OpCalculateScalarFieldValues
Get value at integration points for scalar field.
Definition: UserDataOperators.hpp:82
MoFEM::DMMoFEMGetSnesCtx
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMoFEM.cpp:1094
MoFEM::ProjectionFieldOn10NodeTet::setNodes
bool setNodes
Definition: Projection10NodeCoordsOnField.hpp:61
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
MoFEM::PostProcBrokenMeshInMoab< VolumeElementForcesAndSourcesCore >
Definition: PostProcBrokenMeshInMoabBase.hpp:670
MoFEM::FaceElementForcesAndSourcesCore
Face finite element.
Definition: FaceElementForcesAndSourcesCore.hpp:23
NavierStokesExample::bitLevel
BitRefLevel bitLevel
Definition: navier_stokes.cpp:59
NavierStokesExample::D
SmartPetscObj< Vec > D
Definition: navier_stokes.cpp:66
MoFEM::DMMoFEMSNESSetJacobian
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMoFEM.cpp:759
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
MoFEM::DMoFEMPreProcessFiniteElements
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMoFEM.cpp:546
NavierStokesElement::LoadScale::lambda
static double lambda
Definition: NavierStokesElement.hpp:104
MoFEM::DMMoFEMCreateMoFEM
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMoFEM.cpp:114
NavierStokesExample::neumannForces
boost::ptr_map< std::string, NeumannForcesSurface > neumannForces
Definition: navier_stokes.cpp:77
MoFEM::VolumeElementForcesAndSourcesCore
Volume finite element base.
Definition: VolumeElementForcesAndSourcesCore.hpp:26
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
NavierStokesExample::readInput
MoFEMErrorCode readInput()
[Run problem]
Definition: navier_stokes.cpp:112
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:22
NavierStokesElement::addElement
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.
Definition: NavierStokesElement.hpp:129
MoFEM::OpCalculateVectorFieldGradient
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Definition: UserDataOperators.hpp:1535
NavierStokesElement::LoadScale
Definition: NavierStokesElement.hpp:102
NavierStokesElement::setCalcDragOperators
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.
Definition: NavierStokesElement.cpp:150
NavierStokesExample::snes
SNES snes
Definition: navier_stokes.cpp:62
MoFEM::SnesRhs
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Definition: SnesCtx.cpp:27
Range
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
NavierStokesExample::setupAlgebraicStructures
MoFEMErrorCode setupAlgebraicStructures()
[Setup discrete manager]
Definition: navier_stokes.cpp:381
MoFEM::vectorDuplicate
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
Definition: PetscSmartObj.hpp:221
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::CommInterface
Managing BitRefLevels.
Definition: CommInterface.hpp:21
_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_
#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.
Definition: MeshsetsManager.hpp:71
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
std
Definition: enable_if.hpp:5
NavierStokesExample::setupDiscreteManager
MoFEMErrorCode setupDiscreteManager()
[Define finite elements]
Definition: navier_stokes.cpp:357
BLOCKSET
@ BLOCKSET
Definition: definitions.h:161
NavierStokesExample::solveProblem
MoFEMErrorCode solveProblem()
[Setup SNES]
Definition: navier_stokes.cpp:556
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
lambda
static double lambda
Definition: incompressible_elasticity.cpp:199
NavierStokesExample::postProcess
MoFEMErrorCode postProcess()
[Solve problem]
Definition: navier_stokes.cpp:649
NavierStokesExample::fePostProcPtr
boost::shared_ptr< PostProcVol > fePostProcPtr
Definition: navier_stokes.cpp:79
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
NavierStokesExample::feRhsPtr
boost::shared_ptr< VolumeElementForcesAndSourcesCore > feRhsPtr
Definition: navier_stokes.cpp:71
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
numHoLevels
int numHoLevels
Definition: initial_diffusion.cpp:54
ReactionDiffusionEquation::D
const double D
diffusivity
Definition: reaction_diffusion.cpp:20
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MoFEM::SnesCtx
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:13
MoFEM::VolumeElementForcesAndSourcesCoreOnSide
Base volume element used to integrate on skeleton.
Definition: VolumeElementForcesAndSourcesCoreOnSide.hpp:22
MoFEM::SmartPetscObj< Vec >
MoFEM::SnesMat
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
Definition: SnesCtx.cpp:139
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:586
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::FieldBlas
Basic algebra on fields.
Definition: FieldBlas.hpp:21
EshelbianPlasticity::U
@ U
Definition: EshelbianContact.cpp:197
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
NavierStokesExample::commonData
boost::shared_ptr< NavierStokesElement::CommonData > commonData
Definition: navier_stokes.cpp:64
NavierStokesElement::setNavierStokesOperators
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.
Definition: NavierStokesElement.cpp:22
NavierStokesElement::setStokesOperators
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.
Definition: NavierStokesElement.cpp:79
F
@ F
Definition: free_surface.cpp:394
MoFEM::PostProcBrokenMeshInMoab< FaceElementForcesAndSourcesCore >
Definition: PostProcBrokenMeshInMoabBase.hpp:677
MoFEM::DMMoFEMSetIsPartitioned
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMoFEM.cpp:1123
MoFEM::DMMoFEMSNESSetFunction
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMoFEM.cpp:718
MoFEM::OpPostProcMapInMoab
Post post-proc data at points from hash maps.
Definition: PostProcBrokenMeshInMoabBase.hpp:698
NavierStokesExample::setupSNES
MoFEMErrorCode setupSNES()
[Setup element instances]
Definition: navier_stokes.cpp:520
Generated by Doxygen 1.8.17 and hosted at University of Glasgow