23 {
24
26 "-pc_type lu \n"
27 "-pc_factor_mat_solver_type mumps \n"
28 "-mat_mumps_icntl_20 0 \n"
29 "-ksp_atol 1e-10 \n"
30 "-ksp_rtol 1e-10 \n"
31 "-snes_monitor \n"
32 "-snes_type newtonls \n"
33 "-snes_linesearch_type basic \n"
34 "-snes_max_it 100 \n"
35 "-snes_atol 1e-7 \n"
36 "-snes_rtol 1e-7 \n"
37 "-ts_monitor \n"
38 "-ts_type alpha \n";
39
41 if (!
static_cast<bool>(ifstream(
param_file))) {
42 std::ofstream file(
param_file.c_str(), std::ios::ate);
43 if (file.is_open()) {
45 file.close();
46 }
47 }
48
50
51
52 try {
53
54 moab::Core mb_instance;
55 moab::Interface &moab = mb_instance;
56
57 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab,
MYPCOMM_INDEX);
58 auto moab_comm_wrap =
59 boost::make_shared<WrapMPIComm>(PETSC_COMM_WORLD, false);
60 if (pcomm == NULL)
61 pcomm = new ParallelComm(&moab, moab_comm_wrap->get_comm());
62
63 PetscBool flg = PETSC_TRUE;
67 if (flg != PETSC_TRUE) {
68 SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
69 }
70
73 &flg);
74 if (flg != PETSC_TRUE) {
76 }
77
78
79
80 PetscBool is_partitioned = PETSC_FALSE;
82 &is_partitioned, &flg);
83
84 if (is_partitioned == PETSC_TRUE) {
85
86 const char *option;
87 option = "PARALLEL=BCAST_DELETE;PARALLEL_RESOLVE_SHARED_ENTS;PARTITION="
88 "PARALLEL_PARTITION;";
90 } else {
91 const char *option;
92 option = "";
94 }
95
96
97 Tag th_step_size, th_step;
98 double def_step_size = 1;
99 CHKERR moab.tag_get_handle(
"_STEPSIZE", 1, MB_TYPE_DOUBLE, th_step_size,
100 MB_TAG_CREAT | MB_TAG_MESH, &def_step_size);
101 if (
rval == MB_ALREADY_ALLOCATED)
103
104 int def_step = 1;
105 CHKERR moab.tag_get_handle(
"_STEP", 1, MB_TYPE_INTEGER, th_step,
106 MB_TAG_CREAT | MB_TAG_MESH, &def_step);
107 if (
rval == MB_ALREADY_ALLOCATED)
109
110 const void *tag_data_step_size[1];
112 CHKERR moab.tag_get_by_ptr(th_step_size, &root, 1, tag_data_step_size);
113 double &step_size = *(double *)tag_data_step_size[0];
114 const void *tag_data_step[1];
115 CHKERR moab.tag_get_by_ptr(th_step, &root, 1, tag_data_step);
116 int &step = *(int *)tag_data_step[0];
117
118 CHKERR PetscPrintf(PETSC_COMM_WORLD,
119 "Start step %D and step_size = %6.4e\n", step,
120 step_size);
121
124
125
128 std::vector<BitRefLevel> bit_levels;
131
132 if (step == 1) {
133
134 problem_bit_level = bit_levels.back();
135
136
138 3);
141
143
144
147
148
151
152
153
155 "MESH_NODE_POSITIONS");
156
157
159 "SPATIAL_POSITION");
162 "SPATIAL_POSITION");
164 "ELASTIC", "LAMBDA");
166 "SPATIAL_POSITION");
168 "ELASTIC", "MESH_NODE_POSITIONS");
170
171
173 "LAMBDA");
175 "LAMBDA");
176
178 "LAMBDA");
179
180
182
183
185 "ELASTIC");
187 "ARC_LENGTH");
188
190 "SPRING");
191
192
194 problem_bit_level);
195
196
200
201
202 {
203
205 {
206 const double coords[] = {0, 0, 0};
208 Range range_no_field_vertex;
209 range_no_field_vertex.insert(no_field_vertex);
214 range_no_field_vertex);
215 }
216
218 {
219 CHKERR moab.create_meshset(MESHSET_SET, meshset_fe_arc_length);
220 CHKERR moab.add_entities(meshset_fe_arc_length, &no_field_vertex, 1);
223 }
224
226 meshset_fe_arc_length, "ARC_LENGTH", false);
227 }
228
229
234
239
244
245
248 "SPATIAL_POSITION");
250 "SPATIAL_POSITION");
252 "SPATIAL_POSITION");
254 "NEUMANN_FE", "MESH_NODE_POSITIONS");
256 "NEUMANN_FE");
260 CHKERR moab.get_entities_by_type(it->meshset, MBTRI, tris,
true);
262 "NEUMANN_FE");
263 }
267 CHKERR moab.get_entities_by_type(it->meshset, MBTRI, tris,
true);
269 "NEUMANN_FE");
270 }
271
274 "FORCE_FE");
275 }
276
277
278
279 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr(
281 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr(
283
285 m_field, fe_spring_lhs_ptr, fe_spring_rhs_ptr, "SPATIAL_POSITION",
286 "MESH_NODE_POSITIONS");
287
288 PetscBool linear;
290 &linear);
291
294 CHKERR elastic_materials.setBlocks(elastic.setOfBlocks);
295 CHKERR elastic.addElement(
"ELASTIC",
"SPATIAL_POSITION");
297 false, false);
299 false, false);
301 false, false, false);
302 CHKERR elastic.setOperators(
"SPATIAL_POSITION");
303
304
306 CHKERR post_proc.generateReferenceElementMesh();
308 false);
309 CHKERR post_proc.addFieldValuesPostProc(
"SPATIAL_POSITION");
310 CHKERR post_proc.addFieldValuesPostProc(
"MESH_NODE_POSITIONS");
311 CHKERR post_proc.addFieldValuesGradientPostProc(
"SPATIAL_POSITION");
312 std::map<int, NonlinearElasticElement::BlockData>::iterator sit =
313 elastic.setOfBlocks.begin();
314 for (; sit != elastic.setOfBlocks.end(); sit++) {
316 post_proc.postProcMesh, post_proc.mapGaussPts, "SPATIAL_POSITION",
317 sit->second, post_proc.commonData));
318 }
319
320
322 if (step == 1) {
323
325 "MESH_NODE_POSITIONS");
326 CHKERR m_field.
loop_dofs(
"MESH_NODE_POSITIONS", ent_method_material, 0);
328 "SPATIAL_POSITION");
330 "SPATIAL_POSITION");
332 1., "MESH_NODE_POSITIONS", "SPATIAL_POSITION");
334 "SPATIAL_POSITION");
336 "SPATIAL_POSITION");
337 }
338
339
341
342
344
347
348 if (is_partitioned) {
349 SETERRQ(PETSC_COMM_SELF, 1,
350 "Not implemented, problem with arc-length force multiplayer");
351 } else {
355 }
357
358
363
365
366
369 "ELASTIC_MECHANICS",
COL, &F);
372 Mat Aij;
374 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>("ELASTIC_MECHANICS",
375 &Aij);
376
377 boost::shared_ptr<ArcLengthCtx> arc_ctx = boost::shared_ptr<ArcLengthCtx>(
379
383 CHKERR MatGetLocalSize(Aij, &
m, &
n);
384 boost::scoped_ptr<ArcLengthMatShell> mat_ctx(
386
387 Mat ShellAij;
388 CHKERR MatCreateShell(PETSC_COMM_WORLD,
m,
n,
M,
N, mat_ctx.get(),
389 &ShellAij);
390 CHKERR MatShellSetOperation(ShellAij, MATOP_MULT,
392
393 ArcLengthSnesCtx snes_ctx(m_field, "ELASTIC_MECHANICS", arc_ctx);
394
396
401 CHKERR moab.get_entities_by_type(meshset, MBVERTEX, nodes,
true);
403 node_set.merge(nodes);
404 }
405 PetscPrintf(PETSC_COMM_WORLD, "Nb. nodes in load path: %u\n",
406 node_set.size());
407
409
410 double scaled_reference_load = 1;
411 double *scale_lhs = &(arc_ctx->getFieldData());
412 double *scale_rhs = &(scaled_reference_load);
414 m_field, Aij, arc_ctx->F_lambda, scale_lhs, scale_rhs);
416 neumann_forces.getLoopSpatialFe();
417 if (linear) {
420 }
421 fe_neumann.
uSeF =
true;
423 it)) {
425 }
429 }
430
431 boost::shared_ptr<FEMethod> my_dirichlet_bc =
433 m_field,
"SPATIAL_POSITION", Aij,
D, F));
435 &(my_dirichlet_bc->problemPtr));
437 ->iNitialize();
438
440
441 boost::shared_ptr<ArcLengthCtx> arcPtr;
443
446 : arcPtr(arc_ptr), nodeSet(node_set) {}
447
450
451
455 CHKERR VecGhostUpdateBegin(
snes_f, INSERT_VALUES, SCATTER_FORWARD);
456 CHKERR VecGhostUpdateEnd(
snes_f, INSERT_VALUES, SCATTER_FORWARD);
457 CHKERR VecZeroEntries(arcPtr->F_lambda);
458 CHKERR VecGhostUpdateBegin(arcPtr->F_lambda, INSERT_VALUES,
459 SCATTER_FORWARD);
460 CHKERR VecGhostUpdateEnd(arcPtr->F_lambda, INSERT_VALUES,
461 SCATTER_FORWARD);
462 } break;
463 default:
464 SETERRQ(PETSC_COMM_SELF, 1, "not implemented");
465 }
466
468 }
469
474
475 CHKERR VecGhostUpdateBegin(
snes_f, ADD_VALUES, SCATTER_REVERSE);
476 CHKERR VecGhostUpdateEnd(
snes_f, ADD_VALUES, SCATTER_REVERSE);
479 } break;
480 default:
481 SETERRQ(PETSC_COMM_SELF, 1, "not implemented");
482 }
484 }
485
488 boost::shared_ptr<NumeredDofEntity_multiIndex> numered_dofs_rows =
490 Range::iterator nit = nodeSet.begin();
491 for (; nit != nodeSet.end(); nit++) {
492 NumeredDofEntityByEnt::iterator dit, hi_dit;
493 dit = numered_dofs_rows->get<
Ent_mi_tag>().lower_bound(*nit);
494 hi_dit = numered_dofs_rows->get<
Ent_mi_tag>().upper_bound(*nit);
495 for (; dit != hi_dit; dit++) {
496 PetscPrintf(PETSC_COMM_WORLD, "%s [ %d ] %6.4e -> ", "LAMBDA", 0,
497 arcPtr->getFieldData());
498 PetscPrintf(PETSC_COMM_WORLD, "%s [ %d ] %6.4e\n",
499 dit->get()->getName().c_str(),
500 dit->get()->getDofCoeffIdx(),
501 dit->get()->getFieldData());
502 }
503 }
505 }
506 };
507
508 struct AddLambdaVectorToFInternal :
public FEMethod {
509
510 boost::shared_ptr<ArcLengthCtx> arcPtr;
511 boost::shared_ptr<DirichletSpatialPositionsBc> bC;
512
513 AddLambdaVectorToFInternal(boost::shared_ptr<ArcLengthCtx> &arc_ptr,
514 boost::shared_ptr<FEMethod> &bc)
515 : arcPtr(arc_ptr),
518
522 }
526 }
529 switch (snes_ctx) {
531
532 CHKERR VecGhostUpdateBegin(arcPtr->F_lambda, ADD_VALUES,
533 SCATTER_REVERSE);
534 CHKERR VecGhostUpdateEnd(arcPtr->F_lambda, ADD_VALUES,
535 SCATTER_REVERSE);
536 CHKERR VecAssemblyBegin(arcPtr->F_lambda);
537 CHKERR VecAssemblyEnd(arcPtr->F_lambda);
538 for (std::vector<int>::iterator vit = bC->dofsIndices.begin();
539 vit != bC->dofsIndices.end(); vit++) {
540 CHKERR VecSetValue(arcPtr->F_lambda, *vit, 0, INSERT_VALUES);
541 }
542 CHKERR VecAssemblyBegin(arcPtr->F_lambda);
543 CHKERR VecAssemblyEnd(arcPtr->F_lambda);
544 CHKERR VecDot(arcPtr->F_lambda, arcPtr->F_lambda, &arcPtr->F_lambda2);
545 PetscPrintf(PETSC_COMM_WORLD, "\tFlambda2 = %6.4e\n",
546 arcPtr->F_lambda2);
547
548 CHKERR VecAXPY(snes_f, arcPtr->getFieldData(), arcPtr->F_lambda);
549 PetscPrintf(PETSC_COMM_WORLD, "\tlambda = %6.4e\n",
550 arcPtr->getFieldData());
551 double fnorm;
552 CHKERR VecNorm(snes_f, NORM_2, &fnorm);
553 PetscPrintf(PETSC_COMM_WORLD, "\tfnorm = %6.4e\n", fnorm);
554 } break;
555 default:
556 SETERRQ(PETSC_COMM_SELF, 1, "not implemented");
557 }
559 }
560 };
561
563 AddLambdaVectorToFInternal assemble_F_lambda(arc_ctx, my_dirichlet_bc);
564
565 SNES snes;
566 CHKERR SNESCreate(PETSC_COMM_WORLD, &snes);
567 CHKERR SNESSetApplicationContext(snes, &snes_ctx);
569 CHKERR SNESSetJacobian(snes, ShellAij, Aij,
SnesMat, &snes_ctx);
570 CHKERR SNESSetFromOptions(snes);
571
572 PetscReal my_tol;
574 &flg);
575 if (flg == PETSC_TRUE) {
576 PetscReal atol, rtol, stol;
577 PetscInt maxit, maxf;
578 CHKERR SNESGetTolerances(snes, &atol, &rtol, &stol, &maxit, &maxf);
579 atol = my_tol;
580 rtol = atol * 1e2;
581 CHKERR SNESSetTolerances(snes, atol, rtol, stol, maxit, maxf);
582 }
583
584 KSP ksp;
585 CHKERR SNESGetKSP(snes, &ksp);
586 PC pc;
587 CHKERR KSPGetPC(ksp, &pc);
588 boost::scoped_ptr<PCArcLengthCtx> pc_ctx(
590 CHKERR PCSetType(pc, PCSHELL);
591 CHKERR PCShellSetContext(pc, pc_ctx.get());
594
595 if (flg == PETSC_TRUE) {
596 PetscReal rtol, atol, dtol;
597 PetscInt maxits;
598 CHKERR KSPGetTolerances(ksp, &rtol, &atol, &dtol, &maxits);
599 atol = my_tol * 1e-2;
600 rtol = atol * 1e-2;
601 CHKERR KSPSetTolerances(ksp, rtol, atol, dtol, maxits);
602 }
603
605 snes_ctx.get_loops_to_do_Rhs();
606 snes_ctx.get_preProcess_to_do_Rhs().push_back(my_dirichlet_bc);
607 snes_ctx.get_preProcess_to_do_Rhs().push_back(&pre_post_method);
608 loops_to_do_Rhs.push_back(
610
611 loops_to_do_Rhs.push_back(
613
614
615 loops_to_do_Rhs.push_back(
617
618
619 boost::ptr_map<std::string, EdgeForce> edge_forces;
620 string fe_name_str = "FORCE_FE";
621 edge_forces.insert(fe_name_str,
new EdgeForce(m_field));
623 it)) {
624 CHKERR edge_forces.at(fe_name_str)
625 .addForce("SPATIAL_POSITION", arc_ctx->F_lambda, it->getMeshsetId());
626 }
627 for (boost::ptr_map<std::string, EdgeForce>::iterator eit =
628 edge_forces.begin();
629 eit != edge_forces.end(); eit++) {
630 loops_to_do_Rhs.push_back(
632 }
633
634
635 boost::ptr_map<std::string, NodalForce> nodal_forces;
636
637 nodal_forces.insert(fe_name_str,
new NodalForce(m_field));
639 it)) {
640 CHKERR nodal_forces.at(fe_name_str)
641 .addForce("SPATIAL_POSITION", arc_ctx->F_lambda, it->getMeshsetId());
642 }
643 for (boost::ptr_map<std::string, NodalForce>::iterator fit =
644 nodal_forces.begin();
645 fit != nodal_forces.end(); fit++) {
646 loops_to_do_Rhs.push_back(
648 }
649
650
651 loops_to_do_Rhs.push_back(
653 loops_to_do_Rhs.push_back(
655 snes_ctx.get_postProcess_to_do_Rhs().push_back(&pre_post_method);
656 snes_ctx.get_postProcess_to_do_Rhs().push_back(my_dirichlet_bc);
657
659 snes_ctx.get_loops_to_do_Mat();
660 snes_ctx.get_preProcess_to_do_Mat().push_back(my_dirichlet_bc);
661 loops_to_do_Mat.push_back(
663
664 loops_to_do_Mat.push_back(
666
667 loops_to_do_Mat.push_back(
669 loops_to_do_Mat.push_back(
671 snes_ctx.get_postProcess_to_do_Mat().push_back(my_dirichlet_bc);
672
674 "ELASTIC_MECHANICS",
COL,
D, INSERT_VALUES, SCATTER_FORWARD);
675 CHKERR VecGhostUpdateBegin(
D, INSERT_VALUES, SCATTER_FORWARD);
676 CHKERR VecGhostUpdateEnd(
D, INSERT_VALUES, SCATTER_FORWARD);
677
678 PetscScalar step_size_reduction;
680 &step_size_reduction, &flg);
681 if (flg != PETSC_TRUE) {
682 step_size_reduction = 1.;
683 }
684
685 PetscInt max_steps;
687 &flg);
688 if (flg != PETSC_TRUE) {
689 max_steps = 5;
690 }
691
692 int its_d;
694 &flg);
695 if (flg != PETSC_TRUE) {
696 its_d = 4;
697 }
698 PetscScalar max_reduction = 10, min_reduction = 0.1;
700 &max_reduction, &flg);
702 &min_reduction, &flg);
703
704 double gamma = 0.5, reduction = 1;
705
706 if (step == 1) {
707 step_size = step_size_reduction;
708 } else {
709 reduction = step_size_reduction;
710 step++;
711 }
712 double step_size0 = step_size;
713
714 if (step > 1) {
716 "ELASTIC_MECHANICS",
"SPATIAL_POSITION",
"X0_SPATIAL_POSITION",
COL,
717 arc_ctx->x0, INSERT_VALUES, SCATTER_FORWARD);
718 double x0_nrm;
719 CHKERR VecNorm(arc_ctx->x0, NORM_2, &x0_nrm);
720 CHKERR PetscPrintf(PETSC_COMM_WORLD,
721 "\tRead x0_nrm = %6.4e dlambda = %6.4e\n", x0_nrm,
722 arc_ctx->dLambda);
723 CHKERR arc_ctx->setAlphaBeta(1, 0);
724 } else {
725 CHKERR arc_ctx->setS(step_size);
726 CHKERR arc_ctx->setAlphaBeta(0, 1);
727 }
728
730
733 CHKERR VecDuplicate(arc_ctx->x0, &x00);
734 bool converged_state = false;
735
736 for (int jj = 0; step < max_steps; step++, jj++) {
737
739 CHKERR VecCopy(arc_ctx->x0, x00);
740
741 if (step == 1) {
742
743 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"Load Step %D step_size = %6.4e\n",
744 step, step_size);
745 CHKERR arc_ctx->setS(step_size);
746 CHKERR arc_ctx->setAlphaBeta(0, 1);
747 CHKERR VecCopy(
D, arc_ctx->x0);
748 double dlambda;
749 CHKERR arc_method.calculateInitDlambda(&dlambda);
750 CHKERR arc_method.setDlambdaToX(
D, dlambda);
751
752 } else if (step == 2) {
753
754 CHKERR arc_ctx->setAlphaBeta(1, 0);
755 CHKERR arc_method.calculateDxAndDlambda(
D);
756 step_size = std::sqrt(arc_method.calculateLambdaInt());
757 step_size0 = step_size;
758 CHKERR arc_ctx->setS(step_size);
759 double dlambda = arc_ctx->dLambda;
760 double dx_nrm;
761 CHKERR VecNorm(arc_ctx->dx, NORM_2, &dx_nrm);
762 CHKERR PetscPrintf(PETSC_COMM_WORLD,
763 "Load Step %D step_size = %6.4e dlambda0 = %6.4e "
764 "dx_nrm = %6.4e dx2 = %6.4e\n",
765 step, step_size, dlambda, dx_nrm, arc_ctx->dx2);
766 CHKERR VecCopy(
D, arc_ctx->x0);
767 CHKERR VecAXPY(
D, 1., arc_ctx->dx);
768 CHKERR arc_method.setDlambdaToX(
D, dlambda);
769
770 } else {
771
772 if (jj == 0) {
773 step_size0 = step_size;
774 }
775
776 CHKERR arc_method.calculateDxAndDlambda(
D);
777 step_size *= reduction;
778 if (step_size > max_reduction * step_size0) {
779 step_size = max_reduction * step_size0;
780 } else if (step_size < min_reduction * step_size0) {
781 step_size = min_reduction * step_size0;
782 }
783 CHKERR arc_ctx->setS(step_size);
784 double dlambda = reduction * arc_ctx->dLambda;
785 double dx_nrm;
786 CHKERR VecScale(arc_ctx->dx, reduction);
787 CHKERR VecNorm(arc_ctx->dx, NORM_2, &dx_nrm);
788 CHKERR PetscPrintf(PETSC_COMM_WORLD,
789 "Load Step %D step_size = %6.4e dlambda0 = %6.4e "
790 "dx_nrm = %6.4e dx2 = %6.4e\n",
791 step, step_size, dlambda, dx_nrm, arc_ctx->dx2);
792 CHKERR VecCopy(
D, arc_ctx->x0);
793 CHKERR VecAXPY(
D, 1., arc_ctx->dx);
794 CHKERR arc_method.setDlambdaToX(
D, dlambda);
795 }
796
797 CHKERR SNESSolve(snes, PETSC_NULL,
D);
798 int its;
799 CHKERR SNESGetIterationNumber(snes, &its);
800 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"number of Newton iterations = %D\n",
801 its);
802
803 SNESConvergedReason reason;
804 CHKERR SNESGetConvergedReason(snes, &reason);
805 if (reason < 0) {
806
808 CHKERR VecCopy(x00, arc_ctx->x0);
809
810 double x0_nrm;
811 CHKERR VecNorm(arc_ctx->x0, NORM_2, &x0_nrm);
812 CHKERR PetscPrintf(PETSC_COMM_WORLD,
813 "\tRead x0_nrm = %6.4e dlambda = %6.4e\n", x0_nrm,
814 arc_ctx->dLambda);
815 CHKERR arc_ctx->setAlphaBeta(1, 0);
816
817 reduction = 0.1;
818 converged_state = false;
819
820 continue;
821
822 } else {
823
824 if (step > 1 && converged_state) {
825
826 reduction = pow((double)its_d / (double)(its + 1), gamma);
827 if (step_size >= max_reduction * step_size0 && reduction > 1) {
828 reduction = 1;
829 } else if (step_size <= min_reduction * step_size0 && reduction < 1) {
830 reduction = 1;
831 }
832 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"reduction step_size = %6.4e\n",
833 reduction);
834 }
835
836
838 "ELASTIC_MECHANICS",
COL,
D, INSERT_VALUES, SCATTER_REVERSE);
840 "ELASTIC_MECHANICS",
"SPATIAL_POSITION",
"X0_SPATIAL_POSITION",
COL,
841 arc_ctx->x0, INSERT_VALUES, SCATTER_REVERSE);
842 converged_state = true;
843 }
844
845 if (step % 1 == 0) {
846
847
848
849
850
851
852
853
854
855
857 post_proc);
858 std::ostringstream o1;
859 o1 << "out_" << step << ".h5m";
860 CHKERR post_proc.writeFile(o1.str().c_str());
861 }
862
863 CHKERR pre_post_method.potsProcessLoadPath();
864 }
865
868
869
873 CHKERR MatDestroy(&ShellAij);
874 CHKERR SNESDestroy(&snes);
875 }
877
879
880 return 0;
881}
const std::string default_options
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ NOFIELD
scalar or vector of scalars describe (no true field)
#define MYPCOMM_INDEX
default communicator number PCOMM
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
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 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 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 add_ents_to_finite_element_by_MESHSET(const EntityHandle meshset, const std::string &name, const bool recursive=false)=0
add MESHSET element to finite element database given by name
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 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.
#define _IT_CUBITMESHSETS_BY_NAME_FOR_LOOP_(MESHSET_MANAGER, NAME, IT)
Iterator that loops over Cubit BlockSet having a particular name.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
MoFEMErrorCode partitionProblem(const std::string name, int verb=VERBOSE)
partition problem dofs (collective)
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 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
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
const FTensor::Tensor2< T, Dim, Dim > Vec
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetReal(PetscOptions *, const char pre[], const char name[], PetscReal *dval, PetscBool *set)
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
Store variables for ArcLength analysis.
shell matrix for arc-length method
PetscErrorCode preProcess()
function is run at the beginning of loop
PetscErrorCode postProcess()
function is run at the end of loop
Set Dirichlet boundary conditions on spatial displacements.
Force on edges and lines.
Manage setting parameters and constitutive equations for nonlinear/linear elastic materials.
virtual MoFEMErrorCode postProcess()
function is run at the end of loop
virtual MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
const Problem * problemPtr
raw pointer to problem
virtual moab::Interface & get_moab()=0
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
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.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
structure for User Loop Methods on finite elements
Matrix manager is used to build and partition problems.
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode printForceSet() const
print meshsets with force boundary conditions data structure
MoFEMErrorCode printMaterialsSet() const
print meshsets with material data structure set on it
MoFEMErrorCode printDisplacementSet() const
print meshsets with displacement boundary conditions data structure
auto & getNumeredRowDofsPtr() const
get access to numeredRowDofsPtr storing DOFs on rows
Problem manager is used to build and partition problems.
Projection of edge entities with one mid-node on hierarchical basis.
MoFEM::FEMethodsSequence FEMethodsSequence
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
MoFEMErrorCode addPressure(int ms_id)
MoFEMErrorCode addForce(int ms_id)
NonLinear surface pressure element (obsolete implementation)
structure grouping operators and data used for calculation of nonlinear elastic element
structure for Arc Length pre-conditioner
Implementation of spherical arc-length method.