v0.13.2
Loading...
Searching...
No Matches
Classes | Functions | Variables
PlasticOps Namespace Reference

Classes

struct  CommonData
 [Common data] More...
 
struct  Monitor
 
struct  OpCalculateConstraintsLhs_dEP
 
struct  OpCalculateConstraintsLhs_dTAU
 
struct  OpCalculateConstraintsLhs_dU
 
struct  OpCalculateConstraintsLhs_LogStrain_dU
 
struct  OpCalculateConstraintsRhs
 
struct  OpCalculatePlasticFlowLhs_dEP
 
struct  OpCalculatePlasticFlowLhs_dTAU
 
struct  OpCalculatePlasticFlowLhs_dU
 
struct  OpCalculatePlasticFlowLhs_LogStrain_dU
 
struct  OpCalculatePlasticFlowRhs
 
struct  OpCalculatePlasticInternalForceLhs_dEP
 
struct  OpCalculatePlasticInternalForceLhs_LogStrain_dEP
 
struct  OpCalculatePlasticity
 
struct  OpCalculatePlasticSurface
 [Operators definitions] More...
 
struct  OpPlasticStress
 

Functions

auto diff_tensor ()
 [Operators definitions] More...
 
auto symm_L_tensor ()
 
auto diff_symmetrize ()
 
template<typename T >
double trace (FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress)
 
template<typename T >
auto deviator (FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress, double trace)
 
auto diff_deviator (FTensor::Ddg< double, SPACE_DIM, SPACE_DIM > &&t_diff_stress)
 
double platsic_surface (FTensor::Tensor2_symmetric< double, 3 > &&t_stress_deviator)
 
auto plastic_flow (long double f, FTensor::Tensor2_symmetric< double, 3 > &&t_dev_stress, FTensor::Ddg< double, 3, SPACE_DIM > &&t_diff_deviator)
 
template<typename T >
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)
 
template<typename T >
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)
 
double constrain_diff_sign (double x)
 
double constrian_sign (double x)
 
double constrain_abs (double x)
 
double w (double eqiv, double dot_tau, double f, double sigma_y)
 
double constraint (double eqiv, double dot_tau, double f, double sigma_y, double abs_w)
 
double diff_constrain_ddot_tau (double sign, double eqiv)
 
double diff_constrain_eqiv (double sign, double eqiv, double dot_tau)
 
auto diff_constrain_df (double sign)
 
auto diff_constrain_dsigma_y (double sign)
 
template<typename T >
auto diff_constrain_dstress (double diff_constrain_df, FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_plastic_flow)
 
template<typename T1 , typename T2 >
auto diff_constrain_dstrain (T1 &t_D, T2 &&t_diff_constrain_dstress)
 
template<typename T >
auto equivalent_strain_dot (FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_plastic_strain_dot)
 
template<typename T1 , typename T2 , typename T3 >
auto diff_equivalent_strain_dot (const T1 eqiv, T2 &t_plastic_strain_dot, T3 &t_diff_plastic_strain)
 
static auto get_mat_tensor_sym_dvector (size_t rr, MatrixDouble &mat, FTensor::Number< 2 >)
 [Lambda functions] More...
 
static auto get_mat_tensor_sym_dvector (size_t rr, MatrixDouble &mat, FTensor::Number< 3 >)
 
static auto get_mat_tensor_sym_dtensor_sym (size_t rr, MatrixDouble &mat, FTensor::Number< 2 >)
 
static auto get_mat_tensor_sym_dtensor_sym (size_t rr, MatrixDouble &mat, FTensor::Number< 3 >)
 
static auto get_mat_tensor_sym_dscalar (size_t rr, MatrixDouble &mat, FTensor::Number< 2 >)
 
static auto get_mat_tensor_sym_dscalar (size_t rr, MatrixDouble &mat, FTensor::Number< 3 >)
 
auto get_mat_scalar_dtensor_sym (MatrixDouble &mat, FTensor::Number< 2 >)
 
auto get_mat_scalar_dtensor_sym (MatrixDouble &mat, FTensor::Number< 3 >)
 
static FTensor::Tensor2< FTensor::PackPtr< double *, 3 >, 2, 3 > get_mat_vector_dtensor_sym (size_t rr, MatrixDouble &mat, FTensor::Number< 2 >)
 
static FTensor::Tensor2< FTensor::PackPtr< double *, 6 >, 3, 6 > get_mat_vector_dtensor_sym (size_t rr, MatrixDouble &mat, FTensor::Number< 3 >)
 

Variables

FTensor::Index< 'i', SPACE_DIMi
 [Common data] More...
 
FTensor::Index< 'j', SPACE_DIMj
 
FTensor::Index< 'k', SPACE_DIMk
 
FTensor::Index< 'l', SPACE_DIMl
 
FTensor::Index< 'm', SPACE_DIMm
 
FTensor::Index< 'n', SPACE_DIMn
 
FTensor::Index< 'o', SPACE_DIMo
 
FTensor::Index< 'p', SPACE_DIMp
 
FTensor::Index< 'I', 3 > I
 
FTensor::Index< 'J', 3 > J
 
FTensor::Index< 'M', 3 > M
 
FTensor::Index< 'N', 3 > N
 
FTensor::Index< 'L', size_symmL
 
FTensor::Index< 'O', size_symmO
 

Function Documentation

◆ constrain_abs()

double PlasticOps::constrain_abs ( double  x)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 502 of file PlasticOps.hpp.

502 {
503 const auto y = -x / zeta;
504 if (y > std::numeric_limits<float>::max_exponent10 ||
505 y < std::numeric_limits<float>::min_exponent10) {
506 return std::abs(x);
507 } else {
508 const double e = std::exp(y);
509 return x + 2 * zeta * std::log1p(e);
510 }
511};
double zeta
Definition: plastic.cpp:107

◆ constrain_diff_sign()

double PlasticOps::constrain_diff_sign ( double  x)
inline
Examples
PlasticOps.hpp.

Definition at line 476 of file PlasticOps.hpp.

476 {
477 const auto y = x / zeta;
478 if (y > std::numeric_limits<float>::max_exponent10 ||
479 y < std::numeric_limits<float>::min_exponent10) {
480 return 0;
481 } else {
482 const auto e = std::exp(y);
483 const auto ep1 = e + 1;
484 return (2 / zeta) * (e / (ep1 * ep1));
485 }
486};

◆ constraint()

double PlasticOps::constraint ( double  eqiv,
double  dot_tau,
double  f,
double  sigma_y,
double  abs_w 
)
inline

\[ \dot{\tau} - \frac{1}{2}\left\{\dot{\tau} + (f(\pmb\sigma) - \sigma_y) + \| \dot{\tau} + (f(\pmb\sigma) - \sigma_y) \|\right\} = 0 \\ c_n \sigma_y \dot{\tau} - \frac{1}{2}\left\{c_n\sigma_y \dot{\tau} + (f(\pmb\sigma) - \sigma_y) + \| c_n \sigma_y \dot{\tau} + (f(\pmb\sigma) - \sigma_y) \|\right\} = 0 \]

Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 528 of file PlasticOps.hpp.

529 {
530 return visH * dot_tau + (sigmaY / 2) * ((cn0 * (dot_tau - eqiv) +
531 cn1 * (std::sqrt(eqiv) * dot_tau) -
532 (f - sigma_y) / sigmaY) -
533 abs_w);
534};
double visH
Definition: plastic.cpp:104
double cn0
Definition: plastic.cpp:105
double sigmaY
Definition: plastic.cpp:102
double cn1
Definition: plastic.cpp:106

◆ constrian_sign()

double PlasticOps::constrian_sign ( double  x)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 488 of file PlasticOps.hpp.

488 {
489 const auto y = x / zeta;
490 if (y > std::numeric_limits<float>::max_exponent10 ||
491 y < std::numeric_limits<float>::min_exponent10) {
492 if (x > 0)
493 return 1.;
494 else
495 return -1.;
496 } else {
497 const auto e = std::exp(y);
498 return (e - 1) / (1 + e);
499 }
500};

◆ deviator()

template<typename T >
auto PlasticOps::deviator ( FTensor::Tensor2_symmetric< T, SPACE_DIM > &  t_stress,
double  trace 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 373 of file PlasticOps.hpp.

374 {
376 t_dev(I, J) = 0;
377 for (int ii = 0; ii != SPACE_DIM; ++ii)
378 for (int jj = ii; jj != SPACE_DIM; ++jj)
379 t_dev(ii, jj) = t_stress(ii, jj);
380 t_dev(0, 0) -= trace;
381 t_dev(1, 1) -= trace;
382 t_dev(2, 2) -= trace;
383 return t_dev;
384};
static Index< 'J', 3 > J
constexpr int SPACE_DIM
double trace(FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress)
Definition: PlasticOps.hpp:364
constexpr IntegrationType I

◆ diff_constrain_ddot_tau()

double PlasticOps::diff_constrain_ddot_tau ( double  sign,
double  eqiv 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 536 of file PlasticOps.hpp.

536 {
537 return visH + (sigmaY / 2) * (cn0 + cn1 * std::sqrt(eqiv) * (1 - sign));
538};

◆ diff_constrain_df()

auto PlasticOps::diff_constrain_df ( double  sign)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 545 of file PlasticOps.hpp.

545{ return (-1 - sign) / 2; };

◆ diff_constrain_dsigma_y()

auto PlasticOps::diff_constrain_dsigma_y ( double  sign)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 547 of file PlasticOps.hpp.

547{ return (1 + sign) / 2; }

◆ diff_constrain_dstrain()

template<typename T1 , typename T2 >
auto PlasticOps::diff_constrain_dstrain ( T1 &  t_D,
T2 &&  t_diff_constrain_dstress 
)
inline
Examples
PlasticOps.hpp.

Definition at line 559 of file PlasticOps.hpp.

559 {
560 FTensor::Tensor2_symmetric<double, SPACE_DIM> t_diff_constrain_dstrain;
561 t_diff_constrain_dstrain(k, l) =
562 t_diff_constrain_dstress(i, j) * t_D(i, j, k, l);
563 return t_diff_constrain_dstrain;
564};
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k

◆ diff_constrain_dstress()

template<typename T >
auto PlasticOps::diff_constrain_dstress ( double  diff_constrain_df,
FTensor::Tensor2_symmetric< T, SPACE_DIM > &  t_plastic_flow 
)
inline
Examples
PlasticOps.hpp.

Definition at line 550 of file PlasticOps.hpp.

552 {
553 FTensor::Tensor2_symmetric<double, SPACE_DIM> t_diff_constrain_dstress;
554 t_diff_constrain_dstress(i, j) = diff_constrain_df * t_plastic_flow(i, j);
555 return t_diff_constrain_dstress;
556};

◆ diff_constrain_eqiv()

double PlasticOps::diff_constrain_eqiv ( double  sign,
double  eqiv,
double  dot_tau 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 540 of file PlasticOps.hpp.

540 {
541 return (sigmaY / 2) *
542 (-cn0 + cn1 * dot_tau * (0.5 / std::sqrt(eqiv)) * (1 - sign));
543};

◆ diff_deviator()

auto PlasticOps::diff_deviator ( FTensor::Ddg< double, SPACE_DIM, SPACE_DIM > &&  t_diff_stress)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 387 of file PlasticOps.hpp.

387 {
389 t_diff_deviator(I, J, k, l) = 0;
390 for (int ii = 0; ii != SPACE_DIM; ++ii)
391 for (int jj = ii; jj != SPACE_DIM; ++jj)
392 for (int kk = 0; kk != SPACE_DIM; ++kk)
393 for (int ll = kk; ll != SPACE_DIM; ++ll)
394 t_diff_deviator(ii, jj, kk, ll) = t_diff_stress(ii, jj, kk, ll);
395
396 constexpr double third = boost::math::constants::third<double>();
397
398 t_diff_deviator(0, 0, 0, 0) -= third;
399 t_diff_deviator(0, 0, 1, 1) -= third;
400
401 t_diff_deviator(1, 1, 0, 0) -= third;
402 t_diff_deviator(1, 1, 1, 1) -= third;
403
404 t_diff_deviator(2, 2, 0, 0) -= third;
405 t_diff_deviator(2, 2, 1, 1) -= third;
406
407 if constexpr (SPACE_DIM == 3) {
408 t_diff_deviator(0, 0, 2, 2) -= third;
409 t_diff_deviator(1, 1, 2, 2) -= third;
410 t_diff_deviator(2, 2, 2, 2) -= third;
411 }
412
413 return t_diff_deviator;
414};
constexpr double third

◆ diff_equivalent_strain_dot()

template<typename T1 , typename T2 , typename T3 >
auto PlasticOps::diff_equivalent_strain_dot ( const T1  eqiv,
T2 &  t_plastic_strain_dot,
T3 &  t_diff_plastic_strain 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 576 of file PlasticOps.hpp.

577 {
578 constexpr double A = 2. / 3;
580 t_diff_eqiv(i, j) = A * (t_plastic_strain_dot(k, l) / eqiv) *
581 t_diff_plastic_strain(k, l, i, j);
582 return t_diff_eqiv;
583};
constexpr AssemblyType A

◆ diff_plastic_flow_dstrain()

template<typename T >
auto PlasticOps::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 
)
inline
Examples
PlasticOps.hpp.

Definition at line 467 of file PlasticOps.hpp.

469 {
471 t_diff_flow(i, j, k, l) =
472 t_diff_plastic_flow_dstress(i, j, m, n) * t_D(m, n, k, l);
473 return t_diff_flow;
474};
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m

◆ diff_plastic_flow_dstress()

template<typename T >
auto PlasticOps::diff_plastic_flow_dstress ( long double  f,
FTensor::Tensor2_symmetric< T, SPACE_DIM > &  t_flow,
FTensor::Ddg< double, 3, SPACE_DIM > &&  t_diff_deviator 
)
inline
Examples
PlasticOps.hpp.

Definition at line 455 of file PlasticOps.hpp.

457 {
459 t_diff_flow(i, j, k, l) =
460 (1.5 * (t_diff_deviator(M, N, i, j) * t_diff_deviator(M, N, k, l) -
461 (2. / 3.) * t_flow(i, j) * t_flow(k, l))) /
462 f;
463 return t_diff_flow;
464};
static Index< 'M', 3 > M
const int N
Definition: speed_test.cpp:3

◆ diff_symmetrize()

auto PlasticOps::diff_symmetrize ( )
inline
Examples
PlasticOps.hpp, PlasticOpsLargeStrains.hpp, and PlasticOpsSmallStrains.hpp.

Definition at line 333 of file PlasticOps.hpp.

333 {
335
336 t_diff(i, j, k, l) = 0;
337 t_diff(0, 0, 0, 0) = 1;
338 t_diff(1, 1, 1, 1) = 1;
339
340 t_diff(1, 0, 1, 0) = 0.5;
341 t_diff(1, 0, 0, 1) = 0.5;
342
343 t_diff(0, 1, 0, 1) = 0.5;
344 t_diff(0, 1, 1, 0) = 0.5;
345
346 if constexpr (SPACE_DIM == 3) {
347 t_diff(2, 2, 2, 2) = 1;
348
349 t_diff(2, 0, 2, 0) = 0.5;
350 t_diff(2, 0, 0, 2) = 0.5;
351 t_diff(0, 2, 0, 2) = 0.5;
352 t_diff(0, 2, 2, 0) = 0.5;
353
354 t_diff(2, 1, 2, 1) = 0.5;
355 t_diff(2, 1, 1, 2) = 0.5;
356 t_diff(1, 2, 1, 2) = 0.5;
357 t_diff(1, 2, 2, 1) = 0.5;
358 }
359
360 return t_diff;
361};

◆ diff_tensor()

auto PlasticOps::diff_tensor ( )
inline

[Operators definitions]

[Lambda functions]

Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 308 of file PlasticOps.hpp.

308 {
311 t_diff(i, j, k, l) = (t_kd(i, k) ^ t_kd(j, l)) / 4.;
312 return t_diff;
313};
Kronecker Delta class symmetric.
constexpr auto t_kd

◆ equivalent_strain_dot()

template<typename T >
auto PlasticOps::equivalent_strain_dot ( FTensor::Tensor2_symmetric< T, SPACE_DIM > &  t_plastic_strain_dot)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 567 of file PlasticOps.hpp.

568 {
569 constexpr double A = 2. / 3;
570 return std::sqrt(A * t_plastic_strain_dot(i, j) *
571 t_plastic_strain_dot(i, j)) +
572 std::numeric_limits<double>::epsilon();
573};

◆ get_mat_scalar_dtensor_sym() [1/2]

auto PlasticOps::get_mat_scalar_dtensor_sym ( MatrixDouble &  mat,
FTensor::Number< 2 >   
)
Examples
PlasticOpsGeneric.hpp.

Definition at line 562 of file PlasticOpsGeneric.hpp.

562 {
564 &mat(0, 0), &mat(0, 1), &mat(0, 2)};
565}

◆ get_mat_scalar_dtensor_sym() [2/2]

auto PlasticOps::get_mat_scalar_dtensor_sym ( MatrixDouble &  mat,
FTensor::Number< 3 >   
)

Definition at line 567 of file PlasticOpsGeneric.hpp.

567 {
569 &mat(0, 0), &mat(0, 1), &mat(0, 2), &mat(0, 3), &mat(0, 4), &mat(0, 5)};
570}

◆ get_mat_tensor_sym_dscalar() [1/2]

static auto PlasticOps::get_mat_tensor_sym_dscalar ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 2 >   
)
inlinestatic
Examples
PlasticOpsGeneric.hpp.

Definition at line 492 of file PlasticOpsGeneric.hpp.

493 {
495 &mat(3 * rr + 0, 0), &mat(3 * rr + 1, 0), &mat(3 * rr + 2, 0)};
496}

◆ get_mat_tensor_sym_dscalar() [2/2]

static auto PlasticOps::get_mat_tensor_sym_dscalar ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 3 >   
)
inlinestatic

Definition at line 498 of file PlasticOpsGeneric.hpp.

499 {
501 &mat(6 * rr + 0, 0), &mat(6 * rr + 1, 0), &mat(6 * rr + 2, 0),
502 &mat(6 * rr + 3, 0), &mat(6 * rr + 4, 0), &mat(6 * rr + 5, 0)};
503}

◆ get_mat_tensor_sym_dtensor_sym() [1/2]

static auto PlasticOps::get_mat_tensor_sym_dtensor_sym ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 2 >   
)
inlinestatic
Examples
PlasticOpsGeneric.hpp.

Definition at line 402 of file PlasticOpsGeneric.hpp.

403 {
405 &mat(3 * rr + 0, 0), &mat(3 * rr + 0, 1), &mat(3 * rr + 0, 2),
406 &mat(3 * rr + 1, 0), &mat(3 * rr + 1, 1), &mat(3 * rr + 1, 2),
407 &mat(3 * rr + 2, 0), &mat(3 * rr + 2, 1), &mat(3 * rr + 2, 2)};
408}

◆ get_mat_tensor_sym_dtensor_sym() [2/2]

static auto PlasticOps::get_mat_tensor_sym_dtensor_sym ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 3 >   
)
inlinestatic

Definition at line 410 of file PlasticOpsGeneric.hpp.

411 {
413 &mat(6 * rr + 0, 0), &mat(6 * rr + 0, 1), &mat(6 * rr + 0, 2),
414 &mat(6 * rr + 0, 3), &mat(6 * rr + 0, 4), &mat(6 * rr + 0, 5),
415 &mat(6 * rr + 1, 0), &mat(6 * rr + 1, 1), &mat(6 * rr + 1, 2),
416 &mat(6 * rr + 1, 3), &mat(6 * rr + 1, 4), &mat(6 * rr + 1, 5),
417 &mat(6 * rr + 2, 0), &mat(6 * rr + 2, 1), &mat(6 * rr + 2, 2),
418 &mat(6 * rr + 2, 3), &mat(6 * rr + 2, 4), &mat(6 * rr + 2, 5),
419 &mat(6 * rr + 3, 0), &mat(6 * rr + 3, 1), &mat(6 * rr + 3, 2),
420 &mat(6 * rr + 3, 3), &mat(6 * rr + 3, 4), &mat(6 * rr + 3, 5),
421 &mat(6 * rr + 4, 0), &mat(6 * rr + 4, 1), &mat(6 * rr + 4, 2),
422 &mat(6 * rr + 4, 3), &mat(6 * rr + 4, 4), &mat(6 * rr + 4, 5),
423 &mat(6 * rr + 5, 0), &mat(6 * rr + 5, 1), &mat(6 * rr + 5, 2),
424 &mat(6 * rr + 5, 3), &mat(6 * rr + 5, 4), &mat(6 * rr + 5, 5)};
425}

◆ get_mat_tensor_sym_dvector() [1/2]

static auto PlasticOps::get_mat_tensor_sym_dvector ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 2 >   
)
inlinestatic

[Lambda functions]

[Auxiliary functions functions

Examples
PlasticOps.hpp, PlasticOpsLargeStrains.hpp, and PlasticOpsSmallStrains.hpp.

Definition at line 588 of file PlasticOps.hpp.

589 {
591 &mat(3 * rr + 0, 0), &mat(3 * rr + 0, 1), &mat(3 * rr + 1, 0),
592 &mat(3 * rr + 1, 1), &mat(3 * rr + 2, 0), &mat(3 * rr + 2, 1)};
593}

◆ get_mat_tensor_sym_dvector() [2/2]

static auto PlasticOps::get_mat_tensor_sym_dvector ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 3 >   
)
inlinestatic

Definition at line 595 of file PlasticOps.hpp.

596 {
598 &mat(6 * rr + 0, 0), &mat(6 * rr + 0, 1), &mat(6 * rr + 0, 2),
599 &mat(6 * rr + 1, 0), &mat(6 * rr + 1, 1), &mat(6 * rr + 1, 2),
600 &mat(6 * rr + 2, 0), &mat(6 * rr + 2, 1), &mat(6 * rr + 2, 2),
601 &mat(6 * rr + 3, 0), &mat(6 * rr + 3, 1), &mat(6 * rr + 3, 2),
602 &mat(6 * rr + 4, 0), &mat(6 * rr + 4, 1), &mat(6 * rr + 4, 2),
603 &mat(6 * rr + 5, 0), &mat(6 * rr + 5, 1), &mat(6 * rr + 5, 2)};
604}

◆ get_mat_vector_dtensor_sym() [1/2]

static FTensor::Tensor2< FTensor::PackPtr< double *, 3 >, 2, 3 > PlasticOps::get_mat_vector_dtensor_sym ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 2 >   
)
inlinestatic
Examples
PlasticOpsLargeStrains.hpp, and PlasticOpsSmallStrains.hpp.

Definition at line 20 of file PlasticOpsSmallStrains.hpp.

20 {
22 &mat(2 * rr + 0, 0), &mat(2 * rr + 0, 1), &mat(2 * rr + 0, 2),
23 &mat(2 * rr + 1, 0), &mat(2 * rr + 1, 1), &mat(2 * rr + 1, 2)};
24}

◆ get_mat_vector_dtensor_sym() [2/2]

static FTensor::Tensor2< FTensor::PackPtr< double *, 6 >, 3, 6 > PlasticOps::get_mat_vector_dtensor_sym ( size_t  rr,
MatrixDouble &  mat,
FTensor::Number< 3 >   
)
inlinestatic

Definition at line 27 of file PlasticOpsSmallStrains.hpp.

27 {
29 &mat(3 * rr + 0, 0), &mat(3 * rr + 0, 1), &mat(3 * rr + 0, 2),
30 &mat(3 * rr + 0, 3), &mat(3 * rr + 0, 4), &mat(3 * rr + 0, 5),
31
32 &mat(3 * rr + 1, 0), &mat(3 * rr + 1, 1), &mat(3 * rr + 1, 2),
33 &mat(3 * rr + 1, 3), &mat(3 * rr + 1, 4), &mat(3 * rr + 1, 5),
34
35 &mat(3 * rr + 2, 0), &mat(3 * rr + 2, 1), &mat(3 * rr + 2, 2),
36 &mat(3 * rr + 2, 3), &mat(3 * rr + 2, 4), &mat(3 * rr + 2, 5)};
37}

◆ plastic_flow()

auto PlasticOps::plastic_flow ( long double  f,
FTensor::Tensor2_symmetric< double, 3 > &&  t_dev_stress,
FTensor::Ddg< double, 3, SPACE_DIM > &&  t_diff_deviator 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 445 of file PlasticOps.hpp.

447 {
449 t_diff_f(k, l) =
450 (1.5 * (t_dev_stress(I, J) * t_diff_deviator(I, J, k, l))) / f;
451 return t_diff_f;
452};

◆ platsic_surface()

double PlasticOps::platsic_surface ( FTensor::Tensor2_symmetric< double, 3 > &&  t_stress_deviator)
inline

\[ \begin{split} f&=\sqrt{s_{ij}s_{ij}}\\ A_{ij}&=\frac{\partial f}{\partial \sigma_{ij}}= \frac{1}{f} s_{kl} \frac{\partial s_{kl}}{\partial \sigma_{ij}}\\ \frac{\partial A_{ij}}{\partial \sigma_{kl}}&= \frac{\partial^2 f}{\partial \sigma_{ij}\partial\sigma_{mn}}= \frac{1}{f} \left( \frac{\partial s_{kl}}{\partial \sigma_{mn}}\frac{\partial s_{kl}}{\partial \sigma_{ij}} -A_{mn}A_{ij} \right)\\ \frac{\partial f}{\partial \varepsilon_{ij}}&=A_{mn}D_{mnij} \\ \frac{\partial A_{ij}}{\partial \varepsilon_{kl}}&= \frac{\partial A_{ij}}{\partial \sigma_{mn}} \frac{\partial \sigma_{mn}}{\partial \varepsilon_{kl}}= \frac{\partial A_{ij}}{\partial \sigma_{mn}} D_{mnkl} \end{split} \]

Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 440 of file PlasticOps.hpp.

440 {
441 return std::sqrt(1.5 * t_stress_deviator(I, J) * t_stress_deviator(I, J)) +
442 std::numeric_limits<double>::epsilon();
443};

◆ symm_L_tensor()

auto PlasticOps::symm_L_tensor ( )
inline
Examples
PlasticOps.hpp, PlasticOpsGeneric.hpp, PlasticOpsLargeStrains.hpp, and PlasticOpsSmallStrains.hpp.

Definition at line 315 of file PlasticOps.hpp.

315 {
317 t_L(i, j, L) = 0;
318 if constexpr (SPACE_DIM == 2) {
319 t_L(0, 0, 0) = 1;
320 t_L(1, 0, 1) = 1;
321 t_L(1, 1, 2) = 1;
322 } else {
323 t_L(0, 0, 0) = 1;
324 t_L(1, 0, 1) = 1;
325 t_L(2, 0, 2) = 1;
326 t_L(1, 1, 3) = 1;
327 t_L(2, 1, 4) = 1;
328 t_L(2, 2, 5) = 1;
329 }
330 return t_L;
331}
static Index< 'L', 3 > L

◆ trace()

template<typename T >
double PlasticOps::trace ( FTensor::Tensor2_symmetric< T, SPACE_DIM > &  t_stress)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 364 of file PlasticOps.hpp.

364 {
365 constexpr double third = boost::math::constants::third<double>();
366 if constexpr (SPACE_DIM == 2)
367 return (t_stress(0, 0) + t_stress(1, 1)) * third;
368 else
369 return (t_stress(0, 0) + t_stress(1, 1) + t_stress(2, 2)) * third;
370};

◆ w()

double PlasticOps::w ( double  eqiv,
double  dot_tau,
double  f,
double  sigma_y 
)
inline
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 513 of file PlasticOps.hpp.

513 {
514 return (f - sigma_y) / sigmaY + cn1 * (dot_tau * std::sqrt(eqiv));
515};

Variable Documentation

◆ i

FTensor::Index<'i', SPACE_DIM> PlasticOps::i

◆ I

FTensor::Index<'I', 3> PlasticOps::I
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 101 of file PlasticOps.hpp.

◆ j

FTensor::Index<'j', SPACE_DIM> PlasticOps::j

◆ J

FTensor::Index<'J', 3> PlasticOps::J
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 102 of file PlasticOps.hpp.

◆ k

FTensor::Index<'k', SPACE_DIM> PlasticOps::k

◆ l

FTensor::Index<'l', SPACE_DIM> PlasticOps::l

◆ L

FTensor::Index<'L', size_symm> PlasticOps::L

◆ m

FTensor::Index<'m', SPACE_DIM> PlasticOps::m

◆ M

FTensor::Index<'M', 3> PlasticOps::M
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 103 of file PlasticOps.hpp.

◆ n

FTensor::Index<'n', SPACE_DIM> PlasticOps::n

◆ N

FTensor::Index<'N', 3> PlasticOps::N
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 104 of file PlasticOps.hpp.

◆ o

FTensor::Index<'o', SPACE_DIM> PlasticOps::o
Examples
PlasticOps.hpp.

Definition at line 98 of file PlasticOps.hpp.

◆ O

FTensor::Index<'O', size_symm> PlasticOps::O
Examples
PlasticOps.hpp, and PlasticOpsGeneric.hpp.

Definition at line 107 of file PlasticOps.hpp.

◆ p

FTensor::Index<'p', SPACE_DIM> PlasticOps::p
Examples
PlasticOps.hpp.

Definition at line 99 of file PlasticOps.hpp.