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