v0.10.0
Simple.cpp
Go to the documentation of this file.
1 /** \file Simple.cpp
2  * \brief Implementation of simple interface
3  * \ingroup mofem_simple_interface
4  */
5 
6 /* MoFEM is free software: you can redistribute it and/or modify it under
7  * the terms of the GNU Lesser General Public License as published by the
8  * Free Software Foundation, either version 3 of the License, or (at your
9  * option) any later version.
10  *
11  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14  * License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
18  */
19 
20 namespace MoFEM {
21 
23  UnknownInterface **iface) const {
25  *iface = NULL;
26  if (uuid == IDD_MOFEMSimple) {
27  *iface = const_cast<Simple *>(this);
29  }
30  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
32 }
33 
34 template<int DIM>
36  static_assert(DIM == 2 || DIM == 3, "not implemented");
37  return MOFEM_NOT_IMPLEMENTED;
38 }
39 
40 template<>
41 MoFEMErrorCode Simple::setSkeletonAdjacency<2>() {
42  Interface &m_field = cOre;
44 
45  auto defaultSkeletonEdge = [&](moab::Interface &moab, const Field &field,
46  const EntFiniteElement &fe,
47  Range &adjacency) -> MoFEMErrorCode {
49 
50  CHKERR DefaultElementAdjacency::defaultEdge(moab, field, fe, adjacency);
51 
52  if (std::find(domainFields.begin(), domainFields.end(),
53  field.getName()) != domainFields.end()) {
54 
55  const EntityHandle fe_ent = fe.getEnt();
56  Range bride_adjacency_edge, bride_adjacency;
57  CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, bride_adjacency_edge);
58 
59  switch (field.getSpace()) {
60  case H1:
61  CHKERR moab.get_connectivity(bride_adjacency_edge, bride_adjacency,
62  true);
63  case HCURL:
64  case HDIV:
65  CHKERR moab.get_adjacencies(bride_adjacency_edge, 1, false,
66  bride_adjacency, moab::Interface::UNION);
67  case L2:
68  bride_adjacency.merge(bride_adjacency_edge);
69  break;
70  case NOFIELD:
71  break;
72  default:
73  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
74  "this field is not implemented for TRI finite element");
75  }
76 
77  bride_adjacency = subtract(bride_adjacency, adjacency);
78 
79  for (auto e : bride_adjacency)
80  const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
81  .insert(boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
82 
83  adjacency.merge(bride_adjacency);
84  }
85 
87  };
88 
90  defaultSkeletonEdge);
91 
93 }
94 
95 template<>
96 MoFEMErrorCode Simple::setSkeletonAdjacency<3>() {
97  Interface &m_field = cOre;
99 
100  auto defaultSkeletonEdge = [&](moab::Interface &moab, const Field &field,
101  const EntFiniteElement &fe,
102  Range &adjacency) -> MoFEMErrorCode {
104 
105  CHKERR DefaultElementAdjacency::defaultFace(moab, field, fe, adjacency);
106 
107  if (std::find(domainFields.begin(), domainFields.end(),
108  field.getName()) != domainFields.end()) {
109 
110  const EntityHandle fe_ent = fe.getEnt();
111  Range bride_adjacency_face, bride_adjacency;
112  CHKERR moab.get_adjacencies(&fe_ent, 1, 3, false, bride_adjacency_face);
113 
114  switch (field.getSpace()) {
115  case H1:
116  CHKERR moab.get_connectivity(bride_adjacency_face, bride_adjacency,
117  true);
118  case HCURL:
119  CHKERR moab.get_adjacencies(bride_adjacency_face, 1, false,
120  bride_adjacency, moab::Interface::UNION);
121  case HDIV:
122  CHKERR moab.get_adjacencies(bride_adjacency_face, 2, false,
123  bride_adjacency, moab::Interface::UNION);
124  case L2:
125  bride_adjacency.merge(bride_adjacency_face);
126  break;
127  case NOFIELD:
128  break;
129  default:
130  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
131  "this field is not implemented for TRI finite element");
132  }
133 
134  bride_adjacency = subtract(bride_adjacency, adjacency);
135 
136  for (auto e : bride_adjacency)
137  const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
138  .insert(boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
139 
140  adjacency.merge(bride_adjacency);
141  }
142 
144  };
145 
147  defaultSkeletonEdge);
149  defaultSkeletonEdge);
150 
152 }
153 
154 template<>
157  switch (getDim()) {
158  case 1:
159  THROW_MESSAGE("Not implemented");
160  case 2:
161  return setSkeletonAdjacency<2>();
162  case 3:
163  return setSkeletonAdjacency<3>();
164  default:
165  THROW_MESSAGE("Not implemented");
166  }
168 }
169 
172  if(dim == -1)
173  dim = getDim();
174  switch(dim) {
175  case 2:
176  return setSkeletonAdjacency<2>();
177  case 3:
178  return setSkeletonAdjacency<3>();
179  default:
180  SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED, "Not implemented");
181  }
183 }
184 
185 Simple::Simple(const Core &core)
186  : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
187  meshSet(0), boundaryMeshset(0), nameOfProblem("SimpleProblem"),
188  domainFE("dFE"), boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1) {
189  PetscLogEventRegister("LoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
190  PetscLogEventRegister("buildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
191  PetscLogEventRegister("buildFiniteElements", 0,
193  PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
194  PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
195  strcpy(meshFileName, "mesh.h5m");
196 }
197 
199  PetscBool flg = PETSC_TRUE;
201  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
202  "none");
203  CHKERRG(ierr);
204  ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
205  meshFileName, 255, &flg);
206  CHKERRG(ierr);
207  ierr = PetscOptionsEnd();
208  CHKERRG(ierr);
210 }
211 
212 MoFEMErrorCode Simple::loadFile(const std::string options,
213  const std::string mesh_file_name) {
214  Interface &m_field = cOre;
216  PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
217 
218  if (!mesh_file_name.empty())
219  strcpy(meshFileName, mesh_file_name.c_str());
220 
221  // This is a case of distributed mesh and algebra. In that case each processor
222  // keep only part of the problem.
223  CHKERR m_field.get_moab().load_file(meshFileName, 0, options.c_str());
224  CHKERR m_field.rebuild_database();
225  // determine problem dimension
226  if (dIm == -1) {
227  int nb_ents_3d;
228  CHKERR m_field.get_moab().get_number_entities_by_dimension(
229  meshSet, 3, nb_ents_3d, true);
230  if (nb_ents_3d > 0) {
231  dIm = 3;
232  } else {
233  int nb_ents_2d;
234  CHKERR m_field.get_moab().get_number_entities_by_dimension(
235  meshSet, 2, nb_ents_2d, true);
236  if (nb_ents_2d > 0) {
237  dIm = 2;
238  } else {
239  dIm = 1;
240  }
241  }
242  }
243  Range ents;
244  CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents, true);
245  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
246  false);
247  ParallelComm *pcomm =
248  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
249  if (pcomm == NULL)
250  pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
251  PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
253 }
254 
257  Interface &m_field = cOre;
258  if (m_field.get_comm_size() == 1)
260  else
261  CHKERR loadFile("PARALLEL=READ_PART;"
262  "PARALLEL_RESOLVE_SHARED_ENTS;"
263  "PARTITION=PARALLEL_PARTITION;",
264  meshFileName);
266 }
267 
269 Simple::addDomainField(const std::string &name, const FieldSpace space,
270  const FieldApproximationBase base,
271  const FieldCoefficientsNumber nb_of_cooficients,
272  const TagType tag_type, const enum MoFEMTypes bh,
273  int verb) {
274 
275  Interface &m_field = cOre;
277  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
278  verb);
279  if (space == NOFIELD)
280  noFieldFields.push_back(name);
281  else
282  domainFields.push_back(name);
284 }
285 
287 Simple::addBoundaryField(const std::string &name, const FieldSpace space,
288  const FieldApproximationBase base,
289  const FieldCoefficientsNumber nb_of_cooficients,
290  const TagType tag_type, const enum MoFEMTypes bh,
291  int verb) {
292  Interface &m_field = cOre;
294  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
295  verb);
296  boundaryFields.push_back(name);
297  if (space == NOFIELD)
298  SETERRQ(
299  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
300  "NOFIELD space for boundary filed not implemented in Simple interface");
302 }
303 
305 Simple::addSkeletonField(const std::string &name, const FieldSpace space,
306  const FieldApproximationBase base,
307  const FieldCoefficientsNumber nb_of_cooficients,
308  const TagType tag_type, const enum MoFEMTypes bh,
309  int verb) {
310 
311  Interface &m_field = cOre;
313  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
314  verb);
315  skeletonFields.push_back(name);
316  if (space == NOFIELD)
317  SETERRQ(
318  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
319  "NOFIELD space for boundary filed not implemented in Simple interface");
321 }
322 
324 Simple::addDataField(const std::string &name, const FieldSpace space,
325  const FieldApproximationBase base,
326  const FieldCoefficientsNumber nb_of_cooficients,
327  const TagType tag_type, const enum MoFEMTypes bh,
328  int verb) {
329 
330  Interface &m_field = cOre;
332  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
333  verb);
334  if (space == NOFIELD)
335  noFieldDataFields.push_back(name);
336  else
337  dataFields.push_back(name);
339 }
340 
341 MoFEMErrorCode Simple::removeDomainField(const std::string &name) {
342  Interface &m_field = cOre;
344 
345  auto remove_field_from_list = [&](auto &vec) {
346  auto it = std::find(vec.begin(), vec.end(), name);
347  if (it != vec.end())
348  vec.erase(it);
349  };
350 
351  remove_field_from_list(noFieldFields);
352  remove_field_from_list(domainFields);
353 
355 }
356 
357 MoFEMErrorCode Simple::removeBoundaryField(const std::string &name) {
358  Interface &m_field = cOre;
360 
361  auto remove_field_from_list = [&](auto &vec) {
362  auto it = std::find(vec.begin(), vec.end(), name);
363  if (it != vec.end())
364  vec.erase(it);
365  };
366 
367  remove_field_from_list(boundaryFields);
368 
370 }
371 
372 MoFEMErrorCode Simple::removeSkeletonField(const std::string &name) {
373  Interface &m_field = cOre;
375 
376  auto remove_field_from_list = [&](auto &vec) {
377  auto it = std::find(vec.begin(), vec.end(), name);
378  if (it != vec.end())
379  vec.erase(it);
380  };
381 
382  remove_field_from_list(skeletonFields);
383 
385 }
386 
388  Interface &m_field = cOre;
390 
391  auto clear_rows_and_cols = [&](auto &fe_name) {
393  auto fe_ptr = m_field.get_finite_elements();
394  auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
395  ->get<FiniteElement_name_mi_tag>();
396  auto it_fe = fe_by_name.find(fe_name);
397  if (it_fe != fe_by_name.end()) {
398 
399  if(!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
400  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
401  "modification unsuccessful");
402 
403  if(!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
404  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
405  "modification unsuccessful");
406  }
408  };
409  CHKERR clear_rows_and_cols(domainFE);
410  CHKERR clear_rows_and_cols(boundaryFE);
411  CHKERR clear_rows_and_cols(skeletonFE);
412 
413  // Define finite elements
415 
416  auto add_fields = [&](auto &fe_name, auto &fields) {
418  for (auto &field : fields) {
419  CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
420  CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
421  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
422  }
424  };
425 
426  auto add_data_fields = [&](auto &fe_name, auto &fields) {
428  for (auto &field : fields)
429  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
431  };
432 
433  CHKERR add_fields(domainFE, domainFields);
434  CHKERR add_data_fields(domainFE, dataFields);
435  CHKERR add_fields(domainFE, noFieldFields);
436  CHKERR add_data_fields(domainFE, noFieldDataFields);
437 
438  if (!boundaryFields.empty()) {
440  CHKERR add_fields(boundaryFE, domainFields);
441  CHKERR add_fields(boundaryFE, boundaryFields);
442  CHKERR add_fields(boundaryFE, skeletonFields);
443  CHKERR add_data_fields(boundaryFE, dataFields);
444  CHKERR add_data_fields(boundaryFE, noFieldDataFields);
445  CHKERR add_fields(boundaryFE, noFieldFields);
446  }
447  if (!skeletonFields.empty()) {
449  CHKERR add_fields(skeletonFE, domainFields);
450  CHKERR add_fields(skeletonFE, boundaryFields);
451  CHKERR add_fields(skeletonFE, skeletonFields);
452  CHKERR add_data_fields(skeletonFE, dataFields);
453  CHKERR add_data_fields(skeletonFE, noFieldDataFields);
454  CHKERR add_fields(skeletonFE, noFieldFields);
455  }
457 }
458 
459 MoFEMErrorCode Simple::defineProblem(const PetscBool is_partitioned) {
460  Interface &m_field = cOre;
462  // Create dm instance
463  dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
464  // set dm data structure which created mofem data structures
465  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel);
466  CHKERR DMSetFromOptions(dM);
468  if (!boundaryFields.empty()) {
470  }
471  if (!skeletonFields.empty()) {
473  }
474  for (std::vector<std::string>::iterator fit = otherFEs.begin();
475  fit != otherFEs.end(); ++fit) {
476  CHKERR DMMoFEMAddElement(dM, fit->c_str());
477  }
478  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
480 }
481 
482 MoFEMErrorCode Simple::setFieldOrder(const std::string field_name,
483  const int order, const Range *ents) {
485  fieldsOrder.insert(
486  {field_name, {order, ents == NULL ? Range() : Range(*ents)}});
488 }
489 
491  Interface &m_field = cOre;
493  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
494 
495  auto get_skin = [&](auto meshset) {
496  Range domain_ents;
497  CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
498  domain_ents, true);
499  Skinner skin(&m_field.get_moab());
500  Range domain_skin;
501  CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
502  // filter not owned entities, those are not on boundary
503  ParallelComm *pcomm =
504  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
505  Range proc_domain_skin;
506  CHKERR pcomm->filter_pstatus(domain_skin,
507  PSTATUS_SHARED | PSTATUS_MULTISHARED,
508  PSTATUS_NOT, -1, &proc_domain_skin);
509  return proc_domain_skin;
510  };
511 
512  auto create_boundary_meshset = [&](auto &&proc_domain_skin) {
514  // create boundary meshset
515  if (boundaryMeshset != 0) {
516  CHKERR m_field.get_moab().delete_entities(&boundaryMeshset, 1);
517  }
518  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
519  CHKERR m_field.get_moab().add_entities(boundaryMeshset, proc_domain_skin);
520  for (int dd = 0; dd != dIm - 1; dd++) {
521  Range adj;
522  CHKERR m_field.get_moab().get_adjacencies(proc_domain_skin, dd, false,
523  adj, moab::Interface::UNION);
524  CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
525  }
527  };
528 
529  CHKERR create_boundary_meshset(get_skin(meshSet));
530 
531  auto comm_interface_ptr = m_field.getInterface<CommInterface>();
532  auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
533 
534  // Add entities to the fields
535  auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
537  for (auto &field : fields) {
538  CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
539  CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
540  }
542  };
543 
544  auto make_no_field_ents = [&](auto &fields) {
546  for (auto &field : fields) {
547  std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
548  CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 2);
549  CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
550  CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
551  }
553  };
554 
555  CHKERR add_ents_to_field(meshSet, dIm, domainFields);
556  CHKERR add_ents_to_field(meshSet, dIm, dataFields);
557  CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
558  CHKERR add_ents_to_field(meshSet, dIm - 1, skeletonFields);
559  CHKERR make_no_field_ents(noFieldFields);
560  CHKERR make_no_field_ents(noFieldDataFields);
561 
562  // Set order
563  auto set_order = [&](auto meshset, auto dim, auto &fields) {
565  for (auto &field : fields) {
566  if (fieldsOrder.find(field) == fieldsOrder.end()) {
567  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
568  "Order for field not set %s", field.c_str());
569  }
570  int dds = 0;
571  const Field *field_ptr = m_field.get_field_structure(field);
572  switch (field_ptr->getSpace()) {
573  case L2:
574  dds = dim;
575  break;
576  case HDIV:
577  dds = 2;
578  break;
579  case HCURL:
580  dds = 1;
581  break;
582  case H1:
583  dds = 1;
584  break;
585  default:
586  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
587  "Glasgow we have a problem");
588  }
589 
590  auto set_order = [&](auto field, auto &ents) {
592 
593  auto range = fieldsOrder.equal_range(field);
594  for (auto o = range.first; o != range.second; ++o) {
595  if (!o->second.second.empty())
596  ents = intersect(ents, o->second.second);
597  CHKERR m_field.set_field_order(ents, field, o->second.first);
598  }
599 
601  };
602 
603  if (field_ptr->getSpace() == H1) {
604  if (field_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
605  Range ents;
606  CHKERR m_field.get_field_entities_by_dimension(field, 0, ents);
607  CHKERR set_order(field, ents);
608  } else {
609  CHKERR m_field.set_field_order(meshSet, MBVERTEX, field, 1);
610  }
611  }
612  for (int dd = dds; dd <= dim; dd++) {
613  Range ents;
614  CHKERR m_field.get_field_entities_by_dimension(field, dd, ents);
615  CHKERR set_order(field, ents);
616  }
617  }
619  };
620 
621  CHKERR set_order(meshSet, dIm, domainFields);
622  CHKERR set_order(meshSet, dIm, dataFields);
623  CHKERR set_order(meshSet, dIm - 1, boundaryFields);
624  CHKERR set_order(meshSet, dIm - 1, skeletonFields);
625 
626  // Build fields
627  CHKERR m_field.build_fields();
628  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
630 }
631 
633  Interface &m_field = cOre;
635  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
636  // Add finite elements
638  true);
640  if (!boundaryFields.empty()) {
642  boundaryFE, true);
644  }
645  if (!skeletonFields.empty()) {
647  skeletonFE, true);
649  }
650  for (std::vector<std::string>::iterator fit = otherFEs.begin();
651  fit != otherFEs.end(); ++fit) {
652  CHKERR m_field.build_finite_elements(*fit);
653  }
654  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
656 }
657 
659  Interface &m_field = cOre;
661  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
663  // Set problem by the DOFs on the fields rather that by adding DOFs on the
664  // elements
665  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
666  CHKERR DMSetUp(dM);
667  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
668  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
670 }
671 
672 MoFEMErrorCode Simple::setUp(const PetscBool is_partitioned) {
675  if (!skeletonFields.empty())
677  CHKERR defineProblem(is_partitioned);
682 }
683 
685  Interface &m_field = cOre;
687 
689  if (!skeletonFields.empty())
693 
694  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
696  // Set problem by the DOFs on the fields rather that by adding DOFs on the
697  // elements
698  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
700  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
701  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
702 
704 }
705 
708  CHKERR PetscObjectReference(getPetscObject(dM.get()));
709  *dm = dM.get();
711 }
712 
713 } // namespace MoFEM
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
virtual MoFEMErrorCode rebuild_database(int verb=DEFAULT_VERBOSITY)=0
Clear database and initialize it once again.
std::string domainFE
domain finite element
Definition: Simple.hpp:365
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
MoFEM::Core & cOre
Definition: Simple.hpp:343
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:632
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:672
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:274
Deprecated interface functions.
MoFEM interface unique ID.
Problem manager is used to build and partition problems.
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
virtual MoFEMErrorCode add_ents_to_finite_element_by_dim(const EntityHandle entities, const int dim, const std::string &name, const bool recursive=true)=0
add entities to finite element
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:387
field with continuous normal traction
Definition: definitions.h:179
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
int dIm
dimension of problem
Definition: Simple.hpp:372
Provide data structure for (tensor) field approximation.The Field is intended to provide support for ...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
virtual int get_comm_size() const =0
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:425
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
virtual MoFEMErrorCode get_field_entities_by_dimension(const std::string name, int dim, Range &ents) const =0
get entities in the field by dimension
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1119
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
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.
PetscObject getPetscObject(T obj)
Definition: AuxPETSc.hpp:23
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:628
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.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:369
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEMErrorCode reSetUp()
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:684
char mesh_file_name[255]
FieldSpace getSpace() const
Get field approximation space.
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:29
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.
Definition: Simple.cpp:372
const int dim
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:357
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:360
MoFEMErrorCode addDataField(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_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:324
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
MoFEMErrorCode addBoundaryField(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_ZERO, int verb=-1)
Add field on boundary.
Definition: Simple.cpp:287
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:658
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:362
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
FieldApproximationBase
approximation base
Definition: definitions.h:150
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:358
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:185
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:212
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:490
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
FieldApproximationBase getApproxBase() const
Get approximation base.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: Simple.cpp:22
Managing BitRefLevels.
MoFEMErrorCode addDomainField(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_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:269
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:347
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:348
field with continuous tangents
Definition: definitions.h:178
static Index< 'o', 3 > o
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:105
FieldSpace
approximation spaces
Definition: definitions.h:174
virtual MoFEMErrorCode modify_finite_element_adjacency_table(const std::string &fe_name, const EntityType type, ElementAdjacencyFunct function)=0
modify finite element table, only for advanced user
std::string nameOfProblem
problem name
Definition: Simple.hpp:364
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
Finite element data for entity.
MoFEMErrorCode addSkeletonField(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_ZERO, int verb=-1)
Add field on skeleton.
Definition: Simple.cpp:305
#define CHKERR
Inline error check.
Definition: definitions.h:604
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
constexpr int order
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:292
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
virtual MoFEMErrorCode create_vertices_and_add_to_field(const std::string name, const double coords[], int size, int verb=DEFAULT_VERBOSITY)=0
Create a vertices and add to field object.
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:354
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:349
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:341
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:264
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:350
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:38
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1943
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:198
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
Core (interface) class.
Definition: Core.hpp:77
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
virtual MPI_Comm & get_comm() const =0
continuous field
Definition: definitions.h:177
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:482
auto createSmartDM
Creates smart DM object.
Definition: AuxPETSc.hpp:174
virtual MoFEMErrorCode add_ents_to_field_by_dim(const Range &ents, const int dim, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:982
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:459
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...
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:189
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:351
field with C-1 continuity
Definition: definitions.h:180