v0.14.0
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
10constexpr int SPACE_DIM = 3;
11constexpr AssemblyType A = AssemblyType::SCHUR;
12
13#ifndef __ESHELBIAN_PLASTICITY_HPP__
14#define __ESHELBIAN_PLASTICITY_HPP__
15
16namespace EshelbianPlasticity {
17
18struct ContactTree;
19
22
23using MatrixPtr = boost::shared_ptr<MatrixDouble>;
24using VectorPtr = boost::shared_ptr<VectorDouble>;
25
26using EntData = EntitiesFieldData::EntData;
27using UserDataOperator = ForcesAndSourcesCore::UserDataOperator;
28using VolUserDataOperator = VolumeElementForcesAndSourcesCore::UserDataOperator;
29using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator;
30
33 : public boost::enable_shared_from_this<DataAtIntegrationPts> {
34
35 MatrixDouble approxPAtPts;
36 MatrixDouble approxSigmaAtPts;
37 MatrixDouble divPAtPts;
38 MatrixDouble divSigmaAtPts;
39 MatrixDouble wL2AtPts;
40 MatrixDouble wL2DotAtPts;
41 MatrixDouble wL2DotDotAtPts;
43 MatrixDouble stretchTensorAtPts;
44
49 MatrixDouble rotAxisAtPts;
50 MatrixDouble rotAxisDotAtPts;
51 MatrixDouble WAtPts;
52 MatrixDouble W0AtPts;
53 MatrixDouble GAtPts;
54 MatrixDouble G0AtPts;
55 MatrixDouble wH1AtPts;
56 MatrixDouble contactL2AtPts;
57 MatrixDouble wGradH1AtPts;
59
60 MatrixDouble hAtPts;
61 MatrixDouble hdOmegaAtPts;
62 MatrixDouble hdLogStretchAtPts;
63 MatrixDouble leviKirchhoffAtPts;
64 MatrixDouble leviKirchhoffPAtPts;
67 MatrixDouble adjointPdUAtPts;
69 MatrixDouble adjointPdUdPAtPts;
70 MatrixDouble rotMatAtPts;
71 MatrixDouble diffRotMatAtPts;
72 MatrixDouble diff2RotMatAtPts;
73 MatrixDouble PAtPts;
74 MatrixDouble SigmaAtPts;
75 VectorDouble phiAtPts;
76 MatrixDouble flowL2AtPts;
77
78 MatrixDouble P_du;
79
80 MatrixDouble eigenVals;
81 MatrixDouble eigenVecs;
82 VectorInt nbUniq;
83
84 MatrixDouble matD;
85 MatrixDouble matAxiatorD;
86 MatrixDouble matDeviatorD;
87
88 double mu;
89 double lambda;
90
92 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
94 }
96 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &approxPAtPts);
97 }
98
100 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divPAtPts);
101 }
102
104 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divSigmaAtPts);
105 }
106
108 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2AtPts);
109 }
110
112 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotAtPts);
113 }
114
116 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotDotAtPts);
117 }
118
120 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
122 }
123
125 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
127 }
128
130 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
132 }
133
135 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxisAtPts);
136 }
137
139 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
141 }
142
144 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
145 }
146
148 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
149 }
150
152 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matD);
153 }
154
156 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matAxiatorD);
157 }
158
160 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matDeviatorD);
161 }
162
164 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wH1AtPts);
165 }
166
168 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &contactL2AtPts);
169 }
170
172 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
173 &wGradH1AtPts);
174 }
175
176 boost::shared_ptr<PhysicalEquations> physicsPtr;
177};
178
179struct OpJacobian;
180
182
189
191
194
196 PhysicalEquations(const int size_active, const int size_dependent)
197 : activeVariables(size_active, 0),
198 dependentVariablesPiola(size_dependent, 0),
199 dependentVariablesPiolaDirevatives(size_dependent * size_active, 0) {}
200 virtual ~PhysicalEquations() = default;
201
202 virtual MoFEMErrorCode recordTape(const int tag, DTensor2Ptr *t_h) = 0;
203
204 virtual OpJacobian *
205 returnOpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs,
206 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
207 boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
208
209 std::vector<double> activeVariables;
210 std::vector<double> dependentVariablesPiola;
212
213 /** \name Active variables */
214
215 /**@{*/
216
217 template <int S>
218 inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
219 return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
220 &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
221 }
222
223 template <int S>
224 inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
225 return DTensor0Ptr(&v[S + 0]);
226 }
227
228 template <int S0>
229 inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v,
230 const int nba) {
231
232 const int A00 = nba * 0 + S0;
233 const int A01 = nba * 1 + S0;
234 const int A02 = nba * 2 + S0;
235 const int A10 = nba * 3 + S0;
236 const int A11 = nba * 4 + S0;
237 const int A12 = nba * 5 + S0;
238 const int A20 = nba * 6 + S0;
239 const int A21 = nba * 7 + S0;
240 const int A22 = nba * 8 + S0;
241
242 return DTensor3Ptr(
243
244 &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
245 &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
246
247 &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
248 &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
249
250 &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
251 &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
252
253 );
254 }
255
256 /**@}*/
257
258 /** \name Active variables */
259
260 /**@{*/
261
262 inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
263
264 /**@}*/
265
266 /** \name Dependent variables */
267
268 /**@{*/
269
271 return get_VecTensor2<0>(dependentVariablesPiola);
272 }
273
274 /**@}*/
275
276 /** \name Derivatives of dependent variables */
277
278 /**@{*/
279
281 return get_vecTensor3<0>(dependentVariablesPiolaDirevatives,
282 activeVariables.size());
283 }
285 return get_vecTensor3<3>(dependentVariablesPiolaDirevatives,
286 activeVariables.size());
287 }
289 return get_vecTensor3<6>(dependentVariablesPiolaDirevatives,
290 activeVariables.size());
291 }
292
293 /**@}*/
294};
295
296struct BcDisp {
297 BcDisp(std::string name, std::vector<double> &attr, Range &faces);
298 std::string blockName;
300 VectorDouble3 vals;
301 VectorInt3 flags;
302};
303typedef std::vector<BcDisp> BcDispVec;
304
305struct BcRot {
306 BcRot(std::string name, std::vector<double> &attr, Range &faces);
307 std::string blockName;
309 VectorDouble3 vals;
310 double theta;
311};
312typedef std::vector<BcRot> BcRotVec;
313
314typedef std::vector<Range> TractionFreeBc;
315
317 TractionBc(std::string name, std::vector<double> &attr, Range &faces);
318 std::string blockName;
320 VectorDouble3 vals;
321 VectorInt3 flags;
322};
323typedef std::vector<TractionBc> TractionBcVec;
324
326
327 OpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs,
328 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
329 boost::shared_ptr<PhysicalEquations> &physics_ptr)
330 : UserDataOperator(NOSPACE, OPSPACE), tAg(tag), evalRhs(eval_rhs),
331 evalLhs(eval_lhs), dataAtPts(data_ptr), physicsPtr(physics_ptr) {}
332
333 virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
334 virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
335
336 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
337
338protected:
340
341 int tAg = -1; ///< adol-c tape
342 bool evalRhs = false;
343 bool evalLhs = false;
344
345 boost::shared_ptr<DataAtIntegrationPts>
346 dataAtPts; ///< data at integration pts
347 boost::shared_ptr<PhysicalEquations>
348 physicsPtr; ///< material physical equations
349};
350
351template <typename T> struct OpAssembleBasic : public T {
352
353 using ScaleOff = boost::function<double()>;
355
356 boost::shared_ptr<DataAtIntegrationPts>
357 dataAtPts; ///< data at integration pts
358
359 OpAssembleBasic(const std::string &field_name,
360 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
361 const char type)
362 : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false) {}
363
365 const std::string &row_field, const std::string &col_field,
366 boost::shared_ptr<DataAtIntegrationPts> data_ptr, const char type,
367 const bool assemble_symmetry, ScaleOff scale_off = []() { return 1; })
368 : T(row_field, col_field, type, false), dataAtPts(data_ptr),
369 assembleSymmetry(assemble_symmetry), scaleOff(scale_off) {}
370
371 VectorDouble nF; ///< local right hand side vector
372 MatrixDouble K; ///< local tangent matrix
373 MatrixDouble transposeK;
374
376
377 virtual MoFEMErrorCode integrate(EntData &data) {
379 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
381 }
382
383 virtual MoFEMErrorCode integrate(int row_side, EntityType row_type,
384 EntData &data) {
386 CHKERR integrate(data);
388 }
389
390 virtual MoFEMErrorCode assemble(EntData &data) {
392 double *vec_ptr = &*nF.begin();
393 int nb_dofs = data.getIndices().size();
394 int *ind_ptr = &*data.getIndices().begin();
395 CHKERR VecSetValues(this->getTSf(), nb_dofs, ind_ptr, vec_ptr, ADD_VALUES);
397 }
398
399 virtual MoFEMErrorCode assemble(int row_side, EntityType row_type,
400 EntData &data) {
402 CHKERR assemble(data);
404 }
405
406 virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
408 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
410 }
411
412 virtual MoFEMErrorCode assemble(int row_side, int col_side,
413 EntityType row_type, EntityType col_type,
414 EntData &row_data, EntData &col_data) {
416
417 if (assembleSymmetry) {
418 const auto row_nb_dofs = row_data.getIndices().size();
419 const auto col_nb_dofs = col_data.getIndices().size();
420 transposeK.resize(col_nb_dofs, row_nb_dofs, false);
421 noalias(transposeK) = trans(K);
422 transposeK *= scaleOff();
423 }
424
425 CHKERR MatSetValues<SchurL2Mats>(this->getTSA(), row_data, col_data, K,
426 ADD_VALUES);
427 if (assembleSymmetry) {
428 CHKERR MatSetValues<SchurL2Mats>(this->getTSA(), col_data, row_data,
429 transposeK, ADD_VALUES);
430 }
431
433 }
434
435 MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
437 if (data.getIndices().empty())
439 nF.resize(data.getIndices().size(), false);
440 nF.clear();
441 CHKERR integrate(side, type, data);
442 CHKERR assemble(side, type, data);
444 }
445
446 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
447 EntityType col_type, EntData &row_data,
448 EntData &col_data) {
450 if (row_data.getIndices().empty())
452 if (col_data.getIndices().empty())
454 K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
455 K.clear();
456 CHKERR integrate(row_data, col_data);
457 CHKERR assemble(row_side, col_side, row_type, col_type, row_data, col_data);
459 }
460};
461
462struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
463 OpAssembleVolume(const std::string &field,
464 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
465 const char type)
466 : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
467
468 OpAssembleVolume(const std::string &row_field, const std::string &col_field,
469 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
470 const char type, const bool assemble_symmetry)
471 : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
472 type, assemble_symmetry) {}
473};
474
475struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
476 OpAssembleFace(const std::string &field,
477 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
478 const char type)
479 : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
480
481 OpAssembleFace(const std::string &row_field, const std::string &col_field,
482 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
483 const char type, const bool assemble_symmetry)
484 : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
485 type, assemble_symmetry) {}
486};
487
489 boost::shared_ptr<DataAtIntegrationPts>
490 dataAtPts; ///< data at integration pts
492 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
493 : VolUserDataOperator(NOSPACE, OPSPACE), dataAtPts(data_ptr) {}
494 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
495};
496
498 const double alphaW;
499 const double alphaRho;
501 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
502 const double alpha, const double rho)
503 : OpAssembleVolume(field_name, data_ptr, OPROW), alphaW(alpha),
504 alphaRho(rho) {}
505 MoFEMErrorCode integrate(EntData &data);
506};
507
509 OpSpatialRotation(const std::string &field_name,
510 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
511 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
512 MoFEMErrorCode integrate(EntData &data);
513};
514
516
517 OpSpatialPhysical(const std::string &field_name,
518 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
519 const double alpha_u)
520 : OpAssembleVolume(field_name, data_ptr, OPROW), alphaU(alpha_u) {
521
522 CHK_MOAB_THROW(PetscOptionsGetBool(PETSC_NULL, "", "-poly_convex",
523 &polyConvex, PETSC_NULL),
524 "get ployconvex option failed");
525 }
526
527 MoFEMErrorCode integrate(EntData &data);
528
529 MoFEMErrorCode integratePiola(EntData &data);
530 MoFEMErrorCode integrateHencky(EntData &data);
531 MoFEMErrorCode integratePolyconvexHencky(EntData &data);
532
533private:
534 const double alphaU;
535 PetscBool polyConvex = PETSC_FALSE;
536};
537
540 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
541 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
542 MoFEMErrorCode integrate(EntData &data);
543
544private:
545};
546
549 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
550 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
551 MoFEMErrorCode integrate(EntData &data);
552
553private:
554};
555
558 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
559 : OpAssembleVolume(field_name, data_ptr, OPROW) {}
560 MoFEMErrorCode integrate(EntData &data);
561};
562
563struct OpDispBc : public OpAssembleFace {
564 boost::shared_ptr<BcDispVec> bcDispPtr;
565 OpDispBc(const std::string &field_name,
566 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
567 boost::shared_ptr<BcDispVec> &bc_disp_ptr,
568 std::vector<boost::shared_ptr<ScalingMethod>> smv)
569 : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr),
570 scalingMethodsVec(smv) {}
571 MoFEMErrorCode integrate(EntData &data);
572
573protected:
574 std::vector<boost::shared_ptr<ScalingMethod>> scalingMethodsVec;
575};
576
578 boost::shared_ptr<BcDispVec> bcDispPtr;
579 OpDispBc_dx(const std::string &row_field_name,
580 const std::string &col_field_name,
581 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
582 boost::shared_ptr<BcDispVec> &bc_disp_ptr)
583 : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
584 false),
585 bcDispPtr(bc_disp_ptr) {}
586 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
587};
588
590 boost::shared_ptr<BcRotVec> bcRotPtr;
591 OpRotationBc(const std::string &field_name,
592 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
593 boost::shared_ptr<BcRotVec> &bc_rot_ptr)
594 : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
595 MoFEMErrorCode integrate(EntData &data);
596};
597
599 boost::shared_ptr<BcRotVec> bcRotPtr;
600 OpRotationBc_dx(const std::string &row_field_name,
601 const std::string &col_field_name,
602 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
603 boost::shared_ptr<BcRotVec> &bc_rot_ptr)
604 : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
605 false),
606 bcRotPtr(bc_rot_ptr) {}
607 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
608};
609
610struct FeTractionBc;
611
613 OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
614 : FaceUserDataOperator(row_field, FaceUserDataOperator::OPROW),
615 bcFe(bc_fe) {}
616 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
617
618protected:
620 MatrixDouble matPP;
621 MatrixDouble vecPv;
622};
623
624struct FeTractionBc : public FEMethod {
626 FeTractionBc(MoFEM::Interface &m_field, const std::string field_name,
627 boost::shared_ptr<TractionBcVec> &bc)
628 : FEMethod(), mField(m_field), bcData(bc), fieldName(field_name) {}
629 MoFEMErrorCode preProcess();
630 MoFEMErrorCode postProcess() { return 0; }
631
632 friend MoFEMErrorCode OpTractionBc::doWork(int side, EntityType type,
633 EntData &data);
634
635protected:
637 boost::shared_ptr<TractionBcVec> bcData;
638 std::string fieldName;
639};
640
642 OpSpatialEquilibrium_dw_dP(const std::string &row_field,
643 const std::string &col_field,
644 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
645 const bool assemble_off = false)
646 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
647 assemble_off) {
648 sYmm = false;
649 }
650 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
651};
652
654 const double alphaW;
655 const double alphaRho;
656 OpSpatialEquilibrium_dw_dw(const std::string &row_field,
657 const std::string &col_field,
658 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
659 const double alpha, const double rho)
660 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
661 alphaW(alpha), alphaRho(rho) {
662 sYmm = true;
663 }
664 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
665};
666
668 const double alphaU;
669 OpSpatialPhysical_du_du(const std::string &row_field,
670 const std::string &col_field,
671 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
672 const double alpha)
673 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
674 alphaU(alpha) {
675 sYmm = false;
676
677 CHK_MOAB_THROW(PetscOptionsGetBool(PETSC_NULL, "", "-poly_convex",
678 &polyConvex, PETSC_NULL),
679 "get ployconvex option failed");
680 }
681 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
682 MoFEMErrorCode integratePiola(EntData &row_data, EntData &col_data);
683 MoFEMErrorCode integrateHencky(EntData &row_data, EntData &col_data);
684 MoFEMErrorCode integratePolyconvexHencky(EntData &row_data,
685 EntData &col_data);
686private:
687 PetscBool polyConvex = PETSC_FALSE;
688
689 MatrixDouble dP;
690};
691
693 OpSpatialPhysical_du_dP(const std::string &row_field,
694 const std::string &col_field,
695 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
696 const bool assemble_off = false)
697 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
698 assemble_off) {
699 sYmm = false;
700 }
701
702 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
703};
704
707 const std::string &row_field, const std::string &col_field,
708 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
709 const bool assemble_off = false)
710 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
711 assemble_off) {
712 sYmm = false;
713 }
714
715 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
716};
717
719 OpSpatialPhysical_du_domega(const std::string &row_field,
720 const std::string &col_field,
721 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
722 const bool assemble_off = false)
723 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
724 assemble_off) {
725 sYmm = false;
726 }
727 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
728};
729
731 OpSpatialPhysical_du_dx(const std::string &row_field,
732 const std::string &col_field,
733 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
734 const bool assemble_off = false)
735 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
736 assemble_off) {
737 sYmm = false;
738 }
739 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
740};
741
743 OpSpatialRotation_domega_dP(const std::string &row_field,
744 const std::string &col_field,
745 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
746 const bool assemble_off)
747 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
748 assemble_off) {
749 sYmm = false;
750 }
751 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
752};
753
756 const std::string &row_field, const std::string &col_field,
757 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
758 const bool assemble_off)
759 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
760 assemble_off) {
761 sYmm = false;
762 }
763 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
764};
765
768 const std::string &row_field, const std::string &col_field,
769 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
770 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
771 sYmm = false;
772 }
773 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
774};
775
778 const std::string &row_field, const std::string &col_field,
779 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
780 const bool assemble_off)
781 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
782 assemble_off) {
783 sYmm = false;
784 }
785 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
786
787private:
788};
789
792 const std::string &row_field, const std::string &col_field,
793 boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
794 const bool assemble_off)
795 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
796 assemble_off) {
797 sYmm = false;
798 }
799 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
800
801private:
802};
803
805
806 moab::Interface &postProcMesh;
807 std::vector<EntityHandle> &mapGaussPts;
808 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
809
810 OpPostProcDataStructure(moab::Interface &post_proc_mesh,
811 std::vector<EntityHandle> &map_gauss_pts,
812 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
814 postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
815 dataAtPts(data_ptr) {}
816
817 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
818};
819
821 OpSpatialPrj(const std::string &row_field,
822 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
823 : OpAssembleVolume(row_field, data_ptr, OPROW) {}
824 MoFEMErrorCode integrate(EntData &row_data);
825};
826
828 OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field,
829 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
830 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
831 // FIXME: That is symmetric
832 sYmm = false;
833 }
834 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
835};
836
838 OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field,
839 boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
840 : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
841 sYmm = false;
842 }
843 MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
844};
845
847
849 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
850 boost::shared_ptr<double> &e)
852 dataAtPts(data_at_pts), energy(e) {}
853
854 MoFEMErrorCode doWork(int side, EntityType type,
855 EntitiesFieldData::EntData &data);
856
857private:
858 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
859 boost::shared_ptr<double> energy;
860};
861
863
867
868 static double exponentBase;
869 static boost::function<double(const double)> f;
870 static boost::function<double(const double)> d_f;
871 static boost::function<double(const double)> dd_f;
872 static boost::function<double(const double)> inv_f;
873 static boost::function<double(const double)> inv_d_f;
874 static boost::function<double(const double)> inv_dd_f;
875
876 static double f_log(const double v) {
877 return pow(EshelbianCore::exponentBase, v);
878 }
879 static double d_f_log(const double v) {
880 return pow(exponentBase, v) * log(EshelbianCore::exponentBase);
881 }
882 static double dd_f_log(const double v) {
883 return pow(EshelbianCore::exponentBase, v) *
885 }
886
887 static double inv_f_log(const double v) {
888 return log(v) / log(EshelbianCore::exponentBase);
889 }
890 static double inv_d_f_log(const double v) {
891 return (1. / v) / log(EshelbianCore::exponentBase);
892 }
893 static double inv_dd_f_log(const double v) {
894 return -(1. / (v * v)) / log(EshelbianCore::exponentBase);
895 }
896
897 static double f_linear(const double v) { return v + 1; }
898 static double d_f_linear(const double v) { return 1; }
899 static double dd_f_linear(const double v) { return 0; }
900
901 static double inv_f_linear(const double v) { return v - 1; }
902 static double inv_d_f_linear(const double v) { return 0; }
903 static double inv_dd_f_linear(const double v) { return 0; }
904
905 /**
906 * \brief Getting interface of core database
907 * @param uuid unique ID of interface
908 * @param iface returned pointer to interface
909 * @return error code
910 */
911 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
912 UnknownInterface **iface) const;
913
915
916 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
917 boost::shared_ptr<PhysicalEquations> physicalEquations;
918
919 boost::shared_ptr<VolumeElementForcesAndSourcesCore> elasticFeRhs;
920 boost::shared_ptr<VolumeElementForcesAndSourcesCore> elasticFeLhs;
921 boost::shared_ptr<FaceElementForcesAndSourcesCore> elasticBcLhs;
922 boost::shared_ptr<FaceElementForcesAndSourcesCore> elasticBcRhs;
923 boost::shared_ptr<FaceElementForcesAndSourcesCore> contactRhs;
924 boost::shared_ptr<ContactTree> contactTreeRhs; ///< Make a contact tree
925
926 SmartPetscObj<DM> dM; ///< Coupled problem all fields
927 SmartPetscObj<DM> dmElastic; ///< Elastic problem
928 SmartPetscObj<DM> dmPrjSpatial; ///< Projection spatial displacement
929
930 const std::string piolaStress;
931 const std::string eshelbyStress;
932 const std::string spatialL2Disp;
933 const std::string spatialH1Disp;
934 const std::string contactDisp;
935 const std::string materialL2Disp;
936 const std::string stretchTensor;
937 const std::string rotAxis;
938 const std::string materialGradient;
939 const std::string tauField;
940 const std::string lambdaField;
941 const std::string bubbleField;
942
943 const std::string elementVolumeName;
944 const std::string naturalBcElement;
945 const std::string essentialBcElement;
946 const std::string skinElement;
947 const std::string contactElement;
948
950 virtual ~EshelbianCore();
951
955 double alphaU;
956 double alphaW;
957 double alphaRho;
958 double precEps;
960 double precEpsW;
962
963 MoFEMErrorCode getOptions();
964
965 boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
966 boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
967 boost::shared_ptr<TractionBcVec> bcSpatialTraction;
968 boost::shared_ptr<TractionFreeBc> bcSpatialFreeTraction;
969
970 template <typename BC>
971 MoFEMErrorCode getBc(boost::shared_ptr<BC> &bc_vec_ptr,
972 const std::string block_name, const int nb_attributes) {
974 for (auto it :
975 mField.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
976
977 (boost::format("%s(.*)") % block_name).str()
978
979 ))
980
981 ) {
982 std::vector<double> block_attributes;
983 CHKERR it->getAttributes(block_attributes);
984 if (block_attributes.size() != nb_attributes) {
985 SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
986 "In block %s expected %d attributes, but given %d",
987 it->getName().c_str(), nb_attributes, block_attributes.size());
988 }
989 Range faces;
990 CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
991 true);
992 bc_vec_ptr->emplace_back(it->getName(), block_attributes, faces);
993 }
995 }
996
997 MoFEMErrorCode getSpatialDispBc();
998
999 inline MoFEMErrorCode getSpatialRotationBc() {
1000 bcSpatialRotationVecPtr = boost::make_shared<BcRotVec>();
1001 return getBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC", 4);
1002 }
1003
1004 MoFEMErrorCode getSpatialTractionBc();
1005
1006 /**
1007 * @brief Remove all, but entities where kinematic constrains are applied.
1008 *
1009 * @param meshset
1010 * @param bc_ptr
1011 * @param disp_block_set_name
1012 * @param rot_block_set_name
1013 * @param contact_set_name
1014 * @return MoFEMErrorCode
1015 */
1016 MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset,
1017 boost::shared_ptr<TractionFreeBc> &bc_ptr,
1018 const std::string contact_set_name);
1019
1020 inline MoFEMErrorCode
1023 boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
1024 return getTractionFreeBc(meshset, bcSpatialFreeTraction, "CONTACT");
1025 }
1026
1027 MoFEMErrorCode addFields(const EntityHandle meshset = 0);
1028 MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset = 0);
1029 MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset = 0);
1030 MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
1031
1032 MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape,
1033 const double lambda,
1034 const double mu,
1035 const double sigma_y);
1036
1037 MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha,
1038 const double beta,
1039 const double lambda,
1040 const double sigma_y);
1041
1042 MoFEMErrorCode addMaterial_Hencky(double E, double nu);
1043
1044 MoFEMErrorCode setBaseVolumeElementOps(
1045 const int tag, const bool do_rhs, const bool do_lhs,
1046 boost::shared_ptr<VolumeElementForcesAndSourcesCore> &fe);
1047
1048 MoFEMErrorCode setVolumeElementOps(
1049 const int tag, const bool add_elastic, const bool add_material,
1050 boost::shared_ptr<VolumeElementForcesAndSourcesCore> &fe_rhs,
1051 boost::shared_ptr<VolumeElementForcesAndSourcesCore> &fe_lhs);
1052
1053 MoFEMErrorCode
1054 setFaceElementOps(const bool add_elastic, const bool add_material,
1055 boost::shared_ptr<FaceElementForcesAndSourcesCore> &fe_rhs,
1056 boost::shared_ptr<FaceElementForcesAndSourcesCore> &fe_lhs);
1057
1058 MoFEMErrorCode setContactElementOps(
1059
1060 boost::shared_ptr<ContactTree> &fe_contact_tree,
1061 boost::shared_ptr<FaceElementForcesAndSourcesCore> &fe_rhs,
1062 boost::shared_ptr<FaceElementForcesAndSourcesCore> &fe_lhs
1063
1064 );
1065
1066 MoFEMErrorCode setElasticElementOps(const int tag);
1067 MoFEMErrorCode setElasticElementToTs(DM dm);
1068
1069 MoFEMErrorCode solveElastic(TS ts, Mat m, Vec f, Vec x);
1070
1071 MoFEMErrorCode postProcessResults(const int tag, const std::string file);
1072 MoFEMErrorCode gettingNorms();
1073
1074 struct SetUpSchur {
1075 static boost::shared_ptr<SetUpSchur> createSetUpSchur(
1076
1077 MoFEM::Interface &m_field, SmartPetscObj<Mat> m, SmartPetscObj<Mat> p,
1078 EshelbianCore *ep_core_ptr
1079
1080 );
1081 virtual MoFEMErrorCode setUp(KSP solver) = 0;
1082
1083 virtual MoFEMErrorCode preProc() = 0;
1084 virtual MoFEMErrorCode postProc() = 0;
1085
1086 protected:
1087 SetUpSchur() = default;
1088 };
1089};
1090
1091} // namespace EshelbianPlasticity
1092
1093#endif //__ESHELBIAN_PLASTICITY_HPP__
static Index< 'p', 3 > p
constexpr int SPACE_DIM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ NOSPACE
Definition definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
@ 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()
#define CHKERR
Inline error check.
FTensor::Index< 'm', SPACE_DIM > m
auto bit
set bit
static double lambda
const double v
phase velocity of light in medium (cm/ns)
std::vector< BcDisp > BcDispVec
boost::shared_ptr< VectorDouble > VectorPtr
FaceElementForcesAndSourcesCore::UserDataOperator FaceUserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
std::vector< Range > TractionFreeBc
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
std::vector< BcRot > BcRotVec
enum RotSelector EshelbianCore
std::vector< TractionBc > TractionBcVec
boost::shared_ptr< MatrixDouble > MatrixPtr
double rho
Definition plastic.cpp:191
constexpr auto field_name
BcDisp(std::string name, std::vector< double > &attr, Range &faces)
BcRot(std::string name, std::vector< double > &attr, Range &faces)
boost::shared_ptr< PhysicalEquations > physicsPtr
virtual MoFEMErrorCode setUp(KSP solver)=0
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field, SmartPetscObj< Mat > m, SmartPetscObj< Mat > p, EshelbianCore *ep_core_ptr)
static double inv_d_f_linear(const double v)
MoFEMErrorCode getSpatialDispBc()
[Getting norms]
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Getting interface of core database.
static boost::function< double(const double)> dd_f
boost::shared_ptr< FaceElementForcesAndSourcesCore > elasticBcLhs
boost::shared_ptr< VolumeElementForcesAndSourcesCore > elasticFeLhs
MoFEMErrorCode addMaterial_Hencky(double E, double nu)
static boost::function< double(const double)> inv_f
static boost::function< double(const double)> inv_dd_f
MoFEMErrorCode postProcessResults(const int tag, const std::string file)
MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset=0)
static boost::function< double(const double)> inv_d_f
MoFEMErrorCode addDMs(const BitRefLevel bit=BitRefLevel().set(0))
static double inv_f_log(const double v)
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
boost::shared_ptr< VolumeElementForcesAndSourcesCore > elasticFeRhs
MoFEMErrorCode setElasticElementOps(const int tag)
static double inv_d_f_log(const double v)
static double d_f_linear(const double v)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
static double dd_f_log(const double v)
MoFEMErrorCode setContactElementOps(boost::shared_ptr< ContactTree > &fe_contact_tree, boost::shared_ptr< FaceElementForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< FaceElementForcesAndSourcesCore > &fe_lhs)
MoFEMErrorCode solveElastic(TS ts, Mat m, Vec f, Vec x)
static enum StretchSelector stretchSelector
static double inv_dd_f_log(const double v)
static boost::function< double(const double)> d_f
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_name, const int nb_attributes)
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< VolumeElementForcesAndSourcesCore > &fe)
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
MoFEMErrorCode setFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< FaceElementForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< FaceElementForcesAndSourcesCore > &fe_lhs)
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.
static double d_f_log(const double v)
MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset=0)
SmartPetscObj< DM > dmElastic
Elastic problem.
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
boost::shared_ptr< PhysicalEquations > physicalEquations
static double inv_dd_f_linear(const double v)
static boost::function< double(const double)> f
MoFEMErrorCode setVolumeElementOps(const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< VolumeElementForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< VolumeElementForcesAndSourcesCore > &fe_lhs)
boost::shared_ptr< FaceElementForcesAndSourcesCore > contactRhs
static double dd_f_linear(const double v)
MoFEMErrorCode gettingNorms()
[Getting norms]
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
static double f_linear(const double v)
boost::shared_ptr< TractionBcVec > bcSpatialTraction
MoFEMErrorCode addFields(const EntityHandle meshset=0)
static double f_log(const double v)
boost::shared_ptr< ContactTree > contactTreeRhs
Make a contact tree.
boost::shared_ptr< FaceElementForcesAndSourcesCore > elasticBcRhs
SmartPetscObj< DM > dmPrjSpatial
Projection spatial displacement.
MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset, boost::shared_ptr< TractionFreeBc > &bc_ptr, const std::string contact_set_name)
Remove all, but entities where kinematic constrains are applied.
static double inv_f_linear(const double v)
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)
OpAssembleBasic(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry, ScaleOff scale_off=[]() { return 1;})
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)
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(boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
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, std::vector< boost::shared_ptr< ScalingMethod > > smv)
MoFEMErrorCode integrate(EntData &data)
std::vector< boost::shared_ptr< ScalingMethod > > scalingMethodsVec
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
virtual MoFEMErrorCode evaluateRhs(EntData &data)=0
OpJacobian(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 evaluateLhs(EntData &data)=0
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
moab::Interface & postProcMesh
std::vector< EntityHandle > & mapGaussPts
OpPostProcDataStructure(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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
OpSpatialConsistency_dBubble_domega(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)
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, const bool assemble_off)
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)
MoFEMErrorCode integrateHencky(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 integratePiola(EntData &row_data, EntData &col_data)
MoFEMErrorCode integratePolyconvexHencky(EntData &row_data, EntData &col_data)
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)
MoFEMErrorCode integratePiola(EntData &data)
MoFEMErrorCode integratePolyconvexHencky(EntData &data)
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode integrateHencky(EntData &data)
OpSpatialPhysical(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha_u)
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)
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 doWork(int side, EntityType type, EntData &data)
OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
FTensor::Tensor3< double, 3, 3, 3 > DTensor3
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > DTensor3Ptr
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
static DTensor3Ptr get_vecTensor3(std::vector< double > &v, const int nba)
FTensor::Tensor2< adouble, 3, 3 > ATensor2
FTensor::Tensor2< double, 3, 3 > DTensor2
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 int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)=0
TractionBc(std::string name, std::vector< double > &attr, Range &faces)
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.