v0.14.0
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
94 auto v = t_vals(coeff);
95 if (is_rotation) {
97 coords[0][idx], coords[1][idx], coords[2][idx]);
99 t_angles, t_coords, t_off)(coeff);
100 }
101 if (getCoords) {
102 v += coords[coeff][idx];
103 }
104
105 field_entity_ptr->getEntFieldData()[coeff] = v;
106 ++idx;
107
109 };
110
111 auto zero_lambda =
112 [&](boost::shared_ptr<FieldEntity> field_entity_ptr) {
114 auto size = field_entity_ptr->getEntFieldData().size();
115 for (int i = coeff; i < size; i += nb_field_coeffs)
116 field_entity_ptr->getEntFieldData()[i] = 0;
118 };
119
120 auto verts = bc.second->bcEnts.subset_by_type(MBVERTEX);
121 auto not_verts = subtract(bc.second->bcEnts, verts);
122
123 if (getCoords || is_rotation) {
124 for (auto d : {0, 1, 2})
125 coords[d].resize(verts.size());
126 CHKERR mField.get_moab().get_coords(verts, &*coords[0].begin(),
127 &*coords[1].begin(),
128 &*coords[2].begin());
129 }
130
131 if (disp_bc->data.flag1 || disp_bc->data.flag5 ||
132 disp_bc->data.flag6) {
133 idx = 0;
134 coeff = 0;
135 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
136 CHKERR fb->fieldLambdaOnEntities(zero_lambda, field_name,
137 &not_verts);
138 }
139 if (disp_bc->data.flag2 || disp_bc->data.flag4 ||
140 disp_bc->data.flag6) {
141 idx = 0;
142 coeff = 1;
143 if (nb_field_coeffs > 1) {
144 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
145 CHKERR fb->fieldLambdaOnEntities(zero_lambda, field_name,
146 &not_verts);
147 }
148 }
149 if (disp_bc->data.flag3 || disp_bc->data.flag4 ||
150 disp_bc->data.flag5 || is_rotation) {
151 idx = 0;
152 coeff = 2;
153 if (nb_field_coeffs > 2) {
154 CHKERR fb->fieldLambdaOnEntities(lambda, field_name, &verts);
155 CHKERR fb->fieldLambdaOnEntities(zero_lambda, field_name,
156 &not_verts);
157 }
158 }
159 }
160 }
161 }
162
163 } else {
164 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
165 "Can not lock shared pointer");
166 }
167
169}
170
172 MoFEM::Interface &m_field, boost::shared_ptr<FEMethod> fe_ptr, double diag,
174 : mField(m_field), fePtr(fe_ptr), vDiag(diag), vRhs(rhs) {}
175
177 MOFEM_LOG_CHANNEL("WORLD");
179
180 if (auto fe_method_ptr = fePtr.lock()) {
181
182 auto bc_mng = mField.getInterface<BcManager>();
183 auto vec_mng = mField.getInterface<VecManager>();
184 auto is_mng = mField.getInterface<ISManager>();
185
186 const auto problem_name = fe_method_ptr->problemPtr->getName();
187
188 SmartPetscObj<IS> is_sum;
189
190 for (auto bc : bc_mng->getBcMapByBlockName()) {
191 if (auto disp_bc = bc.second->dispBcPtr) {
192
193 auto &bc_id = bc.first;
194
195 auto regex_str = (boost::format("%s_(.*)") % problem_name).str();
196 if (std::regex_match(bc_id, std::regex(regex_str))) {
197
198 auto [field_name, block_name] =
199 BcManager::extractStringFromBlockId(bc_id, problem_name);
200
201 MOFEM_LOG("WORLD", Sev::noisy)
202 << "Apply EssentialPreProc<DisplacementCubitBcData>: "
203 << problem_name << "_" << field_name << "_" << block_name;
204
205 const bool is_rotation =
206 disp_bc->data.flag4 || disp_bc->data.flag5 || disp_bc->data.flag6;
207
208 auto ents = bc.second->bcEnts;
209
210 std::array<SmartPetscObj<IS>, 3> is_xyz;
211 auto prb_name = fe_method_ptr->problemPtr->getName();
212
213 if (disp_bc->data.flag1 || is_rotation) {
214 CHKERR is_mng->isCreateProblemFieldAndRankLocal(
215 prb_name, ROW, field_name, 0, 0, is_xyz[0], &ents);
216 }
217 if (disp_bc->data.flag2 || is_rotation) {
218 CHKERR is_mng->isCreateProblemFieldAndRankLocal(
219 prb_name, ROW, field_name, 1, 1, is_xyz[1], &ents);
220 }
221 if (disp_bc->data.flag3 || is_rotation) {
222 CHKERR is_mng->isCreateProblemFieldAndRankLocal(
223 prb_name, ROW, field_name, 2, 2, is_xyz[2], &ents);
224 }
225
226 auto get_is_sum = [](auto is1, auto is2) {
227 IS is;
228 CHK_THROW_MESSAGE(ISExpand(is1, is2, &is), "is sum");
229 return SmartPetscObj<IS>(is);
230 };
231
232 for (auto &is : is_xyz) {
233 if (is) {
234 if (!is_sum) {
235 is_sum = is;
236 } else {
237 is_sum = get_is_sum(is_sum, is);
238 }
239 }
240 }
241 }
242 }
243 }
244
245 if (is_sum) {
246 if (auto fe_ptr = fePtr.lock()) {
247
248 auto snes_ctx = fe_ptr->snes_ctx;
249 auto ts_ctx = fe_ptr->ts_ctx;
250
252 vRhs ? vRhs : SmartPetscObj<Vec>(fe_ptr->f, true);
253
254 if (fe_ptr->vecAssembleSwitch) {
255 if ((*fe_ptr->vecAssembleSwitch) && !vRhs) {
256 CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
257 CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
258 CHKERR VecAssemblyBegin(f);
259 CHKERR VecAssemblyEnd(f);
260 *fe_ptr->vecAssembleSwitch = false;
261 }
262 }
263
264 const int *index_ptr;
265 CHKERR ISGetIndices(is_sum, &index_ptr);
266 int size;
267 CHKERR ISGetLocalSize(is_sum, &size);
268 double *a;
269 CHKERR VecGetArray(f, &a);
270
271 auto tmp_x = vectorDuplicate(f);
272 CHKERR vec_mng->setLocalGhostVector(problem_name, ROW, tmp_x,
273 INSERT_VALUES, SCATTER_FORWARD);
274 const double *u;
275 CHKERR VecGetArrayRead(tmp_x, &u);
276
277 if (snes_ctx != FEMethod::CTX_SNESNONE ||
279
280 auto x = fe_ptr->x;
281 CHKERR VecGhostUpdateBegin(x, INSERT_VALUES, SCATTER_FORWARD);
282 CHKERR VecGhostUpdateEnd(x, INSERT_VALUES, SCATTER_FORWARD);
283
284 const double *v;
285 CHKERR VecGetArrayRead(x, &v);
286
287 for (auto i = 0; i != size; ++i) {
288 a[index_ptr[i]] = vDiag * (v[index_ptr[i]] - u[index_ptr[i]]);
289 }
290
291 CHKERR VecRestoreArrayRead(x, &v);
292
293 } else {
294 for (auto i = 0; i != size; ++i) {
295 a[index_ptr[i]] = vDiag * u[index_ptr[i]];
296 }
297 }
298
299 CHKERR VecRestoreArrayRead(tmp_x, &u);
300 CHKERR VecRestoreArray(f, &a);
301 CHKERR ISRestoreIndices(is_sum, &index_ptr);
302 }
303 }
304
305 } else {
306 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
307 "Can not lock shared pointer");
308 }
309
311}
312
314 MoFEM::Interface &m_field, boost::shared_ptr<FEMethod> fe_ptr, double diag,
316 : mField(m_field), fePtr(fe_ptr), vDiag(diag), vLhs(lhs), vAO(ao) {}
317
319 MOFEM_LOG_CHANNEL("WORLD");
321
322 if (auto fe_method_ptr = fePtr.lock()) {
323
324 auto bc_mng = mField.getInterface<BcManager>();
325 auto is_mng = mField.getInterface<ISManager>();
326
327 const auto problem_name = fe_method_ptr->problemPtr->getName();
328
329 SmartPetscObj<IS> is_sum;
330
331 for (auto bc : bc_mng->getBcMapByBlockName()) {
332 if (auto disp_bc = bc.second->dispBcPtr) {
333
334 auto &bc_id = bc.first;
335
336 auto regex_str = (boost::format("%s_(.*)") % problem_name).str();
337 if (std::regex_match(bc_id, std::regex(regex_str))) {
338
339 auto [field_name, block_name] =
340 BcManager::extractStringFromBlockId(bc_id, problem_name);
341
342 MOFEM_LOG("WORLD", Sev::noisy)
343 << "Apply EssentialPreProc<DisplacementCubitBcData>: "
344 << problem_name << "_" << field_name << "_" << block_name;
345
346 const bool is_rotation =
347 disp_bc->data.flag4 || disp_bc->data.flag5 || disp_bc->data.flag6;
348
349 auto ents = bc.second->bcEnts;
350
351 std::array<SmartPetscObj<IS>, 3> is_xyz;
352 auto prb_name = fe_method_ptr->problemPtr->getName();
353
354 if (disp_bc->data.flag1 || is_rotation) {
355 CHKERR is_mng->isCreateProblemFieldAndRank(
356 prb_name, ROW, field_name, 0, 0, is_xyz[0], &ents);
357 }
358 if (disp_bc->data.flag2 || is_rotation) {
359 CHKERR is_mng->isCreateProblemFieldAndRank(
360 prb_name, ROW, field_name, 1, 1, is_xyz[1], &ents);
361 }
362 if (disp_bc->data.flag3 || is_rotation) {
363 CHKERR is_mng->isCreateProblemFieldAndRank(
364 prb_name, ROW, field_name, 2, 2, is_xyz[2], &ents);
365 }
366
367 auto get_is_sum = [](auto is1, auto is2) {
368 IS is;
369 CHK_THROW_MESSAGE(ISExpand(is1, is2, &is), "is sum");
370 return SmartPetscObj<IS>(is);
371 };
372
373 for (auto &is : is_xyz) {
374 if (is) {
375 if (!is_sum) {
376 is_sum = is;
377 } else {
378 is_sum = get_is_sum(is_sum, is);
379 }
380 }
381 }
382 }
383 }
384 }
385
386 if (is_sum) {
387 if (auto fe_ptr = fePtr.lock()) {
389 vLhs ? vLhs : SmartPetscObj<Mat>(fe_ptr->B, true);
390 // User is responsible for assembly if vLhs is provided
391 if ((*fe_ptr->matAssembleSwitch) && !vLhs) {
392 if (*fe_ptr->matAssembleSwitch) {
393 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
394 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
395 *fe_ptr->matAssembleSwitch = false;
396 }
397 }
398 if (vAO) {
399 MOFEM_LOG("WORLD", Sev::noisy) << "Apply AO to IS";
400 CHKERR AOApplicationToPetscIS(vAO, is_sum);
401 }
402 // ISView(is_sum, PETSC_VIEWER_STDOUT_WORLD);
403 CHKERR MatZeroRowsColumnsIS(B, is_sum, vDiag, PETSC_NULL, PETSC_NULL);
404 }
405 }
406
407 } else {
408 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
409 "Can not lock shared pointer");
410 }
411
413}
414
416 MoFEM::Interface &m_field, boost::shared_ptr<FEMethod> fe_ptr,
418 : mField(m_field), fePtr(fe_ptr), vRhs(rhs), sevLevel(sev) {}
419
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}
613
614} // namespace MoFEM
#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
@ ROW
Definition: definitions.h:123
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#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, severity)
Log.
Definition: LogManager.hpp:308
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
FTensor::Index< 'i', SPACE_DIM > i
static double lambda
const double c
speed of light (cm/ns)
const double v
phase velocity of light in medium (cm/ns)
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1884
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
constexpr auto field_name
Simple interface for fast problem set-up.
Definition: BcManager.hpp:25
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
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.
static UId getLoFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
static UId getHiFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
Class (Function) to enforce essential constrains on the left hand side diagonal.
Definition: Essential.hpp:47
Class (Function) to enforce essential constrains on the right hand side diagonal.
Definition: Essential.hpp:55
Class (Function) to enforce essential constrains.
Definition: Essential.hpp:39
Class (Function) to calculate residual side diagonal.
Definition: Essential.hpp:63
Basic algebra on fields.
Definition: FieldBlas.hpp:21
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:23