v0.15.0
Loading...
Searching...
No Matches
ConvectiveMassElement.hpp
Go to the documentation of this file.
1/** \file ConvectiveMassElement.hpp
2 * \brief Operators and data structures for mass and convective mass element
3 * \ingroup convective_mass_elem
4 *
5 */
6
7
8
9#ifndef __CONVECTIVE_MASS_ELEMENT_HPP
10#define __CONVECTIVE_MASS_ELEMENT_HPP
11
12#ifndef WITH_ADOL_C
13#error "MoFEM need to be compiled with ADOL-C"
14#endif
15
16/** \brief structure grouping operators and data used for calculation of mass
17 * (convective) element \ingroup convective_mass_elem \ingroup
18 * nonlinear_elastic_elem
19 *
20 * In order to assemble matrices and right hand vectors, the loops over
21 * elements, entities over that elements and finally loop over integration
22 * points are executed.
23 *
24 * Following implementation separate those three celeries of loops and to each
25 * loop attach operator.
26 *
27 */
29
30 /// \brief definition of volume element
32
33 Mat A;
34 Vec F;
35
37
38 /** \brief it is used to calculate nb. of Gauss integration points
39 *
40 * for more details pleas look
41 * Reference:
42 *
43 * Albert Nijenhuis, Herbert Wilf,
44 * Combinatorial Algorithms for Computers and Calculators,
45 * Second Edition,
46 * Academic Press, 1978,
47 * ISBN: 0-12-519260-6,
48 * LC: QA164.N54.
49 *
50 * More details about algorithm
51 * http://people.sc.fsu.edu/~jburkardt/cpp_src/gm_rule/gm_rule.html
52 **/
53 int getRule(int order);
54
56 double eNergy;
57
60 };
61
62 MyVolumeFE feMassRhs; ///< calculate right hand side for tetrahedral elements
64 return feMassRhs;
65 } ///< get rhs volume element
66 MyVolumeFE feMassLhs; ///< calculate left hand side for tetrahedral
67 ///< elements,i.e. mass element
69 return feMassLhs;
70 } ///< get lhs volume element
71 MyVolumeFE feMassAuxLhs; ///< calculate left hand side for tetrahedral
72 ///< elements for Kuu shell matrix
74 return feMassAuxLhs;
75 } ///< get lhs volume element for Kuu shell matrix
76
77 MyVolumeFE feVelRhs; ///< calculate right hand side for tetrahedral elements
78 MyVolumeFE &getLoopFeVelRhs() { return feVelRhs; } ///< get rhs volume element
79 MyVolumeFE feVelLhs; ///< calculate left hand side for tetrahedral elements
80 MyVolumeFE &getLoopFeVelLhs() { return feVelLhs; } ///< get lhs volume element
81
82 MyVolumeFE feTRhs; ///< calculate right hand side for tetrahedral elements
83 MyVolumeFE &getLoopFeTRhs() { return feTRhs; } ///< get rhs volume element
84 MyVolumeFE feTLhs; ///< calculate left hand side for tetrahedral elements
85 MyVolumeFE &getLoopFeTLhs() { return feTLhs; } ///< get lhs volume element
86
87 MyVolumeFE feEnergy; ///< calculate kinetic energy
89 return feEnergy;
90 } ///< get kinetic energy element
91
94 auto add_ops = [&](auto &fe) {
95 return MoFEM::addHOOpsVol("MESH_NODE_POSITIONS", fe, true, false, false,
96 false);
97 };
98 CHKERR add_ops(feMassRhs);
99 CHKERR add_ops(feMassLhs);
100 CHKERR add_ops(feMassAuxLhs);
101 CHKERR add_ops(feVelRhs);
102 CHKERR add_ops(feTRhs);
103 CHKERR add_ops(feTLhs);
104 CHKERR add_ops(feEnergy);
106 }
107
109 short int tAg;
110
111 ConvectiveMassElement(MoFEM::Interface &m_field, short int tag);
112
113 /** \brief data for calculation inertia forces
114 * \ingroup user_modules
115 */
116 struct BlockData {
117 double rho0; ///< reference density
118 VectorDouble a0; ///< constant acceleration
119 Range tEts; ///< elements in block set
120 };
121 std::map<int, BlockData>
122 setOfBlocks; ///< maps block set id with appropriate BlockData
123
124 /** \brief common data used by volume elements
125 * \ingroup user_modules
126 */
127 struct CommonData {
128
129 bool lInear;
131 CommonData() : lInear(false), staticOnly(false) {}
132
133 std::map<std::string, std::vector<VectorDouble>> dataAtGaussPts;
134 std::map<std::string, std::vector<MatrixDouble>> gradAtGaussPts;
138 std::vector<VectorDouble> valVel;
139 std::vector<std::vector<double *>> jacVelRowPtr;
140 std::vector<MatrixDouble> jacVel;
141 std::vector<VectorDouble> valMass;
142 std::vector<std::vector<double *>> jacMassRowPtr;
143 std::vector<MatrixDouble> jacMass;
144 std::vector<VectorDouble> valT;
145 std::vector<std::vector<double *>> jacTRowPtr;
146 std::vector<MatrixDouble> jacT;
147 };
149
150 boost::ptr_vector<MethodForForceScaling> methodsOp;
151
153 : public VolumeElementForcesAndSourcesCore::UserDataOperator {
154
155 std::vector<VectorDouble> &valuesAtGaussPts;
156 std::vector<MatrixDouble> &gradientAtGaussPts;
157 const EntityType zeroAtType;
158
159 OpGetDataAtGaussPts(const std::string field_name,
160 std::vector<VectorDouble> &values_at_gauss_pts,
161 std::vector<MatrixDouble> &gardient_at_gauss_pts);
162
163 /** \brief operator calculating deformation gradient
164 *
165 */
166 MoFEMErrorCode doWork(int side, EntityType type,
167 EntitiesFieldData::EntData &data);
168 };
169
171 OpGetCommonDataAtGaussPts(const std::string field_name,
172 CommonData &common_data);
173 };
174
176
178 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
180
183 int tAg;
185 bool &lInear;
187
188 boost::ptr_vector<MethodForForceScaling> &methodsOp;
189
190 OpMassJacobian(const std::string field_name, BlockData &data,
191 CommonData &common_data,
192 boost::ptr_vector<MethodForForceScaling> &methods_op,
193 int tag, bool linear = false);
194
198
201 std::vector<double> active;
202
203 MoFEMErrorCode doWork(int row_side, EntityType row_type,
204 EntitiesFieldData::EntData &row_data);
205 };
206
207 struct OpMassRhs : public VolumeElementForcesAndSourcesCore::UserDataOperator,
209
212
213 OpMassRhs(const std::string field_name, BlockData &data,
214 CommonData &common_data);
215
217
218 MoFEMErrorCode doWork(int row_side, EntityType row_type,
219 EntitiesFieldData::EntData &row_data);
220 };
221
223 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
225
229
230 OpMassLhs_dM_dv(const std::string vel_field, const std::string field_name,
231 BlockData &data, CommonData &common_data,
232 Range *forcesonlyonentities_ptr = NULL);
233
235
236 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
237 int gg);
238
239 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
240 EntityType col_type,
241 EntitiesFieldData::EntData &row_data,
242 EntitiesFieldData::EntData &col_data);
243 };
244
246
247 OpMassLhs_dM_dx(const std::string field_name, const std::string col_field,
248 BlockData &data, CommonData &common_data);
249
250 MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg);
251 };
252
254
255 OpMassLhs_dM_dX(const std::string field_name, const std::string col_field,
256 BlockData &data, CommonData &common_data);
257
258 MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg);
259 };
260
261 struct OpEnergy : public VolumeElementForcesAndSourcesCore::UserDataOperator,
263
267 bool &lInear;
268
269 OpEnergy(const std::string field_name, BlockData &data,
270 CommonData &common_data, SmartPetscObj<Vec> v);
271
274
275 MoFEMErrorCode doWork(int row_side, EntityType row_type,
276 EntitiesFieldData::EntData &row_data);
277 };
278
280 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
282
285 int tAg;
287
288 OpVelocityJacobian(const std::string field_name, BlockData &data,
289 CommonData &common_data, int tag, bool jacobian = true);
290
294
295 std::vector<double> active;
296
297 MoFEMErrorCode doWork(int row_side, EntityType row_type,
298 EntitiesFieldData::EntData &row_data);
299 };
300
302 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
304
307
308 OpVelocityRhs(const std::string field_name, BlockData &data,
309 CommonData &common_data);
310
312
313 MoFEMErrorCode doWork(int row_side, EntityType row_type,
314 EntitiesFieldData::EntData &row_data);
315 };
316
318
319 OpVelocityLhs_dV_dv(const std::string vel_field,
320 const std::string field_name, BlockData &data,
321 CommonData &common_data);
322
323 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
324 int gg);
325 };
326
328
329 OpVelocityLhs_dV_dx(const std::string vel_field,
330 const std::string field_name, BlockData &data,
331 CommonData &common_data);
332
333 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
334 int gg);
335 };
336
338
339 OpVelocityLhs_dV_dX(const std::string vel_field,
340 const std::string field_name, BlockData &data,
341 CommonData &common_data);
342
343 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
344 int gg);
345 };
346
348 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
350
353 int tAg;
356
358 BlockData &data,
359 CommonData &common_data, int tag,
360 bool jacobian = true);
361
365
366 MoFEMErrorCode doWork(int row_side, EntityType row_type,
367 EntitiesFieldData::EntData &row_data);
368 };
369
371 : public VolumeElementForcesAndSourcesCore::UserDataOperator,
373
377
379 BlockData &data,
380 CommonData &common_data,
381 Range *forcesonlyonentities_ptr);
382
384
385 MoFEMErrorCode doWork(int row_side, EntityType row_type,
386 EntitiesFieldData::EntData &row_data);
387 };
388
390
391 OpEshelbyDynamicMaterialMomentumLhs_dv(const std::string vel_field,
392 const std::string field_name,
393 BlockData &data,
394 CommonData &common_data,
395 Range *forcesonlyonentities_ptr);
396
397 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
398 int gg);
399 };
400
403
404 OpEshelbyDynamicMaterialMomentumLhs_dx(const std::string vel_field,
405 const std::string field_name,
406 BlockData &data,
407 CommonData &common_data,
408 Range *forcesonlyonentities_ptr);
409
410 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
411 int gg);
412 };
413
416
417 OpEshelbyDynamicMaterialMomentumLhs_dX(const std::string vel_field,
418 const std::string field_name,
419 BlockData &data,
420 CommonData &common_data,
421 Range *forcesonlyonentities_ptr);
422
423 virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data,
424 int gg);
425 };
426
427 /**
428 * @brief Set fields DOT_
429 *
430 * \note This is old solution, to keep rates calculate by TS as a fields. This
431 * is not memory efficient solution.
432 *
433 */
434 struct UpdateAndControl : public FEMethod {
435
437 TS tS;
438 const std::string velocityField;
439 const std::string spatialPositionField;
440
442 UpdateAndControl(MoFEM::Interface &m_field, TS _ts,
443 const std::string velocity_field,
444 const std::string spatial_position_field);
445
446 /**
447 * @brief Scatter values from t_u_dt on the fields
448 *
449 * @return MoFEMErrorCode
450 */
452
453 // This is empty fun cions does nothing
455 };
456
458
459 static MoFEMErrorCode
461 boost::shared_ptr<map<int, BlockData>> &block_sets_ptr);
462
464 string element_name, string velocity_field_name,
465 string spatial_position_field_name,
466 string material_position_field_name = "MESH_NODE_POSITIONS",
467 bool ale = false, BitRefLevel bit = BitRefLevel());
468
470 string element_name, string velocity_field_name,
471 string spatial_position_field_name,
472 string material_position_field_name = "MESH_NODE_POSITIONS",
473 bool ale = false, BitRefLevel bit = BitRefLevel());
474
476 string element_name, string velocity_field_name,
477 string spatial_position_field_name,
478 string material_position_field_name = "MESH_NODE_POSITIONS",
479 bool ale = false, BitRefLevel bit = BitRefLevel(),
480 Range *intersected = NULL);
481
483 string velocity_field_name, string spatial_position_field_name,
484 string material_position_field_name = "MESH_NODE_POSITIONS",
485 bool ale = false, bool linear = false);
486
488 string velocity_field_name, string spatial_position_field_name,
489 string material_position_field_name = "MESH_NODE_POSITIONS",
490 bool ale = false);
491
493 string velocity_field_name, string spatial_position_field_name,
494 string material_position_field_name = "MESH_NODE_POSITIONS",
495 Range *forces_on_entities_ptr = NULL);
496
498 string velocity_field_name, string spatial_position_field_name,
499 string material_position_field_name = "MESH_NODE_POSITIONS",
500 bool linear = false);
501
502 struct MatShellCtx {
503
504 Mat K, M;
505 VecScatter scatterU, scatterV;
506 double ts_a; //,scale;
507
509 MatShellCtx();
510 virtual ~MatShellCtx();
511
512 Mat barK;
513 Vec u, v, Ku, Mv;
515
517
518 friend MoFEMErrorCode MultOpA(Mat A, Vec x, Vec f);
520 };
521
522 /** \brief Mult operator for shell matrix
523 *
524 * \f[
525 \left[
526 \begin{array}{cc}
527 \mathbf{M} & \mathbf{K} \\
528 \mathbf{I} & -\mathbf{I}a
529 \end{array}
530 \right]
531 \left[
532 \begin{array}{c}
533 \mathbf{v} \\
534 \mathbf{u}
535 \end{array}
536 \right] =
537 \left[
538 \begin{array}{c}
539 \mathbf{r}_u \\
540 \mathbf{r}_v
541 \end{array}
542 \right]
543 * \f]
544 *
545 */
546 static MoFEMErrorCode MultOpA(Mat A, Vec x, Vec f) {
548
549 void *void_ctx;
550 CHKERR MatShellGetContext(A, &void_ctx);
551 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
552 if (!ctx->iNitialized) {
553 CHKERR ctx->iNit();
554 }
555 CHKERR VecZeroEntries(f);
556 // Mult Ku
557 CHKERR VecScatterBegin(ctx->scatterU, x, ctx->u, INSERT_VALUES,
558 SCATTER_FORWARD);
559 CHKERR VecScatterEnd(ctx->scatterU, x, ctx->u, INSERT_VALUES,
560 SCATTER_FORWARD);
561 CHKERR MatMult(ctx->K, ctx->u, ctx->Ku);
562 CHKERR VecScatterBegin(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
563 SCATTER_REVERSE);
564 CHKERR VecScatterEnd(ctx->scatterU, ctx->Ku, f, INSERT_VALUES,
565 SCATTER_REVERSE);
566 // Mult Mv
567 CHKERR VecScatterBegin(ctx->scatterV, x, ctx->v, INSERT_VALUES,
568 SCATTER_FORWARD);
569 CHKERR VecScatterEnd(ctx->scatterV, x, ctx->v, INSERT_VALUES,
570 SCATTER_FORWARD);
571 CHKERR MatMult(ctx->M, ctx->v, ctx->Mv);
572 CHKERR VecScatterBegin(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
573 SCATTER_REVERSE);
574 CHKERR VecScatterEnd(ctx->scatterU, ctx->Mv, f, ADD_VALUES,
575 SCATTER_REVERSE);
576 // Velocities
577 CHKERR VecAXPY(ctx->v, -ctx->ts_a, ctx->u);
578 // CHKERR VecScale(ctx->v,ctx->scale);
579 CHKERR VecScatterBegin(ctx->scatterV, ctx->v, f, INSERT_VALUES,
580 SCATTER_REVERSE);
581 CHKERR VecScatterEnd(ctx->scatterV, ctx->v, f, INSERT_VALUES,
582 SCATTER_REVERSE);
583 // Assemble
584 CHKERR VecAssemblyBegin(f);
585 CHKERR VecAssemblyEnd(f);
587 }
588
591
592 void *void_ctx;
593 CHKERR MatShellGetContext(A, &void_ctx);
594 MatShellCtx *ctx = (MatShellCtx *)void_ctx;
595 CHKERR MatZeroEntries(ctx->K);
596 CHKERR MatZeroEntries(ctx->M);
598 }
599
600 struct PCShellCtx {
601
603 bool initPC; ///< check if PC is initialized
604
605 PCShellCtx(Mat shell_mat) : shellMat(shell_mat), initPC(false) {}
606
607 PC pC;
608
610
612
615 friend MoFEMErrorCode PCShellApplyOp(PC pc, Vec f, Vec x);
616 };
617
620
621 void *void_ctx;
622 CHKERR PCShellGetContext(pc, &void_ctx);
623 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
624 CHKERR ctx->iNit();
625 MatShellCtx *shell_mat_ctx;
626 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
627 CHKERR PCSetFromOptions(ctx->pC);
628 CHKERR PCSetOperators(ctx->pC, shell_mat_ctx->barK, shell_mat_ctx->barK);
629 CHKERR PCSetUp(ctx->pC);
631 }
632
635
636 void *void_ctx;
637 CHKERR PCShellGetContext(pc, &void_ctx);
638 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
639 CHKERR ctx->dEstroy();
641 }
642
643 /** \brief apply pre-conditioner for shell matrix
644 *
645 * \f[
646 \left[
647 \begin{array}{cc}
648 \mathbf{M} & \mathbf{K} \\
649 \mathbf{I} & -\mathbf{I}a
650 \end{array}
651 \right]
652 \left[
653 \begin{array}{c}
654 \mathbf{v} \\
655 \mathbf{u}
656 \end{array}
657 \right] =
658 \left[
659 \begin{array}{c}
660 \mathbf{r}_u \\
661 \mathbf{r}_v
662 \end{array}
663 \right]
664 * \f]
665 *
666 * where \f$\mathbf{v} = \mathbf{r}_v + a\mathbf{u}\f$ and
667 \f$\mathbf{u}=(a\mathbf{M}+\mathbf{K})^{-1}(\mathbf{r}_u -
668 \mathbf{M}\mathbf{r}_v\f$.
669 *
670 */
671 static MoFEMErrorCode PCShellApplyOp(PC pc, Vec f, Vec x) {
673
674 void *void_ctx;
675 CHKERR PCShellGetContext(pc, &void_ctx);
676 PCShellCtx *ctx = (PCShellCtx *)void_ctx;
677 MatShellCtx *shell_mat_ctx;
678 CHKERR MatShellGetContext(ctx->shellMat, &shell_mat_ctx);
679 // forward
680 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
681 INSERT_VALUES, SCATTER_FORWARD);
682 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, f, shell_mat_ctx->Ku,
683 INSERT_VALUES, SCATTER_FORWARD);
684 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
685 INSERT_VALUES, SCATTER_FORWARD);
686 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, f, shell_mat_ctx->v,
687 INSERT_VALUES, SCATTER_FORWARD);
688 // CHKERR VecScale(shell_mat_ctx->v,1/shell_mat_ctx->scale);
689 // apply pre-conditioner and calculate u
690 CHKERR MatMult(shell_mat_ctx->M, shell_mat_ctx->v,
691 shell_mat_ctx->Mv); // Mrv
692 CHKERR VecAXPY(shell_mat_ctx->Ku, -1, shell_mat_ctx->Mv); // f-Mrv
693 CHKERR PCApply(ctx->pC, shell_mat_ctx->Ku,
694 shell_mat_ctx->u); // u = (aM+K)^(-1)(ru-Mrv)
695 // VecView(shell_mat_ctx->u,PETSC_VIEWER_STDOUT_WORLD);
696 // calculate velocities
697 CHKERR VecAXPY(shell_mat_ctx->v, shell_mat_ctx->ts_a,
698 shell_mat_ctx->u); // v = v + a*u
699 // VecView(shell_mat_ctx->v,PETSC_VIEWER_STDOUT_WORLD);
700 // reverse
701 CHKERR VecZeroEntries(x);
702 CHKERR VecScatterBegin(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
703 INSERT_VALUES, SCATTER_REVERSE);
704 CHKERR VecScatterEnd(shell_mat_ctx->scatterU, shell_mat_ctx->u, x,
705 INSERT_VALUES, SCATTER_REVERSE);
706 CHKERR VecScatterBegin(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
707 INSERT_VALUES, SCATTER_REVERSE);
708 CHKERR VecScatterEnd(shell_mat_ctx->scatterV, shell_mat_ctx->v, x,
709 INSERT_VALUES, SCATTER_REVERSE);
710 CHKERR VecAssemblyBegin(x);
711 CHKERR VecAssemblyEnd(x);
713 }
714
718
719 // variables bellow need to be set by user
720 MatShellCtx *shellMatCtx; ///< pointer to shell matrix
721
722 /**
723 * @brief Calculate inconsistency between approximation of velocities and
724 * velocities calculated from displacements
725 *
726 * @return MoFEMErrorCode
727 */
729
730 };
731
732#ifdef __DIRICHLET_HPP__
733
734 /** \brief blocked element/problem
735 *
736 * Blocked element run loops for different problem than TS problem. It is
737 * used to calculate matrices of shell matrix.
738 *
739 */
740 struct ShellMatrixElement : public FEMethod {
741
742 MoFEM::Interface &mField;
743 ShellMatrixElement(MoFEM::Interface &m_field);
744
745 typedef std::pair<std::string, FEMethod *> PairNameFEMethodPtr;
746 typedef std::vector<PairNameFEMethodPtr> LoopsToDoType;
747 LoopsToDoType loopK; ///< methods to calculate K shell matrix
748 LoopsToDoType loopM; ///< methods to calculate M shell matrix
749 LoopsToDoType loopAuxM; ///< methods to calculate derivatives of inertia
750 ///< forces over displacements shell matrix
751
752 // variables bellow need to be set by user
753 string problemName; ///< name of shell problem
754 MatShellCtx *shellMatCtx; ///< pointer to shell matrix
755 DirichletDisplacementBc *DirichletBcPtr; ///< boundary conditions
756
757 MoFEMErrorCode preProcess();
758 };
759
760#endif //__DIRICHLET_HPP__
761};
762
763#endif //__CONVECTIVE_MASS_ELEMENT_HPP
764
765/**
766* \defgroup convective_mass_elem Mass Element
767* \ingroup user_modules
768**/
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
constexpr int order
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
ublas::vector< T, ublas::bounded_array< T, N > > VectorBoundedArray
Definition Types.hpp:85
ublas::matrix< T, ublas::row_major, ublas::bounded_array< T, N > > MatrixBoundedArray
Definition Types.hpp:103
VectorBoundedArray< double, 3 > VectorDouble3
Definition Types.hpp:92
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
constexpr auto field_name
data for calculation inertia forces
VectorDouble a0
constant acceleration
common data used by volume elements
std::vector< std::vector< double * > > jacTRowPtr
std::map< std::string, std::vector< VectorDouble > > dataAtGaussPts
std::vector< std::vector< double * > > jacVelRowPtr
std::map< std::string, std::vector< MatrixDouble > > gradAtGaussPts
std::vector< std::vector< double * > > jacMassRowPtr
friend MoFEMErrorCode ZeroEntriesOp(Mat A)
friend MoFEMErrorCode MultOpA(Mat A, Vec x, Vec f)
MoFEMErrorCode postProcess()
function is run at the end of loop
int getRule(int order)
it is used to calculate nb. of Gauss integration points
MoFEMErrorCode preProcess()
function is run at the beginning of loop
OpEnergy(const std::string field_name, BlockData &data, CommonData &common_data, SmartPetscObj< Vec > v)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpEshelbyDynamicMaterialMomentumJacobian(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool jacobian=true)
OpEshelbyDynamicMaterialMomentumLhs_dX(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data, Range *forcesonlyonentities_ptr)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
OpEshelbyDynamicMaterialMomentumLhs_dv(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data, Range *forcesonlyonentities_ptr)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
OpEshelbyDynamicMaterialMomentumLhs_dx(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data, Range *forcesonlyonentities_ptr)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpEshelbyDynamicMaterialMomentumRhs(const std::string field_name, BlockData &data, CommonData &common_data, Range *forcesonlyonentities_ptr)
OpGetCommonDataAtGaussPts(const std::string field_name, CommonData &common_data)
OpGetDataAtGaussPts(const std::string field_name, std::vector< VectorDouble > &values_at_gauss_pts, std::vector< MatrixDouble > &gardient_at_gauss_pts)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
operator calculating deformation gradient
MatrixBoundedArray< adouble, 9 > invH
boost::ptr_vector< MethodForForceScaling > & methodsOp
OpMassJacobian(const std::string field_name, BlockData &data, CommonData &common_data, boost::ptr_vector< MethodForForceScaling > &methods_op, int tag, bool linear=false)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
VectorBoundedArray< adouble, 3 > a_res
VectorBoundedArray< adouble, 3 > dp_dt
VectorBoundedArray< adouble, 3 > dot_W
OpMassLhs_dM_dX(const std::string field_name, const std::string col_field, BlockData &data, CommonData &common_data)
MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpMassLhs_dM_dv(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data, Range *forcesonlyonentities_ptr=NULL)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
OpMassLhs_dM_dx(const std::string field_name, const std::string col_field, BlockData &data, CommonData &common_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpMassRhs(const std::string field_name, BlockData &data, CommonData &common_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpVelocityJacobian(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool jacobian=true)
OpVelocityLhs_dV_dX(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
OpVelocityLhs_dV_dv(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
OpVelocityLhs_dV_dx(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
OpVelocityRhs(const std::string field_name, BlockData &data, CommonData &common_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
friend MoFEMErrorCode PCShellSetUpOp(PC pc)
friend MoFEMErrorCode PCShellDestroy(PC pc)
bool initPC
check if PC is initialized
friend MoFEMErrorCode PCShellApplyOp(PC pc, Vec f, Vec x)
MatShellCtx * shellMatCtx
pointer to shell matrix
MoFEMErrorCode preProcess()
Calculate inconsistency between approximation of velocities and velocities calculated from displaceme...
UpdateAndControl(MoFEM::Interface &m_field, TS _ts, const std::string velocity_field, const std::string spatial_position_field)
MoFEMErrorCode postProcess()
function is run at the end of loop
MoFEMErrorCode preProcess()
Scatter values from t_u_dt on the fields.
structure grouping operators and data used for calculation of mass (convective) element \ nonlinear_e...
static MoFEMErrorCode MultOpA(Mat A, Vec x, Vec f)
Mult operator for shell matrix.
ConvectiveMassElement(MoFEM::Interface &m_field, short int tag)
static MoFEMErrorCode ZeroEntriesOp(Mat A)
MoFEMErrorCode setVelocityOperators(string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false)
MyVolumeFE feEnergy
calculate kinetic energy
MoFEMErrorCode setShellMatrixMassOperators(string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool linear=false)
MyVolumeFE feVelRhs
calculate right hand side for tetrahedral elements
MoFEMErrorCode addVelocityElement(string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel())
MyVolumeFE & getLoopFeVelRhs()
get rhs volume element
static MoFEMErrorCode PCShellDestroy(PC pc)
MoFEMErrorCode addConvectiveMassElement(string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel())
MyVolumeFE & getLoopFeTRhs()
get rhs volume element
static MoFEMErrorCode PCShellApplyOp(PC pc, Vec f, Vec x)
apply pre-conditioner for shell matrix
MoFEMErrorCode setKinematicEshelbyOperators(string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", Range *forces_on_entities_ptr=NULL)
MoFEMErrorCode addEshelbyDynamicMaterialMomentum(string element_name, string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, BitRefLevel bit=BitRefLevel(), Range *intersected=NULL)
MyVolumeFE feTRhs
calculate right hand side for tetrahedral elements
MyVolumeFE feMassRhs
calculate right hand side for tetrahedral elements
MyVolumeFE & getLoopFeMassRhs()
get rhs volume element
MyVolumeFE feTLhs
calculate left hand side for tetrahedral elements
MyVolumeFE & getLoopFeEnergy()
get kinetic energy element
MoFEMErrorCode setConvectiveMassOperators(string velocity_field_name, string spatial_position_field_name, string material_position_field_name="MESH_NODE_POSITIONS", bool ale=false, bool linear=false)
boost::ptr_vector< MethodForForceScaling > methodsOp
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
MyVolumeFE & getLoopFeMassLhs()
get lhs volume element
static MoFEMErrorCode PCShellSetUpOp(PC pc)
MyVolumeFE & getLoopFeVelLhs()
get lhs volume element
MyVolumeFE feVelLhs
calculate left hand side for tetrahedral elements
MyVolumeFE & getLoopFeMassAuxLhs()
get lhs volume element for Kuu shell matrix
MyVolumeFE & getLoopFeTLhs()
get lhs volume element
Set Dirichlet boundary conditions on displacements.
Deprecated interface functions.
structure for User Loop Methods on finite elements
intrusive_ptr for managing petsc objects