v0.14.0
Loading...
Searching...
No Matches
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#include <MoFEM.hpp>
9
10using namespace MoFEM;
11
12static char help[] = "...\n\n";
13
15
16template <int DIM> struct ElementsAndOps {};
17
18template <> struct ElementsAndOps<2> {
20};
21
22template <> struct ElementsAndOps<3> {
24};
25
26constexpr int SPACE_DIM =
27 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
28
33
34struct MyPostProc : public PostProcEle {
35 using PostProcEle::PostProcEle;
36
39
42
43protected:
44 ublas::matrix<int> refEleMap;
46};
47
48struct Example {
49
50 Example(MoFEM::Interface &m_field) : mField(m_field) {}
51
53
54private:
57
67
70};
71
72//! [Run programme]
85}
86//! [Run programme]
87
88//! [Read mesh]
91
92 PetscBool load_file = PETSC_FALSE;
93 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-load_file", &load_file,
94 PETSC_NULL);
95
96 if (load_file == PETSC_FALSE) {
97
98 auto &moab = mField.get_moab();
99
100 if (SPACE_DIM == 3) {
101
102 // create one tet
103 double tet_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
104 EntityHandle nodes[4];
105 for (int nn = 0; nn < 4; nn++) {
106 CHKERR moab.create_vertex(&tet_coords[3 * nn], nodes[nn]);
107 }
108 EntityHandle tet;
109 CHKERR moab.create_element(MBTET, nodes, 4, tet);
110 Range adj;
111 for (auto d : {1, 2})
112 CHKERR moab.get_adjacencies(&tet, 1, d, true, adj);
113 }
114
115 if (SPACE_DIM == 2) {
116
117 // create one triangle
118 double tri_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0};
119 EntityHandle nodes[3];
120 for (int nn = 0; nn < 3; nn++) {
121 CHKERR moab.create_vertex(&tri_coords[3 * nn], nodes[nn]);
122 }
123 EntityHandle tri;
124 CHKERR moab.create_element(MBTRI, nodes, 3, tri);
125 Range adj;
126 CHKERR moab.get_adjacencies(&tri, 1, 1, true, adj);
127 }
128
132
133 // Add all elements to database
134 CHKERR mField.getInterface<BitRefManager>()->setBitRefLevelByDim(
136
137 } else {
138
142 }
143
145}
146//! [Read mesh]
147
148//! [Set up problem]
151 // Add field
152
153 // Declare elements
154 enum bases { AINSWORTH, AINSWORTH_LOBATTO, DEMKOWICZ, BERNSTEIN, LASBASETOP };
155 const char *list_bases[] = {"ainsworth", "ainsworth_lobatto", "demkowicz",
156 "bernstein"};
157
158 PetscBool flg;
159 PetscInt choice_base_value = AINSWORTH;
160 CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, "-base", list_bases, LASBASETOP,
161 &choice_base_value, &flg);
162 if (flg != PETSC_TRUE)
163 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "base not set");
165 if (choice_base_value == AINSWORTH)
167 if (choice_base_value == AINSWORTH_LOBATTO)
169 else if (choice_base_value == DEMKOWICZ)
171 else if (choice_base_value == BERNSTEIN)
173
174 enum spaces { H1SPACE, L2SPACE, HCURLSPACE, HDIVSPACE, LASBASETSPACE };
175 const char *list_spaces[] = {"h1", "l2", "hcurl", "hdiv"};
176 PetscInt choice_space_value = H1SPACE;
177 CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, "-space", list_spaces,
178 LASBASETSPACE, &choice_space_value, &flg);
179 if (flg != PETSC_TRUE)
180 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "space not set");
181 space = H1;
182 if (choice_space_value == H1SPACE)
183 space = H1;
184 else if (choice_space_value == L2SPACE)
185 space = L2;
186 else if (choice_space_value == HCURLSPACE)
187 space = HCURL;
188 else if (choice_space_value == HDIVSPACE)
189 space = HDIV;
190
192
193 int order = 3;
194 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
197
199}
200//! [Set up problem]
201
202//! [Set integration rule]
206}
207//! [Set integration rule]
208
209//! [Create common data]
211//! [Create common data]
212
213//! [Boundary condition]
215//! [Boundary condition]
216
217//! [Push operators to pipeline]
219//! [Push operators to pipeline]
220
221//! [Solve]
223
224//! [Solve]
227
229 pipeline_mng->getDomainLhsFE().reset();
230
231 auto post_proc_fe = boost::make_shared<MyPostProc>(mField);
232 post_proc_fe->generateReferenceElementMesh();
233 pipeline_mng->getDomainRhsFE() = post_proc_fe;
234
235 if (SPACE_DIM == 2) {
236 if (space == HCURL) {
237 auto jac_ptr = boost::make_shared<MatrixDouble>();
238 post_proc_fe->getOpPtrVector().push_back(
239 new OpCalculateHOJacForFace(jac_ptr));
240 post_proc_fe->getOpPtrVector().push_back(new OpMakeHdivFromHcurl());
241 post_proc_fe->getOpPtrVector().push_back(
243 }
244 }
245
246 switch (space) {
247 case H1:
248 case L2:
249
250 {
251
253
254 auto u_ptr = boost::make_shared<VectorDouble>();
255 post_proc_fe->getOpPtrVector().push_back(
256 new OpCalculateScalarFieldValues("U", u_ptr));
257 post_proc_fe->getOpPtrVector().push_back(
258
259 new OpPPMap(
260
261 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
262
263 {{"U", u_ptr}},
264
265 {},
266
267 {},
268
269 {}
270
271 )
272
273 );
274 } break;
275 case HCURL:
276 case HDIV:
277
278 {
280
281 auto u_ptr = boost::make_shared<MatrixDouble>();
282 post_proc_fe->getOpPtrVector().push_back(
283 new OpCalculateHVecVectorField<3>("U", u_ptr));
284
285 post_proc_fe->getOpPtrVector().push_back(
286
287 new OpPPMap(
288
289 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
290
291 {},
292
293 {{"U", u_ptr}},
294
295 {},
296
297 {}
298
299 )
300
301 );
302 } break;
303 default:
304 break;
305 }
306
307 auto scale_tag_val = [&]() {
309 auto &post_proc_mesh = post_proc_fe->getPostProcMesh();
310 Range nodes;
311 CHKERR post_proc_mesh.get_entities_by_type(0, MBVERTEX, nodes);
312 Tag th;
313 CHKERR post_proc_mesh.tag_get_handle("U", th);
314 int length;
315 CHKERR post_proc_mesh.tag_get_length(th, length);
316 std::vector<double> data(nodes.size() * length);
317 CHKERR post_proc_mesh.tag_get_data(th, nodes, &*data.begin());
318 double max_v = 0;
319 for (int i = 0; i != nodes.size(); ++i) {
320 double v = 0;
321 for (int d = 0; d != length; ++d)
322 v += pow(data[length * i + d], 2);
323 v = std::sqrt(v);
324 max_v = std::max(max_v, v);
325 }
326 for (auto &v : data)
327 v /= max_v;
328 CHKERR post_proc_mesh.tag_set_data(th, nodes, &*data.begin());
330 };
331
332 size_t nb = 0;
333 auto dofs_ptr = mField.get_dofs();
334
335 for (auto dof_ptr : (*dofs_ptr)) {
336 MOFEM_LOG("PLOTBASE", Sev::verbose) << *dof_ptr;
337 auto &val = const_cast<double &>(dof_ptr->getFieldData());
338 val = 1;
339 CHKERR pipeline_mng->loopFiniteElements();
340 CHKERR scale_tag_val();
341 CHKERR post_proc_fe->writeFile(
342 "out_base_dof_" + boost::lexical_cast<std::string>(nb) + ".h5m");
343 CHKERR post_proc_fe->getPostProcMesh().delete_mesh();
344 val = 0;
345 ++nb;
346 };
347
349}
350//! [Postprocess results]
351
352//! [Check results]
354//! [Check results]
355
356int main(int argc, char *argv[]) {
357
358 // Initialisation of MoFEM/PETSc and MOAB data structures
359 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
360
361 try {
362
363 //! [Register MoFEM discrete manager in PETSc]
364 DMType dm_name = "DMMOFEM";
365 CHKERR DMRegister_MoFEM(dm_name);
366 //! [Register MoFEM discrete manager in PETSc
367
368 // Add logging channel for example
369 auto core_log = logging::core::get();
370 core_log->add_sink(
372 LogManager::setLog("PLOTBASE");
373 MOFEM_LOG_TAG("PLOTBASE", "plotbase");
374
375 //! [Create MoAB]
376 moab::Core mb_instance; ///< mesh database
377 moab::Interface &moab = mb_instance; ///< mesh database interface
378 //! [Create MoAB]
379
380 //! [Create MoFEM]
381 MoFEM::Core core(moab); ///< finite element database
382 MoFEM::Interface &m_field = core; ///< finite element database insterface
383 //! [Create MoFEM]
384
385 //! [Example]
386 Example ex(m_field);
387 CHKERR ex.runProblem();
388 //! [Example]
389 }
391
393
394 return 0;
395}
396
399 moab::Core core_ref;
400 moab::Interface &moab_ref = core_ref;
401
402 char ref_mesh_file_name[255];
403
404 if (SPACE_DIM == 2)
405 strcpy(ref_mesh_file_name, "ref_mesh2d.h5m");
406 else if (SPACE_DIM == 3)
407 strcpy(ref_mesh_file_name, "ref_mesh3d.h5m");
408 else
409 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
410 "Dimension not implemented");
411
412 CHKERR PetscOptionsGetString(PETSC_NULL, "", "-ref_file", ref_mesh_file_name,
413 255, PETSC_NULL);
414 CHKERR moab_ref.load_file(ref_mesh_file_name, 0, "");
415
416 // Get elements
417 Range elems;
418 CHKERR moab_ref.get_entities_by_dimension(0, SPACE_DIM, elems);
419
420 // Add mid-nodes on edges
421 EntityHandle meshset;
422 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
423 CHKERR moab_ref.add_entities(meshset, elems);
424 CHKERR moab_ref.convert_entities(meshset, true, false, false);
425 CHKERR moab_ref.delete_entities(&meshset, 1);
426
427 // Get nodes on the mesh
428 Range elem_nodes;
429 CHKERR moab_ref.get_connectivity(elems, elem_nodes, false);
430
431 // Map node entity and Gauss pint number
432 std::map<EntityHandle, int> nodes_pts_map;
433
434 // Set gauss points coordinates from the reference mesh
435 gaussPts.resize(SPACE_DIM + 1, elem_nodes.size(), false);
436 gaussPts.clear();
437 Range::iterator nit = elem_nodes.begin();
438 for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
439 double coords[3];
440 CHKERR moab_ref.get_coords(&*nit, 1, coords);
441 for (auto d : {0, 1, 2})
442 gaussPts(d, gg) = coords[d];
443 nodes_pts_map[*nit] = gg;
444 }
445
446 if (SPACE_DIM == 2) {
447 // Set size of adjacency matrix (note ho order nodes 3 nodes and 3 nodes on
448 // edges)
449 refEleMap.resize(elems.size(), 3 + 3);
450 } else if (SPACE_DIM == 3) {
451 refEleMap.resize(elems.size(), 4 + 6);
452 }
453
454 // Set adjacency matrix
455 Range::iterator tit = elems.begin();
456 for (int tt = 0; tit != elems.end(); ++tit, ++tt) {
457 const EntityHandle *conn;
458 int num_nodes;
459 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes, false);
460 for (int nn = 0; nn != num_nodes; ++nn) {
461 refEleMap(tt, nn) = nodes_pts_map[conn[nn]];
462 }
463 }
464
466}
467
470
471 const int num_nodes = gaussPts.size2();
472
473 // Calculate shape functions
474
475 switch (numeredEntFiniteElementPtr->getEntType()) {
476 case MBTRI:
477 shapeFunctions.resize(num_nodes, 3);
479 &gaussPts(0, 0), &gaussPts(1, 0), num_nodes);
480 break;
481 case MBQUAD: {
482 shapeFunctions.resize(num_nodes, 4);
483 for (int gg = 0; gg != num_nodes; gg++) {
484 double ksi = gaussPts(0, gg);
485 double eta = gaussPts(1, gg);
486 shapeFunctions(gg, 0) = N_MBQUAD0(ksi, eta);
487 shapeFunctions(gg, 1) = N_MBQUAD1(ksi, eta);
488 shapeFunctions(gg, 2) = N_MBQUAD2(ksi, eta);
489 shapeFunctions(gg, 3) = N_MBQUAD3(ksi, eta);
490 }
491 } break;
492 case MBTET: {
493 shapeFunctions.resize(num_nodes, 8);
495 &gaussPts(0, 0), &gaussPts(1, 0),
496 &gaussPts(2, 0), num_nodes);
497 } break;
498 case MBHEX: {
499 shapeFunctions.resize(num_nodes, 8);
500 for (int gg = 0; gg != num_nodes; gg++) {
501 double ksi = gaussPts(0, gg);
502 double eta = gaussPts(1, gg);
503 double zeta = gaussPts(2, gg);
504 shapeFunctions(gg, 0) = N_MBHEX0(ksi, eta, zeta);
505 shapeFunctions(gg, 1) = N_MBHEX1(ksi, eta, zeta);
506 shapeFunctions(gg, 2) = N_MBHEX2(ksi, eta, zeta);
507 shapeFunctions(gg, 3) = N_MBHEX3(ksi, eta, zeta);
508 shapeFunctions(gg, 4) = N_MBHEX4(ksi, eta, zeta);
509 shapeFunctions(gg, 5) = N_MBHEX5(ksi, eta, zeta);
510 shapeFunctions(gg, 6) = N_MBHEX6(ksi, eta, zeta);
511 shapeFunctions(gg, 7) = N_MBHEX7(ksi, eta, zeta);
512 }
513 } break;
514 default:
515 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
516 "Not implemented element type");
517 }
518
519 // Create physical nodes
520
521 // MoAB interface allowing for creating nodes and elements in the bulk
522 ReadUtilIface *iface;
523 CHKERR getPostProcMesh().query_interface(iface);
524
525 std::vector<double *> arrays; /// pointers to memory allocated by MoAB for
526 /// storing X, Y, and Z coordinates
527 EntityHandle startv; // Starting handle for vertex
528 // Allocate memory for num_nodes, and return starting handle, and access to
529 // memort.
530 CHKERR iface->get_node_coords(3, num_nodes, 0, startv, arrays);
531
532 mapGaussPts.resize(gaussPts.size2());
533 for (int gg = 0; gg != num_nodes; ++gg)
534 mapGaussPts[gg] = startv + gg;
535
536 Tag th;
537 int def_in_the_loop = -1;
538 CHKERR getPostProcMesh().tag_get_handle("NB_IN_THE_LOOP", 1, MB_TYPE_INTEGER,
539 th, MB_TAG_CREAT | MB_TAG_SPARSE,
540 &def_in_the_loop);
541
542 // Create physical elements
543
544 const int num_el = refEleMap.size1();
545 const int num_nodes_on_ele = refEleMap.size2();
546
547 EntityHandle starte; // Starting handle to first created element
548 EntityHandle *conn; // Access to MOAB memory with connectivity of elements
549
550 // Create tris/tets in the bulk in MoAB database
551 if (SPACE_DIM == 2)
552 CHKERR iface->get_element_connect(num_el, num_nodes_on_ele, MBTRI, 0,
553 starte, conn);
554 else if (SPACE_DIM == 3)
555 CHKERR iface->get_element_connect(num_el, num_nodes_on_ele, MBTET, 0,
556 starte, conn);
557 else
558 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
559 "Dimension not implemented");
560
561 // At this point elements (memory for elements) is allocated, at code bellow
562 // actual connectivity of elements is set.
563 for (unsigned int tt = 0; tt != refEleMap.size1(); ++tt) {
564 for (int nn = 0; nn != num_nodes_on_ele; ++nn)
565 conn[num_nodes_on_ele * tt + nn] = mapGaussPts[refEleMap(tt, nn)];
566 }
567
568 // Finalise elements creation. At that point MOAB updates adjacency tables,
569 // and elements are ready to use.
570 CHKERR iface->update_adjacencies(starte, num_el, num_nodes_on_ele, conn);
571
572 auto physical_elements = Range(starte, starte + num_el - 1);
573 CHKERR getPostProcMesh().tag_clear_data(th, physical_elements, &(nInTheLoop));
574
575 EntityHandle fe_ent = numeredEntFiniteElementPtr->getEnt();
576 int fe_num_nodes;
577 {
578 const EntityHandle *conn;
579 mField.get_moab().get_connectivity(fe_ent, conn, fe_num_nodes, true);
580 coords.resize(3 * fe_num_nodes, false);
581 CHKERR mField.get_moab().get_coords(conn, fe_num_nodes, &coords[0]);
582 }
583
584 // Set physical coordinates to physical nodes
587 &*shapeFunctions.data().begin());
588
590 arrays[0], arrays[1], arrays[2]);
591 const double *t_coords_ele_x = &coords[0];
592 const double *t_coords_ele_y = &coords[1];
593 const double *t_coords_ele_z = &coords[2];
594 for (int gg = 0; gg != num_nodes; ++gg) {
596 t_coords_ele_x, t_coords_ele_y, t_coords_ele_z);
597 t_coords(i) = 0;
598 for (int nn = 0; nn != fe_num_nodes; ++nn) {
599 t_coords(i) += t_n * t_ele_coords(i);
600 for (auto ii : {0, 1, 2})
601 if (std::abs(t_coords(ii)) < std::numeric_limits<float>::epsilon())
602 t_coords(ii) = 0;
603 ++t_ele_coords;
604 ++t_n;
605 }
606 ++t_coords;
607 }
608
610}
611
614 ParallelComm *pcomm_post_proc_mesh =
615 ParallelComm::get_pcomm(coreMeshPtr.get(), MYPCOMM_INDEX);
616 if (pcomm_post_proc_mesh != NULL)
617 delete pcomm_post_proc_mesh;
619};
620
623
624 auto resolve_shared_ents = [&]() {
626
627 ParallelComm *pcomm_post_proc_mesh =
628 ParallelComm::get_pcomm(&(getPostProcMesh()), MYPCOMM_INDEX);
629 if (pcomm_post_proc_mesh == NULL) {
630 // wrapRefMeshComm =
631 // boost::make_shared<WrapMPIComm>(T::mField.get_comm(), false);
632 pcomm_post_proc_mesh = new ParallelComm(
633 &(getPostProcMesh()),
634 PETSC_COMM_WORLD /*(T::wrapRefMeshComm)->get_comm()*/);
635 }
636
637 CHKERR pcomm_post_proc_mesh->resolve_shared_ents(0);
638
640 };
641
642 CHKERR resolve_shared_ents();
643
645}
int main()
Definition: adol-c_atom.cpp:46
constexpr int SPACE_DIM
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
FieldApproximationBase
approximation base
Definition: definitions.h:58
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
FieldSpace
approximation spaces
Definition: definitions.h:82
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_IMPOSSIBLE_CASE
Definition: definitions.h:35
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr int order
#define N_MBQUAD3(x, y)
quad shape function
Definition: fem_tools.h:60
#define N_MBHEX7(x, y, z)
Definition: fem_tools.h:78
#define N_MBHEX3(x, y, z)
Definition: fem_tools.h:74
#define N_MBHEX5(x, y, z)
Definition: fem_tools.h:76
#define N_MBHEX4(x, y, z)
Definition: fem_tools.h:75
#define N_MBHEX0(x, y, z)
Definition: fem_tools.h:71
#define N_MBHEX6(x, y, z)
Definition: fem_tools.h:77
#define N_MBHEX2(x, y, z)
Definition: fem_tools.h:73
#define N_MBQUAD0(x, y)
quad shape function
Definition: fem_tools.h:57
#define N_MBHEX1(x, y, z)
Definition: fem_tools.h:72
#define N_MBQUAD2(x, y)
quad shape function
Definition: fem_tools.h:59
#define N_MBQUAD1(x, y)
quad shape function
Definition: fem_tools.h:58
double eta
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:389
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
FTensor::Index< 'i', SPACE_DIM > i
const 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:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
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
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)
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:13
double zeta
Viscous hardening.
Definition: plastic.cpp:177
static char help[]
Definition: plot_base.cpp:12
constexpr int SPACE_DIM
Definition: plot_base.cpp:26
[Example]
Definition: plastic.cpp:226
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode setIntegrationRules()
[Set up problem]
Definition: plot_base.cpp:203
FieldApproximationBase base
Definition: plot_base.cpp:68
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
Definition: plot_base.cpp:50
Simple * simpleInterface
Definition: helmholtz.cpp:41
MoFEMErrorCode runProblem()
FieldSpace space
Definition: plot_base.cpp:69
MoFEM::Interface & mField
Definition: plastic.cpp:233
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode rebuild_database(int verb=DEFAULT_VERBOSITY)=0
Clear database and initialize it once again.
Core (interface) class.
Definition: Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
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:298
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:344
Get vector field for H-div approximation.
Get value at integration points for scalar field.
Make Hdiv space from Hcurl space in 2d.
Post post-proc data at points from hash maps.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition: Simple.cpp:194
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:264
void setDim(int dim)
Set the problem dimension.
Definition: Simple.hpp:292
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:473
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:611
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition: Simple.hpp:327
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:680
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const int nb)
Calculate shape functions on triangle.
Definition: Tools.hpp:649
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Definition: plot_base.cpp:34
MoFEMErrorCode postProcess()
Definition: plot_base.cpp:621
MoFEMErrorCode generateReferenceElementMesh()
Definition: plot_base.cpp:397
ublas::matrix< int > refEleMap
Definition: plot_base.cpp:44
MoFEMErrorCode setGaussPts(int order)
Definition: plot_base.cpp:468
MoFEMErrorCode preProcess()
Definition: plot_base.cpp:612
MatrixDouble shapeFunctions
Definition: plot_base.cpp:45