v0.13.2
Loading...
Searching...
No Matches
EshelbianPlasticity.hpp
Go to the documentation of this file.
1/**
2 * \file EshelbianPlasticity.hpp
3 * \brief Eshelbian plasticity interface
4 *
5 * \brief Problem implementation for mix element for large-strain elasticity
6 *
7 * \todo Implementation of plasticity
8 */
9
10#ifndef __ESHELBIAN_PLASTICITY_HPP__
11#define __ESHELBIAN_PLASTICITY_HPP__
12
13namespace EshelbianPlasticity {
14
15typedef boost::shared_ptr<MatrixDouble> MatrixPtr;
16typedef boost::shared_ptr<VectorDouble> VectorPtr;
17
19getFTensor3FromMat(MatrixDouble &m) {
21 &m(0, 0), &m(1, 0), &m(2, 0), &m(3, 0), &m(4, 0), &m(5, 0), &m(6, 0),
22 &m(7, 0), &m(8, 0), &m(9, 0), &m(10, 0), &m(11, 0), &m(12, 0), &m(13, 0),
23 &m(14, 0), &m(15, 0), &m(16, 0), &m(17, 0), &m(18, 0), &m(19, 0),
24 &m(20, 0), &m(21, 0), &m(22, 0), &m(23, 0), &m(24, 0), &m(25, 0),
25 &m(26, 0));
26}
27
28using EntData = DataForcesAndSourcesCore::EntData;
29using UserDataOperator = ForcesAndSourcesCore::UserDataOperator;
30using VolUserDataOperator = VolumeElementForcesAndSourcesCore::UserDataOperator;
31using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator;
32
34 SmartPetscObj<Mat> Suu;
35 SmartPetscObj<AO> aoSuu;
36 SmartPetscObj<Mat> SBubble;
37 SmartPetscObj<AO> aoSBubble;
38 SmartPetscObj<Mat> SOmega;
39 SmartPetscObj<AO> aoSOmega;
40 SmartPetscObj<Mat> Sw;
41 SmartPetscObj<AO> aoSw;
42 EpElementBase() = default;
43 virtual ~EpElementBase() = default;
44
45 MoFEMErrorCode addStreachSchurMatrix(SmartPetscObj<Mat> &Suu,
46 SmartPetscObj<AO> &aoSuu) {
48 this->Suu = Suu;
49 this->aoSuu = aoSuu;
51 }
52
53 MoFEMErrorCode addBubbleSchurMatrix(SmartPetscObj<Mat> &SBubble,
54 SmartPetscObj<AO> &aoSBubble) {
56 this->SBubble = SBubble;
57 this->aoSBubble = aoSBubble;
59 }
60
61 MoFEMErrorCode addSpatialDispStressSchurMatrix(SmartPetscObj<Mat> &Sw,
62 SmartPetscObj<AO> &aoSw) {
64 this->Sw = Sw;
65 this->aoSw = aoSw;
67 }
68
69 MoFEMErrorCode addOmegaSchurMatrix(SmartPetscObj<Mat> &SOmega,
70 SmartPetscObj<AO> &aoSOmega) {
72 this->SOmega = SOmega;
73 this->aoSOmega = aoSOmega;
75 }
76};
77
78template <typename E> struct EpElement : public E, EpElementBase {
79 EpElement(MoFEM::Interface &m_field) : E(m_field), EpElementBase() {}
80};
81
82template <> struct EpElement<BasicMethod> : public FEMethod, EpElementBase {
84};
85
86template <> struct EpElement<FEMethod> : public FEMethod, EpElementBase {
88};
89
90struct EpFEMethod : EpElement<FEMethod> {
92 MoFEMErrorCode preProcess() {
94 if (Suu)
95 CHKERR MatZeroEntries(Suu);
96 if (SBubble)
97 CHKERR MatZeroEntries(SBubble);
98 if (Sw)
99 CHKERR MatZeroEntries(Sw);
100 if (SOmega)
101 CHKERR MatZeroEntries(SOmega);
103 }
104
105 MoFEMErrorCode postProcess() {
107 auto assemble = [](Mat a) {
109 if (a) {
110 CHKERR MatAssemblyBegin(a, MAT_FINAL_ASSEMBLY);
111 CHKERR MatAssemblyEnd(a, MAT_FINAL_ASSEMBLY);
112 }
114 };
115 CHKERR assemble(Suu);
116 CHKERR assemble(SBubble);
117 CHKERR assemble(SOmega);
118 CHKERR assemble(Sw);
119 // std::string wait;
120 // CHKERR MatView(SOmega, PETSC_VIEWER_DRAW_WORLD);
121 // std::cin >> wait;
122 // CHKERR MatView(Sw, PETSC_VIEWER_DRAW_WORLD);
123 // std::cin >> wait;
124 // CHKERR MatView(SBubble, PETSC_VIEWER_DRAW_WORLD);
125 // std::cin >> wait;
126 // CHKERR MatView(Suu, PETSC_VIEWER_DRAW_WORLD);
127 // std::cin >> wait;
128 // CHKERR MatView(getTSB(), PETSC_VIEWER_DRAW_WORLD);
129 // std::cin >> wait;
131 }
132};
133
134struct PhysicalEquations;
136 : public boost::enable_shared_from_this<DataAtIntegrationPts> {
137
138 MatrixDouble approxPAtPts;
139 MatrixDouble approxSigmaAtPts;
140 MatrixDouble divPAtPts;
141 MatrixDouble divSigmaAtPts;
142 MatrixDouble wAtPts;
143 MatrixDouble wDotAtPts;
144 MatrixDouble wDotDotAtPts;
146 MatrixDouble streachTensorAtPts;
147
152 MatrixDouble rotAxisAtPts;
153 MatrixDouble rotAxisDotAtPts;
154 MatrixDouble WAtPts;
155 MatrixDouble W0AtPts;
156 MatrixDouble GAtPts;
157 MatrixDouble G0AtPts;
158
159 MatrixDouble hAtPts;
160 MatrixDouble rotMatAtPts;
161 MatrixDouble diffRotMatAtPts;
162 MatrixDouble PAtPts;
163 MatrixDouble SigmaAtPts;
164 VectorDouble phiAtPts;
165 MatrixDouble flowAtPts;
166
167 MatrixDouble P_dh0;
168 MatrixDouble P_dh1;
169 MatrixDouble P_dh2;
170 MatrixDouble P_dH0;
171 MatrixDouble P_dH1;
172 MatrixDouble P_dH2;
173 MatrixDouble Sigma_dh0;
174 MatrixDouble Sigma_dh1;
175 MatrixDouble Sigma_dh2;
176 MatrixDouble Sigma_dH0;
177 MatrixDouble Sigma_dH1;
178 MatrixDouble Sigma_dH2;
179 MatrixDouble phi_dh;
180 MatrixDouble phi_dH;
181 MatrixDouble Flow_dh0;
182 MatrixDouble Flow_dh1;
183 MatrixDouble Flow_dh2;
184 MatrixDouble Flow_dH0;
185 MatrixDouble Flow_dH1;
186 MatrixDouble Flow_dH2;
187
188 MatrixDouble eigenVals;
189 MatrixDouble eigenVecs;
190 VectorInt nbUniq;
191
193 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
195 }
197 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &approxPAtPts);
198 }
199
201 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divPAtPts);
202 }
203
205 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divSigmaAtPts);
206 }
207
209 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wAtPts);
210 }
211
213 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wDotAtPts);
214 }
215
217 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wDotDotAtPts);
218 }
219
221 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
223 }
224
226 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
228 }
229
231 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
233 }
234
236 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxisAtPts);
237 }
238
240 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
242 }
243
245 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
246 }
247
249 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
250 }
251
252 // Not really data at integration points, used to calculate Schur complement
253
255
256 std::string rowField;
257 std::string colField;
262 MatrixDouble M;
263 VectorInt rowInd;
264 VectorInt colInd;
265
267
268 BlockMatData(const std::string row_field, const std::string col_field,
269 EntityType row_type, EntityType col_type, int row_side,
270 int col_side, const MatrixDouble &m, const VectorInt row_ind,
271 VectorInt col_ind)
272 : rowField(row_field), colField(col_field), rowType(row_type),
273 colType(col_type), rowSide(row_side), colSide(col_side),
274 setAtElement(true) {
275
276 M.resize(m.size1(), m.size2(), false);
277 noalias(M) = m;
278 rowInd.resize(row_ind.size(), false);
279 noalias(rowInd) = row_ind;
280 colInd.resize(col_ind.size(), false);
281 noalias(colInd) = col_ind;
282 }
283
284 void setInd(const VectorInt &row_ind, const VectorInt &col_ind) const {
285 auto &const_row_ind = const_cast<VectorInt &>(rowInd);
286 auto &const_col_ind = const_cast<VectorInt &>(colInd);
287 const_row_ind.resize(row_ind.size(), false);
288 noalias(const_row_ind) = row_ind;
289 const_col_ind.resize(col_ind.size(), false);
290 noalias(const_col_ind) = col_ind;
291 }
292
293 void setMat(const MatrixDouble &m) const {
294 auto &const_m = const_cast<MatrixDouble &>(M);
295 const_m.resize(m.size1(), m.size2(), false);
296 noalias(const_m) = m;
297 }
298
299 void addMat(const MatrixDouble &m) const {
300 auto &const_m = const_cast<MatrixDouble &>(M);
301 const_m += m;
302 }
303
304 void clearMat() const {
305 auto &const_m = const_cast<MatrixDouble &>(M);
306 const_m.clear();
307 }
308
309 void setSetAtElement() const {
310 bool &set = const_cast<bool &>(setAtElement);
311 set = true;
312 }
313
314 void unSetAtElement() const {
315 bool &set = const_cast<bool &>(setAtElement);
316 set = false;
317 }
318 };
319
320 typedef multi_index_container<
321 BlockMatData,
322 indexed_by<
323
324 ordered_unique<
325
326 composite_key<
327 BlockMatData,
328
329 member<BlockMatData, std::string, &BlockMatData::rowField>,
330 member<BlockMatData, std::string, &BlockMatData::colField>,
331 member<BlockMatData, EntityType, &BlockMatData::rowType>,
332 member<BlockMatData, EntityType, &BlockMatData::colType>,
333 member<BlockMatData, int, &BlockMatData::rowSide>,
334 member<BlockMatData, int, &BlockMatData::colSide>
335
336 >>,
337 ordered_non_unique<
338
339 composite_key<
340 BlockMatData,
341
342 member<BlockMatData, std::string, &BlockMatData::rowField>,
343 member<BlockMatData, std::string, &BlockMatData::colField>,
344 member<BlockMatData, EntityType, &BlockMatData::rowType>,
345 member<BlockMatData, EntityType, &BlockMatData::colType>
346
347 >>,
348 ordered_non_unique<
349
350 composite_key<
351 BlockMatData,
352
353 member<BlockMatData, std::string, &BlockMatData::rowField>,
354 member<BlockMatData, std::string, &BlockMatData::colField>
355
356 >>,
357 ordered_non_unique<
358 member<BlockMatData, std::string, &BlockMatData::rowField>>,
359 ordered_non_unique<
360 member<BlockMatData, std::string, &BlockMatData::colField>>
361
362 >>
364
368
369 boost::shared_ptr<PhysicalEquations> physicsPtr;
370};
371
372struct OpJacobian;
373
375
380
383
385 PhysicalEquations(const int size_active, const int size_dependent)
386 : activeVariables(size_active, 0), dependentVariables(size_dependent, 0),
387 dependentVariablesDirevatives(size_dependent * size_active, 0) {}
389
390 virtual MoFEMErrorCode recordTape(const int tag, DTensor2Ptr *t_h) = 0;
391
392 virtual OpJacobian *
393 returnOpJacobian(const std::string &field_name, const int tag,
394 const bool eval_rhs, const bool eval_lhs,
395 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
396 boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
397
398 std::vector<double> activeVariables;
399 std::vector<double> dependentVariables;
400 std::vector<double> dependentVariablesDirevatives;
401
402 /** \name Active variables */
403
404 /**@{*/
405
406 template <int S>
407 inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
408 return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
409 &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
410 }
411
412 template <int S>
413 inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
414 return DTensor0Ptr(&v[S + 0]);
415 }
416
417 template <int S0, int S1, int S2>
418 inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v) {
419 constexpr int A00 = 18 * S0 + 18 * 0 + 9 * S1 + 3 * S2;
420 constexpr int A01 = 18 * S0 + 18 * 1 + 9 * S1 + 3 * S2;
421 constexpr int A02 = 18 * S0 + 18 * 2 + 9 * S1 + 3 * S2;
422 constexpr int A10 = 18 * S0 + 18 * 3 + 9 * S1 + 3 * S2;
423 constexpr int A11 = 18 * S0 + 18 * 4 + 9 * S1 + 3 * S2;
424 constexpr int A12 = 18 * S0 + 18 * 5 + 9 * S1 + 3 * S2;
425 constexpr int A20 = 18 * S0 + 18 * 6 + 9 * S1 + 3 * S2;
426 constexpr int A21 = 18 * S0 + 18 * 7 + 9 * S1 + 3 * S2;
427 constexpr int A22 = 18 * S0 + 18 * 8 + 9 * S1 + 3 * S2;
428 return DTensor3Ptr(
429
430 &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
431 &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
432
433 &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
434 &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
435
436 &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
437 &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
438
439 );
440 }
441
442 /**@}*/
443
444 /** \name Active variables */
445
446 /**@{*/
447
448 inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
449 inline DTensor2Ptr get_H() { return get_VecTensor2<9>(activeVariables); }
450
451 /**@}*/
452
453 /** \name Dependent variables */
454
455 /**@{*/
456
457 inline DTensor2Ptr get_P() { return get_VecTensor2<0>(dependentVariables); }
459 return get_VecTensor2<9>(dependentVariables);
460 }
462 return get_VecTensor0<18>(dependentVariables);
463 }
464 inline double &get_PhiRef() { return dependentVariables[18]; }
466 return get_VecTensor2<9 + 9 + 1>(dependentVariables);
467 }
468
469 /**@}*/
470
471 /** \name Derivatives of dependent variables */
472
473 /**@{*/
474
476 return get_vecTensor3<0, 0, 0>(dependentVariablesDirevatives);
477 }
479 return get_vecTensor3<0, 1, 0>(dependentVariablesDirevatives);
480 }
482 return get_vecTensor3<0, 0, 1>(dependentVariablesDirevatives);
483 }
485 return get_vecTensor3<0, 1, 1>(dependentVariablesDirevatives);
486 }
488 return get_vecTensor3<0, 0, 2>(dependentVariablesDirevatives);
489 }
491 return get_vecTensor3<0, 1, 2>(dependentVariablesDirevatives);
492 }
493
495 return get_vecTensor3<9, 0, 0>(dependentVariablesDirevatives);
496 }
498 return get_vecTensor3<9, 1, 0>(dependentVariablesDirevatives);
499 }
501 return get_vecTensor3<9, 0, 1>(dependentVariablesDirevatives);
502 }
504 return get_vecTensor3<9, 1, 1>(dependentVariablesDirevatives);
505 }
507 return get_vecTensor3<9, 0, 2>(dependentVariablesDirevatives);
508 }
510 return get_vecTensor3<9, 1, 2>(dependentVariablesDirevatives);
511 }
512
514 return get_VecTensor2<18 * (9 + 9) + 0>(dependentVariablesDirevatives);
515 }
517 return get_VecTensor2<18 * (9 + 9) + 9>(dependentVariablesDirevatives);
518 }
519
521 return get_vecTensor3<9 + 9 + 1, 0, 0>(dependentVariablesDirevatives);
522 }
524 return get_vecTensor3<9 + 9 + 1, 1, 0>(dependentVariablesDirevatives);
525 }
527 return get_vecTensor3<9 + 9 + 1, 0, 1>(dependentVariablesDirevatives);
528 }
530 return get_vecTensor3<9 + 9 + 1, 1, 1>(dependentVariablesDirevatives);
531 }
533 return get_vecTensor3<9 + 9 + 1, 0, 2>(dependentVariablesDirevatives);
534 }
536 return get_vecTensor3<9 + 9 + 1, 1, 2>(dependentVariablesDirevatives);
537 }
538
539 /**@}*/
540};
541
542struct BcDisp {
543 BcDisp(std::string name, std::vector<double> &attr, Range &faces);
544 std::string blockName;
546 VectorDouble3 vals;
547 VectorInt3 flags;
548};
549typedef std::vector<BcDisp> BcDispVec;
550
551struct BcRot {
552 BcRot(std::string name, std::vector<double> &attr, Range &faces);
553 std::string blockName;
555 VectorDouble3 vals;
556 double theta;
557};
558typedef std::vector<BcRot> BcRotVec;
559
560typedef std::vector<Range> TractionFreeBc;
561
563 TractionBc(std::string name, std::vector<double> &attr, Range &faces);
564 std::string blockName;
566 VectorDouble3 vals;
567 VectorInt3 flags;
568};
569typedef std::vector<TractionBc> TractionBcVec;
570
572 const int tAg; ///< adol-c tape
573 const bool evalRhs;
574 const bool evalLhs;
575 boost::shared_ptr<DataAtIntegrationPts>
576 dataAtPts; ///< data at integration pts
577 boost::shared_ptr<PhysicalEquations>
578 physicsPtr; ///< material physical equations
579
580 OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs,
581 const bool eval_lhs,
582 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
583 boost::shared_ptr<PhysicalEquations> &physics_ptr)
584 : UserDataOperator(field_name, OPROW), tAg(tag), evalRhs(eval_rhs),
585 evalLhs(eval_lhs), dataAtPts(data_ptr), physicsPtr(physics_ptr) {}
586
587 virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
588 virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
589
590 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
591};
592
593template <typename T> struct OpAssembleBasic : public T {
594
596
597 boost::shared_ptr<DataAtIntegrationPts>
598 dataAtPts; ///< data at integration pts
599
600 OpAssembleBasic(const std::string &field_name,
601 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
602 const char type)
603 : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false) {}
604
605 OpAssembleBasic(const std::string &row_field, const std::string &col_field,
606 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
607 const char type, const bool assemble_symmetry)
608 : T(row_field, col_field, type, false), dataAtPts(data_ptr),
609 assembleSymmetry(assemble_symmetry) {}
610
611 VectorDouble nF; ///< local right hand side vector
612 MatrixDouble K; ///< local tangent matrix
613 MatrixDouble transposeK;
614
615 virtual MoFEMErrorCode integrate(EntData &data) {
617 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
619 }
620
621 virtual MoFEMErrorCode integrate(int row_side, EntityType row_type,
622 EntData &data) {
624 CHKERR integrate(data);
626 }
627
628 virtual MoFEMErrorCode assemble(EntData &data) {
630 double *vec_ptr = &*nF.begin();
631 int nb_dofs = data.getIndices().size();
632 int *ind_ptr = &*data.getIndices().begin();
633 CHKERR VecSetValues(this->getTSf(), nb_dofs, ind_ptr, vec_ptr, ADD_VALUES);
635 }
636
637 virtual MoFEMErrorCode assemble(int row_side, EntityType row_type,
638 EntData &data) {
640 CHKERR assemble(data);
642 }
643
644 virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
646 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
648 }
649
650 virtual MoFEMErrorCode assemble(int row_side, int col_side,
651 EntityType row_type, EntityType col_type,
652 EntData &row_data, EntData &col_data) {
654
655 auto &bmc = dataAtPts->blockMatContainor;
656 int *row_ind_ptr = &*row_data.getIndices().begin();
657 int *col_ind_ptr = &*col_data.getIndices().begin();
658 int row_nb_dofs = row_data.getIndices().size();
659 int col_nb_dofs = col_data.getIndices().size();
660
661 auto add_block = [&](auto &row_name, auto &col_name, auto row_side,
662 auto col_side, auto row_type, auto col_type,
663 const auto &m, const auto &row_ind,
664 const auto &col_ind) {
665 auto it = bmc.get<0>().find(boost::make_tuple(
666 row_name, col_name, row_type, col_type, row_side, col_side));
667 if (it != bmc.get<0>().end()) {
668 it->setMat(m);
669 it->setInd(row_ind, col_ind);
670 it->setSetAtElement();
671 } else
673 row_name, col_name, row_type, col_type, row_side, col_side, m,
674 row_ind, col_ind));
675 };
676
677 add_block(this->rowFieldName, this->colFieldName, row_side, col_side,
678 row_type, col_type, K, row_data.getIndices(),
679 col_data.getIndices());
680 if (assembleSymmetry) {
681 transposeK.resize(col_nb_dofs, row_nb_dofs, false);
682 noalias(transposeK) = trans(K);
683 add_block(this->colFieldName, this->rowFieldName, col_side, row_side,
684 col_type, row_type, transposeK, col_data.getIndices(),
685 row_data.getIndices());
686 }
687
688 double *mat_ptr = &*K.data().begin();
689 CHKERR MatSetValues(this->getTSB(), row_nb_dofs, row_ind_ptr, col_nb_dofs,
690 col_ind_ptr, mat_ptr, ADD_VALUES);
691 if (assembleSymmetry) {
692 double *mat_ptr = &*transposeK.data().begin();
693 CHKERR MatSetValues(this->getTSB(), col_nb_dofs, col_ind_ptr, row_nb_dofs,
694 row_ind_ptr, mat_ptr, ADD_VALUES);
695 }
696
698 }
699
700 MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
702 if (data.getIndices().empty())
704 nF.resize(data.getIndices().size(), false);
705 nF.clear();
706 CHKERR integrate(side, type, data);
707 CHKERR assemble(side, type, data);
709 }
710
711 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
712 EntityType col_type, EntData &row_data,
713 EntData &col_data) {
715 if (row_data.getIndices().empty())
717 if (col_data.getIndices().empty())
719 K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
720 K.clear();
721 CHKERR integrate(row_data, col_data);
722 CHKERR assemble(row_side, col_side, row_type, col_type, row_data, col_data);
724 }
725};
726
727struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
728 OpAssembleVolume(const std::string &field,
729 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
730 const char type)
731 : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
732
733 OpAssembleVolume(const std::string &row_field, const std::string &col_field,
734 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
735 const char type, const bool assemble_symmetry)
736 : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
737 type, assemble_symmetry) {}
738};
739
740struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
741 OpAssembleFace(const std::string &field,
742 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
743 const char type)
744 : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
745
746 OpAssembleFace(const std::string &row_field, const std::string &col_field,
747 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
748 const char type, const bool assemble_symmetry)
749 : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
750 type, assemble_symmetry) {}
751};
752
754 boost::shared_ptr<DataAtIntegrationPts>
755 dataAtPts; ///< data at integration pts
757 const std::string &field_name,
758 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
759 : VolUserDataOperator(field_name, OPROW), dataAtPts(data_ptr) {}
760 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
761};
762
764 const double alphaW;
765 const double alphaRho;
767 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
768 const double alpha, const double rho)
769 : OpAssembleVolume(field_name, data_ptr, OPROW), alphaW(alpha),
770 alphaRho(rho) {}
771 MoFEMErrorCode integrate(EntData &data);
772};
773
775 OpSpatialRotation(const std::string &field_name,
776 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
777 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
778 MoFEMErrorCode integrate(EntData &data);
779};
780
782 const double alphaU;
783 OpSpatialPhysical(const std::string &field_name,
784 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
785 const double alpha)
786 : OpAssembleVolume(field_name, data_ptr, OPROW), alphaU(alpha) {}
787 MoFEMErrorCode integrate(EntData &data);
788};
789
792 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
793 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
794 MoFEMErrorCode integrate(EntData &data);
795};
796
799 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
800 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
801 MoFEMErrorCode integrate(EntData &data);
802};
803
806 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
807 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
808 MoFEMErrorCode integrate(EntData &data);
809};
810
811struct OpDispBc : public OpAssembleFace {
812 boost::shared_ptr<BcDispVec> bcDispPtr;
813 OpDispBc(const std::string &field_name,
814 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
815 boost::shared_ptr<BcDispVec> &bc_disp_ptr)
816 : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr) {}
817 MoFEMErrorCode integrate(EntData &data);
818};
819
821 boost::shared_ptr<BcDispVec> bcDispPtr;
822 OpDispBc_dx(const std::string &row_field_name,
823 const std::string &col_field_name,
824 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
825 boost::shared_ptr<BcDispVec> &bc_disp_ptr)
826 : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
827 false),
828 bcDispPtr(bc_disp_ptr) {}
829 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
830};
831
833 boost::shared_ptr<BcRotVec> bcRotPtr;
834 OpRotationBc(const std::string &field_name,
835 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
836 boost::shared_ptr<BcRotVec> &bc_rot_ptr)
837 : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
838 MoFEMErrorCode integrate(EntData &data);
839};
840
842 boost::shared_ptr<BcRotVec> bcRotPtr;
843 OpRotationBc_dx(const std::string &row_field_name,
844 const std::string &col_field_name,
845 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
846 boost::shared_ptr<BcRotVec> &bc_rot_ptr)
847 : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
848 false),
849 bcRotPtr(bc_rot_ptr) {}
850 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
851};
852
853struct FeTractionBc;
854
856 OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
857 : FaceUserDataOperator(row_field, FaceUserDataOperator::OPROW),
858 bcFe(bc_fe) {}
859 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
860
861protected:
863 MatrixDouble matPP;
864 MatrixDouble vecPv;
865};
866
867struct FeTractionBc : public FEMethod {
869 FeTractionBc(MoFEM::Interface &m_field, const std::string field_name,
870 boost::shared_ptr<TractionBcVec> &bc)
871 : FEMethod(), mField(m_field), bcData(bc), fieldName(field_name) {}
872 MoFEMErrorCode preProcess();
873 MoFEMErrorCode postProcess() { return 0; }
874
875 friend MoFEMErrorCode OpTractionBc::doWork(int side, EntityType type,
876 EntData &data);
877
878protected:
880 boost::shared_ptr<TractionBcVec> bcData;
881 std::string fieldName;
882};
883
884template <>
886 EpElement(MoFEM::Interface &m_field, const std::string field_name,
887 boost::shared_ptr<TractionBcVec> &bc)
888 : FeTractionBc(m_field, field_name, bc), EpElementBase() {}
889 MoFEMErrorCode postProcess() {
893 }
894};
895
897 OpSpatialEquilibrium_dw_dP(const std::string &row_field,
898 const std::string &col_field,
899 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
900 const bool assemble_off = false)
901 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
902 assemble_off) {
903 sYmm = false;
904 }
905 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
906};
907
909 const double alphaW;
910 const double alphaRho;
911 OpSpatialEquilibrium_dw_dw(const std::string &row_field,
912 const std::string &col_field,
913 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
914 const double alpha, const double rho)
915 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
916 alphaW(alpha), alphaRho(rho) {
917 sYmm = true;
918 }
919 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
920};
921
923 const double alphaU;
924 OpSpatialPhysical_du_du(const std::string &row_field,
925 const std::string &col_field,
926 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
927 const double alpha)
928 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
929 alphaU(alpha) {
930 sYmm = false;
931 }
932 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
933};
934
936 OpSpatialPhysical_du_dP(const std::string &row_field,
937 const std::string &col_field,
938 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
939 const bool assemble_off = false)
940 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
941 assemble_off) {
942 sYmm = false;
943 }
944
945 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
946};
947
950 const std::string &row_field, const std::string &col_field,
951 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
952 const bool assemble_off = false)
953 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
954 assemble_off) {
955 sYmm = false;
956 }
957
958 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
959};
960
962 OpSpatialPhysical_du_domega(const std::string &row_field,
963 const std::string &col_field,
964 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
965 const bool assemble_off = false)
966 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
967 assemble_off) {
968 sYmm = false;
969 }
970 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
971};
972
974 OpSpatialPhysical_du_dx(const std::string &row_field,
975 const std::string &col_field,
976 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
977 const bool assemble_off = false)
978 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
979 assemble_off) {
980 sYmm = false;
981 }
982 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
983};
984
986 OpSpatialRotation_domega_dP(const std::string &row_field,
987 const std::string &col_field,
988 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
989 const bool assemble_off)
990 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
991 assemble_off) {
992 sYmm = false;
993 }
994 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
995};
996
999 const std::string &row_field, const std::string &col_field,
1000 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1001 const bool assemble_off)
1002 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
1003 assemble_off) {
1004 sYmm = false;
1005 }
1006 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1007};
1008
1011 const std::string &row_field, const std::string &col_field,
1012 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1013 const bool assemble_off = false)
1014 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
1015 assemble_off) {
1016 sYmm = false;
1017 }
1018 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1019};
1020
1023 const std::string &row_field, const std::string &col_field,
1024 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1025 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1026 sYmm = false;
1027 }
1028 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1029};
1030
1033 const std::string &row_field, const std::string &col_field,
1034 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1035 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1036 sYmm = false;
1037 }
1038 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1039};
1040
1042
1043 moab::Interface &postProcMesh;
1044 std::vector<EntityHandle> &mapGaussPts;
1045 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1046
1047 OpPostProcDataStructure(moab::Interface &post_proc_mesh,
1048 std::vector<EntityHandle> &map_gauss_pts,
1049 const std::string field_name,
1050 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1052 postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
1053 dataAtPts(data_ptr) {}
1054
1055 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
1056};
1057
1059 OpSpatialPrj(const std::string &row_field,
1060 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1061 : OpAssembleVolume(row_field, data_ptr, OPROW) {}
1062 MoFEMErrorCode integrate(EntData &row_data);
1063};
1064
1066 OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field,
1067 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1068 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1069 // FIXME: That is symmetric
1070 sYmm = false;
1071 }
1072 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1073};
1074
1076 OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field,
1077 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1078 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1079 sYmm = false;
1080 }
1081 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1082};
1083
1085 OpSpatialSchurBegin(const std::string &row_field,
1086 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1087 : OpAssembleVolume(row_field, data_ptr, OPROW) {
1088 sYmm = false;
1089 }
1090 MoFEMErrorCode integrate(EntData &row_data) { return 0; }
1091 MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data);
1092};
1093
1095 OpSpatialPreconditionMass(const std::string &row_field, MatrixPtr m_ptr)
1096 : OpAssembleVolume(row_field, nullptr, OPROW), mPtr(m_ptr) {
1097 sYmm = false;
1098 }
1099 MoFEMErrorCode integrate(EntData &row_data);
1100 MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data) {
1101 return 0;
1102 }
1103
1104private:
1106};
1107
1109 OpSpatialSchurEnd(const std::string &row_field,
1110 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1111 const double eps)
1112 : OpAssembleVolume(row_field, data_ptr, OPROW), eps(eps) {
1113 sYmm = false;
1114 }
1115
1116 MoFEMErrorCode integrate(EntData &row_data) { return 0; }
1117 MoFEMErrorCode assemble(int side, EntityType type, EntData &data);
1118 MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
1119 return assemble(side, type, data);
1120 }
1121
1122private:
1123 const double eps;
1124
1125 MatrixDouble invMat;
1126 VectorInt iPIV;
1127 VectorDouble lapackWork;
1128
1129 map<std::string, MatrixDouble> invBlockMat;
1130 map<std::string, DataAtIntegrationPts::BlockMatContainor> blockMat;
1131};
1132
1134
1136 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
1137 boost::shared_ptr<double> &e)
1139 dataAtPts(data_at_pts), energy(e) {}
1140
1141 MoFEMErrorCode doWork(int side, EntityType type,
1142 DataForcesAndSourcesCore::EntData &data);
1143
1144private:
1145 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1146 boost::shared_ptr<double> energy;
1147};
1148
1150
1152
1153 MoFEMErrorCode doWork(int side, EntityType type,
1154 DataForcesAndSourcesCore::EntData &data);
1155};
1156
1158
1159 /**
1160 * \brief Getting interface of core database
1161 * @param uuid unique ID of interface
1162 * @param iface returned pointer to interface
1163 * @return error code
1164 */
1165 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
1166 UnknownInterface **iface) const;
1167
1169
1170 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1171 boost::shared_ptr<PhysicalEquations> physicalEquations;
1172
1173 boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeRhs;
1174 boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeLhs;
1175 boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcLhs;
1176 boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcRhs;
1177 boost::shared_ptr<EpFEMethod> schurAssembly;
1178
1179 SmartPetscObj<DM> dM; ///< Coupled problem all fields
1180 SmartPetscObj<DM> dmElastic; ///< Elastic problem
1181 SmartPetscObj<DM> dmElasticSchurStreach; ///< Sub problem of dmElastic Schur
1182 SmartPetscObj<DM> dmElasticSchurBubble; ///< Sub problem of dmElastic Schur
1183 SmartPetscObj<DM>
1184 dmElasticSchurSpatialDisp; ///< Sub problem of dmElastic Schur
1185 SmartPetscObj<DM> dmElasticSchurOmega; ///< Sub problem of dmElastic Schur
1186 SmartPetscObj<DM> dmMaterial; ///< Material problem
1187
1188 const std::string piolaStress;
1189 const std::string eshelbyStress;
1190 const std::string spatialDisp;
1191 const std::string materialDisp;
1192 const std::string streachTensor;
1193 const std::string rotAxis;
1194 const std::string materialGradient;
1195 const std::string tauField;
1196 const std::string lambdaField;
1197 const std::string bubbleField;
1198
1199 const std::string elementVolumeName;
1200 const std::string naturalBcElement;
1201 const std::string essentialBcElement;
1202
1204 virtual ~EshelbianCore();
1205
1208 double alphaU;
1209 double alphaW;
1210 double alphaRho;
1211 double precEps;
1212
1213 MoFEMErrorCode getOptions();
1214
1215 boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
1216 boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
1217 boost::shared_ptr<TractionBcVec> bcSpatialTraction;
1218 boost::shared_ptr<TractionFreeBc> bcSpatialFreeTraction;
1219
1220 template <typename BC>
1221 MoFEMErrorCode getBc(boost::shared_ptr<BC> &bc_vec_ptr,
1222 const std::string block_set_name,
1223 const int nb_attributes) {
1226 auto block_name = it->getName();
1227 if (block_name.compare(0, block_set_name.length(), block_set_name) == 0) {
1228 std::vector<double> block_attributes;
1229 CHKERR it->getAttributes(block_attributes);
1230 if (block_attributes.size() != nb_attributes) {
1231 SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1232 "In block %s expected %d attributes, but given %d",
1233 it->getName().c_str(), nb_attributes,
1234 block_attributes.size());
1235 }
1236 Range faces;
1237 CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
1238 true);
1239 bc_vec_ptr->emplace_back(block_name, block_attributes, faces);
1240 }
1241 }
1243 }
1244
1245 inline MoFEMErrorCode getSpatialDispBc() {
1246 bcSpatialDispVecPtr = boost::make_shared<BcDispVec>();
1247 return getBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC", 6);
1248 }
1249
1250 inline MoFEMErrorCode getSpatialRotationBc() {
1251 bcSpatialRotationVecPtr = boost::make_shared<BcRotVec>();
1252 return getBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC", 4);
1253 }
1254
1255 inline MoFEMErrorCode getSpatialTractionBc() {
1256 bcSpatialTraction = boost::make_shared<TractionBcVec>();
1257 return getBc(bcSpatialTraction, "SPATIAL_TRACTION_BC", 6);
1258 }
1259
1260 MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset,
1261 boost::shared_ptr<TractionFreeBc> &bc_ptr,
1262 const std::string disp_block_set_name,
1263 const std::string rot_block_set_name);
1264 inline MoFEMErrorCode
1267 boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
1268 return getTractionFreeBc(meshset, bcSpatialFreeTraction, "SPATIAL_DISP_BC",
1269 "SPATIAL_ROTATION_BC");
1270 }
1271
1272 MoFEMErrorCode addFields(const EntityHandle meshset = 0);
1273 MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset = 0);
1274 MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset = 0);
1275 MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
1276
1277 MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape,
1278 const double lambda,
1279 const double mu,
1280 const double sigma_y);
1281
1282 MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha,
1283 const double beta,
1284 const double lambda,
1285 const double sigma_y);
1286
1287 MoFEMErrorCode setBaseVolumeElementOps(
1288 const int tag, const bool do_rhs, const bool do_lhs,
1289 boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe);
1290
1291 MoFEMErrorCode setGenericVolumeElementOps(
1292 const int tag, const bool add_elastic, const bool add_material,
1293 boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_rhs,
1294 boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_lhs);
1295
1296 MoFEMErrorCode setGenericFaceElementOps(
1297 const bool add_elastic, const bool add_material,
1298 boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_rhs,
1299 boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_lhs);
1300
1301 MoFEMErrorCode setElasticElementOps(const int tag);
1302 MoFEMErrorCode setElasticElementToTs(DM dm);
1303
1304 MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f);
1305 MoFEMErrorCode solveElastic(TS ts, Vec x);
1306
1307 MoFEMErrorCode postProcessResults(const int tag, const std::string file);
1308};
1309
1310} // namespace EshelbianPlasticity
1311
1312#endif //__ESHELBIAN_PLASTICITY_HPP__
constexpr double a
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ L2
field with C-1 continuity
Definition: definitions.h:88
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ BLOCKSET
Definition: definitions.h:148
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#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 MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FTensor::Index< 'm', SPACE_DIM > m
constexpr double lambda
#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.
auto bit
set bit
const double v
phase velocity of light in medium (cm/ns)
const double T
std::vector< BcDisp > BcDispVec
DataForcesAndSourcesCore::EntData EntData
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
boost::shared_ptr< VectorDouble > VectorPtr
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat(MatrixDouble &m)
std::vector< Range > TractionFreeBc
double f(const double v)
ForcesAndSourcesCore::UserDataOperator UserDataOperator
std::vector< BcRot > BcRotVec
std::vector< TractionBc > TractionBcVec
boost::shared_ptr< MatrixDouble > MatrixPtr
FaceElementForcesAndSourcesCore::UserDataOperator FaceUserDataOperator
double rho
Definition: plastic.cpp:101
constexpr auto field_name
constexpr double mu
BlockMatData(const std::string row_field, const std::string col_field, EntityType row_type, EntityType col_type, int row_side, int col_side, const MatrixDouble &m, const VectorInt row_ind, VectorInt col_ind)
void setInd(const VectorInt &row_ind, const VectorInt &col_ind) const
multi_index_container< BlockMatData, indexed_by< ordered_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField >, member< BlockMatData, EntityType, &BlockMatData::rowType >, member< BlockMatData, EntityType, &BlockMatData::colType >, member< BlockMatData, int, &BlockMatData::rowSide >, member< BlockMatData, int, &BlockMatData::colSide > > >, ordered_non_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField >, member< BlockMatData, EntityType, &BlockMatData::rowType >, member< BlockMatData, EntityType, &BlockMatData::colType > > >, ordered_non_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField > > >, ordered_non_unique< member< BlockMatData, std::string, &BlockMatData::rowField > >, ordered_non_unique< member< BlockMatData, std::string, &BlockMatData::colField > > > > BlockMatContainor
boost::shared_ptr< PhysicalEquations > physicsPtr
EpElement(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
MoFEMErrorCode addBubbleSchurMatrix(SmartPetscObj< Mat > &SBubble, SmartPetscObj< AO > &aoSBubble)
MoFEMErrorCode addSpatialDispStressSchurMatrix(SmartPetscObj< Mat > &Sw, SmartPetscObj< AO > &aoSw)
MoFEMErrorCode addOmegaSchurMatrix(SmartPetscObj< Mat > &SOmega, SmartPetscObj< AO > &aoSOmega)
MoFEMErrorCode addStreachSchurMatrix(SmartPetscObj< Mat > &Suu, SmartPetscObj< AO > &aoSuu)
EpElement(MoFEM::Interface &m_field)
MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f)
SmartPetscObj< DM > dmMaterial
Material problem.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Getting interface of core database.
SmartPetscObj< DM > dmElasticSchurBubble
Sub problem of dmElastic Schur.
MoFEMErrorCode postProcessResults(const int tag, const std::string file)
MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset=0)
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
SmartPetscObj< DM > dmElasticSchurStreach
Sub problem of dmElastic Schur.
MoFEMErrorCode addDMs(const BitRefLevel bit=BitRefLevel().set(0))
MoFEMErrorCode setGenericFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > &fe_rhs, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > &fe_lhs)
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
MoFEMErrorCode setElasticElementOps(const int tag)
MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset, boost::shared_ptr< TractionFreeBc > &bc_ptr, const std::string disp_block_set_name, const std::string rot_block_set_name)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode setGenericVolumeElementOps(const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > &fe_rhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > &fe_lhs)
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
SmartPetscObj< DM > dmElasticSchurOmega
Sub problem of dmElastic Schur.
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
SmartPetscObj< DM > dmElasticSchurSpatialDisp
Sub problem of dmElastic Schur.
MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha, const double beta, const double lambda, const double sigma_y)
SmartPetscObj< DM > dM
Coupled problem all fields.
MoFEMErrorCode solveElastic(TS ts, Vec x)
MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset=0)
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > &fe)
SmartPetscObj< DM > dmElastic
Elastic problem.
boost::shared_ptr< EpFEMethod > schurAssembly
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
boost::shared_ptr< PhysicalEquations > physicalEquations
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
boost::shared_ptr< TractionBcVec > bcSpatialTraction
MoFEMErrorCode addFields(const EntityHandle meshset=0)
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
boost::shared_ptr< TractionBcVec > bcData
FeTractionBc(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
virtual MoFEMErrorCode assemble(EntData &data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode integrate(int row_side, EntityType row_type, EntData &data)
virtual MoFEMErrorCode integrate(EntData &data)
virtual MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MatrixDouble K
local tangent matrix
VectorDouble nF
local right hand side vector
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
virtual MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpAssembleBasic(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
OpAssembleBasic(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
OpAssembleFace(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry)
OpAssembleFace(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
OpAssembleVolume(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
OpAssembleVolume(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
OpCalculateRotationAndSpatialGradient(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
OpCalculateStrainEnergy(const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< double > &e)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpDispBc_dx(const std::string &row_field_name, const std::string &col_field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcDispVec > &bc_disp_ptr)
boost::shared_ptr< BcDispVec > bcDispPtr
boost::shared_ptr< BcDispVec > bcDispPtr
OpDispBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcDispVec > &bc_disp_ptr)
MoFEMErrorCode integrate(EntData &data)
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)
virtual MoFEMErrorCode evaluateRhs(EntData &data)=0
virtual MoFEMErrorCode evaluateLhs(EntData &data)=0
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
moab::Interface & postProcMesh
OpPostProcDataStructure(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
std::vector< EntityHandle > & mapGaussPts
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpRotationBc_dx(const std::string &row_field_name, const std::string &col_field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< BcRotVec > bcRotPtr
MoFEMErrorCode integrate(EntData &data)
OpRotationBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
boost::shared_ptr< BcRotVec > bcRotPtr
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialConsistency_dBubble_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialConsistency_dP_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyBubble(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyDivTerm(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyP(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialEquilibrium_dw_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialEquilibrium_dw_dw(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha, const double rho)
OpSpatialEquilibrium(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha, const double rho)
OpSpatialPhysical_du_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPhysical_du_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialPhysical_du_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPhysical_du_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPhysical_du_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialPhysical(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha)
MoFEMErrorCode integrate(EntData &data)
OpSpatialPreconditionMass(const std::string &row_field, MatrixPtr m_ptr)
MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data)
OpSpatialPrj(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialRotation_domega_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialRotation_domega_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialRotation_domega_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &data)
OpSpatialRotation(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MoFEMErrorCode integrate(EntData &row_data)
OpSpatialSchurBegin(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int side, EntityType type, EntData &data)
map< std::string, MatrixDouble > invBlockMat
MoFEMErrorCode integrate(EntData &row_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpSpatialSchurEnd(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double eps)
map< std::string, DataAtIntegrationPts::BlockMatContainor > blockMat
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > DTensor3Ptr
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
FTensor::Tensor2< adouble, 3, 3 > ATensor2
static DTensor3Ptr get_vecTensor3(std::vector< double > &v)
PhysicalEquations(const int size_active, const int size_dependent)
static DTensor2Ptr get_VecTensor2(std::vector< double > &v)
FTensor::Tensor1< adouble, 3 > ATensor1
FTensor::Tensor3< adouble, 3, 3, 3 > ATensor3
virtual MoFEMErrorCode recordTape(const int tag, DTensor2Ptr *t_h)=0
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > DTensor0Ptr
static DTensor0Ptr get_VecTensor0(std::vector< double > &v)
virtual OpJacobian * returnOpJacobian(const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)=0
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
base class for all interface classes