v0.9.1
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 L2:
64  case HCURL:
65  case HDIV:
66  CHKERR moab.get_adjacencies(bride_adjacency_edge, 1, false,
67  bride_adjacency, moab::Interface::UNION);
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 L2:
122  case HDIV:
123  CHKERR moab.get_adjacencies(bride_adjacency_face, 2, false,
124  bride_adjacency, moab::Interface::UNION);
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 
256 Simple::addDomainField(const std::string &name, const FieldSpace space,
257  const FieldApproximationBase base,
258  const FieldCoefficientsNumber nb_of_cooficients,
259  const TagType tag_type, const enum MoFEMTypes bh,
260  int verb) {
261 
262  Interface &m_field = cOre;
264  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
265  verb);
266  if (space == NOFIELD)
267  noFieldFields.push_back(name);
268  else
269  domainFields.push_back(name);
271 }
272 
274 Simple::addBoundaryField(const std::string &name, const FieldSpace space,
275  const FieldApproximationBase base,
276  const FieldCoefficientsNumber nb_of_cooficients,
277  const TagType tag_type, const enum MoFEMTypes bh,
278  int verb) {
279  Interface &m_field = cOre;
281  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
282  verb);
283  boundaryFields.push_back(name);
284  if (space == NOFIELD)
285  SETERRQ(
286  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
287  "NOFIELD space for boundary filed not implemented in Simple interface");
289 }
290 
292 Simple::addSkeletonField(const std::string &name, const FieldSpace space,
293  const FieldApproximationBase base,
294  const FieldCoefficientsNumber nb_of_cooficients,
295  const TagType tag_type, const enum MoFEMTypes bh,
296  int verb) {
297 
298  Interface &m_field = cOre;
300  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
301  verb);
302  skeletonFields.push_back(name);
303  if (space == NOFIELD)
304  SETERRQ(
305  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
306  "NOFIELD space for boundary filed not implemented in Simple interface");
308 }
309 
311 Simple::addDataField(const std::string &name, const FieldSpace space,
312  const FieldApproximationBase base,
313  const FieldCoefficientsNumber nb_of_cooficients,
314  const TagType tag_type, const enum MoFEMTypes bh,
315  int verb) {
316 
317  Interface &m_field = cOre;
319  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
320  verb);
321  if (space == NOFIELD)
322  noFieldDataFields.push_back(name);
323  else
324  dataFields.push_back(name);
326 }
327 
328 MoFEMErrorCode Simple::removeDomainField(const std::string &name) {
329  Interface &m_field = cOre;
331 
332  auto remove_field_from_list = [&](auto &vec) {
333  auto it = std::find(vec.begin(), vec.end(), name);
334  if (it != vec.end())
335  vec.erase(it);
336  };
337 
338  remove_field_from_list(noFieldFields);
339  remove_field_from_list(domainFields);
340 
342 }
343 
344 MoFEMErrorCode Simple::removeBoundaryField(const std::string &name) {
345  Interface &m_field = cOre;
347 
348  auto remove_field_from_list = [&](auto &vec) {
349  auto it = std::find(vec.begin(), vec.end(), name);
350  if (it != vec.end())
351  vec.erase(it);
352  };
353 
354  remove_field_from_list(boundaryFields);
355 
357 }
358 
359 MoFEMErrorCode Simple::removeSkeletonField(const std::string &name) {
360  Interface &m_field = cOre;
362 
363  auto remove_field_from_list = [&](auto &vec) {
364  auto it = std::find(vec.begin(), vec.end(), name);
365  if (it != vec.end())
366  vec.erase(it);
367  };
368 
369  remove_field_from_list(skeletonFields);
370 
372 }
373 
375  Interface &m_field = cOre;
377 
378  auto clear_rows_and_cols = [&](auto &fe_name) {
380  const FiniteElement_multiIndex *fe_ptr;
381  CHKERR m_field.get_finite_elements(&fe_ptr);
382  auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
383  ->get<FiniteElement_name_mi_tag>();
384  auto it_fe = fe_by_name.find(fe_name);
385  if (it_fe != fe_by_name.end()) {
386 
387  if(!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
388  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
389  "modification unsuccessful");
390 
391  if(!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
392  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
393  "modification unsuccessful");
394  }
396  };
397  CHKERR clear_rows_and_cols(domainFE);
398  CHKERR clear_rows_and_cols(boundaryFE);
399  CHKERR clear_rows_and_cols(skeletonFE);
400 
401  // Define finite elements
403 
404  auto add_fields = [&](auto &fe_name, auto &fields) {
406  for (auto &field : fields) {
407  CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
408  CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
409  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
410  }
412  };
413 
414  auto add_data_fields = [&](auto &fe_name, auto &fields) {
416  for (auto &field : fields)
417  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
419  };
420 
421  CHKERR add_fields(domainFE, domainFields);
422  CHKERR add_data_fields(domainFE, dataFields);
423  CHKERR add_fields(domainFE, noFieldFields);
424  CHKERR add_data_fields(domainFE, noFieldDataFields);
425 
426  if (!boundaryFields.empty()) {
428  CHKERR add_fields(boundaryFE, domainFields);
429  CHKERR add_fields(boundaryFE, boundaryFields);
430  CHKERR add_fields(boundaryFE, skeletonFields);
431  CHKERR add_data_fields(boundaryFE, dataFields);
432  CHKERR add_data_fields(boundaryFE, noFieldDataFields);
433  CHKERR add_fields(boundaryFE, noFieldFields);
434  }
435  if (!skeletonFields.empty()) {
437  CHKERR add_fields(skeletonFE, domainFields);
438  CHKERR add_fields(skeletonFE, boundaryFields);
439  CHKERR add_fields(skeletonFE, skeletonFields);
440  CHKERR add_data_fields(skeletonFE, dataFields);
441  CHKERR add_data_fields(skeletonFE, noFieldDataFields);
442  CHKERR add_fields(skeletonFE, noFieldFields);
443  }
445 }
446 
447 MoFEMErrorCode Simple::defineProblem(const PetscBool is_partitioned) {
448  Interface &m_field = cOre;
450  // Create dm instance
451  dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
452  // set dm data structure which created mofem data structures
453  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel);
454  CHKERR DMSetFromOptions(dM);
456  if (!boundaryFields.empty()) {
458  }
459  if (!skeletonFields.empty()) {
461  }
462  for (std::vector<std::string>::iterator fit = otherFEs.begin();
463  fit != otherFEs.end(); ++fit) {
464  CHKERR DMMoFEMAddElement(dM, fit->c_str());
465  }
466  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
468 }
469 
470 MoFEMErrorCode Simple::setFieldOrder(const std::string field_name,
471  const int order, const Range *ents) {
473  fieldsOrder.insert(
474  {field_name, {order, ents == NULL ? Range() : Range(*ents)}});
476 }
477 
479  Interface &m_field = cOre;
481  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
482 
483  auto get_skin = [&](auto meshset) {
484  Range domain_ents;
485  CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
486  domain_ents, true);
487  Skinner skin(&m_field.get_moab());
488  Range domain_skin;
489  CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
490  // filter not owned entities, those are not on boundary
491  ParallelComm *pcomm =
492  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
493  Range proc_domain_skin;
494  CHKERR pcomm->filter_pstatus(domain_skin,
495  PSTATUS_SHARED | PSTATUS_MULTISHARED,
496  PSTATUS_NOT, -1, &proc_domain_skin);
497  return proc_domain_skin;
498  };
499 
500  auto create_boundary_meshset = [&](auto &&proc_domain_skin) {
502  // create boundary meshset
503  if (boundaryMeshset != 0) {
504  CHKERR m_field.get_moab().delete_entities(&boundaryMeshset, 1);
505  }
506  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
507  CHKERR m_field.get_moab().add_entities(boundaryMeshset, proc_domain_skin);
508  for (int dd = 0; dd != dIm - 1; dd++) {
509  Range adj;
510  CHKERR m_field.get_moab().get_adjacencies(proc_domain_skin, dd, false,
511  adj, moab::Interface::UNION);
512  CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
513  }
515  };
516 
517  CHKERR create_boundary_meshset(get_skin(meshSet));
518 
519  auto comm_interface_ptr = m_field.getInterface<CommInterface>();
520  auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
521 
522  // Add entities to the fields
523  auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
525  for (auto &field : fields) {
526  CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
527  CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
528  }
530  };
531 
532  auto make_no_field_ents = [&](auto &fields) {
534  for (auto &field : fields) {
535  std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
536  CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 2);
537  CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
538  CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
539  }
541  };
542 
543  CHKERR add_ents_to_field(meshSet, dIm, domainFields);
544  CHKERR add_ents_to_field(meshSet, dIm, dataFields);
545  CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
546  CHKERR add_ents_to_field(meshSet, dIm - 1, skeletonFields);
547  CHKERR make_no_field_ents(noFieldFields);
548  CHKERR make_no_field_ents(noFieldDataFields);
549 
550  // Set order
551  auto set_order = [&](auto meshset, auto dim, auto &fields) {
553  for (auto &field : fields) {
554  if (fieldsOrder.find(field) == fieldsOrder.end()) {
555  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
556  "Order for field not set %s", field.c_str());
557  }
558  int dds = 0;
559  const Field *field_ptr = m_field.get_field_structure(field);
560  switch (field_ptr->getSpace()) {
561  case L2:
562  dds = dim;
563  break;
564  case HDIV:
565  dds = 2;
566  break;
567  case HCURL:
568  dds = 1;
569  break;
570  case H1:
571  dds = 1;
572  break;
573  default:
574  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
575  "Glasgow we have a problem");
576  }
577 
578  auto set_order = [&](auto field, auto &ents) {
580 
581  auto range = fieldsOrder.equal_range(field);
582  for (auto o = range.first; o != range.second; ++o) {
583  if (!o->second.second.empty())
584  ents = intersect(ents, o->second.second);
585  CHKERR m_field.set_field_order(ents, field, o->second.first);
586  }
587 
589  };
590 
591  if (field_ptr->getSpace() == H1) {
592  if (field_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
593  Range ents;
594  CHKERR m_field.get_field_entities_by_dimension(field, 0, ents);
595  CHKERR set_order(field, ents);
596  } else {
597  CHKERR m_field.set_field_order(meshSet, MBVERTEX, field, 1);
598  }
599  }
600  for (int dd = dds; dd <= dim; dd++) {
601  Range ents;
602  CHKERR m_field.get_field_entities_by_dimension(field, dd, ents);
603  CHKERR set_order(field, ents);
604  }
605  }
607  };
608 
609  CHKERR set_order(meshSet, dIm, domainFields);
610  CHKERR set_order(meshSet, dIm, dataFields);
611  CHKERR set_order(meshSet, dIm - 1, boundaryFields);
612  CHKERR set_order(meshSet, dIm - 1, skeletonFields);
613 
614  // Build fields
615  CHKERR m_field.build_fields();
616  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
618 }
619 
621  Interface &m_field = cOre;
623  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
624  // Add finite elements
626  true);
628  if (!boundaryFields.empty()) {
630  boundaryFE, true);
632  }
633  if (!skeletonFields.empty()) {
635  skeletonFE, true);
637  }
638  for (std::vector<std::string>::iterator fit = otherFEs.begin();
639  fit != otherFEs.end(); ++fit) {
640  CHKERR m_field.build_finite_elements(*fit);
641  }
642  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
644 }
645 
647  Interface &m_field = cOre;
649  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
651  // Set problem by the DOFs on the fields rather that by adding DOFs on the
652  // elements
653  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
654  CHKERR DMSetUp(dM);
655  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
656  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
658 }
659 
660 MoFEMErrorCode Simple::setUp(const PetscBool is_partitioned) {
663  if (!skeletonFields.empty())
665  CHKERR defineProblem(is_partitioned);
670 }
671 
673  Interface &m_field = cOre;
675 
677  if (!skeletonFields.empty())
681 
682  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
684  // Set problem by the DOFs on the fields rather that by adding DOFs on the
685  // elements
686  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
688  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
689  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
690 
692 }
693 
696  CHKERR PetscObjectReference(getPetscObject(dM.get()));
697  *dm = dM.get();
699 }
700 
701 } // namespace MoFEM
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363
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:362
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:340
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:620
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:350
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:660
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:271
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:374
field with continuous normal traction
Definition: definitions.h:178
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
char meshFileName[255]
mesh file name
Definition: Simple.hpp:368
int dIm
dimension of problem
Definition: Simple.hpp:369
Provide data structure for (tensor) field approximation.The Field is intended to provide support for ...
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:411
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:372
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:482
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:1105
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
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:353
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:625
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.
Core (interface) class.
Definition: Core.hpp:50
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:356
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:366
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEMErrorCode reSetUp()
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:672
char mesh_file_name[255]
FieldSpace getSpace() const
Get field approximation space.
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:29
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.
Definition: Simple.cpp:359
virtual MoFEMErrorCode get_finite_elements(const FiniteElement_multiIndex **fe_ptr) const =0
Get finite elements multi-index.
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:344
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:357
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:311
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:274
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:646
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:359
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:149
MoFEMErrorCode loadFile(const std::string options="PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;", const std::string mesh_file_name="")
Load mesh file.
Definition: Simple.cpp:212
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:355
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:185
constexpr int order
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:478
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
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:256
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:344
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:345
field with continuous tangents
Definition: definitions.h:177
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:173
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:361
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:292
#define CHKERR
Inline error check.
Definition: definitions.h:601
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:290
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:351
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:346
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:328
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:261
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:347
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:38
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
virtual MPI_Comm & get_comm() const =0
continuous field
Definition: definitions.h:176
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:470
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:342
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:968
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:447
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:188
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:348
field with C-1 continuity
Definition: definitions.h:179