30 return A * std::pow(1.5, 0.5 * (1 + n)) / std::pow(sigma0, n);
33 template <
int DisplacementDim>
34 std::optional<std::tuple<typename CreepBGRa<DisplacementDim>::KelvinVector,
36 DisplacementDim>::MaterialStateVariables>,
45 auto const& eps_m = std::get<MPL::SymmetricTensor<DisplacementDim>>(
46 variable_array[
static_cast<int>(MPL::Variable::mechanical_strain)]);
47 auto const& eps_m_prev = std::get<MPL::SymmetricTensor<DisplacementDim>>(
48 variable_array_prev[
static_cast<int>(
49 MPL::Variable::mechanical_strain)]);
50 auto const& sigma_prev = std::get<MPL::SymmetricTensor<DisplacementDim>>(
51 variable_array_prev[
static_cast<int>(MPL::Variable::stress)]);
52 auto const T = std::get<double>(
57 Eigen::FullPivLU<Eigen::Matrix<double, KelvinVectorSize, KelvinVectorSize,
61 const auto C = this->getElasticTensor(t, x, T);
62 KelvinVector sigma_try = sigma_prev + C * (eps_m - eps_m_prev);
64 auto const& deviatoric_matrix = Invariants::deviatoric_projection;
66 double const norm_s_try =
67 Invariants::FrobeniusNorm(deviatoric_matrix * sigma_try);
69 if (norm_s_try < std::numeric_limits<double>::epsilon() * C(0, 0))
71 return {std::make_tuple(sigma_try, createMaterialStateVariables(), C)};
76 const double A = _a(t, x)[0];
77 const double n = _n(t, x)[0];
78 const double sigma0 = _sigma_f(t, x)[0];
79 const double Q = _q(t, x)[0];
81 const double constant_coefficient =
85 dt * constant_coefficient *
88 double const G2b = 2.0 * b * this->_mp.mu(t, x);
90 auto const update_jacobian = [&solution, &G2b, &n](
JacobianMatrix& jacobian)
92 auto const& D = Invariants::deviatoric_projection;
94 double const norm_s_n1 = Invariants::FrobeniusNorm(s_n1);
95 double const pow_norm_s_n1_n_minus_one_2b_G =
96 G2b * std::pow(norm_s_n1, n - 1);
97 jacobian = KelvinMatrix::Identity() +
98 (pow_norm_s_n1_n_minus_one_2b_G * D +
99 (n - 1) * G2b * std::pow(norm_s_n1, n - 3) * s_n1 *
103 auto const update_residual =
106 KelvinVector const s_n1 = Invariants::deviatoric_projection * solution;
107 double const norm_s_n1 = Invariants::FrobeniusNorm(s_n1);
108 double const pow_norm_s_n1_n_minus_one_2b_G =
109 G2b * std::pow(norm_s_n1, n - 1);
110 r = solution - sigma_try + pow_norm_s_n1_n_minus_one_2b_G * s_n1;
114 { solution += increment; };
119 decltype(update_residual),
120 decltype(update_solution)>(
121 linear_solver, update_jacobian, update_residual, update_solution,
122 _nonlinear_solver_parameters);
125 auto const success_iterations = newton_solver.
solve(jacobian);
127 if (!success_iterations)
136 (*success_iterations == 0) ? C : linear_solver.solve(C);
138 return {std::make_tuple(solution, createMaterialStateVariables(),
142 template <
int DisplacementDim>
145 double const T,
double const deviatoric_stress_norm)
const
147 const double A = _a(t, x)[0];
148 const double n = _n(t, x)[0];
149 const double sigma0 = _sigma_f(t, x)[0];
150 const double Q = _q(t, x)[0];
152 const double constant_coefficient =
155 return 2.0 * constant_coefficient *
158 this->_mp.mu(t, x) * std::pow(deviatoric_stress_norm, n - 1) * dt *
A class for computing the BGRa creep model, which reads.
Eigen::Matrix< double, KelvinVectorSize, 1 > ResidualVectorType
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > KelvinMatrix
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > KelvinVector
double getTemperatureRelatedCoefficient(double const t, double const dt, ParameterLib::SpatialPosition const &x, double const T, double const deviatoric_stress_norm) const override
std::optional< std::tuple< KelvinVector, std::unique_ptr< typename MechanicsBase< DisplacementDim >::MaterialStateVariables >, KelvinMatrix > > integrateStress(MaterialPropertyLib::VariableArray const &variable_array_prev, MaterialPropertyLib::VariableArray const &variable_array, double const t, ParameterLib::SpatialPosition const &x, double const dt, typename MechanicsBase< DisplacementDim >::MaterialStateVariables const &material_state_variables) const override
Eigen::Matrix< double, KelvinVectorSize, KelvinVectorSize, Eigen::RowMajor > JacobianMatrix
std::optional< int > solve(JacobianMatrix &jacobian) const
constexpr double IdealGasConstant
double getCreepConstantCoefficient(const double A, const double n, const double sigma0)
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray