v0.13.0
PlasticOps.hpp
Go to the documentation of this file.
1 /* This file is part of MoFEM.
2  * MoFEM is free software: you can redistribute it and/or modify it under
3  * the terms of the GNU Lesser General Public License as published by the
4  * Free Software Foundation, either version 3 of the License, or (at your
5  * option) any later version.
6  *
7  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10  * License for more details.
11  *
12  * You should have received a copy of the GNU Lesser General Public
13  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14 
15 /** \file PlasticOps.hpp
16  * \example PlasticOps.hpp
17 
18 \f[
19 \left\{
20 \begin{array}{ll}
21 \frac{\partial \sigma_{ij}}{\partial x_j} - b_i = 0 & \forall x \in \Omega \\
22 \varepsilon_{ij} = \frac{1}{2}\left( \frac{\partial u_i}{\partial x_j} +
23 \frac{\partial u_j}{\partial x_i} \right)\\
24 \sigma_{ij} = D_{ijkl}\left(\varepsilon_{kl}-\varepsilon^p_{kl}\right) \\
25 \dot{\varepsilon}^p_{kl} - \dot{\tau} \left( \left. \frac{\partial f}{\partial
26 \sigma_{kl}} \right|_{(\sigma,\tau) } \right) = 0 \\
27 f(\sigma, \tau) \leq 0,\; \dot{\tau} \geq 0,\;\dot{\tau}f(\sigma, \tau)=0\\
28 u_i = \overline{u}_i & \forall x \in \partial\Omega_u \\
29 \sigma_{ij}n_j = \overline{t}_i & \forall x \in \partial\Omega_\sigma \\
30 \Omega_u \cup \Omega_\sigma = \Omega \\
31 \Omega_u \cap \Omega_\sigma = \emptyset
32 \end{array}
33 \right.
34 \f]
35 
36 \f[
37 \left\{
38 \begin{array}{ll}
39 \left(\frac{\partial \delta u_i}{\partial x_j},\sigma_{ij}\right)_\Omega-(\delta
40 u_i,b_i)_\Omega -(\delta u_i,\overline{t}_i)_{\partial\Omega_\sigma}=0 & \forall
41 \delta u_i \in H^1(\Omega)\\ \left(\delta\varepsilon^p_{kl} ,D_{ijkl}\left(
42 \dot{\varepsilon}^p_{kl} - \dot{\tau} A_{kl} \right)\right) = 0
43 & \forall \delta\varepsilon^p_{ij} \in L^2(\Omega) \cap \mathcal{S} \\
44 \left(\delta\tau,c_n\dot{\tau} - \frac{1}{2}\left\{c_n \dot{\tau} +
45 (f(\pmb\sigma,\tau) - \sigma_y) +
46 \| c_n \dot{\tau} + (f(\pmb\sigma,\tau) - \sigma_y) \|\right\}\right) = 0 &
47 \forall \delta\tau \in L^2(\Omega) \end{array} \right.
48 \f]
49 
50 */
51 
52 namespace PlasticOps {
53 
54 //! [Common data]
55 struct CommonData : public boost::enable_shared_from_this<CommonData> {
56  boost::shared_ptr<MatrixDouble> mDPtr;
57  boost::shared_ptr<MatrixDouble> mDPtr_Axiator;
58  boost::shared_ptr<MatrixDouble> mDPtr_Deviator;
59  boost::shared_ptr<MatrixDouble> mGradPtr;
60  boost::shared_ptr<MatrixDouble> mStrainPtr;
61  boost::shared_ptr<MatrixDouble> mStressPtr;
62 
70 
71  inline auto getPlasticTauPtr() {
72  return boost::shared_ptr<VectorDouble>(shared_from_this(), &plasticTau);
73  }
74  inline auto getPlasticTauDotPtr() {
75  return boost::shared_ptr<VectorDouble>(shared_from_this(), &plasticTauDot);
76  }
77  inline auto getPlasticStrainPtr() {
78  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &plasticStrain);
79  }
80  inline auto getPlasticStrainDotPtr() {
81  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
83  }
84  inline auto getTempValPtr() {
85  return boost::shared_ptr<VectorDouble>(shared_from_this(), &tempVal);
86  }
87 
88 };
89 //! [Common data]
90 
97 
102 
103 FTensor::Index<'L', (SPACE_DIM * (SPACE_DIM + 1)) / 2> L;
104 
105 //! [Operators definitions]
107  OpCalculatePlasticSurface(const std::string field_name,
108  boost::shared_ptr<CommonData> common_data_ptr);
109  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
110 
111 protected:
112  boost::shared_ptr<CommonData> commonDataPtr;
113 };
114 
115 struct OpPlasticStress : public DomainEleOp {
116  OpPlasticStress(const std::string field_name,
117  boost::shared_ptr<CommonData> common_data_ptr,
118  boost::shared_ptr<MatrixDouble> mDPtr,
119  const double scale = 1);
120  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
121 
122 private:
123  boost::shared_ptr<MatrixDouble> mDPtr;
124  const double scaleStress;
125  boost::shared_ptr<CommonData> commonDataPtr;
126 };
127 
129  OpCalculatePlasticFlowRhs(const std::string field_name,
130  boost::shared_ptr<CommonData> common_data_ptr);
132 
133 private:
134  boost::shared_ptr<CommonData> commonDataPtr;
135 };
136 
138  OpCalculateContrainsRhs(const std::string field_name,
139  boost::shared_ptr<CommonData> common_data_ptr);
141 
142 private:
143  boost::shared_ptr<CommonData> commonDataPtr;
144 };
145 
148  const std::string row_field_name, const std::string col_field_name,
149  boost::shared_ptr<CommonData> common_data_ptr,
150  boost::shared_ptr<MatrixDouble> m_D_ptr);
152  EntitiesFieldData::EntData &col_data);
153 
154 private:
155  boost::shared_ptr<CommonData> commonDataPtr;
156  boost::shared_ptr<MatrixDouble> mDPtr;
157 };
158 
160  : public AssemblyDomainEleOp {
162  const std::string row_field_name, const std::string col_field_name,
163  boost::shared_ptr<CommonData> common_data_ptr,
164  boost::shared_ptr<HenckyOps::CommonData> common_henky_data_ptr,
165  boost::shared_ptr<MatrixDouble> m_D_ptr);
167  EntitiesFieldData::EntData &col_data);
168 
169 private:
170  boost::shared_ptr<CommonData> commonDataPtr;
171  boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
172  boost::shared_ptr<MatrixDouble> mDPtr;
174 };
175 
177  OpCalculatePlasticFlowLhs_dU(const std::string row_field_name,
178  const std::string col_field_name,
179  boost::shared_ptr<CommonData> common_data_ptr,
180  boost::shared_ptr<MatrixDouble> m_D_ptr);
182  EntitiesFieldData::EntData &col_data);
183 
184 private:
185  boost::shared_ptr<CommonData> commonDataPtr;
186  boost::shared_ptr<MatrixDouble> mDPtr;
187 };
188 
191  const std::string row_field_name, const std::string col_field_name,
192  boost::shared_ptr<CommonData> common_data_ptr,
193  boost::shared_ptr<HenckyOps::CommonData> comman_henky_data_ptr,
194  boost::shared_ptr<MatrixDouble> m_D_ptr);
196  EntitiesFieldData::EntData &col_data);
197 
198 private:
199  boost::shared_ptr<CommonData> commonDataPtr;
200  boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
201  boost::shared_ptr<MatrixDouble> mDPtr;
202 };
203 
205  OpCalculatePlasticFlowLhs_dEP(const std::string row_field_name,
206  const std::string col_field_name,
207  boost::shared_ptr<CommonData> common_data_ptr,
208  boost::shared_ptr<MatrixDouble> m_D_ptr);
210  EntitiesFieldData::EntData &col_data);
211 
212 private:
213  boost::shared_ptr<CommonData> commonDataPtr;
214  boost::shared_ptr<MatrixDouble> mDPtr;
215 };
216 
218  OpCalculatePlasticFlowLhs_dTAU(const std::string row_field_name,
219  const std::string col_field_name,
220  boost::shared_ptr<CommonData> common_data_ptr);
222  EntitiesFieldData::EntData &col_data);
223 
224 private:
225  boost::shared_ptr<CommonData> commonDataPtr;
226 };
227 
229  OpCalculateContrainsLhs_dU(const std::string row_field_name,
230  const std::string col_field_name,
231  boost::shared_ptr<CommonData> common_data_ptr,
232  boost::shared_ptr<MatrixDouble> m_D_ptr);
234  EntitiesFieldData::EntData &col_data);
235 
236 private:
237  boost::shared_ptr<CommonData> commonDataPtr;
238  boost::shared_ptr<MatrixDouble> mDPtr;
239 };
240 
243  const std::string row_field_name, const std::string col_field_name,
244  boost::shared_ptr<CommonData> common_data_ptr,
245  boost::shared_ptr<HenckyOps::CommonData> comman_henky_data_ptr,
246  boost::shared_ptr<MatrixDouble> m_D_ptr);
248  EntitiesFieldData::EntData &col_data);
249 
250 private:
251  boost::shared_ptr<CommonData> commonDataPtr;
252  boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
253  boost::shared_ptr<MatrixDouble> mDPtr;
254 };
255 
257  OpCalculateContrainsLhs_dEP(const std::string row_field_name,
258  const std::string col_field_name,
259  boost::shared_ptr<CommonData> common_data_ptr,
260  boost::shared_ptr<MatrixDouble> mat_D_ptr);
262  EntitiesFieldData::EntData &col_data);
263 
264 private:
265  boost::shared_ptr<CommonData> commonDataPtr;
266  boost::shared_ptr<MatrixDouble> mDPtr;
267 };
268 
270  OpCalculateContrainsLhs_dTAU(const std::string row_field_name,
271  const std::string col_field_name,
272  boost::shared_ptr<CommonData> common_data_ptr);
274  EntitiesFieldData::EntData &col_data);
275 private:
276  boost::shared_ptr<CommonData> commonDataPtr;
277 };
278 
280  OpPostProcPlastic(const std::string field_name,
281  moab::Interface &post_proc_mesh,
282  std::vector<EntityHandle> &map_gauss_pts,
283  boost::shared_ptr<CommonData> common_data_ptr);
284  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
285 
286 private:
288  std::vector<EntityHandle> &mapGaussPts;
289  boost::shared_ptr<CommonData> commonDataPtr;
290 };
291 //! [Operators definitions]
292 
293 //! [Lambda functions]
294 inline auto diff_tensor() {
297  t_diff(i, j, k, l) = (t_kd(i, k) ^ t_kd(j, l)) / 4.;
298  return t_diff;
299 };
300 
301 inline auto symm_L_tensor() {
302  constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
304  t_L(i, j, L) = 0;
305  if (SPACE_DIM == 2) {
306  t_L(0, 0, 0) = 1;
307  t_L(1, 0, 1) = 1;
308  t_L(1, 1, 2) = 1;
309  } else if (SPACE_DIM == 3) {
310  t_L(0, 0, 0) = 1;
311  t_L(1, 0, 1) = 1;
312  t_L(2, 0, 2) = 1;
313  t_L(1, 1, 3) = 1;
314  t_L(2, 1, 4) = 1;
315  t_L(2, 2, 5) = 1;
316  }
317  return t_L;
318 }
319 
320 inline auto diff_symmetrize() {
322 
323  t_diff(i, j, k, l) = 0;
324  t_diff(0, 0, 0, 0) = 1;
325  t_diff(1, 1, 1, 1) = 1;
326 
327  t_diff(1, 0, 1, 0) = 0.5;
328  t_diff(1, 0, 0, 1) = 0.5;
329 
330  t_diff(0, 1, 0, 1) = 0.5;
331  t_diff(0, 1, 1, 0) = 0.5;
332 
333  if (SPACE_DIM == 3) {
334  t_diff(2, 2, 2, 2) = 1;
335 
336  t_diff(2, 0, 2, 0) = 0.5;
337  t_diff(2, 0, 0, 2) = 0.5;
338  t_diff(0, 2, 0, 2) = 0.5;
339  t_diff(0, 2, 2, 0) = 0.5;
340 
341  t_diff(2, 1, 2, 1) = 0.5;
342  t_diff(2, 1, 1, 2) = 0.5;
343  t_diff(1, 2, 1, 2) = 0.5;
344  t_diff(1, 2, 2, 1) = 0.5;
345  }
346 
347  return t_diff;
348 };
349 
350 template <typename T>
352  constexpr double third = boost::math::constants::third<double>();
353  if (SPACE_DIM == 2)
354  return (t_stress(0, 0) + t_stress(1, 1)) * third;
355  else
356  return (t_stress(0, 0) + t_stress(1, 1) + t_stress(2, 2)) * third;
357 };
358 
359 template <typename T>
361  double trace) {
363  t_dev(I, J) = 0;
364  for (int ii = 0; ii != SPACE_DIM; ++ii)
365  for (int jj = ii; jj != SPACE_DIM; ++jj)
366  t_dev(ii, jj) = t_stress(ii, jj);
367  t_dev(0, 0) -= trace;
368  t_dev(1, 1) -= trace;
369  t_dev(2, 2) -= trace;
370  return t_dev;
371 };
372 
373 inline auto
375  FTensor::Ddg<double, 3, SPACE_DIM> t_diff_deviator;
376  t_diff_deviator(I, J, k, l) = 0;
377  for (int ii = 0; ii != SPACE_DIM; ++ii)
378  for (int jj = ii; jj != SPACE_DIM; ++jj)
379  for (int kk = 0; kk != SPACE_DIM; ++kk)
380  for (int ll = kk; ll != SPACE_DIM; ++ll)
381  t_diff_deviator(ii, jj, kk, ll) = t_diff_stress(ii, jj, kk, ll);
382 
383  constexpr double third = boost::math::constants::third<double>();
384 
385  t_diff_deviator(0, 0, 0, 0) -= third;
386  t_diff_deviator(0, 0, 1, 1) -= third;
387 
388  t_diff_deviator(1, 1, 0, 0) -= third;
389  t_diff_deviator(1, 1, 1, 1) -= third;
390 
391  t_diff_deviator(2, 2, 0, 0) -= third;
392  t_diff_deviator(2, 2, 1, 1) -= third;
393 
394  if (SPACE_DIM == 3) {
395  t_diff_deviator(0, 0, 2, 2) -= third;
396  t_diff_deviator(1, 1, 2, 2) -= third;
397  t_diff_deviator(2, 2, 2, 2) -= third;
398  }
399 
400  return t_diff_deviator;
401 };
402 
403 /**
404  *
405 
406 \f[
407 \begin{split}
408 f&=\sqrt{s_{ij}s_{ij}}\\
409 A_{ij}&=\frac{\partial f}{\partial \sigma_{ij}}=
410 \frac{1}{f} s_{kl} \frac{\partial s_{kl}}{\partial \sigma_{ij}}\\
411 \frac{\partial A_{ij}}{\partial \sigma_{kl}}&= \frac{\partial^2 f}{\partial
412 \sigma_{ij}\partial\sigma_{mn}}= \frac{1}{f} \left( \frac{\partial
413 s_{kl}}{\partial \sigma_{mn}}\frac{\partial s_{kl}}{\partial \sigma_{ij}}
414 -A_{mn}A_{ij}
415 \right)\\
416 \frac{\partial f}{\partial \varepsilon_{ij}}&=A_{mn}D_{mnij}
417 \\
418 \frac{\partial A_{ij}}{\partial \varepsilon_{kl}}&=
419 \frac{\partial A_{ij}}{\partial \sigma_{mn}} \frac{\partial
420 \sigma_{mn}}{\partial \varepsilon_{kl}}= \frac{\partial A_{ij}}{\partial
421 \sigma_{mn}} D_{mnkl}
422 \end{split}
423 \f]
424 
425  */
426 inline double
428  return std::sqrt(1.5 * t_stress_deviator(I, J) * t_stress_deviator(I, J)) +
429  std::numeric_limits<double>::epsilon();
430 };
431 
432 inline auto plastic_flow(long double f,
434  FTensor::Ddg<double, 3, SPACE_DIM> &&t_diff_deviator) {
436  t_diff_f(k, l) =
437  (1.5 * (t_dev_stress(I, J) * t_diff_deviator(I, J, k, l))) / f;
438  return t_diff_f;
439 };
440 
441 template <typename T>
443  long double f, FTensor::Tensor2_symmetric<T, SPACE_DIM> &t_flow,
444  FTensor::Ddg<double, 3, SPACE_DIM> &&t_diff_deviator) {
446  t_diff_flow(i, j, k, l) =
447  (1.5 * (t_diff_deviator(M, N, i, j) * t_diff_deviator(M, N, k, l) -
448  (2. / 3.) * t_flow(i, j) * t_flow(k, l))) /
449  f;
450  return t_diff_flow;
451 };
452 
453 template <typename T>
456  FTensor::Ddg<double, SPACE_DIM, SPACE_DIM> &&t_diff_plastic_flow_dstress) {
458  t_diff_flow(i, j, k, l) =
459  t_diff_plastic_flow_dstress(i, j, m, n) * t_D(m, n, k, l);
460  return t_diff_flow;
461 };
462 
463 inline double constrain_abs(long double x) {
464  return std::abs(x);
465 };
466 
467 inline double constrian_sign(long double x) {
468  if (x > 0)
469  return 1;
470  else if (x < 0)
471  return -1;
472  else
473  return 0;
474 };
475 
476 inline double w(long double dot_tau, long double f, long double sigma_y) {
477  return (f - sigma_y) / sigmaY + cn * dot_tau;
478 };
479 
480 /**
481 
482 \f[
483 \dot{\tau} - \frac{1}{2}\left\{\dot{\tau} + (f(\pmb\sigma) - \sigma_y) +
484 \| \dot{\tau} + (f(\pmb\sigma) - \sigma_y) \|\right\} = 0 \\
485 c_n \sigma_y \dot{\tau} - \frac{1}{2}\left\{c_n\sigma_y \dot{\tau} +
486 (f(\pmb\sigma) - \sigma_y) +
487 \| c_n \sigma_y \dot{\tau} + (f(\pmb\sigma) - \sigma_y) \|\right\} = 0
488 \f]
489 
490  */
491 inline double constrain(long double dot_tau, long double f,
492  long double sigma_y) {
493  return visH * dot_tau +
494  (sigmaY / 2) * ((cn * dot_tau - (f - sigma_y) / sigmaY) -
495  constrain_abs(w(dot_tau, f, sigma_y)));
496 };
497 
498 inline double diff_constrain_ddot_tau(long double dot_tau, long double f,
499  long double sigma_y) {
500  return visH +
501  (sigmaY / 2) * (cn - cn * constrian_sign(w(dot_tau, f, sigma_y)));
502 };
503 
504 inline auto diff_constrain_df(long double dot_tau, long double f,
505  long double sigma_y) {
506  return (-1 - constrian_sign(w(dot_tau, f, sigma_y))) / 2;
507 };
508 
509 inline auto diff_constrain_dsigma_y(long double dot_tau, long double f,
510  long double sigma_y) {
511  return (1 + constrian_sign(w(dot_tau, f, sigma_y))) / 2;
512 }
513 
514 template <typename T>
516  double &&diff_constrain_df,
518  FTensor::Tensor2_symmetric<double, SPACE_DIM> t_diff_constrain_dstress;
519  t_diff_constrain_dstress(i, j) = diff_constrain_df * t_plastic_flow(i, j);
520  return t_diff_constrain_dstress;
521 };
522 
523 template <typename T1, typename T2>
524 inline auto diff_constrain_dstrain(T1 &t_D, T2 &&t_diff_constrain_dstress) {
525  FTensor::Tensor2_symmetric<double, SPACE_DIM> t_diff_constrain_dstrain;
526  t_diff_constrain_dstrain(k, l) =
527  t_diff_constrain_dstress(i, j) * t_D(i, j, k, l);
528  return t_diff_constrain_dstrain;
529 };
530 //! [Lambda functions]
531 
532 //! [Auxiliary functions functions
536  &mat(3 * rr + 0, 0), &mat(3 * rr + 0, 1), &mat(3 * rr + 1, 0),
537  &mat(3 * rr + 1, 1), &mat(3 * rr + 1, 0), &mat(3 * rr + 1, 1),
538  &mat(3 * rr + 2, 0), &mat(3 * rr + 2, 1)};
539 }
540 
544  &mat(6 * rr + 0, 0), &mat(6 * rr + 0, 1), &mat(6 * rr + 0, 2),
545  &mat(6 * rr + 1, 0), &mat(6 * rr + 1, 1), &mat(6 * rr + 1, 2),
546  &mat(6 * rr + 2, 0), &mat(6 * rr + 2, 1), &mat(6 * rr + 2, 2),
547  &mat(6 * rr + 1, 0), &mat(6 * rr + 1, 1), &mat(6 * rr + 1, 2),
548  &mat(6 * rr + 3, 0), &mat(6 * rr + 3, 1), &mat(6 * rr + 3, 2),
549  &mat(6 * rr + 4, 0), &mat(6 * rr + 4, 1), &mat(6 * rr + 4, 2),
550  &mat(6 * rr + 2, 0), &mat(6 * rr + 2, 1), &mat(6 * rr + 2, 2),
551  &mat(6 * rr + 4, 0), &mat(6 * rr + 4, 1), &mat(6 * rr + 4, 2),
552  &mat(6 * rr + 5, 0), &mat(6 * rr + 5, 1), &mat(6 * rr + 5, 2)};
553 }
554 
555 static inline auto get_mat_scalar_dvector(MatrixDouble &mat,
558  &mat(0, 1)};
559 }
560 
561 static inline auto get_mat_scalar_dvector(MatrixDouble &mat,
564  &mat(0, 0), &mat(0, 1), &mat(0, 2)};
565 }
566 //! [Auxiliary functions functions
567 
568 
569 }; // namespace PlasticOps
570 
571 #include <PlasticOpsGeneric.hpp>
574 #include <PlasticOpsMonitor.hpp>
constexpr double third
EntitiesFieldData::EntData EntData
constexpr int SPACE_DIM
Kronecker Delta class symmetric.
constexpr auto t_kd
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
auto symm_L_tensor()
Definition: PlasticOps.hpp:301
auto deviator(FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress, double trace)
Definition: PlasticOps.hpp:360
auto diff_constrain_dstress(double &&diff_constrain_df, FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_plastic_flow)
Definition: PlasticOps.hpp:515
double platsic_surface(FTensor::Tensor2_symmetric< double, 3 > &&t_stress_deviator)
Definition: PlasticOps.hpp:427
FTensor::Index< 'j', SPACE_DIM > j
[Common data]
Definition: PlasticOps.hpp:91
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:100
auto diff_constrain_dstrain(T1 &t_D, T2 &&t_diff_constrain_dstress)
Definition: PlasticOps.hpp:524
auto plastic_flow(long double f, FTensor::Tensor2_symmetric< double, 3 > &&t_dev_stress, FTensor::Ddg< double, 3, SPACE_DIM > &&t_diff_deviator)
Definition: PlasticOps.hpp:432
static auto get_mat_scalar_dvector(MatrixDouble &mat, FTensor::Number< 2 >)
Definition: PlasticOps.hpp:555
FTensor::Index< 'l', SPACE_DIM > l
Definition: PlasticOps.hpp:93
FTensor::Index< 'k', SPACE_DIM > k
Definition: PlasticOps.hpp:92
auto diff_tensor()
[Operators definitions]
Definition: PlasticOps.hpp:294
FTensor::Index< 'N', 3 > N
Definition: PlasticOps.hpp:101
FTensor::Index< 'I', 3 > I
Definition: PlasticOps.hpp:98
double constrain(long double dot_tau, long double f, long double sigma_y)
Definition: PlasticOps.hpp:491
auto diff_deviator(FTensor::Ddg< double, SPACE_DIM, SPACE_DIM > &&t_diff_stress)
Definition: PlasticOps.hpp:374
FTensor::Index< 'i', SPACE_DIM > i
Definition: PlasticOps.hpp:95
double w(long double dot_tau, long double f, long double sigma_y)
Definition: PlasticOps.hpp:476
FTensor::Index< 'm', SPACE_DIM > m
Definition: PlasticOps.hpp:94
auto diff_constrain_dsigma_y(long double dot_tau, long double f, long double sigma_y)
Definition: PlasticOps.hpp:509
FTensor::Index< 'J', 3 > J
Definition: PlasticOps.hpp:99
double trace(FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress)
Definition: PlasticOps.hpp:351
auto diff_symmetrize()
Definition: PlasticOps.hpp:320
double constrian_sign(long double x)
Definition: PlasticOps.hpp:467
auto diff_plastic_flow_dstress(long double f, FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_flow, FTensor::Ddg< double, 3, SPACE_DIM > &&t_diff_deviator)
Definition: PlasticOps.hpp:442
auto diff_plastic_flow_dstrain(FTensor::Ddg< T, SPACE_DIM, SPACE_DIM > &t_D, FTensor::Ddg< double, SPACE_DIM, SPACE_DIM > &&t_diff_plastic_flow_dstress)
Definition: PlasticOps.hpp:454
double diff_constrain_ddot_tau(long double dot_tau, long double f, long double sigma_y)
Definition: PlasticOps.hpp:498
double constrain_abs(long double x)
Definition: PlasticOps.hpp:463
auto diff_constrain_df(long double dot_tau, long double f, long double sigma_y)
Definition: PlasticOps.hpp:504
static FTensor::Tensor3< FTensor::PackPtr< double *, 2 >, 2, 2, 2 > get_mat_tensor_sym_dvector(size_t rr, MatrixDouble &mat, FTensor::Number< 2 >)
[Lambda functions]
Definition: PlasticOps.hpp:534
FTensor::Index< 'L',(SPACE_DIM *(SPACE_DIM+1))/2 > L
Definition: PlasticOps.hpp:103
FTensor::Index< 'n', SPACE_DIM > n
Definition: PlasticOps.hpp:96
double visH
Definition: plastic.cpp:110
double scale
Definition: plastic.cpp:103
double sigmaY
Definition: plastic.cpp:108
double cn
Definition: plastic.cpp:111
Data on single entity (This is passed as argument to DataOperator::doWork)
boost::shared_ptr< MatrixDouble > mStrainPtr
Definition: PlasticOps.hpp:60
VectorDouble tempVal
Definition: PlasticOps.hpp:69
boost::shared_ptr< MatrixDouble > mGradPtr
Definition: PlasticOps.hpp:59
VectorDouble plasticTauDot
Definition: PlasticOps.hpp:66
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:56
VectorDouble plasticSurface
Definition: PlasticOps.hpp:63
boost::shared_ptr< MatrixDouble > mDPtr_Axiator
Definition: PlasticOps.hpp:57
MatrixDouble plasticStrain
Definition: PlasticOps.hpp:67
boost::shared_ptr< MatrixDouble > mStressPtr
Definition: PlasticOps.hpp:61
MatrixDouble plasticFlow
Definition: PlasticOps.hpp:64
MatrixDouble plasticStrainDot
Definition: PlasticOps.hpp:68
VectorDouble plasticTau
Definition: PlasticOps.hpp:65
boost::shared_ptr< MatrixDouble > mDPtr_Deviator
Definition: PlasticOps.hpp:58
OpCalculateContrainsLhs_LogStrain_dU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< HenckyOps::CommonData > comman_henky_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:251
boost::shared_ptr< HenckyOps::CommonData > commonHenckyDataPtr
Definition: PlasticOps.hpp:252
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:253
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:266
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpCalculateContrainsLhs_dEP(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > mat_D_ptr)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:265
OpCalculateContrainsLhs_dTAU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:276
OpCalculateContrainsLhs_dU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:238
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:237
OpCalculateContrainsRhs(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:143
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpCalculatePlasticFlowLhs_LogStrain_dU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< HenckyOps::CommonData > comman_henky_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:201
boost::shared_ptr< HenckyOps::CommonData > commonHenckyDataPtr
Definition: PlasticOps.hpp:200
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:199
OpCalculatePlasticFlowLhs_dEP(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:214
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:213
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:225
OpCalculatePlasticFlowLhs_dTAU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpCalculatePlasticFlowLhs_dU(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:186
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:185
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:134
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
OpCalculatePlasticFlowRhs(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
[Calculate stress]
boost::shared_ptr< HenckyOps::CommonData > commonHenckyDataPtr
Definition: PlasticOps.hpp:171
OpCalculatePlasticInternalForceLhs_LogStrain_dEP(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< HenckyOps::CommonData > common_henky_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpCalculatePlasticInternalForceLhs_dEP(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:155
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:156
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculatePlasticSurface(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:112
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:125
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[Calculate stress]
OpPlasticStress(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > mDPtr, const double scale=1)
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOps.hpp:123
Definition: PlasticOps.hpp:279
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOps.hpp:289
std::vector< EntityHandle > & mapGaussPts
Definition: PlasticOps.hpp:288
OpPostProcPlastic(const std::string field_name, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > common_data_ptr)
moab::Interface & postProcMesh
Definition: PlasticOps.hpp:287
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[Postprocessing]