v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Protected Attributes | List of all members
MoFEM::EssentialPreProcReaction< DisplacementCubitBcData > Struct Reference

Specialization for DisplacementCubitBcData. More...

#include <src/boundary_conditions/EssentialDisplacementCubitBcData.hpp>

Collaboration diagram for MoFEM::EssentialPreProcReaction< DisplacementCubitBcData >:
[legend]

Public Member Functions

 EssentialPreProcReaction (MoFEM::Interface &m_field, boost::shared_ptr< FEMethod > fe_ptr, SmartPetscObj< Vec > rhs=nullptr, LogManager::SeverityLevel sev=Sev::inform)
 
MoFEMErrorCode operator() ()
 

Protected Attributes

MoFEM::InterfacemField
 
boost::weak_ptr< FEMethodfePtr
 
SmartPetscObj< Vec > vRhs
 
LogManager::SeverityLevel sevLevel
 

Detailed Description

Specialization for DisplacementCubitBcData.

Template Parameters

Definition at line 151 of file EssentialDisplacementCubitBcData.hpp.

Constructor & Destructor Documentation

◆ EssentialPreProcReaction()

MoFEM::EssentialPreProcReaction< DisplacementCubitBcData >::EssentialPreProcReaction ( MoFEM::Interface m_field,
boost::shared_ptr< FEMethod fe_ptr,
SmartPetscObj< Vec >  rhs = nullptr,
LogManager::SeverityLevel  sev = Sev::inform 
)

Member Function Documentation

◆ operator()()

Definition at line 420 of file EssentialDisplacementCubitBcData.cpp.

420 {
421 MOFEM_LOG_CHANNEL("WORLD");
423
424 enum { X = 0, Y, Z, MX, MY, MZ, LAST };
425
426 if (auto fe_ptr = fePtr.lock()) {
427
428 SmartPetscObj<Vec> f = vRhs ? vRhs : SmartPetscObj<Vec>(fe_ptr->f, true);
429
430 if (fe_ptr->vecAssembleSwitch) {
431 if ((*fe_ptr->vecAssembleSwitch) && !vRhs) {
432 CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
433 CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
434 CHKERR VecAssemblyBegin(f);
435 CHKERR VecAssemblyEnd(f);
436 *fe_ptr->vecAssembleSwitch = false;
437 }
438 }
439
440 auto get_low_hi_uid_by_entities = [](auto bit_number, auto f, auto s) {
441 return std::make_pair(DofEntity::getLoFieldEntityUId(bit_number, f),
442 DofEntity::getHiFieldEntityUId(bit_number, s));
443 };
444
445 auto get_low_hi = [fe_ptr](auto lo_uid, auto hi_uid) {
446 auto it = fe_ptr->problemPtr->numeredRowDofsPtr->get<Unique_mi_tag>()
447 .lower_bound(lo_uid);
448 auto hi_it = fe_ptr->problemPtr->numeredRowDofsPtr->get<Unique_mi_tag>()
449 .upper_bound(hi_uid);
450 return std::make_pair(it, hi_it);
451 };
452
453 auto mpi_array_reduce = [this](auto &array) {
454 std::array<double, LAST> array_sum{0, 0, 0, 0, 0, 0};
455 MPI_Allreduce(&array[0], &array_sum[0], LAST, MPI_DOUBLE, MPI_SUM,
456 mField.get_comm());
457 return array_sum;
458 };
459
460 const double *a;
461 CHKERR VecGetArrayRead(f, &a);
462
463 auto bc_mng = mField.getInterface<BcManager>();
464 const auto problem_name = fe_ptr->problemPtr->getName();
465 const auto nb_local_dofs = fe_ptr->problemPtr->nbLocDofsRow;
466
467 std::array<double, LAST> total_reactions{0, 0, 0, 0, 0, 0};
468
469 for (auto bc : bc_mng->getBcMapByBlockName()) {
470 if (auto disp_bc = bc.second->dispBcPtr) {
471
472 auto &bc_id = bc.first;
473
474 auto regex_str = (boost::format("%s_(.*)") % problem_name).str();
475 if (std::regex_match(bc_id, std::regex(regex_str))) {
476
477 auto [field_name, block_name] =
478 BcManager::extractStringFromBlockId(bc_id, problem_name);
479
480 MOFEM_TAG_AND_LOG("WORLD", sevLevel, "Essential")
481 << "EssentialPreProc<DisplacementCubitBcData>: " << problem_name
482 << "_" << field_name << "_" << block_name;
483 auto bit_number = mField.get_field_bit_number(field_name);
484
485 FTensor::Tensor1<double, 3> t_off{0., 0., 0.};
486 if (auto ext_disp_bc =
487 dynamic_cast<DisplacementCubitBcDataWithRotation const *>(
488 disp_bc.get())) {
489 for (int a = 0; a != 3; ++a)
490 t_off(a) = ext_disp_bc->rotOffset[a];
491 }
492
493 auto verts = bc.second->bcEnts.subset_by_type(MBVERTEX);
494
498
499 auto get_coords_vec = [&]() {
500 VectorDouble coords_vec(3 * verts.size());
501 if (verts.size()) {
502 CHKERR mField.get_moab().get_coords(verts, &*coords_vec.begin());
503 }
504 return coords_vec;
505 };
506
507 auto coords_vec = get_coords_vec();
508 std::array<double, LAST> reactions{0, 0, 0, 0, 0, 0};
509
510 for (auto pit = verts.const_pair_begin();
511 pit != verts.const_pair_end(); ++pit) {
512 auto [lo_uid, hi_uid] =
513 get_low_hi_uid_by_entities(bit_number, pit->first, pit->second);
514 auto [lo, hi] = get_low_hi(lo_uid, hi_uid);
515
516 for (; lo != hi; ++lo) {
517 const auto loc_dof = (*lo)->getPetscLocalDofIdx();
518 if (loc_dof < nb_local_dofs) {
519
520 const auto coeff = (*lo)->getDofCoeffIdx();
521
522 if (
523
524 ((disp_bc->data.flag1 || disp_bc->data.flag4) &&
525 coeff == 0) ||
526 ((disp_bc->data.flag2 || disp_bc->data.flag5) &&
527 coeff == 1) ||
528 ((disp_bc->data.flag3 || disp_bc->data.flag6) &&
529 coeff == 2)) {
530
531 const auto ent = (*lo)->getEnt();
532 reactions[coeff] += a[loc_dof];
533
534 auto force = [&]() {
535 FTensor::Tensor1<double, 3> t_force{0., 0., 0.};
536 t_force(coeff) = a[loc_dof];
537 return t_force;
538 };
539
540 auto coord = [&]() {
541 const auto idx = verts.index(ent);
543 coords_vec[3 * idx + X], coords_vec[3 * idx + Y],
544 coords_vec[3 * idx + Z]};
545 t_coords(i) -= t_off(i);
546 return t_coords;
547 };
548
549 auto moment = [&](auto &&t_force, auto &&t_coords) {
551 t_moment(i) =
552 (FTensor::levi_civita<double>(i, j, k) * t_coords(k)) *
553 t_force(j);
554 return t_moment;
555 };
556
557 auto t_moment = moment(force(), coord());
558 reactions[MX] += t_moment(X);
559 reactions[MY] += t_moment(Y);
560 reactions[MZ] += t_moment(Z);
561 }
562 }
563 }
564 }
565
566 FTensor::Tensor1<double, 3> t_force{reactions[X], reactions[Y],
567 reactions[Z]};
568 FTensor::Tensor1<double, 3> t_moment{reactions[MX], reactions[MY],
569 reactions[MZ]};
571 &total_reactions[X], &total_reactions[Y], &total_reactions[Z]};
573 &total_reactions[MX], &total_reactions[MY], &total_reactions[MZ]};
574 t_total_force(i) += t_force(i);
575 t_total_moment(i) +=
576 t_moment(i) +
577 (FTensor::levi_civita<double>(i, j, k) * t_off(k)) * t_force(j);
578
579 auto mpi_reactions = mpi_array_reduce(reactions);
580 MOFEM_TAG_AND_LOG_C("WORLD", sevLevel, "Essential",
581 "Offset: %4.3e %4.3e %4.3e", t_off(X), t_off(Y),
582 t_off(Z));
583 MOFEM_TAG_AND_LOG_C("WORLD", sevLevel, "Essential",
584 "Force: %4.3e %4.3e %4.3e", mpi_reactions[X],
585 mpi_reactions[Y], mpi_reactions[Z]);
586 MOFEM_TAG_AND_LOG_C("WORLD", sevLevel, "Essential",
587 "Moment: %4.3e %4.3e %4.3e", mpi_reactions[MX],
588 mpi_reactions[MY], mpi_reactions[MZ]);
589
590
591 }
592 }
593 }
594
595 CHKERR VecRestoreArrayRead(f, &a);
596
597 auto mpi_total_reactions = mpi_array_reduce(total_reactions);
599 "WORLD", sevLevel, "Essential", "Total force: %4.3e %4.3e %4.3e",
600 mpi_total_reactions[X], mpi_total_reactions[Y], mpi_total_reactions[Z]);
601 MOFEM_TAG_AND_LOG_C("WORLD", sevLevel, "Essential",
602 "Total moment: %4.3e %4.3e %4.3e",
603 mpi_total_reactions[MX], mpi_total_reactions[MY],
604 mpi_total_reactions[MZ]);
605
606 } else {
607 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
608 "Can not lock shared pointer");
609 }
610
612}
#define MOFEM_TAG_AND_LOG_C(channel, severity, tag, format,...)
Tag and log in channel.
Definition: LogManager.hpp:370
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
constexpr double a
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#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 MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
constexpr auto field_name
static std::pair< std::string, std::string > extractStringFromBlockId(const std::string block_id, const std::string prb_name)
Extract block name and block name form block id.
Definition: BcManager.cpp:1218
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
static UId getLoFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
static UId getHiFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

Member Data Documentation

◆ fePtr

Definition at line 161 of file EssentialDisplacementCubitBcData.hpp.

◆ mField

Definition at line 160 of file EssentialDisplacementCubitBcData.hpp.

◆ sevLevel

Definition at line 163 of file EssentialDisplacementCubitBcData.hpp.

◆ vRhs

Definition at line 162 of file EssentialDisplacementCubitBcData.hpp.


The documentation for this struct was generated from the following files: