v0.13.0
plot_base.cpp
Go to the documentation of this file.
1 /**
2  * \file plot_base.cpp
3  * \example plot_base.cpp
4  *
5  * Utility for plotting base functions for different spaces, polynomial bases
6  */
7 
8 /* This file is part of MoFEM.
9  * MoFEM is free software: you can redistribute it and/or modify it under
10  * the terms of the GNU Lesser General Public License as published by the
11  * Free Software Foundation, either version 3 of the License, or (at your
12  * option) any later version.
13  *
14  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
15  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
16  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17  * License for more details.
18  *
19  * You should have received a copy of the GNU Lesser General Public
20  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
21 
22 #include <MoFEM.hpp>
23 
24 using namespace MoFEM;
25 
26 static char help[] = "...\n\n";
27 
28 #include <BasicFiniteElements.hpp>
29 
30 template <int DIM> struct ElementsAndOps {};
31 
32 template <> struct ElementsAndOps<2> {
36 };
37 
38 template <> struct ElementsAndOps<3> {
42 };
43 
44 constexpr int SPACE_DIM =
45  EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
46 
51 
52 struct MyPostProc : public PostProcEle {
54 
55  MoFEMErrorCode generateReferenceElementMesh();
56  MoFEMErrorCode setGaussPts(int order);
57 
58 protected:
59  ublas::matrix<int> refEleMap;
61 };
62 
63 struct Example {
64 
65  Example(MoFEM::Interface &m_field) : mField(m_field) {}
66 
68 
69 private:
70  MoFEM::Interface &mField;
71  Simple *simpleInterface;
72 
75  MoFEMErrorCode setIntegrationRules();
82 
85 };
86 
87 //! [Run programme]
90  CHKERR readMesh();
91  CHKERR setupProblem();
92  CHKERR setIntegrationRules();
93  CHKERR createCommonData();
94  CHKERR boundaryCondition();
95  CHKERR assembleSystem();
96  CHKERR solveSystem();
97  CHKERR outputResults();
98  CHKERR checkResults();
100 }
101 //! [Run programme]
102 
103 //! [Read mesh]
106 
107  PetscBool load_file = PETSC_FALSE;
108  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-load_file", &load_file,
109  PETSC_NULL);
110 
111  if (load_file == PETSC_FALSE) {
112 
113  auto &moab = mField.get_moab();
114 
115  if (SPACE_DIM == 3) {
116 
117  // create one tet
118  double tet_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
119  EntityHandle nodes[4];
120  for (int nn = 0; nn < 4; nn++) {
121  CHKERR moab.create_vertex(&tet_coords[3 * nn], nodes[nn]);
122  }
123  EntityHandle tet;
124  CHKERR moab.create_element(MBTET, nodes, 4, tet);
125  Range adj;
126  for (auto d : {1, 2})
127  CHKERR moab.get_adjacencies(&tet, 1, d, true, adj);
128  }
129 
130  if (SPACE_DIM == 2) {
131 
132  // create one triangle
133  double tri_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0};
134  EntityHandle nodes[3];
135  for (int nn = 0; nn < 3; nn++) {
136  CHKERR moab.create_vertex(&tri_coords[3 * nn], nodes[nn]);
137  }
138  EntityHandle tri;
139  CHKERR moab.create_element(MBTRI, nodes, 3, tri);
140  Range adj;
141  CHKERR moab.get_adjacencies(&tri, 1, 1, true, adj);
142  }
143 
144  CHKERR mField.rebuild_database();
145  CHKERR mField.getInterface(simpleInterface);
146  simpleInterface->setDim(SPACE_DIM);
147 
148  // Add all elements to database
149  CHKERR mField.getInterface<BitRefManager>()->setBitRefLevelByDim(
150  0, SPACE_DIM, simpleInterface->getBitRefLevel());
151 
152  } else {
153 
154  CHKERR mField.getInterface(simpleInterface);
155  CHKERR simpleInterface->getOptions();
156  CHKERR simpleInterface->loadFile();
157 
158  }
159 
161 }
162 //! [Read mesh]
163 
164 //! [Set up problem]
167  // Add field
168 
169  // Declare elements
170  enum bases { AINSWORTH, AINSWORTH_LOBATTO, DEMKOWICZ, BERNSTEIN, LASBASETOP };
171  const char *list_bases[] = {"ainsworth", "ainsworth_lobatto", "demkowicz",
172  "bernstein"};
173 
174  PetscBool flg;
175  PetscInt choice_base_value = AINSWORTH;
176  CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, "-base", list_bases, LASBASETOP,
177  &choice_base_value, &flg);
178  if (flg != PETSC_TRUE)
179  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "base not set");
181  if (choice_base_value == AINSWORTH)
183  if (choice_base_value == AINSWORTH_LOBATTO)
184  base = AINSWORTH_LOBATTO_BASE;
185  else if (choice_base_value == DEMKOWICZ)
186  base = DEMKOWICZ_JACOBI_BASE;
187  else if (choice_base_value == BERNSTEIN)
189 
190  enum spaces { H1SPACE, L2SPACE, HCURLSPACE, HDIVSPACE, LASBASETSPACE };
191  const char *list_spaces[] = {"h1", "l2", "hcurl", "hdiv"};
192  PetscInt choice_space_value = H1SPACE;
193  CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, "-space", list_spaces,
194  LASBASETSPACE, &choice_space_value, &flg);
195  if (flg != PETSC_TRUE)
196  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "space not set");
197  space = H1;
198  if (choice_space_value == H1SPACE)
199  space = H1;
200  else if (choice_space_value == L2SPACE)
201  space = L2;
202  else if (choice_space_value == HCURLSPACE)
203  space = HCURL;
204  else if (choice_space_value == HDIVSPACE)
205  space = HDIV;
206 
207  CHKERR simpleInterface->addDomainField("U", space, base, 1);
208 
209  int order = 3;
210  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
211  CHKERR simpleInterface->setFieldOrder("U", order);
212  CHKERR simpleInterface->setUp();
213 
215 }
216 //! [Set up problem]
217 
218 //! [Set integration rule]
222 }
223 //! [Set integration rule]
224 
225 //! [Create common data]
227 //! [Create common data]
228 
229 //! [Boundary condition]
231 //! [Boundary condition]
232 
233 //! [Push operators to pipeline]
235 //! [Push operators to pipeline]
236 
237 //! [Solve]
238 MoFEMErrorCode Example::solveSystem() { return 0; }
239 
240 //! [Solve]
243 
244  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
245  pipeline_mng->getDomainLhsFE().reset();
246 
247  auto post_proc_fe = boost::make_shared<MyPostProc>(mField);
248  post_proc_fe->generateReferenceElementMesh();
249  pipeline_mng->getDomainRhsFE() = post_proc_fe;
250 
251  if (SPACE_DIM == 2) {
252  if (space == HCURL) {
253  auto jac_ptr = boost::make_shared<MatrixDouble>();
254  post_proc_fe->getOpPtrVector().push_back(
255  new OpCalculateHOJacForFace(jac_ptr));
256  post_proc_fe->getOpPtrVector().push_back(new OpMakeHdivFromHcurl());
257  post_proc_fe->getOpPtrVector().push_back(
259  }
260  }
261 
262  post_proc_fe->addFieldValuesPostProc("U");
263 
264  auto scale_tag_val = [&]() {
266  auto &post_proc_mesh = post_proc_fe->postProcMesh;
267  Range nodes;
268  CHKERR post_proc_mesh.get_entities_by_type(0, MBVERTEX, nodes);
269  Tag th;
270  CHKERR post_proc_mesh.tag_get_handle("U", th);
271  int length;
272  CHKERR post_proc_mesh.tag_get_length(th, length);
273  std::vector<double> data(nodes.size() * length);
274  CHKERR post_proc_mesh.tag_get_data(th, nodes, &*data.begin());
275  double max_v = 0;
276  for (int i = 0; i != nodes.size(); ++i) {
277  double v = 0;
278  for (int d = 0; d != length; ++d)
279  v += pow(data[length * i + d], 2);
280  v = std::sqrt(v);
281  max_v = std::max(max_v, v);
282  }
283  for (auto &v : data)
284  v /= max_v;
285  CHKERR post_proc_mesh.tag_set_data(th, nodes, &*data.begin());
287  };
288 
289  size_t nb = 0;
290  auto dofs_ptr = mField.get_dofs();
291 
292  for (auto dof_ptr : (*dofs_ptr)) {
293  MOFEM_LOG("PLOTBASE", Sev::verbose) << *dof_ptr;
294  auto &val = const_cast<double &>(dof_ptr->getFieldData());
295  val = 1;
296  CHKERR pipeline_mng->loopFiniteElements();
297  CHKERR scale_tag_val();
298  CHKERR post_proc_fe->writeFile(
299  "out_base_dof_" + boost::lexical_cast<std::string>(nb) + ".h5m");
300  CHKERR post_proc_fe->postProcMesh.delete_mesh();
301  val = 0;
302  ++nb;
303  };
304 
306 }
307 //! [Postprocess results]
308 
309 //! [Check results]
311 //! [Check results]
312 
313 int main(int argc, char *argv[]) {
314 
315  // Initialisation of MoFEM/PETSc and MOAB data structures
316  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
317 
318  try {
319 
320  //! [Register MoFEM discrete manager in PETSc]
321  DMType dm_name = "DMMOFEM";
322  CHKERR DMRegister_MoFEM(dm_name);
323  //! [Register MoFEM discrete manager in PETSc
324 
325  // Add logging channel for example
326  auto core_log = logging::core::get();
327  core_log->add_sink(
329  LogManager::setLog("PLOTBASE");
330  MOFEM_LOG_TAG("PLOTBASE", "plotbase");
331 
332  //! [Create MoAB]
333  moab::Core mb_instance; ///< mesh database
334  moab::Interface &moab = mb_instance; ///< mesh database interface
335  //! [Create MoAB]
336 
337  //! [Create MoFEM]
338  MoFEM::Core core(moab); ///< finite element database
339  MoFEM::Interface &m_field = core; ///< finite element database insterface
340  //! [Create MoFEM]
341 
342  //! [Example]
343  Example ex(m_field);
344  CHKERR ex.runProblem();
345  //! [Example]
346  }
347  CATCH_ERRORS;
348 
350 
351  return 0;
352 }
353 
356  moab::Core core_ref;
357  moab::Interface &moab_ref = core_ref;
358 
359  char ref_mesh_file_name[255];
360 
361  if (SPACE_DIM == 2)
362  strcpy(ref_mesh_file_name, "ref_mesh2d.h5m");
363  else if (SPACE_DIM == 3)
364  strcpy(ref_mesh_file_name, "ref_mesh3d.h5m");
365  else
366  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
367  "Dimension not implemented");
368 
369  CHKERR PetscOptionsGetString(PETSC_NULL, "", "-ref_file", ref_mesh_file_name,
370  255, PETSC_NULL);
371  CHKERR moab_ref.load_file(ref_mesh_file_name, 0, "");
372 
373  // Get elements
374  Range elems;
375  CHKERR moab_ref.get_entities_by_dimension(0, SPACE_DIM, elems);
376 
377  // Add mid-nodes on edges
378  EntityHandle meshset;
379  CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
380  CHKERR moab_ref.add_entities(meshset, elems);
381  CHKERR moab_ref.convert_entities(meshset, true, false, false);
382  CHKERR moab_ref.delete_entities(&meshset, 1);
383 
384  // Get nodes on the mesh
385  Range elem_nodes;
386  CHKERR moab_ref.get_connectivity(elems, elem_nodes, false);
387 
388  // Map node entity and Gauss pint number
389  std::map<EntityHandle, int> nodes_pts_map;
390 
391  // Set gauss points coordinates from the reference mesh
392  gaussPts.resize(SPACE_DIM + 1, elem_nodes.size(), false);
393  gaussPts.clear();
394  Range::iterator nit = elem_nodes.begin();
395  for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
396  double coords[3];
397  CHKERR moab_ref.get_coords(&*nit, 1, coords);
398  for (auto d : {0, 1, 2})
399  gaussPts(d, gg) = coords[d];
400  nodes_pts_map[*nit] = gg;
401  }
402 
403  if (SPACE_DIM == 2) {
404  // Set size of adjacency matrix (note ho order nodes 3 nodes and 3 nodes on
405  // edges)
406  refEleMap.resize(elems.size(), 3 + 3);
407  } else if (SPACE_DIM == 3) {
408  refEleMap.resize(elems.size(), 4 + 6);
409  }
410 
411  // Set adjacency matrix
412  Range::iterator tit = elems.begin();
413  for (int tt = 0; tit != elems.end(); ++tit, ++tt) {
414  const EntityHandle *conn;
415  int num_nodes;
416  CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes, false);
417  for (int nn = 0; nn != num_nodes; ++nn) {
418  refEleMap(tt, nn) = nodes_pts_map[conn[nn]];
419  }
420  }
421 
423 }
424 
427 
428  const int num_nodes = gaussPts.size2();
429 
430  // Calculate sheape functions
431 
432  switch (numeredEntFiniteElementPtr->getEntType()) {
433  case MBTRI:
434  shapeFunctions.resize(num_nodes, 3);
435  CHKERR Tools::shapeFunMBTRI(&*shapeFunctions.data().begin(),
436  &gaussPts(0, 0), &gaussPts(1, 0), num_nodes);
437  break;
438  case MBQUAD: {
439  shapeFunctions.resize(num_nodes, 4);
440  for (int gg = 0; gg != num_nodes; gg++) {
441  double ksi = gaussPts(0, gg);
442  double eta = gaussPts(1, gg);
443  shapeFunctions(gg, 0) = N_MBQUAD0(ksi, eta);
444  shapeFunctions(gg, 1) = N_MBQUAD1(ksi, eta);
445  shapeFunctions(gg, 2) = N_MBQUAD2(ksi, eta);
446  shapeFunctions(gg, 3) = N_MBQUAD3(ksi, eta);
447  }
448  } break;
449  case MBTET: {
450  shapeFunctions.resize(num_nodes, 8);
451  CHKERR Tools::shapeFunMBTET(&*shapeFunctions.data().begin(),
452  &gaussPts(0, 0), &gaussPts(1, 0),
453  &gaussPts(2, 0), num_nodes);
454  } break;
455  case MBHEX: {
456  shapeFunctions.resize(num_nodes, 8);
457  for (int gg = 0; gg != num_nodes; gg++) {
458  double ksi = gaussPts(0, gg);
459  double eta = gaussPts(1, gg);
460  double zeta = gaussPts(2, gg);
461  shapeFunctions(gg, 0) = N_MBHEX0(ksi, eta, zeta);
462  shapeFunctions(gg, 1) = N_MBHEX1(ksi, eta, zeta);
463  shapeFunctions(gg, 2) = N_MBHEX2(ksi, eta, zeta);
464  shapeFunctions(gg, 3) = N_MBHEX3(ksi, eta, zeta);
465  shapeFunctions(gg, 4) = N_MBHEX4(ksi, eta, zeta);
466  shapeFunctions(gg, 5) = N_MBHEX5(ksi, eta, zeta);
467  shapeFunctions(gg, 6) = N_MBHEX6(ksi, eta, zeta);
468  shapeFunctions(gg, 7) = N_MBHEX7(ksi, eta, zeta);
469  }
470  } break;
471  default:
472  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
473  "Not implemented element type");
474  }
475 
476 
477  // Create physical nodes
478 
479  // MoAB interface allowing for creating nodes and elements in the bulk
480  ReadUtilIface *iface;
481  CHKERR postProcMesh.query_interface(iface);
482 
483  std::vector<double *> arrays; /// pointers to memory allocated by MoAB for
484  /// storing X, Y, and Z coordinates
485  EntityHandle startv; // Starting handle for vertex
486  // Allocate memory for num_nodes, and return starting handle, and access to
487  // memort.
488  CHKERR iface->get_node_coords(3, num_nodes, 0, startv, arrays);
489 
490  mapGaussPts.resize(gaussPts.size2());
491  for (int gg = 0; gg != num_nodes; ++gg)
492  mapGaussPts[gg] = startv + gg;
493 
494  Tag th;
495  int def_in_the_loop = -1;
496  CHKERR postProcMesh.tag_get_handle("NB_IN_THE_LOOP", 1, MB_TYPE_INTEGER, th,
497  MB_TAG_CREAT | MB_TAG_SPARSE,
498  &def_in_the_loop);
499 
500  // Create physical elements
501 
502  const int num_el = refEleMap.size1();
503  const int num_nodes_on_ele = refEleMap.size2();
504 
505  EntityHandle starte; // Starting handle to first created element
506  EntityHandle *conn; // Access to MOAB memory with connectivity of elements
507 
508  // Create tris/tets in the bulk in MoAB database
509  if (SPACE_DIM == 2)
510  CHKERR iface->get_element_connect(num_el, num_nodes_on_ele, MBTRI, 0,
511  starte, conn);
512  else if (SPACE_DIM == 3)
513  CHKERR iface->get_element_connect(num_el, num_nodes_on_ele, MBTET, 0,
514  starte, conn);
515  else
516  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
517  "Dimension not implemented");
518 
519  // At this point elements (memory for elements) is allocated, at code bellow
520  // actual connectivity of elements is set.
521  for (unsigned int tt = 0; tt != refEleMap.size1(); ++tt) {
522  for (int nn = 0; nn != num_nodes_on_ele; ++nn)
523  conn[num_nodes_on_ele * tt + nn] = mapGaussPts[refEleMap(tt, nn)];
524  }
525 
526  // Finalise elements creation. At that point MOAB updates adjacency tables,
527  // and elements are ready to use.
528  CHKERR iface->update_adjacencies(starte, num_el, num_nodes_on_ele, conn);
529 
530  auto physical_elements = Range(starte, starte + num_el - 1);
531  CHKERR postProcMesh.tag_clear_data(th, physical_elements, &(nInTheLoop));
532 
533  EntityHandle fe_ent = numeredEntFiniteElementPtr->getEnt();
534  int fe_num_nodes;
535  {
536  const EntityHandle *conn;
537  mField.get_moab().get_connectivity(fe_ent, conn, fe_num_nodes, true);
538  coords.resize(3 * fe_num_nodes, false);
539  CHKERR mField.get_moab().get_coords(conn, fe_num_nodes, &coords[0]);
540  }
541 
542  // Set physical coordinates to physical nodes
545  &*shapeFunctions.data().begin());
546 
548  arrays[0], arrays[1], arrays[2]);
549  const double *t_coords_ele_x = &coords[0];
550  const double *t_coords_ele_y = &coords[1];
551  const double *t_coords_ele_z = &coords[2];
552  for (int gg = 0; gg != num_nodes; ++gg) {
554  t_coords_ele_x, t_coords_ele_y, t_coords_ele_z);
555  t_coords(i) = 0;
556  for (int nn = 0; nn != fe_num_nodes; ++nn) {
557  t_coords(i) += t_n * t_ele_coords(i);
558  for (auto ii : {0, 1, 2})
559  if (std::abs(t_coords(ii)) < std::numeric_limits<float>::epsilon())
560  t_coords(ii) = 0;
561  ++t_ele_coords;
562  ++t_n;
563  }
564  ++t_coords;
565  }
566 
568 }
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
FieldApproximationBase
approximation base
Definition: definitions.h:71
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:75
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:79
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:77
FieldSpace
approximation spaces
Definition: definitions.h:95
@ L2
field with C-1 continuity
Definition: definitions.h:101
@ H1
continuous field
Definition: definitions.h:98
@ HCURL
field with continuous tangents
Definition: definitions.h:99
@ HDIV
field with continuous normal traction
Definition: definitions.h:100
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_IMPOSIBLE_CASE
Definition: definitions.h:48
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define N_MBQUAD3(x, y)
quad shape function
Definition: fem_tools.h:70
#define N_MBHEX7(x, y, z)
Definition: fem_tools.h:88
#define N_MBHEX3(x, y, z)
Definition: fem_tools.h:84
#define N_MBHEX5(x, y, z)
Definition: fem_tools.h:86
#define N_MBHEX4(x, y, z)
Definition: fem_tools.h:85
#define N_MBHEX0(x, y, z)
Definition: fem_tools.h:81
#define N_MBHEX6(x, y, z)
Definition: fem_tools.h:87
#define N_MBHEX2(x, y, z)
Definition: fem_tools.h:83
#define N_MBQUAD0(x, y)
quad shape function
Definition: fem_tools.h:67
#define N_MBHEX1(x, y, z)
Definition: fem_tools.h:82
#define N_MBQUAD2(x, y)
quad shape function
Definition: fem_tools.h:69
#define N_MBQUAD1(x, y)
quad shape function
Definition: fem_tools.h:68
constexpr double eta
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
FTensor::Index< 'i', SPACE_DIM > i
double v
phase velocity of light in medium (cm/ns)
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
CoreTmp< 0 > Core
Definition: Core.hpp:1096
OpSetContravariantPiolaTransformOnFace2DImpl< 2 > OpSetContravariantPiolaTransformOnFace2D
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:24
int main(int argc, char *argv[])
[Check results]
Definition: plot_base.cpp:313
EntitiesFieldData::EntData EntData
Definition: plot_base.cpp:47
static char help[]
Definition: plot_base.cpp:26
constexpr int SPACE_DIM
Definition: plot_base.cpp:44
ElementsAndOps< SPACE_DIM >::PostProcEle PostProcEle
Definition: plot_base.cpp:50
PipelineManager::FaceEle DomainEle
[Example]
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode setIntegrationRules()
[Set up problem]
Definition: plot_base.cpp:219
FieldApproximationBase base
Definition: plot_base.cpp:83
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
Definition: plot_base.cpp:65
MoFEMErrorCode runProblem()
FieldSpace space
Definition: plot_base.cpp:84
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
Managing BitRefLevels.
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
Make Hdiv space from Hcurl space in 2d.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi, const double *eta, const double *zeta, const double nb)
Calculate shape functions on tetrahedron.
Definition: Tools.hpp:607
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const int nb)
Calculate shape functions on triangle.
Definition: Tools.hpp:576
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
VolumeElementForcesAndSourcesCoreBase::UserDataOperator UserDataOperator
Definition: plot_base.cpp:52
MoFEMErrorCode generateReferenceElementMesh()
Definition: plot_base.cpp:354
ublas::matrix< int > refEleMap
Definition: plot_base.cpp:59
MoFEMErrorCode setGaussPts(int order)
Definition: plot_base.cpp:425
MatrixDouble shapeFunctions
Definition: plot_base.cpp:60
Postprocess on face.
Post processing.