v0.13.1
Functions | Variables
forces_and_sources_h1_continuity_check.cpp File Reference
#include <MoFEM.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 
static const double face_coords [4][9]
 
static const double edge_coords [6][6]
 

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 34 of file forces_and_sources_h1_continuity_check.cpp.

34 {
35
36 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
37
38 try {
39
40 moab::Core mb_instance;
41 moab::Interface &moab = mb_instance;
42 int rank;
43 MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
44
45 PetscBool flg = PETSC_TRUE;
46 char mesh_file_name[255];
47#if PETSC_VERSION_GE(3, 6, 4)
48 ierr = PetscOptionsGetString(PETSC_NULL, "", "-my_file", mesh_file_name,
49 255, &flg);
51#else
52 ierr = PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
53 mesh_file_name, 255, &flg);
55#endif
56 if (flg != PETSC_TRUE) {
57 SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
58 }
59 const char *option;
60 option = "";
61 CHKERR moab.load_file(mesh_file_name, 0, option);
62
63 // Create MoFEM (Joseph) database
64 MoFEM::Core core(moab);
65 MoFEM::Interface &m_field = core;
66
67 // set entitities bit level
68 BitRefLevel bit_level0;
69 bit_level0.set(0);
70 EntityHandle meshset_level0;
71 CHKERR moab.create_meshset(MESHSET_SET, meshset_level0);
72
73 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
74 0, 3, bit_level0);
75
76
77 // Fields
78 CHKERR m_field.add_field("MESH_NODE_POSITIONS", H1, AINSWORTH_LEGENDRE_BASE,
79 3);
80
81 CHKERR m_field.add_field("H1", H1, AINSWORTH_LEGENDRE_BASE, 1);
82
83
84 // FE
85 CHKERR m_field.add_finite_element("TET_FE");
86 CHKERR m_field.add_finite_element("TRI_FE");
87 CHKERR m_field.add_finite_element("SKIN_FE");
88 CHKERR m_field.add_finite_element("EDGE_FE");
89
90
91 // Define rows/cols and element data
92 CHKERR m_field.modify_finite_element_add_field_row("TET_FE", "H1");
93 CHKERR m_field.modify_finite_element_add_field_col("TET_FE", "H1");
94 CHKERR m_field.modify_finite_element_add_field_data("TET_FE", "H1");
96 "MESH_NODE_POSITIONS");
97
98
99 CHKERR m_field.modify_finite_element_add_field_row("SKIN_FE", "H1");
100 CHKERR m_field.modify_finite_element_add_field_col("SKIN_FE", "H1");
101 CHKERR m_field.modify_finite_element_add_field_data("SKIN_FE", "H1");
103 "MESH_NODE_POSITIONS");
104
105
106 CHKERR m_field.modify_finite_element_add_field_row("TRI_FE", "H1");
107 CHKERR m_field.modify_finite_element_add_field_col("TRI_FE", "H1");
108 CHKERR m_field.modify_finite_element_add_field_data("TRI_FE", "H1");
110 "MESH_NODE_POSITIONS");
111
112
113 CHKERR m_field.modify_finite_element_add_field_row("EDGE_FE", "H1");
114 CHKERR m_field.modify_finite_element_add_field_col("EDGE_FE", "H1");
115 CHKERR m_field.modify_finite_element_add_field_data("EDGE_FE", "H1");
117 "MESH_NODE_POSITIONS");
118
119 // Problem
120 CHKERR m_field.add_problem("TEST_PROBLEM");
121
122
123 // set finite elements for problem
124 CHKERR m_field.modify_problem_add_finite_element("TEST_PROBLEM", "TET_FE");
125 CHKERR m_field.modify_problem_add_finite_element("TEST_PROBLEM", "SKIN_FE");
126 CHKERR m_field.modify_problem_add_finite_element("TEST_PROBLEM", "TRI_FE");
127 CHKERR m_field.modify_problem_add_finite_element("TEST_PROBLEM", "EDGE_FE");
128
129 // set refinement level for problem
130 CHKERR m_field.modify_problem_ref_level_add_bit("TEST_PROBLEM", bit_level0);
131
132
133 // meshset consisting all entities in mesh
134 EntityHandle root_set = moab.get_root_set();
135 // add entities to field
136 CHKERR m_field.add_ents_to_field_by_type(root_set, MBTET, "H1");
137
138
139 // add entities to finite element
140 CHKERR
141 m_field.add_ents_to_finite_element_by_type(root_set, MBTET, "TET_FE");
142
143
144 Range tets;
145 CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
146 BitRefLevel().set(0), BitRefLevel().set(), MBTET, tets);
147
148 Skinner skin(&moab);
149 Range skin_faces; // skin faces from 3d ents
150 CHKERR skin.find_skin(0, tets, false, skin_faces);
151 CHKERR m_field.add_ents_to_finite_element_by_type(skin_faces, MBTRI,
152 "SKIN_FE");
153
154 Range faces;
155 CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
156 BitRefLevel().set(0), BitRefLevel().set(), MBTRI, faces);
157
158 faces = subtract(faces, skin_faces);
159 CHKERR m_field.add_ents_to_finite_element_by_type(faces, MBTRI, "TRI_FE");
160
161 Range edges;
162 CHKERR moab.get_adjacencies(faces, 1, false, edges, moab::Interface::UNION);
163 CHKERR m_field.add_ents_to_finite_element_by_type(edges, MBEDGE, "EDGE_FE");
164
165 Range skin_edges;
166 CHKERR moab.get_adjacencies(skin_faces, 1, false, skin_edges,
167 moab::Interface::UNION);
168 CHKERR m_field.add_ents_to_finite_element_by_type(edges, MBEDGE, "EDGE_FE");
169
170
171 // set app. order
172 int order = 4;
173 CHKERR m_field.set_field_order(root_set, MBTET, "H1", order);
174 CHKERR m_field.set_field_order(root_set, MBTRI, "H1", order);
175 CHKERR m_field.set_field_order(root_set, MBEDGE, "H1", order);
176
177 CHKERR m_field.add_ents_to_field_by_type(0, MBTET, "MESH_NODE_POSITIONS");
178 CHKERR m_field.set_field_order(0, MBVERTEX, "MESH_NODE_POSITIONS", 1);
179 CHKERR m_field.set_field_order(0, MBEDGE, "MESH_NODE_POSITIONS", 2);
180 CHKERR m_field.set_field_order(0, MBTRI, "MESH_NODE_POSITIONS", 2);
181 CHKERR m_field.set_field_order(0, MBTET, "MESH_NODE_POSITIONS", 2);
182
183 /****/
184 // build database
185 // build field
186 CHKERR m_field.build_fields();
187
188 // build finite elemnts
190
191 // build adjacencies
192 CHKERR m_field.build_adjacencies(bit_level0);
193
194
195 ProblemsManager *prb_mng_ptr;
196 CHKERR m_field.getInterface(prb_mng_ptr);
197
198 // build problem
199 CHKERR prb_mng_ptr->buildProblem("TEST_PROBLEM", true);
200
201
202 // project geometry form 10 node tets on higher order approx. functions
203 Projection10NodeCoordsOnField ent_method(m_field, "MESH_NODE_POSITIONS");
204 CHKERR m_field.loop_dofs("MESH_NODE_POSITIONS", ent_method);
205
206 /****/
207 // mesh partitioning
208 // partition
209 CHKERR prb_mng_ptr->partitionSimpleProblem("TEST_PROBLEM");
210 CHKERR prb_mng_ptr->partitionFiniteElements("TEST_PROBLEM");
211
212 // what are ghost nodes, see Petsc Manual
213 CHKERR prb_mng_ptr->partitionGhostDofs("TEST_PROBLEM");
214
215
216 Vec v;
217 CHKERR m_field.getInterface<VecManager>()->vecCreateGhost("TEST_PROBLEM",
218 ROW, &v);
219 CHKERR VecSetRandom(v, PETSC_NULL);
220
221 CHKERR m_field.getInterface<VecManager>()->setLocalGhostVector(
222 "TEST_PROBLEM", ROW, v, INSERT_VALUES, SCATTER_REVERSE);
223
224 CHKERR VecDestroy(&v);
225
226
227 typedef tee_device<std::ostream, std::ofstream> TeeDevice;
228 typedef stream<TeeDevice> TeeStream;
229 std::ofstream ofs("forces_and_sources_hdiv_continuity_check.txt");
230 TeeDevice my_tee(std::cout, ofs);
231 TeeStream my_split(my_tee);
232
233 struct OpTetFluxes
235
236 MoFEM::Interface &mField;
237 Tag tH;
238
239 OpTetFluxes(MoFEM::Interface &m_field, Tag th)
241 "H1", UserDataOperator::OPROW),
242 mField(m_field), tH(th) {}
243
244 MoFEMErrorCode doWork(int side, EntityType type,
247
248 if (data.getFieldData().size() == 0)
250
251 if (type == MBTRI) {
252
253 boost::shared_ptr<const NumeredEntFiniteElement> mofem_fe =
254 getNumeredEntFiniteElementPtr();
255 SideNumber_multiIndex &side_table = mofem_fe->getSideNumberTable();
256 EntityHandle face = side_table.get<1>()
257 .find(boost::make_tuple(type, side))
258 ->get()
259 ->ent;
260 int sense = side_table.get<1>()
261 .find(boost::make_tuple(type, side))
262 ->get()
263 ->sense;
264
265 double t = 0;
266 int nb_dofs = data.getN().size2();
267 for (int dd = 0; dd < nb_dofs; dd++) {
268 t += data.getN(side)(dd) * data.getFieldData()[dd];
269 }
270
271 double *t_ptr;
272 CHKERR mField.get_moab().tag_get_by_ptr(tH, &face, 1,
273 (const void **)&t_ptr);
274
275 *t_ptr += sense * t;
276 }
277
278 if (type == MBEDGE) {
279
280 boost::shared_ptr<const NumeredEntFiniteElement> mofem_fe =
281 getNumeredEntFiniteElementPtr();
282 SideNumber_multiIndex &side_table = mofem_fe->getSideNumberTable();
283 EntityHandle edge = side_table.get<1>()
284 .find(boost::make_tuple(type, side))
285 ->get()
286 ->ent;
287 Range adj_tets;
288 CHKERR mField.get_moab().get_adjacencies(&edge, 1, 3, false, adj_tets,
289 moab::Interface::UNION);
290 const int nb_adj_tets = adj_tets.size();
291
292 double t = 0;
293 int nb_dofs = data.getN().size2();
294 for (int dd = 0; dd < nb_dofs; dd++) {
295 t += data.getN(4 + side)(dd) * data.getFieldData()[dd];
296 }
297
298 double *t_ptr;
299 CHKERR mField.get_moab().tag_get_by_ptr(tH, &edge, 1,
300 (const void **)&t_ptr);
301
302 *t_ptr += t / nb_adj_tets;
303 }
304
306 }
307 };
308
309 struct MyTetFE : public VolumeElementForcesAndSourcesCore {
310
311 MyTetFE(MoFEM::Interface &m_field)
313 int getRule(int order) { return -1; };
314
315 MatrixDouble N_tri;
317
319
320 try {
321
322 N_tri.resize(1, 3);
323 CHKERR ShapeMBTRI(&N_tri(0, 0), G_TRI_X1, G_TRI_Y1, 1);
324
325
326 gaussPts.resize(4, 4 + 6);
327 int ff = 0;
328 for (; ff < 4; ff++) {
329 int dd = 0;
330 for (; dd < 3; dd++) {
331 gaussPts(dd, ff) =
332 cblas_ddot(3, &N_tri(0, 0), 1, &face_coords[ff][dd], 3);
333 }
334 gaussPts(3, ff) = G_TRI_W1[0];
335 }
336
337 int ee = 0;
338 for (; ee < 6; ee++) {
339 int dd = 0;
340 for (; dd < 3; dd++) {
341 gaussPts(dd, 4 + ee) =
342 (edge_coords[ee][0 + dd] + edge_coords[ee][3 + dd]) * 0.5;
343 }
344 gaussPts(3, 4 + ee) = 1;
345 }
346
347 // std::cerr << gaussPts << std::endl;
348
349 } catch (std::exception &ex) {
350 std::ostringstream ss;
351 ss << "thorw in method: " << ex.what() << " at line " << __LINE__
352 << " in file " << __FILE__;
353 SETERRQ(PETSC_COMM_SELF, MOFEM_STD_EXCEPTION_THROW, ss.str().c_str());
354 }
355
357 }
358 };
359
360 struct OpFacesSkinFluxes
362
363 MoFEM::Interface &mField;
364 Tag tH1, tH2;
365 TeeStream &mySplit;
366
367 OpFacesSkinFluxes(MoFEM::Interface &m_field, Tag th1, Tag th2,
368 TeeStream &my_split)
370 "H1", UserDataOperator::OPROW),
371 mField(m_field), tH1(th1), tH2(th2), mySplit(my_split) {}
372
373 MoFEMErrorCode doWork(int side, EntityType type,
376
377 if (type != MBTRI)
379 EntityHandle face = getNumeredEntFiniteElementPtr()->getEnt();
380
381 double *t_ptr;
382 CHKERR mField.get_moab().tag_get_by_ptr(tH1, &face, 1,
383 (const void **)&t_ptr);
384
385 double *tn_ptr;
386 CHKERR mField.get_moab().tag_get_by_ptr(tH2, &face, 1,
387 (const void **)&tn_ptr);
388
389
390 *tn_ptr = *t_ptr;
391
392 int dd = 0;
393 int nb_dofs = data.getN().size2();
394 for (; dd < nb_dofs; dd++) {
395 double val = data.getFieldData()[dd];
396 *tn_ptr += -data.getN()(0, dd) * val;
397 }
398
399 const double eps = 1e-8;
400 if (fabs(*tn_ptr) > eps) {
401 SETERRQ1(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
402 "H1 continuity failed %6.4e", *tn_ptr);
403 }
404
405 mySplit.precision(5);
406 mySplit << face << " " << /*std::fixed <<*/ fabs(*tn_ptr) << std::endl;
407
409 }
410 };
411
412 struct OpFacesFluxes
414
415 MoFEM::Interface &mField;
416 Tag tH1, tH2;
417 TeeStream &mySplit;
418
419 OpFacesFluxes(MoFEM::Interface &m_field, Tag _th1, Tag _th2,
420 TeeStream &my_split)
422 "H1", UserDataOperator::OPROW),
423 mField(m_field), tH1(_th1), tH2(_th2), mySplit(my_split) {}
424
428
429 if (type != MBTRI)
431
432 EntityHandle face = getNumeredEntFiniteElementPtr()->getEnt();
433
434 double *t_ptr;
435 CHKERR mField.get_moab().tag_get_by_ptr(tH1, &face, 1,
436 (const void **)&t_ptr);
437
438 double *tn_ptr;
439 CHKERR mField.get_moab().tag_get_by_ptr(tH2, &face, 1,
440 (const void **)&tn_ptr);
441
442
443 *tn_ptr = *t_ptr;
444
445 const double eps = 1e-8;
446 if (fabs(*tn_ptr) > eps) {
447 SETERRQ1(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
448 "H1 continuity failed %6.4e", *tn_ptr);
449 }
450
451 mySplit.precision(5);
452
453 mySplit << face << " " << /*std::fixed <<*/ fabs(*tn_ptr) << std::endl;
454
456 }
457 };
458
459 struct MyTriFE : public FaceElementForcesAndSourcesCore {
460
461 MyTriFE(MoFEM::Interface &m_field)
463 int getRule(int order) { return -1; };
464
467
468 gaussPts.resize(3, 1);
469 gaussPts(0, 0) = G_TRI_X1[0];
470 gaussPts(1, 0) = G_TRI_Y1[0];
471 gaussPts(2, 0) = G_TRI_W1[0];
472
474 }
475 };
476
477 struct OpEdgesFluxes
479
480 MoFEM::Interface &mField;
481 Tag tH1, tH2;
482 TeeStream &mySplit;
483
484 OpEdgesFluxes(MoFEM::Interface &m_field, Tag _th1, Tag _th2,
485 TeeStream &my_split)
487 "H1", UserDataOperator::OPROW),
488 mField(m_field), tH1(_th1), tH2(_th2), mySplit(my_split) {}
489
490 MoFEMErrorCode doWork(int side, EntityType type,
493
494 if (type != MBEDGE)
496
497 EntityHandle edge = getNumeredEntFiniteElementPtr()->getEnt();
498
499 double *t_ptr;
500 CHKERR mField.get_moab().tag_get_by_ptr(tH1, &edge, 1,
501 (const void **)&t_ptr);
502
503 double *tn_ptr;
504 CHKERR mField.get_moab().tag_get_by_ptr(tH2, &edge, 1,
505 (const void **)&tn_ptr);
506
507
508 *tn_ptr = *t_ptr;
509
510 double tn = 0;
511 int nb_dofs = data.getN().size2();
512 int dd = 0;
513 for (; dd < nb_dofs; dd++) {
514 double val = data.getFieldData()[dd];
515 tn += data.getN()(0, dd) * val;
516 }
517
518 // mySplit << *tn_ptr << " " << tn << " " << getLength() << endl;
519 *tn_ptr -= tn;
520
521 // mySplit << getTangetAtGaussPts() << " " << getDirection() << endl;
522
523 // cerr << t_ptr[0] << " " << t_ptr[1] << " " << t_ptr[2] << endl;
524
525 // const double eps = 1e-8;
526 // if(fabs(*tn_ptr)>eps) {
527 // SETERRQ1(
528 // PETSC_COMM_SELF,
529 // MOFEM_ATOM_TEST_INVALID,
530 // "H1 continuity failed %6.4e",
531 // *tn_ptr
532 // );
533 // }
534
535 mySplit.precision(5);
536
537 mySplit << edge << " " << /*std::fixed <<*/ fabs(*tn_ptr) << std::endl;
538
540 }
541 };
542
543 struct MyEdgeFE : public EdgeElementForcesAndSourcesCore {
544
545 MyEdgeFE(MoFEM::Interface &m_field)
547 int getRule(int order) { return -1; };
548
551
552 gaussPts.resize(2, 1);
553 gaussPts(0, 0) = 0.5;
554 gaussPts(1, 0) = 1;
555
557 }
558 };
559
560 MyTetFE tet_fe(m_field);
561 MyTriFE tri_fe(m_field);
562 MyTriFE skin_fe(m_field);
563 MyEdgeFE edge_fe(m_field);
564
565 Tag th1;
566 double def_val[] = {0, 0, 0};
567 CHKERR moab.tag_get_handle("T", 3, MB_TYPE_DOUBLE, th1,
568 MB_TAG_CREAT | MB_TAG_SPARSE, &def_val);
569
570 tet_fe.getOpPtrVector().push_back(
571 new OpCalculateHOCoords("MESH_NODE_POSITIONS"));
572 tet_fe.getOpPtrVector().push_back(new OpTetFluxes(m_field, th1));
573
574 Tag th2;
575 CHKERR moab.tag_get_handle("TN", 1, MB_TYPE_DOUBLE, th2,
576 MB_TAG_CREAT | MB_TAG_SPARSE, &def_val);
577
578 tri_fe.getOpPtrVector().push_back(
579 new OpCalculateHOCoords("MESH_NODE_POSITIONS"));
580 tri_fe.getOpPtrVector().push_back(
581 new OpGetHONormalsOnFace("MESH_NODE_POSITIONS"));
582 tri_fe.getOpPtrVector().push_back(
583 new OpFacesFluxes(m_field, th1, th2, my_split));
584 skin_fe.getOpPtrVector().push_back(
585 new OpFacesSkinFluxes(m_field, th1, th2, my_split));
586
587 edge_fe.getOpPtrVector().push_back(
588 new OpGetHOTangentsOnEdge("MESH_NODE_POSITIONS"));
589 edge_fe.getOpPtrVector().push_back(
590 new OpCalculateHOCoords("MESH_NODE_POSITIONS"));
591 edge_fe.getOpPtrVector().push_back(
592 new OpEdgesFluxes(m_field, th1, th2, my_split));
593
594 for (Range::iterator fit = faces.begin(); fit != faces.end(); fit++) {
595 CHKERR moab.tag_set_data(th1, &*fit, 1, &def_val);
596
597 CHKERR moab.tag_set_data(th2, &*fit, 1, &def_val);
598
599 }
600
601 CHKERR m_field.loop_finite_elements("TEST_PROBLEM", "TET_FE", tet_fe);
602 my_split << "internal\n";
603 CHKERR m_field.loop_finite_elements("TEST_PROBLEM", "TRI_FE", tri_fe);
604 my_split << "skin\n";
605 CHKERR m_field.loop_finite_elements("TEST_PROBLEM", "SKIN_FE", skin_fe);
606 my_split << "edges\n";
607 CHKERR m_field.loop_finite_elements("TEST_PROBLEM", "EDGE_FE", edge_fe);
608
609
610 EntityHandle meshset;
611 CHKERR moab.create_meshset(MESHSET_SET, meshset);
612 CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
613 BitRefLevel().set(0), BitRefLevel().set(), MBTRI, meshset);
614
615 CHKERR moab.write_file("out.vtk", "VTK", "", &meshset, 1);
616
617 }
619
621}
ForcesAndSourcesCore::UserDataOperator UserDataOperator
static const double eps
@ ROW
Definition: definitions.h:136
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
@ MOFEM_STD_EXCEPTION_THROW
Definition: definitions.h:52
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:53
#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
static const double G_TRI_Y1[]
Definition: fem_tools.h:323
static const double G_TRI_W1[]
Definition: fem_tools.h:324
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:194
static const double G_TRI_X1[]
Definition: fem_tools.h:322
static const double face_coords[4][9]
static const double edge_coords[6][6]
tee_device< std::ostream, std::ofstream > TeeDevice
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
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
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
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
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
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.
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
virtual MoFEMErrorCode loop_finite_elements(const std::string &problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
MoFEMErrorCode partitionSimpleProblem(const std::string name, int verb=VERBOSE)
partition problem dofs
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
MoFEMErrorCode partitionFiniteElements(const std::string name, bool part_from_moab=false, int low_proc=-1, int hi_proc=-1, int verb=VERBOSE)
partition finite elements
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string &name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
double v
phase velocity of light in medium (cm/ns)
char mesh_file_name[255]
const FTensor::Tensor2< T, Dim, Dim > Vec
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
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
CoreTmp< 0 > Core
Definition: Core.hpp:1096
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
constexpr double t
plate stiffness
Definition: plate.cpp:76
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
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.
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)
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorDouble & getFieldData() const
get dofs values
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
MatrixDouble gaussPts
Matrix of integration points.
Calculate HO coordinates at gauss points.
Calculate normals at Gauss points of triangle element.
Calculate tangent vector on edge form HO geometry approximation.
Problem manager is used to build and partition problems.
Projection of edge entities with one mid-node on hierarchical basis.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:33
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)

Variable Documentation

◆ edge_coords

const double edge_coords[6][6]
static
Initial value:
= {
{0, 0, 0, 1, 0, 0}, {1, 0, 0, 0, 1, 0}, {0, 1, 0, 0, 0, 0},
{0, 0, 0, 0, 0, 1}, {1, 0, 0, 0, 0, 1}, {0, 1, 0, 0, 0, 1}}
Examples
prism_elements_from_surface.cpp.

Definition at line 30 of file forces_and_sources_h1_continuity_check.cpp.

◆ face_coords

const double face_coords[4][9]
static
Initial value:
= {{0, 0, 0, 1, 0, 0, 0, 0, 1},
{1, 0, 0, 0, 1, 0, 0, 0, 1},
{0, 0, 0, 0, 1, 0, 0, 0, 1},
{0, 0, 0, 1, 0, 0, 0, 1, 0}}

Definition at line 25 of file forces_and_sources_h1_continuity_check.cpp.

◆ help

char help[] = "...\n\n"
static

Definition at line 23 of file forces_and_sources_h1_continuity_check.cpp.