It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.
Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }
317 {
319
320 auto make_vtk = [&]() {
326 boost::lexical_cast<std::string>(
ts_step) +
327 ".h5m");
329 };
330
331 auto save_skeleton = [&]() {
333
336 "out_skeleton_" + boost::lexical_cast<std::string>(
ts_step) +
337 ".h5m");
339 };
340
341 auto save_skin = [&]() {
345
348 "out_skin_" + boost::lexical_cast<std::string>(
ts_step) +
".h5m");
350 };
351
352 double max_disp, min_disp;
353
354 auto print_max_min = [&](auto &tuple, const std::string msg) {
356 CHKERR VecScatterBegin(std::get<1>(tuple),
ts_u, std::get<0>(tuple),
357 INSERT_VALUES, SCATTER_FORWARD);
358 CHKERR VecScatterEnd(std::get<1>(tuple),
ts_u, std::get<0>(tuple),
359 INSERT_VALUES, SCATTER_FORWARD);
360 CHKERR VecMax(std::get<0>(tuple), PETSC_NULL, &max_disp);
361 CHKERR VecMin(std::get<0>(tuple), PETSC_NULL, &min_disp);
362 PetscPrintf(PETSC_COMM_WORLD, "%s time %3.4e min %3.4e max %3.4e\n",
363 msg.c_str(),
ts_t, min_disp, max_disp);
365 };
366
367 double gauss_area;
368
372
373 int roller_nb = (*cache).rollerDataVec.size();
374 l_reactions.resize(roller_nb * 3);
375 std::fill(l_reactions.begin(), l_reactions.end(), 0);
376
378
379 if (true) {
380 std::vector<double> gauss_pts(2, 0);
382 CHKERR MPIU_Allreduce(lgauss_pts.data(), gauss_pts.data(), 2,
383 MPIU_REAL, MPIU_SUM,
384 PetscObjectComm((PetscObject)dm));
386 "\t TS %d Total contact area (ref): %g / %g",
ts_step,
387 gauss_pts[0], gauss_pts[1]);
388 gauss_area = gauss_pts[0];
389 lgauss_pts[0] = lgauss_pts[1] = 0;
390
391 vector<double> g_reactions(l_reactions.size() * 3, 0);
392 CHKERR MPIU_Allreduce(l_reactions.data(), g_reactions.data(),
393 roller_nb * 3, MPIU_REAL, MPIU_SUM,
394 PetscObjectComm((PetscObject)dm));
395 for (
int dd = 0;
dd != roller_nb; ++
dd)
397 "\t Body %d reactions Time %g Force X: %3.4e Y: "
398 "%3.4e Z: %3.4e",
399 (*cache).rollerDataVec[
dd].iD,
ts_t,
400 g_reactions[3 *
dd + 0], g_reactions[3 *
dd + 1],
401 g_reactions[3 *
dd + 2]);
402 }
403 }
404
406 std::ostringstream ostrm, ostrm2;
407 ostrm <<
"out_debug_" <<
ts_step <<
".vtk";
408 ostrm2 <<
"out_debug_" <<
ts_step <<
".h5m";
412 else
414 "PARALLEL=WRITE_PART");
415
417 }
418 }
419
422 CHKERR VecZeroEntries(vec);
423 CHKERR VecGhostUpdateBegin(vec, INSERT_VALUES, SCATTER_FORWARD);
424 CHKERR VecGhostUpdateEnd(vec, INSERT_VALUES, SCATTER_FORWARD);
429
430 ParallelComm *pcomm =
433
436
437
438 double def_value = 0.;
440 ents, &def_value);
441
442
443
444 CHKERR VecAssemblyBegin(vec);
445 CHKERR VecAssemblyEnd(vec);
446 CHKERR VecGhostUpdateBegin(vec, ADD_VALUES, SCATTER_REVERSE);
447 CHKERR VecGhostUpdateEnd(vec, ADD_VALUES, SCATTER_REVERSE);
448 CHKERR VecGhostUpdateBegin(vec, INSERT_VALUES, SCATTER_FORWARD);
449 CHKERR VecGhostUpdateEnd(vec, INSERT_VALUES, SCATTER_FORWARD);
450
451 const double *react_ptr;
452 CHKERR VecGetArrayRead(vec, &react_ptr);
453
454 PetscPrintf(PETSC_COMM_WORLD,
455 "Reactions time %3.4e X: %3.4e Y: %3.4e Z: %3.4e \n",
ts_t,
456 react_ptr[0], react_ptr[1], react_ptr[2]);
457
458 CHKERR VecRestoreArrayRead(vec, &react_ptr);
459 } else {
462 }
463
466
471 if (m_field.get_comm_rank() == 0) {
472 std::ostringstream ost;
473 for (auto &rol : (*cache).rollerDataVec) {
476 roller_coords += rol.originCoords;
477
480 }
481
482 ost <<
"out_roller_" <<
ts_step <<
".vtk";
485 }
486 }
487
489
490 CHKERR VecGhostUpdateBegin(
ts_u, INSERT_VALUES, SCATTER_FORWARD);
491 CHKERR VecGhostUpdateEnd(
ts_u, INSERT_VALUES, SCATTER_FORWARD);
492
493 string rest_name =
494 "restart_" + to_string(
ts_step) +
"_" + to_string(
ts_t) +
".dat";
495 PetscViewer viewer;
496 CHKERR PetscViewerBinaryOpen(PETSC_COMM_WORLD, rest_name.c_str(),
497 FILE_MODE_WRITE, &viewer);
498 CHKERR PetscViewerBinarySetSkipInfo(viewer, PETSC_TRUE);
499 VecView(
ts_u, viewer);
500 PetscViewerDestroy(&viewer);
501 }
502
506
508 case 1: {
510 const double *react_ptr;
512 if (fabs(react_ptr[0] + 0.12519621572) > 1e-6)
514 "atom test diverged!");
516 }
517 break;
518 }
519 case 2: {
521 const double *react_ptr;
523
524 if ((0.03663003663 - react_ptr[2]) / 0.03663003663 > 1e-2)
526 "atom test diverged!");
528 }
529 break;
530 }
531 case 3: {
533 double min_disp;
535 if (fabs((min_disp + 2.70) / 2.70) > 5e-2)
537 "atom test diverged!");
538 }
539 break;
540 }
541 case 4: {
543 const double *react_ptr;
545
546
547
548 if (abs(0.056090728 - (react_ptr[2] * 4.0)) / 0.056090728 > 5.0e-2)
550 "atom test diverged!");
551 }
552 case 5: {
554 const double *react_ptr;
555
557 if (abs(19030.7 + react_ptr[2]) / 19030.7 > 1e-3)
559 "atom test diverged!");
560 }
561 break;
562 }
563 break;
564 }
565
566 default:
567 break;
568 }
569
571 }
#define MOFEM_LOG_C(channel, severity, format,...)
void simple(double P1[], double P2[], double P3[], double c[], const int N)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MYPCOMM_INDEX
default communicator number PCOMM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_ATOM_TEST_INVALID
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
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)
virtual int get_comm_size() const =0
virtual moab::Interface & get_moab()=0
Simple interface for fast problem set-up.
PetscInt ts_step
time step number
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.