136 {
137
139
140 try {
141
144 int rank;
145 MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
146
147
148 PetscBool flg = PETSC_TRUE;
151 if(flg != PETSC_TRUE) {
152 SETERRQ(PETSC_COMM_SELF,
MOFEM_NOT_FOUND,
"*** ERROR -my_file (MESH FILE NEEDED)");
153 }
154
157 if(flg != PETSC_TRUE) {
159 }
160
165 if(flg != PETSC_TRUE) {
167 }
168
169
170 const char *option;
171 option = "";
173 ParallelComm* pcomm = ParallelComm::get_pcomm(&moab,
MYPCOMM_INDEX);
174 if(pcomm == NULL) pcomm = new ParallelComm(&moab,PETSC_COMM_WORLD);
175
176
177 PetscLogDouble t1,t2;
178 PetscLogDouble v1,v2;
179 ierr = PetscTime(&v1); CHKERRQ(
ierr);
180 ierr = PetscGetCPUTime(&t1); CHKERRQ(
ierr);
181
182
185
186 vector<BitRefLevel> bit_levels;
187 {
188 Tag th_meshset_info;
189 int def_meshset_info[2] = {0,0};
190 rval = moab.tag_get_handle(
191 "MESHSET_INFO",2,MB_TYPE_INTEGER,th_meshset_info,MB_TAG_CREAT|MB_TAG_SPARSE,&def_meshset_info
192 );
193 int meshset_data[2];
194 EntityHandle root = moab.get_root_set();
196 if(meshset_data[0]==0) {
197 meshset_data[0] = 1;
199
200 }
201 bit_levels.push_back(
BitRefLevel().set(meshset_data[0]-1));
202 }
205
206
209
210 Range preriodic_prisms;
212
213
214
215
216
217 Range SurTrisNeg;
218 ierr = m_field.get_cubit_msId_entities_by_dimension(101,
SIDESET,2,SurTrisNeg,
true); CHKERRQ(
ierr);
219 ierr = PetscPrintf(PETSC_COMM_WORLD,
"number of SideSet 101 = %d\n",SurTrisNeg.size()); CHKERRQ(
ierr);
221 double TriCen[3], coords_Tri[9];
222
223
224 double roundfact=1e3;
225
226 for(Range::iterator it = SurTrisNeg.begin(); it!=SurTrisNeg.end(); it++) {
227 const EntityHandle* conn_face; int num_nodes_Tri;
228
229
231
233
234
235 TriCen[0]= (coords_Tri[0]+coords_Tri[3]+coords_Tri[6])/3.0;
236 TriCen[1]= (coords_Tri[1]+coords_Tri[4]+coords_Tri[7])/3.0;
237 TriCen[2]= (coords_Tri[2]+coords_Tri[5]+coords_Tri[8])/3.0;
238
239
240 TriCen[0]=round(TriCen[0]*roundfact)/roundfact;
241 TriCen[1]=round(TriCen[1]*roundfact)/roundfact;
242 TriCen[2]=round(TriCen[2]*roundfact)/roundfact;
243
244
245
246 Face_CenPos_Handle_varNeg.insert(
Face_CenPos_Handle(TriCen[0], TriCen[1], TriCen[2], *it));
247
248
249 }
250
251
252 Range SurTrisPos;
253 ierr = m_field.get_cubit_msId_entities_by_dimension(102,
SIDESET,2,SurTrisPos,
true); CHKERRQ(
ierr);
254 ierr = PetscPrintf(PETSC_COMM_WORLD,
"number of SideSet 102 = %d\n",SurTrisPos.size()); CHKERRQ(
ierr);
255 for(Range::iterator it = SurTrisPos.begin(); it!=SurTrisPos.end(); it++) {
256 const EntityHandle* conn_face; int num_nodes_Tri;
257
258
260
262
263
264 TriCen[0]= (coords_Tri[0]+coords_Tri[3]+coords_Tri[6])/3.0;
265 TriCen[1]= (coords_Tri[1]+coords_Tri[4]+coords_Tri[7])/3.0;
266 TriCen[2]= (coords_Tri[2]+coords_Tri[5]+coords_Tri[8])/3.0;
267
268
269 TriCen[0]=round(TriCen[0]*roundfact)/roundfact;
270 TriCen[1]=round(TriCen[1]*roundfact)/roundfact;
271 TriCen[2]=round(TriCen[2]*roundfact)/roundfact;
272
273
274
275 Face_CenPos_Handle_varPos.insert(
Face_CenPos_Handle(TriCen[0], TriCen[1], TriCen[2], *it));
276 }
277
278
279 double XcoordMin, YcoordMin, ZcoordMin, XcoordMax, YcoordMax, ZcoordMax;
280 typedef Face_CenPos_Handle_multiIndex::index<xcoord_tag>::type::iterator Tri_Xcoord_iterator;
281 typedef Face_CenPos_Handle_multiIndex::index<ycoord_tag>::type::iterator Tri_Ycoord_iterator;
282 typedef Face_CenPos_Handle_multiIndex::index<zcoord_tag>::type::iterator Tri_Zcoord_iterator;
283 Tri_Xcoord_iterator XcoordMin_it, XcoordMax_it;
284 Tri_Ycoord_iterator YcoordMin_it, YcoordMax_it;
285 Tri_Zcoord_iterator ZcoordMin_it, ZcoordMax_it;
286
287
288 XcoordMin_it=Face_CenPos_Handle_varNeg.get<
xcoord_tag>().begin(); XcoordMin=XcoordMin_it->xcoord;
289 XcoordMax_it=Face_CenPos_Handle_varPos.get<
xcoord_tag>().end(); XcoordMax_it--; XcoordMax=XcoordMax_it->xcoord;
290 YcoordMin_it=Face_CenPos_Handle_varNeg.get<
ycoord_tag>().begin(); YcoordMin=YcoordMin_it->ycoord;
291 YcoordMax_it=Face_CenPos_Handle_varPos.get<
ycoord_tag>().end(); YcoordMax_it--; YcoordMax=YcoordMax_it->ycoord;
292 ZcoordMin_it=Face_CenPos_Handle_varNeg.get<
zcoord_tag>().begin(); ZcoordMin=ZcoordMin_it->zcoord;
293 ZcoordMax_it=Face_CenPos_Handle_varPos.get<
zcoord_tag>().end(); ZcoordMax_it--; ZcoordMax=ZcoordMax_it->zcoord;
294
295
296
297
298
299
300
301 typedef Face_CenPos_Handle_multiIndex::index<Tri_Hand_tag>::type::iterator Tri_Hand_iterator;
302 Tri_Hand_iterator Tri_Neg;
303 typedef Face_CenPos_Handle_multiIndex::index<Composite_xyzcoord>::type::iterator xyzcoord_iterator;
304 xyzcoord_iterator Tri_Pos;
305 double XPos, YPos, ZPos;
306
307
308
309 for(Range::iterator it = SurTrisNeg.begin(); it!=SurTrisNeg.end(); it++) {
310
311 Tri_Neg=Face_CenPos_Handle_varNeg.get<
Tri_Hand_tag>().find(*it);
312
313
314
315 if(Tri_Neg->xcoord==XcoordMin){XPos=XcoordMax; YPos=Tri_Neg->ycoord; ZPos=Tri_Neg->zcoord;};
316 if(Tri_Neg->ycoord==YcoordMin){XPos=YPos=Tri_Neg->xcoord; YPos=YcoordMax; ZPos=Tri_Neg->zcoord;};
317 if(Tri_Neg->zcoord==ZcoordMin){XPos=YPos=Tri_Neg->xcoord; YPos=Tri_Neg->ycoord; ZPos=ZcoordMax; };
318 Tri_Pos=Face_CenPos_Handle_varPos.get<
Composite_xyzcoord>().find(boost::make_tuple(XPos, YPos, ZPos));
319
320
321
322
323
324 EntityHandle PrismNodes[6];
325 vector<EntityHandle> TriNodesNeg, TriNodesPos;
326 double CoordNodeNeg[9], CoordNodePos[9];
331 for(int ii=0; ii<3; ii++){
332 PrismNodes[ii]=TriNodesNeg[ii];
333 }
334
335
336
337
338
339
340
341
342 double XNodeNeg, YNodeNeg, ZNodeNeg, XNodePos, YNodePos, ZNodePos;
343 for(int ii=0; ii<3; ii++){
344 if(Tri_Neg->xcoord==XcoordMin){XNodeNeg=XcoordMax; YNodeNeg=CoordNodeNeg[3*ii+1]; ZNodeNeg=CoordNodeNeg[3*ii+2];};
345 if(Tri_Neg->ycoord==YcoordMin){XNodeNeg=CoordNodeNeg[3*ii]; YNodeNeg=YcoordMax; ZNodeNeg=CoordNodeNeg[3*ii+2];};
346 if(Tri_Neg->zcoord==ZcoordMin){XNodeNeg=CoordNodeNeg[3*ii]; YNodeNeg=CoordNodeNeg[3*ii+1]; ZNodeNeg=ZcoordMax;};
347 for(int jj=0; jj<3; jj++){
348
349 XNodeNeg=round(XNodeNeg*roundfact)/roundfact;
350 YNodeNeg=round(YNodeNeg*roundfact)/roundfact;
351 ZNodeNeg=round(ZNodeNeg*roundfact)/roundfact;
352
353
354 XNodePos=CoordNodePos[3*jj]; YNodePos=CoordNodePos[3*jj+1]; ZNodePos=CoordNodePos[3*jj+2];
355 XNodePos=round(XNodePos*roundfact)/roundfact;
356 YNodePos=round(YNodePos*roundfact)/roundfact;
357 ZNodePos=round(ZNodePos*roundfact)/roundfact;
358
359 if(XNodeNeg==XNodePos && YNodeNeg==YNodePos && ZNodeNeg==ZNodePos){
360 PrismNodes[3+ii]=TriNodesPos[jj];
361 break;
362 }
363 }
364 }
365
366 double CoordNodesPrisms[18];
368
369
370
371
372
373 EntityHandle PeriodicPrism;
375 preriodic_prisms.insert(PeriodicPrism);
376
377 }
378
379
381
382 }
383
384 EntityHandle out_meshset;
386
387 ierr = m_field.get_entities_by_ref_level(bit_levels.back(),
BitRefLevel().set(),out_meshset); CHKERRQ(
ierr);
388 Range LatestRefinedTets;
389 rval = moab.get_entities_by_type(out_meshset, MBTET,LatestRefinedTets,
true);
CHKERRQ_MOAB(
rval);
390 Range LatestRefinedPrisms;
391 rval = moab.get_entities_by_type(out_meshset, MBPRISM,LatestRefinedPrisms,
true);
CHKERRQ_MOAB(
rval);
392
393 Range prims_on_problem_bit_level;
394 ierr = m_field.get_entities_by_type_and_ref_level(
395 problem_bit_level,
BitRefLevel().set(),MBPRISM,prims_on_problem_bit_level
397
398 EntityHandle meshset_prims_on_problem_bit_level;
400 rval = moab.add_entities(meshset_prims_on_problem_bit_level,prims_on_problem_bit_level);
CHKERRQ_MOAB(
rval);
401 ierr = m_field.seed_ref_level_MESHSET(meshset_prims_on_problem_bit_level,
BitRefLevel().set()); CHKERRQ(
ierr);
402
403
404 int field_rank=3;
409
410
414 if(it->getName().compare(0,12,"AllBoundSurf") == 0 || it->getMeshsetId() == 103) {
415 Range tris;
418 }
419 }
420 }
421
424
426
427 EntityHandle no_field_vertex;
428 {
429 const double coords[] = {0,0,0};
431 Range range_no_field_vertex;
432 range_no_field_vertex.insert(no_field_vertex);
434 EntityHandle meshset;
437 }
438 Range surface_negative;
439 ierr = m_field.get_cubit_msId_entities_by_dimension(101,
SIDESET,2,surface_negative,
true); CHKERRQ(
ierr);
440 ierr = PetscPrintf(PETSC_COMM_WORLD,
"number of SideSet 101 = %d\n",surface_negative.size()); CHKERRQ(
ierr);
442 }
443
446
448
450 EntityHandle no_field_vertex;
451 {
452 const double coords[] = {0,0,0};
454 Range range_no_field_vertex;
455 range_no_field_vertex.insert(no_field_vertex);
457 EntityHandle meshset;
464 }
465 }
466
467
468
473
478
479 PetscBool fo_boundary = PETSC_FALSE;
481 if(fo_boundary) {
483 if(it->getName().compare(0,12,"AllBoundSurf") == 0 || it->getMeshsetId() == 103) {
484 Range tris;
486 Range tris_edges;
487 rval = moab.get_adjacencies(tris,1,
false,tris_edges,moab::Interface::UNION);
CHKERRQ_MOAB(
rval);
490 }
491 }
492 }
493
498 }
499
504 }
505
506
508
510 ierr = m_field.
loop_dofs(
"MESH_NODE_POSITIONS",ent_method_material); CHKERRQ(
ierr);
511
512 boost::shared_ptr<Hooke<adouble> > hooke_adouble(
new Hooke<adouble>());
513 boost::shared_ptr<Hooke<double> > hooke_double(
new Hooke<double>());
514
517 if(it->getName() != "MAT_ELASTIC_1") continue;
519 ierr = it->getAttributeDataStructure(mydata); CHKERRQ(
ierr);
520 int id = it->getMeshsetId();
521 EntityHandle meshset = it->getMeshset();
523 meshset,MBTET,iso_elastic.setOfBlocks[id].tEts,true
525 iso_elastic.setOfBlocks[id].iD = id;
526 iso_elastic.setOfBlocks[id].E = mydata.
data.Young;
527 iso_elastic.setOfBlocks[id].PoissonRatio = mydata.
data.Poisson;
528 iso_elastic.setOfBlocks[id].materialDoublePtr = hooke_double;
529 iso_elastic.setOfBlocks[id].materialAdoublePtr = hooke_adouble;
530 ierr = m_field.seed_finite_elements(iso_elastic.setOfBlocks[
id].tEts); CHKERRQ(
ierr);
531 }
532 ierr = iso_elastic.addElement(
"ELASTIC",
"DISPLACEMENT"); CHKERRQ(
ierr);
533 ierr = iso_elastic.setOperators(
"DISPLACEMENT",
"MESH_NODE_POSITIONS",
false,
true); CHKERRQ(
ierr);
536 }
537
539 trans_elastic.commonData.spatialPositions = "DISPLACEMENT";
540 trans_elastic.commonData.meshPositions = "MESH_NODE_POSITIONS";
541 std::map<int,boost::shared_ptr<SmallStrainTranverslyIsotropicADouble> > tranversly_isotropic_adouble_ptr_map;
542 std::map<int,boost::shared_ptr<SmallStrainTranverslyIsotropicDouble> > tranversly_isotropic_double_ptr_map;
543 bool trans_iso_blocks = false;
545
546 string name = it->getName();
547 if (name.compare(0,20,"MAT_ELASTIC_TRANSISO") == 0) {
548 trans_iso_blocks = true;
549 int id = it->getMeshsetId();
551 ierr = it->getAttributeDataStructure(mydata); CHKERRQ(
ierr);
552 tranversly_isotropic_adouble_ptr_map[id] = boost::make_shared<SmallStrainTranverslyIsotropicADouble>();
553 tranversly_isotropic_double_ptr_map[id] = boost::make_shared<SmallStrainTranverslyIsotropicDouble>();
554
555 tranversly_isotropic_adouble_ptr_map.at(
id)->E_p = mydata.
data.Youngp;
556 tranversly_isotropic_double_ptr_map.at(
id)->E_p = mydata.
data.Youngp;
557 tranversly_isotropic_adouble_ptr_map.at(
id)->E_z = mydata.
data.Youngz;
558 tranversly_isotropic_double_ptr_map.at(
id)->E_z = mydata.
data.Youngz;
559 tranversly_isotropic_adouble_ptr_map.at(
id)->nu_p = mydata.
data.Poissonp;
560 tranversly_isotropic_double_ptr_map.at(
id)->nu_p = mydata.
data.Poissonp;
561 tranversly_isotropic_adouble_ptr_map.at(
id)->nu_pz = mydata.
data.Poissonpz;
562 tranversly_isotropic_double_ptr_map.at(
id)->nu_pz = mydata.
data.Poissonpz;
563 double shear_zp;
564 if(mydata.
data.Shearzp!=0) {
565 shear_zp = mydata.
data.Shearzp;
566 } else {
567 shear_zp = mydata.
data.Youngz/(2*(1+mydata.
data.Poissonpz));
568 }
569 tranversly_isotropic_adouble_ptr_map.at(it->getMeshsetId())->G_zp = shear_zp;
570 tranversly_isotropic_double_ptr_map.at(it->getMeshsetId())->G_zp = shear_zp;
571
572 EntityHandle meshset = it->getMeshset();
574 meshset,MBTET,trans_elastic.setOfBlocks[id].tEts,true
576
577 trans_elastic.setOfBlocks[id].iD = id;
578
579 trans_elastic.setOfBlocks[id].E = 0;
580 trans_elastic.setOfBlocks[id].PoissonRatio = 0;
581 trans_elastic.setOfBlocks[id].materialDoublePtr = tranversly_isotropic_double_ptr_map.at(id);
582 trans_elastic.setOfBlocks[id].materialAdoublePtr = tranversly_isotropic_adouble_ptr_map.at(id);
583 ierr = m_field.seed_finite_elements(trans_elastic.setOfBlocks[
id].tEts); CHKERRQ(
ierr);
584 }
585 }
586 if(trans_iso_blocks) {
595 }
596 for(
597 map<int,NonlinearElasticElement::BlockData>::iterator sit = trans_elastic.setOfBlocks.begin();
598 sit!=trans_elastic.setOfBlocks.end();sit++
599 ) {
601 }
602 }
603 if(trans_iso_blocks) {
604
605 trans_elastic.feRhs.getOpPtrVector().push_back(
607 );
608 trans_elastic.feRhs.getOpPtrVector().push_back(
610 );
612 trans_elastic.feRhs.getOpPtrVector().push_back(
614 );
615 }
616 map<int,NonlinearElasticElement::BlockData>::iterator sit = trans_elastic.setOfBlocks.begin();
617 for(;sit!=trans_elastic.setOfBlocks.end();sit++) {
618 trans_elastic.feRhs.getOpPtrVector().push_back(
620 "DISPLACEMENT",sit->second,trans_elastic.commonData,2,false,false,true
621 )
622 );
623 trans_elastic.feRhs.getOpPtrVector().push_back(
625 "DISPLACEMENT",sit->second,trans_elastic.commonData
626 )
627 );
628 }
629
630
631 trans_elastic.feLhs.getOpPtrVector().push_back(
633 );
634 trans_elastic.feLhs.getOpPtrVector().push_back(
636 );
638 trans_elastic.feLhs.getOpPtrVector().push_back(
640 );
641 }
642 sit = trans_elastic.setOfBlocks.begin();
643 for(;sit!=trans_elastic.setOfBlocks.end();sit++) {
644 trans_elastic.feLhs.getOpPtrVector().push_back(
646 "DISPLACEMENT",sit->second,trans_elastic.commonData,2,true,false,true
647 )
648 );
649 trans_elastic.feLhs.getOpPtrVector().push_back(
651 "DISPLACEMENT","DISPLACEMENT",sit->second,trans_elastic.commonData
652 )
653 );
654 }
655 }
656
659 lagrangian_element_disp.addLagrangiangElement(
660 "LAGRANGE_ELEM","DISPLACEMENT","LAGRANGE_MUL_DISP","MESH_NODE_POSITIONS"
661 );
662 }
663
668 lagrangian_element_trac.addLagrangiangElement(
669 "LAGRANGE_ELEM","DISPLACEMENT","LAGRANGE_MUL_TRAC","MESH_NODE_POSITIONS"
670 );
671 lagrangian_element_trac.addLagrangiangElement(
672 "LAGRANGE_ELEM_TRANS","DISPLACEMENT","LAGRANGE_MUL_RIGID_TRANS","MESH_NODE_POSITIONS"
673 );
674 lagrangian_element_trac.addLagrangiangElement(
675 "LAGRANGE_ELEM_ROT","DISPLACEMENT","LAGRANGE_MUL_RIGID_ROT","MESH_NODE_POSITIONS"
676 );
677 }
678
681 lagrangian_element_periodic.addLagrangiangElement(
682 "LAGRANGE_ELEM","DISPLACEMENT","LAGRANGE_MUL_PERIODIC","MESH_NODE_POSITIONS",preriodic_prisms
683 );
684 lagrangian_element_trac.addLagrangiangElement(
685 "LAGRANGE_ELEM_TRANS","DISPLACEMENT","LAGRANGE_MUL_RIGID_TRANS","MESH_NODE_POSITIONS"
686 );
687 }
688
689 struct MinMaxNodes {
690 enum MINAMX { C0,MAXLAST };
691 EntityHandle entMinMax[MAXLAST];
692 ublas::vector<int> rowIndices;
694 MinMaxNodes() {
695 rowIndices.resize(3*MAXLAST);
696 for(int ii = 0;ii<6;ii++) {
697 rHs[ii].resize(3*MAXLAST);
698 }
699 }
700 };
701 MinMaxNodes minMaxNodes;
702
709 EntityHandle condensed_traction_element_meshset;
710 rval = moab.create_meshset(MESHSET_TRACK_OWNER,condensed_traction_element_meshset); CHKERRQ(
ierr);
711 Range nodes;
713 if(it->getName().compare(0,12,"AllBoundSurf") == 0 || it->getMeshsetId() == 103) {
714 Range tris;
717 Range tris_nodes;
719 nodes.merge(tris_nodes);
720 }
721 }
722
723 {
725 x.resize(nodes.size(),false);
726 y.resize(nodes.size(),false);
727 z.resize(nodes.size(),false);
729 int bc_nb = 0;
730 for(int sx = -1; sx<=+1; sx+=2) {
731 for(int sy = -1; sy<=+1; sy+=2) {
732 for(int sz = -1; sz<=+1; sz+=2) {
733 if(bc_nb == MinMaxNodes::MAXLAST) break;
735 dist_up_right.resize(x.size(),false);
736 dist_up_right.clear();
737 for(unsigned int nn = 0;nn<x.size();nn++) {
738 if(
739 ((sx*x[nn])>0)&&
740 ((sy*y[nn])>0)&&
741 ((sz*z[nn])>0)
742 ) {
743 dist_up_right[nn] = sx*x[nn]+sy*y[nn]+sz*z[nn];
744 }
745 }
746 VectorDouble::iterator dist_it;
747 dist_it = max_element(dist_up_right.begin(),dist_up_right.end());
748 minMaxNodes.entMinMax[bc_nb++] = nodes[std::distance(dist_up_right.begin(),dist_it)];
749 }
750 }
751 }
752 }
753
754 }
755
756
758
760
761
763
765 if(trans_iso_blocks) {
767 "ELASTIC_MECHANICS","TRAN_ISOTROPIC_ELASTIC"
769 }
772 }
777 }
781 }
784 }
785
786
788
789 ierr = m_field.build_problems(); CHKERRQ(
ierr);
790
791
792
793
794
795 ierr = m_field.partition_problem(
"ELASTIC_MECHANICS"); CHKERRQ(
ierr);
796 ierr = m_field.partition_finite_elements(
799
800 ierr = m_field.partition_ghost_dofs(
"ELASTIC_MECHANICS"); CHKERRQ(
ierr);
801
802
804 vector<Vec> F(6);
805 ierr = m_field.VecCreateGhost(
"ELASTIC_MECHANICS",
ROW,&F[0]); CHKERRQ(
ierr);
806 for(int ii = 1;ii<6;ii++) {
807 ierr = VecDuplicate(F[0],&F[ii]); CHKERRQ(
ierr);
808 }
809 ierr = m_field.VecCreateGhost(
"ELASTIC_MECHANICS",
COL,&
D); CHKERRQ(
ierr);
810
811 Mat Aij;
812 ierr = m_field.MatCreateMPIAIJWithArrays(
"ELASTIC_MECHANICS",&Aij); CHKERRQ(
ierr);
813 ierr = MatSetOption(Aij,MAT_STRUCTURALLY_SYMMETRIC,PETSC_TRUE); CHKERRQ(
ierr);
814 ierr = MatSetOption(Aij,MAT_NEW_NONZERO_LOCATIONS,PETSC_TRUE); CHKERRQ(
ierr);
815 ierr = MatSetOption(Aij,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE); CHKERRQ(
ierr);
816 ierr = MatSetOption(Aij,MAT_USE_INODES,PETSC_TRUE); CHKERRQ(
ierr);
817 ierr = MatSetOption(Aij,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE); CHKERRQ(
ierr);
818 ierr = MatSetOption(Aij,MAT_KEEP_NONZERO_PATTERN,PETSC_FALSE); CHKERRQ(
ierr);
819
820
821
822
823
824
825
827 ierr = VecGhostUpdateBegin(
D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
828 ierr = VecGhostUpdateEnd(
D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
829 ierr = m_field.set_global_ghost_vector(
830 "ELASTIC_MECHANICS",
ROW,
D,INSERT_VALUES,SCATTER_REVERSE
832 for(int ii = 0;ii<6;ii++) {
833 ierr = VecZeroEntries(F[ii]); CHKERRQ(
ierr);
834 ierr = VecGhostUpdateBegin(F[ii],INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
835 ierr = VecGhostUpdateEnd(F[ii],INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
836 }
837 ierr = MatZeroEntries(Aij); CHKERRQ(
ierr);
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1053 int volume_vec_ghost[] = { 0 };
1054 ierr = VecCreateGhost(
1055 PETSC_COMM_WORLD,(!m_field.
get_comm_rank())?1:0,1,1,volume_vec_ghost,&volume_vec
1057 ierr = VecZeroEntries(volume_vec); CHKERRQ(
ierr);
1058
1059 iso_elastic.getLoopFeLhs().getOpPtrVector().push_back(
new VolumeCalculation(
"DISPLACEMENT",volume_vec));
1060 trans_elastic.getLoopFeLhs().getOpPtrVector().push_back(
new VolumeCalculation(
"DISPLACEMENT",volume_vec));
1061
1062
1063 iso_elastic.getLoopFeLhs().snes_x =
D;
1064 iso_elastic.getLoopFeLhs().snes_B = Aij;
1065 trans_elastic.getLoopFeLhs().snes_x =
D;
1066 trans_elastic.getLoopFeLhs().snes_B = Aij;
1070
1072 lagrangian_element_disp.setRVEBCsOperators("DISPLACEMENT","LAGRANGE_MUL_DISP","MESH_NODE_POSITIONS",Aij,F);
1075 }
1077 lagrangian_element_trac.setRVEBCsOperators("DISPLACEMENT","LAGRANGE_MUL_TRAC","MESH_NODE_POSITIONS",Aij,F);
1080 lagrangian_element_rigid_body_trans.setRVEBCsRigidBodyTranOperators(
1081 "DISPLACEMENT","LAGRANGE_MUL_RIGID_TRANS",Aij,lagrangian_element_trac.setOfRVEBC
1082 );
1084 "ELASTIC_MECHANICS","LAGRANGE_ELEM_TRANS",lagrangian_element_rigid_body_trans.getLoopFeRVEBCLhs()
1086 lagrangian_element_rigid_body_rot.setRVEBCsRigidBodyRotOperators(
1087 "DISPLACEMENT","LAGRANGE_MUL_RIGID_ROT",Aij,lagrangian_element_trac.setOfRVEBC
1088 );
1090 "ELASTIC_MECHANICS","LAGRANGE_ELEM_ROT",lagrangian_element_rigid_body_rot.getLoopFeRVEBCLhs()
1092 }
1094 lagrangian_element_periodic.setRVEBCsOperators("DISPLACEMENT","LAGRANGE_MUL_PERIODIC","MESH_NODE_POSITIONS",Aij,F);
1096 "ELASTIC_MECHANICS","LAGRANGE_ELEM",lagrangian_element_periodic.getLoopFeRVEBCLhs()
1099 "ELASTIC_MECHANICS","LAGRANGE_ELEM",lagrangian_element_periodic.getLoopFeRVEBCRhs()
1101 lagrangian_element_rigid_body_trans.setRVEBCsRigidBodyTranOperators(
1102 "DISPLACEMENT","LAGRANGE_MUL_RIGID_TRANS",Aij,lagrangian_element_periodic.setOfRVEBC
1103 );
1105 "ELASTIC_MECHANICS","LAGRANGE_ELEM_TRANS",lagrangian_element_rigid_body_trans.getLoopFeRVEBCLhs()
1107 }
1108
1109 for(int ii = 0;ii<6;ii++) {
1110 ierr = VecGhostUpdateBegin(F[ii],ADD_VALUES,SCATTER_REVERSE); CHKERRQ(
ierr);
1111 ierr = VecGhostUpdateEnd(F[ii],ADD_VALUES,SCATTER_REVERSE); CHKERRQ(
ierr);
1112 ierr = VecAssemblyBegin(F[ii]); CHKERRQ(
ierr);
1113 ierr = VecAssemblyEnd(F[ii]); CHKERRQ(
ierr);
1114 }
1115 ierr = MatAssemblyBegin(Aij,MAT_FINAL_ASSEMBLY); CHKERRQ(
ierr);
1116 ierr = MatAssemblyEnd(Aij,MAT_FINAL_ASSEMBLY); CHKERRQ(
ierr);
1117
1118 {
1119
1120
1121
1122 }
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171 ierr = VecAssemblyBegin(volume_vec); CHKERRQ(
ierr);
1172 ierr = VecAssemblyEnd(volume_vec); CHKERRQ(
ierr);
1173 double rve_volume;
1174 ierr = VecSum(volume_vec,&rve_volume); CHKERRQ(
ierr);
1175 ierr = PetscPrintf(PETSC_COMM_WORLD,
"RVE Volume %3.2g\n",rve_volume); CHKERRQ(
ierr);
1176
1177 ierr = VecDestroy(&volume_vec);
1178
1179
1180 KSP solver;
1181 ierr = KSPCreate(PETSC_COMM_WORLD,&solver); CHKERRQ(
ierr);
1182 ierr = KSPSetOperators(solver,Aij,Aij); CHKERRQ(
ierr);
1183 ierr = KSPSetFromOptions(solver); CHKERRQ(
ierr);
1184 ierr = KSPSetUp(solver); CHKERRQ(
ierr);
1185
1187 Dmat.resize(6,6);
1188 Dmat.clear();
1189
1190
1192 int stress_homo_ghost[] = { 0,1,2,3,4,5,6 };
1195 PetscBool stress_by_boundary_integral = PETSC_FALSE;
1196 ierr = VecCreateGhost(
1197 PETSC_COMM_WORLD,(!m_field.
get_comm_rank())?6:0,6,6,stress_homo_ghost,&stress_homo
1199 switch(choise_value) {
1201 lagrangian_element_disp.setRVEBCsHomoStressOperators(
1202 "DISPLACEMENT","LAGRANGE_MUL_DISP","MESH_NODE_POSITIONS",stress_homo
1203 );
1204 break;
1206 lagrangian_element_trac.setRVEBCsHomoStressOperators(
1207 "DISPLACEMENT","LAGRANGE_MUL_TRAC","MESH_NODE_POSITIONS",stress_homo
1208 );
1209 break;
1211 lagrangian_element_periodic.setRVEBCsHomoStressOperators(
1212 "DISPLACEMENT","LAGRANGE_MUL_PERIODIC","MESH_NODE_POSITIONS",stress_homo
1213 );
1214 break;
1216
1217
1218
1219
1220
1221
1222
1223
1224
1225
1226
1227
1228
1229
1230
1231
1232
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243
1244
1245
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259
1260
1261
1262
1263
1264
1265
1266
1267
1268
1269
1270
1271
1272
1273
1274
1275
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325
1326
1327
1328
1329
1330
1331
1332
1333
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
1347 break;
1350 }
1351
1353
1354 bool doPreProcess;
1355 bool doPostProcess;
1356
1359 doPreProcess(true),
1360 doPostProcess(true)
1361 {}
1362
1363 void setDoPreProcess() { doPreProcess = true; }
1364 void unSetDoPreProcess() { doPreProcess = false; }
1365 void setDoPostProcess() { doPostProcess = true; }
1366 void unSetDoPostProcess() { doPostProcess = false; }
1367
1368
1369
1370 PetscErrorCode preProcess() {
1371 PetscFunctionBegin;
1372 if(doPreProcess) {
1374 }
1375 PetscFunctionReturn(0);
1376 }
1377 PetscErrorCode postProcess() {
1378 PetscFunctionBegin;
1379 if(doPostProcess) {
1381 }
1382 PetscFunctionReturn(0);
1383 }
1384 };
1385
1387
1388 ierr = post_proc.generateReferenceElementMesh(); CHKERRQ(
ierr);
1389 ierr = post_proc.addFieldValuesPostProc(
"DISPLACEMENT"); CHKERRQ(
ierr);
1390 ierr = post_proc.addFieldValuesGradientPostProc(
"DISPLACEMENT"); CHKERRQ(
ierr);
1391 ierr = post_proc.addFieldValuesPostProc(
"MESH_NODE_POSITIONS"); CHKERRQ(
ierr);
1392 ierr = post_proc.addFieldValuesGradientPostProc(
"MESH_NODE_POSITIONS"); CHKERRQ(
ierr);
1393 if(trans_iso_blocks) {
1394 ierr = post_proc.addFieldValuesGradientPostProc(
"POTENTIAL_FIELD"); CHKERRQ(
ierr);
1395 }
1396 for(
1397 map<int,NonlinearElasticElement::BlockData>::iterator sit = iso_elastic.setOfBlocks.begin();
1398 sit != iso_elastic.setOfBlocks.end(); sit++
1399 ) {
1400 post_proc.getOpPtrVector().push_back(
1402 post_proc.postProcMesh,
1403 post_proc.mapGaussPts,
1404 "DISPLACEMENT",
1405 sit->second,
1406 post_proc.commonData,
1407 false
1408 )
1409 );
1410 }
1411 for(
1412 map<int,NonlinearElasticElement::BlockData>::iterator sit = trans_elastic.setOfBlocks.begin();
1413 sit != trans_elastic.setOfBlocks.end(); sit++
1414 ) {
1415 post_proc.getOpPtrVector().push_back(
1417 post_proc.postProcMesh,
1418 post_proc.mapGaussPts,
1419 "DISPLACEMENT",
1420 sit->second,
1421 post_proc.commonData
1422 )
1423 );
1424 }
1425
1426 PetscScalar *avec;
1427 ierr = VecGetArray(stress_homo,&avec); CHKERRQ(
ierr);
1428 for(int ii = 0;ii<6;ii++) {
1429 ierr = VecZeroEntries(
D); CHKERRQ(
ierr);
1430 ierr = KSPSolve(solver,F[ii],
D); CHKERRQ(
ierr);
1431 ierr = VecGhostUpdateBegin(
D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
1432 ierr = VecGhostUpdateEnd(
D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
1433 ierr = m_field.set_global_ghost_vector(
1434 "ELASTIC_MECHANICS",
ROW,
D,INSERT_VALUES,SCATTER_REVERSE
1436 post_proc.setDoPreProcess();
1437 post_proc.unSetDoPostProcess();
1439 "ELASTIC_MECHANICS","ELASTIC",post_proc
1441 post_proc.unSetDoPreProcess();
1442 post_proc.setDoPostProcess();
1443
1446 "ELASTIC_MECHANICS","TRAN_ISOTROPIC_ELASTIC",post_proc
1448 {
1449 ostringstream sss;
1450 sss <<
"mode_" <<
homo_bc_names[choise_value] <<
"_" << ii <<
".h5m";
1451 rval = post_proc.postProcMesh.write_file(
1452 sss.str().c_str(),"MOAB","PARALLEL=WRITE_PART"
1454 }
1455 ierr = VecZeroEntries(stress_homo); CHKERRQ(
ierr);
1458 "ELASTIC_MECHANICS","LAGRANGE_ELEM",lagrangian_element_disp.getLoopFeRVEBCStress()
1460 }
1463 "ELASTIC_MECHANICS","LAGRANGE_ELEM",lagrangian_element_trac.getLoopFeRVEBCStress()
1465 }
1468 "ELASTIC_MECHANICS","LAGRANGE_ELEM",lagrangian_element_periodic.getLoopFeRVEBCStress()
1470 }
1471
1472
1473
1474
1475
1476
1477
1478
1479
1480
1481
1482
1483
1484
1485
1486
1487
1488
1489
1490
1491
1492
1493
1495 PETSC_NULL,"-my_rve_volume",&rve_volume,PETSC_NULL
1497 ierr = VecAssemblyBegin(stress_homo); CHKERRQ(
ierr);
1498 ierr = VecAssemblyEnd(stress_homo); CHKERRQ(
ierr);
1499 ierr = VecScale(stress_homo,1.0/rve_volume); CHKERRQ(
ierr);
1500 ierr = VecGhostUpdateBegin(stress_homo,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
1501 ierr = VecGhostUpdateEnd(stress_homo,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(
ierr);
1502 for(int jj=0; jj<6; jj++) {
1503 Dmat(jj,ii) = avec[jj];
1504 }
1505 }
1506 ierr = VecRestoreArray(stress_homo,&avec); CHKERRQ(
ierr);
1507
1508 PetscPrintf(
1509 PETSC_COMM_WORLD,"\nHomogenised Stiffens Matrix = \n\n"
1510 );
1511
1512 for(int ii=0; ii<6; ii++) {
1513 PetscPrintf(
1514 PETSC_COMM_WORLD,
1515 "stress %d\t\t%8.5e\t\t%8.5e\t\t%8.5e\t\t%8.5e\t\t%8.5e\t\t%8.5e\n",
1516 ii,Dmat(ii,0),Dmat(ii,1),Dmat(ii,2),Dmat(ii,3),Dmat(ii,4),Dmat(ii,5)
1517 );
1518 }
1519
1520
1521
1522
1523 char output_file_Dmat[255];
1525 if(flg) {
1526
1527
1528 if(pcomm->rank()==0){
1529 int fd;
1530 PetscViewer view_out;
1531 PetscViewerBinaryOpen(PETSC_COMM_SELF,output_file_Dmat,FILE_MODE_WRITE,&view_out);
1532 PetscViewerBinaryGetDescriptor(view_out,&fd);
1533 PetscBinaryWrite(fd,&Dmat(0,0),36,PETSC_DOUBLE,PETSC_FALSE);
1534 PetscViewerDestroy(&view_out);
1535 }
1536
1537
1538
1539
1540
1541
1542
1543
1544
1545
1546
1547
1548
1549 }
1550
1551
1552 for(int ii = 0;ii<6;ii++) {
1553 ierr = VecDestroy(&F[ii]); CHKERRQ(
ierr);
1554 }
1556 ierr = MatDestroy(&Aij); CHKERRQ(
ierr);
1557 ierr = KSPDestroy(&solver); CHKERRQ(
ierr);
1558 ierr = VecDestroy(&stress_homo); CHKERRQ(
ierr);
1559
1560 ierr = PetscTime(&v2);CHKERRQ(
ierr);
1561 ierr = PetscGetCPUTime(&t2);CHKERRQ(
ierr);
1562
1563 PetscSynchronizedPrintf(PETSC_COMM_WORLD,"Total Rank %d Time = %f CPU Time = %f\n",pcomm->rank(),v2-v1,t2-t1);
1564
1565 }
1567
1568
1570}
DEPRECATED typedef PostProcStress PostPorcStress
#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.
@ NOFIELD
scalar or vector of scalars describe (no true field)
#define MYPCOMM_INDEX
default communicator number PCOMM
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
#define CHKERRQ_MOAB(a)
check error code of MoAB function
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 bool check_field(const std::string &name) const =0
check if field is in database
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_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
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
multi_index_container< Face_CenPos_Handle, indexed_by< ordered_non_unique< tag< xcoord_tag >, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::xcoord > >, ordered_non_unique< tag< ycoord_tag >, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::ycoord > >, ordered_non_unique< tag< zcoord_tag >, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::zcoord > >, ordered_unique< tag< Tri_Hand_tag >, member< Face_CenPos_Handle, const EntityHandle,&Face_CenPos_Handle::Tri_Hand > >, ordered_unique< tag< Composite_xyzcoord >, composite_key< Face_CenPos_Handle, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::xcoord >, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::ycoord >, member< Face_CenPos_Handle, double,&Face_CenPos_Handle::zcoord > > > > > Face_CenPos_Handle_multiIndex
const FTensor::Tensor2< T, Dim, Dim > Vec
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
UBlasMatrix< double > MatrixDouble
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
UBlasVector< double > VectorDouble
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 PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
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
const double D
diffusivity
const char * homo_bc_names[]
virtual int get_comm_size() const =0
virtual moab::Interface & get_moab()=0
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
virtual bool check_finite_element(const std::string &name) const =0
Check if finite element is in database.
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.
virtual int get_comm_rank() const =0
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.
DEPRECATED MoFEMErrorCode seed_ref_level(const Range &ents, const BitRefLevel &bit, const bool only_tets=true, int verb=-1)
seed entities in the range and their adjacencies in a particular BitRefLevel
DEPRECATED MoFEMErrorCode seed_ref_level_3D(const EntityHandle meshset, const BitRefLevel &bit, int verb=-1)
seed 2D entities in the meshset and their adjacencies (only TETs adjacencies) in a particular BitRefL...
Transverse Isotropic material data structure.
Elastic material data structure.
Projection of edge entities with one mid-node on hierarchical basis.
definition of volume element
Operator performs automatic differentiation.
structure grouping operators and data used for calculation of nonlinear elastic element
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode postProcess()
function is run at the end of loop