v0.13.2
Loading...
Searching...
No Matches
EssentialDisplacementCubitBcData.cpp
Go to the documentation of this file.
1/**
2 * @file EssentialDisplacementCubitBcData.cpp
3 * @brief Essential boundary conditions
4 * @version 13.1
5 * @date 2022-09-03
6 *
7 * @copyright Copyright (c) 2022
8 *
9 */
10
11#include <MoFEM.hpp>
12
13namespace MoFEM {
14
16 MoFEM::Interface &m_field, boost::shared_ptr<FEMethod> fe_ptr,
17 std::vector<boost::shared_ptr<ScalingMethod>> smv, bool get_coords)
18 : mField(m_field), fePtr(fe_ptr), vecOfTimeScalingMethods(smv),
19 getCoords(get_coords) {}
20
22 MOFEM_LOG_CHANNEL("WORLD");
24
25 if (auto fe_method_ptr = fePtr.lock()) {
26
27 auto bc_mng = mField.getInterface<BcManager>();
28 auto fb = mField.getInterface<FieldBlas>();
29 const auto problem_name = fe_method_ptr->problemPtr->getName();
30
31 for (auto bc : bc_mng->getBcMapByBlockName()) {
32 if (auto disp_bc = bc.second->dispBcPtr) {
33
34 auto &bc_id = bc.first;
35
36 auto regex_str = (boost::format("%s_(.*)") % problem_name).str();
37 if (std::regex_match(bc_id, std::regex(regex_str))) {
38
39 auto [field_name, block_name] =
40 BcManager::extractStringFromBlockId(bc_id, problem_name);
41
42 auto get_field_coeffs = [&](auto field_name) {
43 auto field_ptr = mField.get_field_structure(field_name);
44 return field_ptr->getNbOfCoeffs();
45 };
46 const auto nb_field_coeffs = get_field_coeffs(field_name);
47
48 MOFEM_LOG("WORLD", Sev::noisy)
49 << "Apply EssentialPreProc<DisplacementCubitBcData>: "
50 << problem_name << "_" << field_name << "_" << block_name;
51
52 FTensor::Tensor1<double, 3> t_angles(0., 0., 0.);
53 FTensor::Tensor1<double, 3> t_vals(0., 0., 0.);
54 FTensor::Tensor1<double, 3> t_off(0., 0., 0.);
55
56 if (auto ext_disp_bc =
57 dynamic_cast<DisplacementCubitBcDataWithRotation const *>(
58 disp_bc.get())) {
59 for (int a = 0; a != 3; ++a)
60 t_off(a) = ext_disp_bc->rotOffset[a];
61 }
62
63 auto scale_value = [&](const double &c) {
64 double val = c;
65 for (auto s : vecOfTimeScalingMethods) {
66 val *= s->getScale(fe_method_ptr->ts_t);
67 }
68 return val;
69 };
70
71 if (disp_bc->data.flag1 == 1)
72 t_vals(0) = scale_value(-disp_bc->data.value1);
73 if (disp_bc->data.flag2 == 1)
74 t_vals(1) = scale_value(-disp_bc->data.value2);
75 if (disp_bc->data.flag3 == 1)
76 t_vals(2) = scale_value(-disp_bc->data.value3);
77 if (disp_bc->data.flag4 == 1)
78 t_angles(0) = scale_value(-disp_bc->data.value4);
79 if (disp_bc->data.flag5 == 1)
80 t_angles(1) = scale_value(-disp_bc->data.value5);
81 if (disp_bc->data.flag6 == 1)
82 t_angles(2) = scale_value(-disp_bc->data.value6);
83
84 int coeff;
85 std::array<std::vector<double>, 3> coords;
86 int idx;
87
88 const bool is_rotation =
89 disp_bc->data.flag4 || disp_bc->data.flag5 || disp_bc->data.flag6;
90
91 auto lambda = [&](boost::shared_ptr<FieldEntity> field_entity_ptr) {
93 auto v = t_vals(coeff);
94 if (is_rotation) {
96 coords[0][idx], coords[1][idx], coords[2][idx]);
98 t_angles, t_coords, t_off)(coeff);
99 }
100 if (getCoords) {
101 v += coords[coeff][idx];
102 }
103
104 field_entity_ptr->getEntFieldData()[coeff] = v;
105 ++idx;
106
108 };
109
110 auto verts = bc.second->bcEnts.subset_by_type(MBVERTEX);
111 if (getCoords || is_rotation) {
112 for (auto d : {0, 1, 2})
113 coords[d].resize(verts.size());
114 CHKERR mField.get_moab().get_coords(verts, &*coords[0].begin(),
115 &*coords[1].begin(),
116 &*coords[2].begin());
117 }
118
119 if (disp_bc->data.flag1 || disp_bc->data.flag5 ||
120 disp_bc->data.flag6) {
121 idx = 0;
122 coeff = 0;
123 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
124 }
125 if (disp_bc->data.flag2 || disp_bc->data.flag4 ||
126 disp_bc->data.flag6) {
127 idx = 0;
128 coeff = 1;
129 if (nb_field_coeffs > 1)
130 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
131 }
132 if (disp_bc->data.flag3 || disp_bc->data.flag4 ||
133 disp_bc->data.flag5 || is_rotation) {
134 idx = 0;
135 coeff = 2;
136 if (nb_field_coeffs > 2)
137 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
138 }
139 }
140 }
141 }
142
143 } else {
144 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
145 "Can not lock shared pointer");
146 }
147
149}
150
151} // namespace MoFEM
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
constexpr double lambda
surface tension
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
const double c
speed of light (cm/ns)
const double v
phase velocity of light in medium (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
constexpr auto field_name
Simple interface for fast problem set-up.
Definition: BcManager.hpp:24
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:1231
Deprecated interface functions.
A specialized version of DisplacementCubitBcData that includes an additional rotation offset.
static FTensor::Tensor1< double, 3 > GetRotDisp(const FTensor::Tensor1< double, 3 > &angles, FTensor::Tensor1< double, 3 > coordinates, FTensor::Tensor1< double, 3 > offset=FTensor::Tensor1< double, 3 >{ 0., 0., 0.})
Calculates the rotated displacement given the rotation angles, coordinates, and an optional offset.
Class (Function) to enforce essential constrains.
Definition: Essential.hpp:39
Basic algebra on fields.
Definition: FieldBlas.hpp:21
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.