v0.14.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
31 struct MyVolumeFE : public VolumeElementForcesAndSourcesCore {
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
55 SmartPetscObj<Vec> V;
56 double eNergy;
57
58 MoFEMErrorCode preProcess();
59 MoFEMErrorCode postProcess();
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;
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
199 VectorBoundedArray<adouble, 3> a, dot_W, dp_dt, a_res;
200 MatrixBoundedArray<adouble, 9> h, H, invH, F, g, G;
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
216 VectorDouble nf;
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
234 MatrixDouble k, jac;
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
266 SmartPetscObj<Vec> V;
267 bool &lInear;
268
269 OpEnergy(const std::string field_name, BlockData &data,
270 CommonData &common_data, SmartPetscObj<Vec> v);
271
272 MatrixDouble3by3 h, H, invH, F;
273 VectorDouble3 v;
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
291 VectorBoundedArray<adouble, 3> a_res, v, dot_w, dot_W, dot_u;
292 MatrixBoundedArray<adouble, 9> h, H, invH, F;
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
311 VectorDouble nf;
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
362 VectorBoundedArray<adouble, 3> a, v, a_T;
363 MatrixBoundedArray<adouble, 9> g, H, invH, h, F, G;
364 VectorDouble active;
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
383 VectorDouble nf;
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 */
451 MoFEMErrorCode preProcess();
452
453 // This is empty fun cions does nothing
454 MoFEMErrorCode postProcess();
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;
514 MoFEMErrorCode iNit();
515
516 MoFEMErrorCode dEstroy();
517
518 friend MoFEMErrorCode MultOpA(Mat A, Vec x, Vec f);
519 friend MoFEMErrorCode ZeroEntriesOp(Mat A);
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
609 MoFEMErrorCode iNit();
610
611 MoFEMErrorCode dEstroy();
612
613 friend MoFEMErrorCode PCShellSetUpOp(PC pc);
614 friend MoFEMErrorCode PCShellDestroy(PC pc);
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 */
728 MoFEMErrorCode preProcess();
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()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#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
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
constexpr AssemblyType A
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)
int getRule(int order)
it is used to calculate nb. of Gauss integration points
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
operator calculating deformation gradient
MatrixBoundedArray< adouble, 9 > invH
MatrixBoundedArray< adouble, 9 > G
boost::ptr_vector< MethodForForceScaling > & methodsOp
MatrixBoundedArray< adouble, 9 > H
MatrixBoundedArray< adouble, 9 > g
MatrixBoundedArray< adouble, 9 > h
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
VectorBoundedArray< adouble, 3 > a_res
VectorBoundedArray< adouble, 3 > dp_dt
MatrixBoundedArray< adouble, 9 > F
VectorBoundedArray< adouble, 3 > dot_W
VectorBoundedArray< adouble, 3 > a
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)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
virtual MoFEMErrorCode getJac(EntitiesFieldData::EntData &col_data, int gg)
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...
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.
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.
Definition: DirichletBC.hpp:57
Deprecated interface functions.