v0.8.20
Classes | Typedefs | Functions | Variables
#include <BasicFiniteElements.hpp>
#include <boost/program_options.hpp>
#include <Hooke.hpp>

Go to the source code of this file.

Classes

struct  BlockOptionData
 stores some optional data More...
 
struct  VolRule
 Set integration rule. More...
 

Typedefs

using BlockData = NonlinearElasticElement::BlockData
 
using VolUserDataOperator = VolumeElementForcesAndSourcesCore::UserDataOperator
 

Functions

int main (int argc, char *argv[])
 

Variables

static char help []
 

Typedef Documentation

◆ BlockData

◆ VolUserDataOperator

Definition at line 88 of file elasticity.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)
Examples:
elasticity.cpp.

Definition at line 95 of file elasticity.cpp.

95  {
96 
97  // Initialize PETSCc
98  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
99 
100  try {
101 
102  PetscBool flg_block_config, flg_file;
103  char mesh_file_name[255];
104  char block_config_file[255];
105  PetscInt order = 2;
106  PetscBool is_partitioned = PETSC_FALSE;
107 
108  // Read options from command line
109  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Elastic Config", "none");
110  CHKERR(ierr);
111  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
112  mesh_file_name, 255, &flg_file);
113 
114  CHKERR PetscOptionsInt("-my_order", "default approximation order", "",
115  order, &order, PETSC_NULL);
116 
117  CHKERR PetscOptionsBool("-my_is_partitioned",
118  "set if mesh is partitioned (this result that each "
119  "process keeps only one part of the mesh)",
120  "", is_partitioned, &is_partitioned, PETSC_NULL);
121 
122  CHKERR PetscOptionsString("-my_block_config", "elastic configure file name",
123  "", "block_conf.in", block_config_file, 255,
124  &flg_block_config);
125 
126  ierr = PetscOptionsEnd();
127  CHKERRG(ierr);
128 
129  // Throw error if file with mesh is not provided
130  if (flg_file != PETSC_TRUE) {
131  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
132  }
133 
134  // Create mesh database
135  moab::Core mb_instance;
136  moab::Interface &moab = mb_instance;
137 
138  // Create moab communicator
139  // Create separate MOAB communicator, it is duplicate of PETSc communicator.
140  // NOTE That this should eliminate potential communication problems between
141  // MOAB and PETSC functions.
142  MPI_Comm moab_comm_world;
143  MPI_Comm_dup(PETSC_COMM_WORLD, &moab_comm_world);
144  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
145  if (pcomm == NULL)
146  pcomm = new ParallelComm(&moab, moab_comm_world);
147 
148  // Read whole mesh or part of it if partitioned
149  if (is_partitioned == PETSC_TRUE) {
150  // This is a case of distributed mesh and algebra. In this case each
151  // processor keeps only one part of the problem.
152  const char *option;
153  option = "PARALLEL=READ_PART;"
154  "PARALLEL_RESOLVE_SHARED_ENTS;"
155  "PARTITION=PARALLEL_PARTITION;";
156  CHKERR moab.load_file(mesh_file_name, 0, option);
157  } else {
158  // If that case we have distributed algebra, i.e. assembly of vectors and
159  // matrices is in parallel, but whole mesh is stored on all processors.
160  // Solver and matrix scales well, however problem set-up of problem is
161  // not fully parallel.
162  const char *option;
163  option = "";
164  CHKERR moab.load_file(mesh_file_name, 0, option);
165  }
166 
167  // Create MoFEM database and link it to MoAB
168  MoFEM::Core core(moab);
169  MoFEM::Interface &m_field = core;
170 
171  // Print boundary conditions and material parameters
172  MeshsetsManager *meshsets_mng_ptr;
173  CHKERR m_field.getInterface(meshsets_mng_ptr);
174  CHKERR meshsets_mng_ptr->printDisplacementSet();
175  CHKERR meshsets_mng_ptr->printForceSet();
176  CHKERR meshsets_mng_ptr->printMaterialsSet();
177 
178  // Set bit refinement level to all entities (we do not refine mesh in this
179  // example so all entities have the same bit refinement level)
180  BitRefLevel bit_level0;
181  bit_level0.set(0);
182  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
183  0, 3, bit_level0);
184 
185  // Declare approximation fields
186  CHKERR m_field.add_field("DISPLACEMENT", H1, AINSWORTH_LOBATTO_BASE, 3,
187  MB_TAG_DENSE, MF_ZERO);
188 
189  // We can use higher oder geometry to define body
190  CHKERR m_field.add_field("MESH_NODE_POSITIONS", H1, AINSWORTH_LEGENDRE_BASE,
191  3, MB_TAG_DENSE, MF_ZERO);
192 
193  // Declare problem
194 
195  // Add entities (by tets) to the field (all entities in the mesh, root_set
196  // = 0 )
197  CHKERR m_field.add_ents_to_field_by_type(0, MBTET, "DISPLACEMENT");
198  CHKERR m_field.add_ents_to_field_by_type(0, MBTET, "MESH_NODE_POSITIONS");
199 
200  // Set approximation order.
201  // See Hierarchical Finite Element Bases on Unstructured Tetrahedral Meshes.
202  CHKERR m_field.set_field_order(0, MBTET, "DISPLACEMENT", order);
203  CHKERR m_field.set_field_order(0, MBTRI, "DISPLACEMENT", order);
204  CHKERR m_field.set_field_order(0, MBEDGE, "DISPLACEMENT", order);
205  CHKERR m_field.set_field_order(0, MBVERTEX, "DISPLACEMENT", 1);
206 
207  // Set order of approximation of geometry.
208  // Apply 2nd order only on skin (or in whole body)
209  auto setting_second_order_geometry = [&m_field]() {
211  // Setting geometry order everywhere
212  Range tets, edges;
213  CHKERR m_field.get_moab().get_entities_by_type(0, MBTET, tets);
214  CHKERR m_field.get_moab().get_adjacencies(tets, 1, false, edges,
215  moab::Interface::UNION);
216 
217  // Setting 2nd geometry order only on skin
218  // Range tets, faces, edges;
219  // Skinner skin(&m_field.get_moab());
220  // CHKERR skin.find_skin(0,tets,false,faces);
221  // CHKERR m_field.get_moab().get_adjacencies(
222  // faces,1,false,edges,moab::Interface::UNION
223  // );
224  // CHKERR m_field.synchronise_entities(edges);
225 
226  CHKERR m_field.set_field_order(edges, "MESH_NODE_POSITIONS", 2);
227  CHKERR m_field.set_field_order(0, MBVERTEX, "MESH_NODE_POSITIONS", 1);
228 
230  };
231  CHKERR setting_second_order_geometry();
232 
233  // Configure blocks by parsing config file. It allows setting approximation
234  // order for each block independently.
235  std::map<int, BlockOptionData> block_data;
236  auto setting_blocks_data_and_order_from_config_file =
237  [&m_field, &moab, &block_data, flg_block_config, block_config_file,
238  order](boost::shared_ptr<std::map<int, BlockData>> &block_sets_ptr) {
240  if (flg_block_config) {
241  ifstream ini_file(block_config_file);
242  po::variables_map vm;
243  po::options_description config_file_options;
245  it)) {
246  std::ostringstream str_order;
247  str_order << "block_" << it->getMeshsetId()
248  << ".displacement_order";
249  config_file_options.add_options()(
250  str_order.str().c_str(),
251  po::value<int>(&block_data[it->getMeshsetId()].oRder)
252  ->default_value(order));
253 
254  std::ostringstream str_cond;
255  str_cond << "block_" << it->getMeshsetId() << ".young_modulus";
256  config_file_options.add_options()(
257  str_cond.str().c_str(),
258  po::value<double>(&block_data[it->getMeshsetId()].yOung)
259  ->default_value(-1));
260 
261  std::ostringstream str_capa;
262  str_capa << "block_" << it->getMeshsetId() << ".poisson_ratio";
263  config_file_options.add_options()(
264  str_capa.str().c_str(),
265  po::value<double>(&block_data[it->getMeshsetId()].pOisson)
266  ->default_value(-2));
267 
268  std::ostringstream str_init_temp;
269  str_init_temp << "block_" << it->getMeshsetId()
270  << ".initial_temperature";
271  config_file_options.add_options()(
272  str_init_temp.str().c_str(),
273  po::value<double>(&block_data[it->getMeshsetId()].initTemp)
274  ->default_value(0));
275  }
276  po::parsed_options parsed =
277  parse_config_file(ini_file, config_file_options, true);
278  store(parsed, vm);
279  po::notify(vm);
281  it)) {
282  if (block_data[it->getMeshsetId()].oRder == -1)
283  continue;
284  if (block_data[it->getMeshsetId()].oRder == order)
285  continue;
286  PetscPrintf(PETSC_COMM_WORLD, "Set block %d order to %d\n",
287  it->getMeshsetId(),
288  block_data[it->getMeshsetId()].oRder);
289  Range block_ents;
290  CHKERR moab.get_entities_by_handle(it->getMeshset(), block_ents,
291  true);
292  Range ents_to_set_order;
293  CHKERR moab.get_adjacencies(block_ents, 3, false,
294  ents_to_set_order,
295  moab::Interface::UNION);
296  ents_to_set_order = ents_to_set_order.subset_by_type(MBTET);
297  CHKERR moab.get_adjacencies(block_ents, 2, false,
298  ents_to_set_order,
299  moab::Interface::UNION);
300  CHKERR moab.get_adjacencies(block_ents, 1, false,
301  ents_to_set_order,
302  moab::Interface::UNION);
303  CHKERR m_field.synchronise_entities(ents_to_set_order);
304 
305  CHKERR m_field.set_field_order(
306  ents_to_set_order, "DISPLACEMENT",
307  block_data[it->getMeshsetId()].oRder);
308  }
309  std::vector<std::string> additional_parameters;
310  additional_parameters =
311  collect_unrecognized(parsed.options, po::include_positional);
312  for (std::vector<std::string>::iterator vit =
313  additional_parameters.begin();
314  vit != additional_parameters.end(); vit++) {
315  CHKERR PetscPrintf(PETSC_COMM_WORLD,
316  "** WARNING Unrecognized option %s\n",
317  vit->c_str());
318  }
319  }
320 
321  // Update material parameters. Set material parameters block by block.
323  m_field, BLOCKSET | MAT_ELASTICSET, it)) {
324  const int id = it->getMeshsetId();
325  auto &bd = (*block_sets_ptr)[id];
326  if (block_data[id].yOung > 0)
327  bd.E = block_data[id].yOung;
328  if (block_data[id].pOisson >= -1)
329  bd.PoissonRatio = block_data[id].pOisson;
330  CHKERR PetscPrintf(PETSC_COMM_WORLD, "Block %d\n", id);
331  CHKERR PetscPrintf(PETSC_COMM_WORLD, "\tYoung modulus %3.4g\n",
332  bd.E);
333  CHKERR PetscPrintf(PETSC_COMM_WORLD, "\tPoisson ratio %3.4g\n",
334  bd.PoissonRatio);
335  }
336 
338  };
339 
340  // Add elastic element
341 
342  boost::shared_ptr<std::map<int, BlockData>> block_sets_ptr =
343  boost::make_shared<std::map<int, BlockData>>();
344  CHKERR HookeElement::setBlocks(m_field, block_sets_ptr);
345  CHKERR setting_blocks_data_and_order_from_config_file(block_sets_ptr);
346 
347  boost::shared_ptr<ForcesAndSourcesCore> fe_lhs_ptr(
348  new VolumeElementForcesAndSourcesCore(m_field));
349  boost::shared_ptr<ForcesAndSourcesCore> fe_rhs_ptr(
350  new VolumeElementForcesAndSourcesCore(m_field));
351  fe_lhs_ptr->getRuleHook = VolRule();
352  fe_rhs_ptr->getRuleHook = VolRule();
353 
354  CHKERR HookeElement::addElasticElement(m_field, block_sets_ptr, "ELASTIC",
355  "DISPLACEMENT",
356  "MESH_NODE_POSITIONS", false);
357  CHKERR HookeElement::setOperators(fe_lhs_ptr, fe_rhs_ptr, block_sets_ptr,
358  "DISPLACEMENT", "MESH_NODE_POSITIONS",
359  false, true);
360 
361  // Add spring boundary condition applied on surfaces.
362  // This is only declaration not implementation.
363  CHKERR MetaSpringBC::addSpringElements(m_field, "DISPLACEMENT",
364  "MESH_NODE_POSITIONS");
365 
366  // Implementation of spring element
367  // Create new instances of face elements for springs
368  boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr(
369  new FaceElementForcesAndSourcesCore(m_field));
370  boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr(
371  new FaceElementForcesAndSourcesCore(m_field));
372 
373  CHKERR MetaSpringBC::setSpringOperators(m_field, fe_spring_lhs_ptr,
374  fe_spring_rhs_ptr, "DISPLACEMENT",
375  "MESH_NODE_POSITIONS");
376 
377  // Add body force element. This is only declaration of element. not its
378  // implementation.
379  CHKERR m_field.add_finite_element("BODY_FORCE");
380  CHKERR m_field.modify_finite_element_add_field_row("BODY_FORCE",
381  "DISPLACEMENT");
382  CHKERR m_field.modify_finite_element_add_field_col("BODY_FORCE",
383  "DISPLACEMENT");
384  CHKERR m_field.modify_finite_element_add_field_data("BODY_FORCE",
385  "DISPLACEMENT");
386  CHKERR m_field.modify_finite_element_add_field_data("BODY_FORCE",
387  "MESH_NODE_POSITIONS");
388 
390  m_field, BLOCKSET | BODYFORCESSET, it)) {
391  Range tets;
392  CHKERR m_field.get_moab().get_entities_by_type(it->meshset, MBTET, tets,
393  true);
394  CHKERR m_field.add_ents_to_finite_element_by_type(tets, MBTET, "BODY_FORCE");
395  }
396 
397  // Add Neumann forces, i.e. pressure or traction forces applied on body
398  // surface. This is only declaration not implementation.
399  CHKERR MetaNeummanForces::addNeumannBCElements(m_field, "DISPLACEMENT");
400  CHKERR MetaNodalForces::addElement(m_field, "DISPLACEMENT");
401  CHKERR MetaEdgeForces::addElement(m_field, "DISPLACEMENT");
402 
403  // Add fluid pressure finite elements. This is special pressure on the
404  // surface from fluid, i.e. pressure which linearly change with the depth.
405  FluidPressure fluid_pressure_fe(m_field);
406  // This function only declare element. Element is implemented by operators
407  // in class FluidPressure.
408  fluid_pressure_fe.addNeumannFluidPressureBCElements("DISPLACEMENT");
409 
410  // Add elements for thermo elasticity if temperature field is defined.
411  ThermalStressElement thermal_stress_elem(m_field);
412  // Check if TEMP field exist, and then add element.
413  if (!m_field.check_field("TEMP")) {
414  bool add_temp_field = false;
416  if (block_data[it->getMeshsetId()].initTemp != 0) {
417  add_temp_field = true;
418  break;
419  }
420  }
421  if (add_temp_field) {
422  CHKERR m_field.add_field("TEMP", H1, AINSWORTH_LEGENDRE_BASE, 1,
423  MB_TAG_SPARSE, MF_ZERO);
424 
425  CHKERR m_field.add_ents_to_field_by_type(0, MBTET, "TEMP");
426 
427  CHKERR m_field.set_field_order(0, MBVERTEX, "TEMP", 1);
428  }
429  }
430  if (m_field.check_field("TEMP")) {
431  CHKERR thermal_stress_elem.addThermalStressElement(
432  "ELASTIC", "DISPLACEMENT", "TEMP");
433  }
434 
435  // All is declared, at this point build fields first,
436  CHKERR m_field.build_fields();
437  // If 10-node test are on the mesh, use mid nodes to set HO-geometry. Class
438  // Projection10NodeCoordsOnField
439  // do the trick.
440  Projection10NodeCoordsOnField ent_method_material(m_field,
441  "MESH_NODE_POSITIONS");
442  CHKERR m_field.loop_dofs("MESH_NODE_POSITIONS", ent_method_material);
443  if (m_field.check_field("TEMP")) {
445  if (block_data[it->getMeshsetId()].initTemp != 0) {
446  PetscPrintf(PETSC_COMM_WORLD, "Set block %d temperature to %3.2g\n",
447  it->getMeshsetId(),
448  block_data[it->getMeshsetId()].initTemp);
449  Range block_ents;
450  CHKERR moab.get_entities_by_handle(it->meshset, block_ents, true);
451  Range vertices;
452  CHKERR moab.get_connectivity(block_ents, vertices, true);
453  CHKERR m_field.getInterface<FieldBlas>()->setField(
454  block_data[it->getMeshsetId()].initTemp, MBVERTEX, vertices,
455  "TEMP");
456  }
457  }
458  }
459 
460  // Build database for elements. Actual implementation of element is not need
461  // here, only elements has to be declared.
462  CHKERR m_field.build_finite_elements();
463  // Build adjacencies between elements and field entities
464  CHKERR m_field.build_adjacencies(bit_level0);
465 
466  // Register MOFEM DM implementation in PETSc
468 
469  // Create DM manager
470  DM dm;
471  CHKERR DMCreate(PETSC_COMM_WORLD, &dm);
472  CHKERR DMSetType(dm, "MOFEM");
473  CHKERR DMMoFEMCreateMoFEM(dm, &m_field, "ELASTIC_PROB", bit_level0);
474  CHKERR DMSetFromOptions(dm);
475  CHKERR DMMoFEMSetIsPartitioned(dm, is_partitioned);
476  // Add elements to DM manager
477  CHKERR DMMoFEMAddElement(dm, "ELASTIC");
478  CHKERR DMMoFEMAddElement(dm, "SPRING");
479  CHKERR DMMoFEMAddElement(dm, "BODY_FORCE");
480  CHKERR DMMoFEMAddElement(dm, "FLUID_PRESSURE_FE");
481  CHKERR DMMoFEMAddElement(dm, "FORCE_FE");
482  CHKERR DMMoFEMAddElement(dm, "PRESSURE_FE");
483  CHKERR DMSetUp(dm);
484 
485  // Create matrices & vectors. Note that native PETSc DM interface is used,
486  // but under the PETSc interface MOFEM implementation is running.
487  Vec F, D, D0;
488  CHKERR DMCreateGlobalVector(dm, &F);
489  CHKERR VecDuplicate(F, &D);
490  CHKERR VecDuplicate(F, &D0);
491  Mat Aij;
492  CHKERR DMCreateMatrix(dm, &Aij);
493  CHKERR MatSetOption(Aij, MAT_SPD, PETSC_TRUE);
494 
495  // Assign global matrix/vector contributed by springs
496  fe_spring_lhs_ptr->ksp_B = Aij;
497  fe_spring_rhs_ptr->ksp_f = F;
498 
499  // Zero vectors and matrices
500  CHKERR VecZeroEntries(F);
501  CHKERR VecGhostUpdateBegin(F, INSERT_VALUES, SCATTER_FORWARD);
502  CHKERR VecGhostUpdateEnd(F, INSERT_VALUES, SCATTER_FORWARD);
503  CHKERR VecZeroEntries(D);
504  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
505  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
506  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
507  CHKERR MatZeroEntries(Aij);
508 
509  // This controls how kinematic constrains are set, by blockset or nodeset.
510  // Cubit sets kinetic boundary conditions by blockset.
511  bool flag_cubit_disp = false;
513  m_field, NODESET | DISPLACEMENTSET, it)) {
514  flag_cubit_disp = true;
515  }
516 
517  // Below particular implementations of finite elements are used to assemble
518  // problem matrixes and vectors. Implementation of element does not change
519  // how element is declared.
520 
521  // Assemble Aij and F. Define Dirichlet bc element, which stets constrains
522  // to MatrixDouble and the right hand side vector.
523  boost::shared_ptr<FEMethod> dirichlet_bc_ptr;
524 
525  // if normally defined boundary conditions are not found, try to use
526  // DISPLACEMENT blockset. To implementations available here, depending how
527  // BC is defined on mesh file.
528  if (!flag_cubit_disp) {
529  dirichlet_bc_ptr =
530  boost::shared_ptr<FEMethod>(new DirichletSetFieldFromBlockWithFlags(
531  m_field, "DISPLACEMENT", "DISPLACEMENT", Aij, D0, F));
532  } else {
533  dirichlet_bc_ptr = boost::shared_ptr<FEMethod>(
534  new DirichletDisplacementBc(m_field, "DISPLACEMENT", Aij, D0, F));
535  }
536  // That sets Dirichlet bc objects that problem is linear, i.e. no newton
537  // (SNES) solver is run for this problem.
538  dirichlet_bc_ptr->snes_ctx = FEMethod::CTX_SNESNONE;
539  dirichlet_bc_ptr->ts_ctx = FEMethod::CTX_TSNONE;
540 
541  // D0 vector will store initial displacements
542  CHKERR VecZeroEntries(D0);
543  CHKERR VecGhostUpdateBegin(D0, INSERT_VALUES, SCATTER_FORWARD);
544  CHKERR VecGhostUpdateEnd(D0, INSERT_VALUES, SCATTER_FORWARD);
545  CHKERR DMoFEMMeshToLocalVector(dm, D0, INSERT_VALUES, SCATTER_REVERSE);
546  // Run dirichlet_bc, from that on the mesh set values in vector D0. Run
547  // implementation
548  // of particular dirichlet_bc.
549  CHKERR DMoFEMPreProcessFiniteElements(dm, dirichlet_bc_ptr.get());
550  // Set values from D0 on the field (on the mesh)
551 
552  CHKERR VecGhostUpdateBegin(D0, INSERT_VALUES, SCATTER_FORWARD);
553  CHKERR VecGhostUpdateEnd(D0, INSERT_VALUES, SCATTER_FORWARD);
554  CHKERR DMoFEMMeshToLocalVector(dm, D0, INSERT_VALUES, SCATTER_REVERSE);
555 
556  // Calculate residual forces as result of applied kinematic constrains. Run
557  // implementation
558  // of particular finite element implementation. Look how
559  // NonlinearElasticElement is implemented,
560  // in that case. We run NonlinearElasticElement with hook material.
561  // Calculate right hand side vector
562  fe_rhs_ptr->snes_f = F;
563  PetscPrintf(PETSC_COMM_WORLD, "Assemble external force vector ...");
564  CHKERR DMoFEMLoopFiniteElements(dm, "ELASTIC", fe_rhs_ptr);
565  PetscPrintf(PETSC_COMM_WORLD, " done\n");
566  // Assemble matrix
567  fe_lhs_ptr->snes_B = Aij;
568  PetscPrintf(PETSC_COMM_WORLD, "Calculate stiffness matrix ...");
569  CHKERR DMoFEMLoopFiniteElements(dm, "ELASTIC", fe_lhs_ptr);
570  PetscPrintf(PETSC_COMM_WORLD, " done\n");
571 
572  // Assemble springs
573  CHKERR DMoFEMLoopFiniteElements(dm, "SPRING", fe_spring_lhs_ptr);
574  CHKERR DMoFEMLoopFiniteElements(dm, "SPRING", fe_spring_rhs_ptr);
575 
576  // Assemble pressure and traction forces. Run particular implemented for do
577  // this, see
578  // MetaNeummanForces how this is implemented.
579  boost::ptr_map<std::string, NeummanForcesSurface> neumann_forces;
580  CHKERR MetaNeummanForces::setMomentumFluxOperators(m_field, neumann_forces,
581  F, "DISPLACEMENT");
582 
583  {
584  boost::ptr_map<std::string, NeummanForcesSurface>::iterator mit =
585  neumann_forces.begin();
586  for (; mit != neumann_forces.end(); mit++) {
587  CHKERR DMoFEMLoopFiniteElements(dm, mit->first.c_str(),
588  &mit->second->getLoopFe());
589  }
590  }
591  // Assemble forces applied to nodes, see implementation in NodalForce
592  boost::ptr_map<std::string, NodalForce> nodal_forces;
593  CHKERR MetaNodalForces::setOperators(m_field, nodal_forces, F,
594  "DISPLACEMENT");
595 
596  {
597  boost::ptr_map<std::string, NodalForce>::iterator fit =
598  nodal_forces.begin();
599  for (; fit != nodal_forces.end(); fit++) {
600  CHKERR DMoFEMLoopFiniteElements(dm, fit->first.c_str(),
601  &fit->second->getLoopFe());
602  }
603  }
604  // Assemble edge forces
605  boost::ptr_map<std::string, EdgeForce> edge_forces;
606  CHKERR MetaEdgeForces::setOperators(m_field, edge_forces, F,
607  "DISPLACEMENT");
608  {
609  boost::ptr_map<std::string, EdgeForce>::iterator fit =
610  edge_forces.begin();
611  for (; fit != edge_forces.end(); fit++) {
612  CHKERR DMoFEMLoopFiniteElements(dm, fit->first.c_str(),
613  &fit->second->getLoopFe());
614  }
615  }
616  // Assemble body forces, implemented in BodyForceConstantField
617  BodyForceConstantField body_forces_methods(m_field);
619  m_field, BLOCKSET | BODYFORCESSET, it)) {
620  CHKERR body_forces_methods.addBlock("DISPLACEMENT", F,
621  it->getMeshsetId());
622  }
623  CHKERR DMoFEMLoopFiniteElements(dm, "BODY_FORCE",
624  &body_forces_methods.getLoopFe());
625  // Assemble fluid pressure forces
626  CHKERR fluid_pressure_fe.setNeumannFluidPressureFiniteElementOperators(
627  "DISPLACEMENT", F, false, true);
628 
629  CHKERR DMoFEMLoopFiniteElements(dm, "FLUID_PRESSURE_FE",
630  &fluid_pressure_fe.getLoopFe());
631  // Apply kinematic constrain to right hand side vector and matrix
632  CHKERR DMoFEMPostProcessFiniteElements(dm, dirichlet_bc_ptr.get());
633 
634  // Matrix View
635  // MatView(Aij,PETSC_VIEWER_STDOUT_WORLD);
636  // MatView(Aij,PETSC_VIEWER_DRAW_WORLD);//PETSC_VIEWER_STDOUT_WORLD);
637  // std::string wait;
638  // std::cin >> wait;
639 
640  // Set matrix positive defined and symmetric for Cholesky and icc
641  // pre-conditioner
642  CHKERR MatSetOption(Aij, MAT_SPD, PETSC_TRUE);
643  CHKERR VecGhostUpdateBegin(F, ADD_VALUES, SCATTER_REVERSE);
644  CHKERR VecGhostUpdateEnd(F, ADD_VALUES, SCATTER_REVERSE);
645  CHKERR VecAssemblyBegin(F);
646  CHKERR VecAssemblyEnd(F);
647  CHKERR VecScale(F, -1);
648 
649  // Create solver
650  KSP solver;
651  CHKERR KSPCreate(PETSC_COMM_WORLD, &solver);
652  CHKERR KSPSetDM(solver, dm);
653  CHKERR KSPSetFromOptions(solver);
654  CHKERR KSPSetOperators(solver, Aij, Aij);
655  // Setup multi-grid pre-conditioner if set from command line
656  {
657  // from PETSc example ex42.c
658  PetscBool same = PETSC_FALSE;
659  PC pc;
660  CHKERR KSPGetPC(solver, &pc);
661  PetscObjectTypeCompare((PetscObject)pc, PCMG, &same);
662  if (same) {
663  PCMGSetUpViaApproxOrdersCtx pc_ctx(dm, Aij, true);
664  CHKERR PCMGSetUpViaApproxOrders(pc, &pc_ctx);
665  CHKERR PCSetFromOptions(pc);
666  } else {
667  // Operators are already set, do not use DM for doing that
668  CHKERR KSPSetDMActive(solver, PETSC_FALSE);
669  }
670  }
671  CHKERR KSPSetInitialGuessKnoll(solver, PETSC_FALSE);
672  CHKERR KSPSetInitialGuessNonzero(solver, PETSC_TRUE);
673  // Set up solver
674  CHKERR KSPSetUp(solver);
675 
676  // Set up post-processor. It is some generic implementation of finite
677  // element.
678  PostProcVolumeOnRefinedMesh post_proc(m_field);
679  // Add operators to the elements, starting with some generic operators
680  CHKERR post_proc.generateReferenceElementMesh();
681  CHKERR post_proc.addFieldValuesPostProc("DISPLACEMENT");
682  CHKERR post_proc.addFieldValuesPostProc("MESH_NODE_POSITIONS");
683  CHKERR post_proc.addFieldValuesGradientPostProc("DISPLACEMENT");
684  // Add problem specific operator on element to post-process stresses
685  post_proc.getOpPtrVector().push_back(new PostProcHookStress(
686  m_field, post_proc.postProcMesh, post_proc.mapGaussPts, "DISPLACEMENT",
687  post_proc.commonData, block_sets_ptr.get()));
688 
689  // Temperature field is defined on the mesh
690  if (m_field.check_field("TEMP")) {
691 
692  // Create thermal vector
693  Vec F_thermal;
694  CHKERR VecDuplicate(F, &F_thermal);
695 
696  // Set up implementation for calculation of thermal stress vector. Look
697  // how thermal stresses and vector is assembled in ThermalStressElement.
698  CHKERR thermal_stress_elem.setThermalStressRhsOperators(
699  "DISPLACEMENT", "TEMP", F_thermal);
700 
701  SeriesRecorder *recorder_ptr;
702  CHKERR m_field.getInterface(recorder_ptr);
703 
704  // Read time series and do thermo-elastic analysis, this is when time
705  // dependent
706  // temperature problem was run before on the mesh. It means that before
707  // non-stationary
708  // problem was solved for temperature and filed "TEMP" is stored for
709  // subsequent time
710  // steps in the recorder.
711  if (recorder_ptr->check_series("THEMP_SERIES")) {
712  // This is time dependent case, so loop of data series stored by tape
713  // recorder.
714  // Loop over time steps
715  for (_IT_SERIES_STEPS_BY_NAME_FOR_LOOP_(recorder_ptr, "THEMP_SERIES",
716  sit)) {
717  PetscPrintf(PETSC_COMM_WORLD, "Process step %d\n",
718  sit->get_step_number());
719  // Load field data for this time step
720  CHKERR recorder_ptr->load_series_data("THEMP_SERIES",
721  sit->get_step_number());
722 
723  CHKERR VecZeroEntries(F_thermal);
724  CHKERR VecGhostUpdateBegin(F_thermal, INSERT_VALUES, SCATTER_FORWARD);
725  CHKERR VecGhostUpdateEnd(F_thermal, INSERT_VALUES, SCATTER_FORWARD);
726 
727  // Calculate the right-hand side vector as result of thermal stresses.
728  // It MetaNodalForces
729  // that on "ELASTIC" element data structure the element implementation
730  // in thermal_stress_elem
731  // is executed.
733  dm, "ELASTIC", &thermal_stress_elem.getLoopThermalStressRhs());
734 
735  // Assemble vector
736  CHKERR VecAssemblyBegin(F_thermal);
737  CHKERR VecAssemblyEnd(F_thermal);
738  // Accumulate ghost dofs
739  CHKERR VecGhostUpdateBegin(F_thermal, ADD_VALUES, SCATTER_REVERSE);
740  CHKERR VecGhostUpdateEnd(F_thermal, ADD_VALUES, SCATTER_REVERSE);
741 
742  // Calculate norm of vector and print values
743  PetscReal nrm_F;
744  CHKERR VecNorm(F, NORM_2, &nrm_F);
745 
746  PetscPrintf(PETSC_COMM_WORLD, "norm2 F = %6.4e\n", nrm_F);
747  PetscReal nrm_F_thermal;
748  CHKERR VecNorm(F_thermal, NORM_2, &nrm_F_thermal);
749  PetscPrintf(PETSC_COMM_WORLD, "norm2 F_thermal = %6.4e\n",
750  nrm_F_thermal);
751 
752  CHKERR VecScale(F_thermal, -1);
753  // check this !!!
754  CHKERR VecAXPY(F_thermal, 1, F);
755 
756  // Set dirichlet boundary to thermal stresses vector
757  dirichlet_bc_ptr->snes_x = D;
758  dirichlet_bc_ptr->snes_f = F_thermal;
759  CHKERR DMoFEMPostProcessFiniteElements(dm, dirichlet_bc_ptr.get());
760 
761  // Solve problem
762  CHKERR KSPSolve(solver, F_thermal, D);
763  // Add boundary conditions vector
764  CHKERR VecAXPY(D, 1., D0);
765  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
766  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
767 
768  // Save data on the mesh
769  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
770 
771  // Save data on mesh
772  CHKERR DMoFEMPreProcessFiniteElements(dm, dirichlet_bc_ptr.get());
773  // Post-process results
774  CHKERR DMoFEMLoopFiniteElements(dm, "ELASTIC", &post_proc);
775 
776  std::ostringstream o1;
777  o1 << "out_" << sit->step_number << ".h5m";
778  CHKERR post_proc.writeFile(o1.str().c_str());
779  }
780  } else {
781 
782  // This is a case when stationary problem for temperature was solved.
783  CHKERR VecZeroEntries(F_thermal);
784  CHKERR VecGhostUpdateBegin(F_thermal, INSERT_VALUES, SCATTER_FORWARD);
785  CHKERR VecGhostUpdateEnd(F_thermal, INSERT_VALUES, SCATTER_FORWARD);
786 
787  // Calculate the right-hand side vector with thermal stresses
789  dm, "ELASTIC", &thermal_stress_elem.getLoopThermalStressRhs());
790 
791  // Assemble vector
792  CHKERR VecAssemblyBegin(F_thermal);
793  CHKERR VecAssemblyEnd(F_thermal);
794 
795  // Accumulate ghost dofs
796  CHKERR VecGhostUpdateBegin(F_thermal, ADD_VALUES, SCATTER_REVERSE);
797  CHKERR VecGhostUpdateEnd(F_thermal, ADD_VALUES, SCATTER_REVERSE);
798 
799  // Calculate norm
800  PetscReal nrm_F;
801  CHKERR VecNorm(F, NORM_2, &nrm_F);
802 
803  PetscPrintf(PETSC_COMM_WORLD, "norm2 F = %6.4e\n", nrm_F);
804  PetscReal nrm_F_thermal;
805  CHKERR VecNorm(F_thermal, NORM_2, &nrm_F_thermal);
806 
807  PetscPrintf(PETSC_COMM_WORLD, "norm2 F_thermal = %6.4e\n",
808  nrm_F_thermal);
809 
810  // Add thermal stress vector and other forces vector
811  CHKERR VecScale(F_thermal, -1);
812  CHKERR VecAXPY(F_thermal, 1, F);
813 
814  // Apply kinetic boundary conditions
815  dirichlet_bc_ptr->snes_x = D;
816  dirichlet_bc_ptr->snes_f = F_thermal;
817  CHKERR DMoFEMPostProcessFiniteElements(dm, dirichlet_bc_ptr.get());
818 
819  // Solve problem
820  CHKERR KSPSolve(solver, F_thermal, D);
821  CHKERR VecAXPY(D, 1., D0);
822 
823  // Update ghost values for solution vector
824  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
825  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
826 
827  // Save data on mesh
828  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
829  CHKERR DMoFEMLoopFiniteElements(dm, "ELASTIC", &post_proc);
830  // Save results to file
831  PetscPrintf(PETSC_COMM_WORLD, "Write output file ...");
832  CHKERR post_proc.writeFile("out.h5m");
833  PetscPrintf(PETSC_COMM_WORLD, " done\n");
834  }
835 
836  // Destroy vector, no needed any more
837  CHKERR VecDestroy(&F_thermal);
838 
839  } else {
840  // Elastic analysis no temperature field
841  // Solve for vector D
842  CHKERR KSPSolve(solver, F, D);
843  // Add kinetic boundary conditions
844  CHKERR VecAXPY(D, 1., D0);
845  // Update ghost values
846  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
847  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
848  // Save data from vector on mesh
849  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
850  // Post-process results
851  CHKERR DMoFEMLoopFiniteElements(dm, "ELASTIC", &post_proc);
852  // CHKERR DMoFEMLoopFiniteElements(dm, "SPRING", &post_proc);
853  // Write mesh in parallel (using h5m MOAB format, writing is in parallel)
854  PetscPrintf(PETSC_COMM_WORLD, "Write output file ..,");
855  CHKERR post_proc.writeFile("out.h5m");
856  PetscPrintf(PETSC_COMM_WORLD, " done\n");
857  }
858 
859  // Calculate elastic energy
860  auto calculate_strain_energy = [dm, &block_sets_ptr]() {
862 
863  Vec v_energy;
864  CHKERR HookeElement::calculateEnergy(dm, block_sets_ptr, "DISPLACEMENT",
865  "MESH_NODE_POSITIONS", false, true,
866  &v_energy);
867 
868  // Print elastic energy
869  const double *eng_ptr;
870  CHKERR VecGetArrayRead(v_energy, &eng_ptr);
871  PetscPrintf(PETSC_COMM_WORLD, "Elastic energy %6.4e\n", *eng_ptr);
872  CHKERR VecRestoreArrayRead(v_energy, &eng_ptr);
873 
874  CHKERR VecDestroy(&v_energy);
875 
877  };
878  CHKERR calculate_strain_energy();
879 
880  // Destroy matrices, vecors, solver and DM
881  CHKERR VecDestroy(&F);
882  CHKERR VecDestroy(&D);
883  CHKERR VecDestroy(&D0);
884  CHKERR MatDestroy(&Aij);
885  CHKERR KSPDestroy(&solver);
886  CHKERR DMDestroy(&dm);
887 
888  MPI_Comm_free(&moab_comm_world);
889  }
890  CATCH_ERRORS;
891 
893 
894  return 0;
895 }
Set integration rule.
static MoFEMErrorCode setMomentumFluxOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, NeummanForcesSurface > &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.
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
#define _IT_SERIES_STEPS_BY_NAME_FOR_LOOP_(RECORDER, NAME, IT)
loop over recorded series step
Fluid pressure forces.
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
Implentation of thermal stress element.
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
static MoFEMErrorCode addElasticElement(MoFEM::Interface &m_field, boost::shared_ptr< map< int, BlockData >> &block_sets_ptr, const std::string element_name, const std::string x_field, const std::string X_field, const bool ale)
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.
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode load_series_data(const std::string &serie_name, const int step_number)
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:414
static MoFEMErrorCode setBlocks(MoFEM::Interface &m_field, boost::shared_ptr< map< int, BlockData >> &block_sets_ptr)
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:394
static MoFEMErrorCode calculateEnergy(DM dm, boost::shared_ptr< map< int, BlockData >> &block_sets_ptr, const std::string x_field, const std::string X_field, const bool ale, const bool field_disp, Vec *v_energy_ptr)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
static MoFEMErrorCode setOperators(boost::shared_ptr< ForcesAndSourcesCore > &fe_lhs_ptr, boost::shared_ptr< ForcesAndSourcesCore > &fe_rhs_ptr, boost::shared_ptr< map< int, BlockData >> &block_sets_ptr, const std::string x_field, const std::string X_field, const bool ale, const bool field_disp)
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
virtual MoFEMErrorCode add_field(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
Interface for managing meshsets containing materials and boundary conditions.
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
static MoFEMErrorCode addElement(MoFEM::Interface &m_field, const std::string field_name, Range *intersect_ptr=NULL)
Add element taking information from NODESET.
Definition: EdgeForce.hpp:74
Core (interface) class.
Definition: Core.hpp:50
Volume finite elementUser is implementing own operator at Gauss point level, by class derived from Vo...
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:476
static char help[]
Definition: elasticity.cpp:72
char mesh_file_name[255]
static MoFEMErrorCode setSpringOperators(MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Implementation of spring element. Set operators to calculate LHS and RHS.
Projection of edge entities with one mid-node on hierarchical basis.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
Operator post-procesing stresses for Hook isotropic material.
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
static MoFEMErrorCode setOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, NodalForce > &nodal_forces, Vec F, const std::string field_name)
Set integration point operators.
Definition: NodalForce.hpp:140
virtual bool check_series(const std::string &name) const
check if series is in database
Body forces elements.
Definition: BodyForce.hpp:24
MoFEMErrorCode printForceSet() const
print meshsets with force boundary conditions data structure
Basic algebra on fields.
Definition: FieldBlas.hpp:34
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
Add boundary conditions form block set having 6 attributes.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
static MoFEMErrorCode addSpringElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Declare spring element.
Set Dirichlet boundary conditions on displacements.
Definition: DirichletBC.hpp:43
Managing BitRefLevels.
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:143
MoFEMErrorCode printDisplacementSet() const
print meshsets with displacement boundary conditions data structure
static MoFEMErrorCode addElement(MoFEM::Interface &m_field, const std::string field_name, Range *intersect_ptr=NULL)
Add element taking information from NODESET.
Definition: NodalForce.hpp:104
Post processing.
block name is "BODY_FORCES"
Definition: definitions.h:223
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
Set data structures of MG pre-conditioner via approximation orders.
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:437
MoFEMErrorCode printMaterialsSet() const
print meshsets with material data structure set on it
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:447
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: DMMMoFEM.cpp:109
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
Face finite elementUser is implementing own operator at Gauss point level, by own object derived from...
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#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...
#define CHKERR
Inline error check.
Definition: definitions.h:594
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
block name is "MAT_ELASTIC"
Definition: definitions.h:220
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:283
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
virtual MoFEMErrorCode synchronise_entities(Range &ent, int verb=DEFAULT_VERBOSITY)=0
static MoFEMErrorCode setOperators(MoFEM::Interface &m_field, boost::ptr_map< std::string, EdgeForce > &edge_forces, Vec F, const std::string field_name)
Set integration point operators.
Definition: EdgeForce.hpp:109
virtual bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:431
continuous field
Definition: definitions.h:169
MoFEMErrorCode PCMGSetUpViaApproxOrders(PC pc, PCMGSetUpViaApproxOrdersCtx *ctx, int verb)
Function build MG structure.
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:807
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elementsBuild finite element data structures. Have to be run before problem and adjacenc...
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61

Variable Documentation

◆ help

char help[]
static
Initial value:
= "-my_block_config set block data\n"
"\n"
Examples:
elasticity.cpp.

Definition at line 72 of file elasticity.cpp.