OGS
MaterialLib::Solids::Ehlers Namespace Reference

Classes

class  Damage
struct  DamageProperties
struct  DamagePropertiesParameters
struct  MaterialProperties
struct  MaterialPropertiesParameters
struct  OnePlusGamma_pTheta
 Holds powers of 1 + gamma_p*theta to base 0, m_p, and m_p-1. More...
struct  PhysicalStressWithInvariants
struct  PlasticStrain
class  SolidEhlers
struct  StateVariables

Enumerations

enum class  TangentType { Elastic , PlasticDamageSecant , Plastic }

Functions

std::unique_ptr< DamagePropertiesParameterscreateDamageProperties (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config)
template<int DisplacementDim>
std::unique_ptr< SolidEhlers< DisplacementDim > > createEhlers (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config)
template<int DisplacementDim>
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > sOdotS (MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &v)
template<int DisplacementDim>
double plasticFlowVolumetricPart (PhysicalStressWithInvariants< DisplacementDim > const &s, double const sqrtPhi, double const alpha_p, double const beta_p, double const delta_p, double const epsilon_p)
template<int DisplacementDim>
SolidEhlers< DisplacementDim >::KelvinVector plasticFlowDeviatoricPart (PhysicalStressWithInvariants< DisplacementDim > const &s, OnePlusGamma_pTheta const &one_gt, double const sqrtPhi, typename SolidEhlers< DisplacementDim >::KelvinVector const &dtheta_dsigma, double const gamma_p, double const m_p)
template<int DisplacementDim>
double yieldFunction (MaterialProperties const &mp, PhysicalStressWithInvariants< DisplacementDim > const &s, double const k)
template<int DisplacementDim>
std::pair< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > thetaSigmaDerivatives (double theta, PhysicalStressWithInvariants< DisplacementDim > const &s)
template<int DisplacementDim>
SolidEhlers< DisplacementDim >::ResidualVectorType calculatePlasticResidual (MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps_D, double const eps_V, PhysicalStressWithInvariants< DisplacementDim > const &s, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps_p_D, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps_p_D_dot, double const eps_p_V, double const eps_p_V_dot, double const eps_p_eff_dot, double const lambda, double const k, MaterialProperties const &mp)
template<int DisplacementDim>
SolidEhlers< DisplacementDim >::JacobianMatrix calculatePlasticJacobian (double const dt, PhysicalStressWithInvariants< DisplacementDim > const &s, double const lambda, MaterialProperties const &mp)
template<int DisplacementDim>
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > calculateDResidualDEps (double const K, double const G)
double calculateIsotropicHardening (double const kappa, double const hardening_coefficient, double const eps_p_eff)
template<int DisplacementDim>
SolidEhlers< DisplacementDim >::KelvinVector predict_sigma (double const G, double const K, typename SolidEhlers< DisplacementDim >::KelvinVector const &sigma_prev, typename SolidEhlers< DisplacementDim >::KelvinVector const &eps, typename SolidEhlers< DisplacementDim >::KelvinVector const &eps_prev, double const eps_V)
template<typename ResidualVector, typename KelvinVector>
std::tuple< KelvinVector, PlasticStrain< KelvinVector >, double > splitSolutionVector (ResidualVector const &solution)
template<>
MathLib::KelvinVector::KelvinMatrixType< 3 > sOdotS< 3 > (MathLib::KelvinVector::KelvinVectorType< 3 > const &v)
template<>
MathLib::KelvinVector::KelvinMatrixType< 2 > sOdotS< 2 > (MathLib::KelvinVector::KelvinVectorType< 2 > const &v)
TangentType makeTangentType (std::string const &s)

Enumeration Type Documentation

◆ TangentType

Enumerator
Elastic 
PlasticDamageSecant 
Plastic 

Definition at line 33 of file Ehlers.h.

Function Documentation

◆ calculateDResidualDEps()

template<int DisplacementDim>
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > MaterialLib::Solids::Ehlers::calculateDResidualDEps ( double const K,
double const G )

Calculates the derivative of the residuals with respect to total strain. Implementation fully implicit only.

Definition at line 432 of file Ehlers.cpp.

434{
435 static int const KelvinVectorSize =
438
439 auto const& P_dev = Invariants::deviatoric_projection;
440 auto const& P_sph = Invariants::spherical_projection;
441 auto const& I =
443
444 return -2. * I * P_dev - 3. * K / G * I * P_sph;
445}
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType

References MathLib::KelvinVector::kelvin_vector_dimensions().

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ calculateIsotropicHardening()

double MaterialLib::Solids::Ehlers::calculateIsotropicHardening ( double const kappa,
double const hardening_coefficient,
double const eps_p_eff )
inline

Definition at line 447 of file Ehlers.cpp.

450{
451 return kappa * (1. + eps_p_eff * hardening_coefficient);
452}

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ calculatePlasticJacobian()

template<int DisplacementDim>
SolidEhlers< DisplacementDim >::JacobianMatrix MaterialLib::Solids::Ehlers::calculatePlasticJacobian ( double const dt,
PhysicalStressWithInvariants< DisplacementDim > const & s,
double const lambda,
MaterialProperties const & mp )

Definition at line 251 of file Ehlers.cpp.

256{
257 static int const KelvinVectorSize =
260 using KelvinVector =
262 using KelvinMatrix =
264
265 auto const& P_dev = Invariants::deviatoric_projection;
266 auto const& identity2 = Invariants::identity2;
267
268 double const theta = s.J_3 / (s.J_2 * std::sqrt(s.J_2));
269 OnePlusGamma_pTheta const one_gt{mp.gamma_p, theta, mp.m_p};
270
271 auto const [dtheta_dsigma, d2theta_dsigma2] =
272 thetaSigmaDerivatives(theta, s);
273
274 // deviatoric flow
275 double const sqrtPhi = std::sqrt(
276 s.J_2 * one_gt.pow_m_p + mp.alpha_p / 2. * boost::math::pow<2>(s.I_1) +
277 boost::math::pow<2>(mp.delta_p) * boost::math::pow<4>(s.I_1));
278 KelvinVector const flow_D = plasticFlowDeviatoricPart(
279 s, one_gt, sqrtPhi, dtheta_dsigma, mp.gamma_p, mp.m_p);
280 KelvinVector const lambda_flow_D = lambda * flow_D;
281
283 SolidEhlers<DisplacementDim>::JacobianMatrix::Zero();
284
285 // G_11
286 jacobian.template block<KelvinVectorSize, KelvinVectorSize>(0, 0)
287 .noalias() = KelvinMatrix::Identity();
288
289 // G_12
290 jacobian
291 .template block<KelvinVectorSize, KelvinVectorSize>(0, KelvinVectorSize)
292 .noalias() = 2 * KelvinMatrix::Identity();
293
294 // G_13
295 jacobian.template block<KelvinVectorSize, 1>(0, 2 * KelvinVectorSize)
296 .noalias() = mp.K / mp.G * identity2;
297
298 // G_14 and G_15 are zero
299
300 // G_21 -- derivative of deviatoric flow
301
302 double const gm_p = mp.gamma_p * mp.m_p;
303 // intermediate variable for derivative of deviatoric flow
304 KelvinVector const M0 = s.J_2 / one_gt.value * dtheta_dsigma;
305 // derivative of Phi w.r.t. sigma
306 KelvinVector const dPhi_dsigma =
307 one_gt.pow_m_p * (s.D + gm_p * M0) +
308 (mp.alpha_p * s.I_1 +
309 4 * boost::math::pow<2>(mp.delta_p) * boost::math::pow<3>(s.I_1)) *
310 identity2;
311
312 // intermediate variable for derivative of deviatoric flow
313 KelvinMatrix const M1 =
314 one_gt.pow_m_p *
315 (s.D * dPhi_dsigma.transpose() + gm_p * M0 * dPhi_dsigma.transpose());
316 // intermediate variable for derivative of deviatoric flow
317 KelvinMatrix const M2 =
318 one_gt.pow_m_p * (P_dev + s.D * gm_p * M0.transpose());
319
320 // intermediate variable for derivative of deviatoric flow
321 KelvinMatrix const M3 =
322 gm_p * one_gt.pow_m_p1 *
323 ((s.D + (gm_p - mp.gamma_p) * M0) * dtheta_dsigma.transpose() +
324 s.J_2 * d2theta_dsigma2);
325
326 // derivative of flow_D w.r.t. sigma
327 KelvinMatrix const dflow_D_dsigma =
328 (-M1 / (4 * boost::math::pow<3>(sqrtPhi)) + (M2 + M3) / (2 * sqrtPhi)) *
329 mp.G;
330 jacobian
331 .template block<KelvinVectorSize, KelvinVectorSize>(KelvinVectorSize, 0)
332 .noalias() = -lambda * dflow_D_dsigma;
333
334 // G_22
335 jacobian
336 .template block<KelvinVectorSize, KelvinVectorSize>(KelvinVectorSize,
337 KelvinVectorSize)
338 .noalias() = KelvinMatrix::Identity() / dt;
339
340 // G_23 and G_24 are zero
341
342 // G_25
343 jacobian
344 .template block<KelvinVectorSize, 1>(KelvinVectorSize,
345 2 * KelvinVectorSize + 2)
346 .noalias() = -flow_D;
347
348 // G_31
349 {
350 // derivative of flow_V w.r.t. sigma
351 KelvinVector const dflow_V_dsigma =
352 3 * mp.G *
353 (-(mp.alpha_p * s.I_1 + 4 * boost::math::pow<2>(mp.delta_p) *
354 boost::math::pow<3>(s.I_1)) /
355 (4 * boost::math::pow<3>(sqrtPhi)) * dPhi_dsigma +
356 (mp.alpha_p * identity2 +
357 12 * boost::math::pow<2>(mp.delta_p * s.I_1) * identity2) /
358 (2 * sqrtPhi) +
359 2 * mp.epsilon_p * identity2);
360
361 jacobian.template block<1, KelvinVectorSize>(2 * KelvinVectorSize, 0)
362 .noalias() = -lambda * dflow_V_dsigma.transpose();
363 }
364
365 // G_32 is zero
366
367 // G_33
368 jacobian(2 * KelvinVectorSize, 2 * KelvinVectorSize) = 1. / dt;
369
370 // G_34 is zero
371
372 // G_35
373 {
374 double const flow_V = plasticFlowVolumetricPart<DisplacementDim>(
375 s, sqrtPhi, mp.alpha_p, mp.beta_p, mp.delta_p, mp.epsilon_p);
376 jacobian(2 * KelvinVectorSize, 2 * KelvinVectorSize + 2) = -flow_V;
377 }
378
379 // increment of effectiv plastic strain
380 double const eff_flow =
381 std::sqrt(2. / 3. * lambda_flow_D.transpose() * lambda_flow_D);
382
383 if (eff_flow > 0)
384 {
385 // intermediate variable for derivative of plastic jacobian
386 KelvinVector const eff_flow23_lambda_flow_D =
387 -2 / 3. / eff_flow * lambda_flow_D;
388 // G_41
389 jacobian
390 .template block<1, KelvinVectorSize>(2 * KelvinVectorSize + 1, 0)
391 .noalias() = lambda * dflow_D_dsigma * eff_flow23_lambda_flow_D;
392 // G_45
393 jacobian(2 * KelvinVectorSize + 1, 2 * KelvinVectorSize + 2) =
394 eff_flow23_lambda_flow_D.transpose() * flow_D;
395 }
396
397 // G_42 and G_43 are zero
398
399 // G_44
400 jacobian(2 * KelvinVectorSize + 1, 2 * KelvinVectorSize + 1) = 1. / dt;
401
402 // G_51
403 {
404 double const one_gt_pow_m = std::pow(one_gt.value, mp.m);
405 double const gm = mp.gamma * mp.m;
406 // derivative of yield function w.r.t. sigma
407 KelvinVector const dF_dsigma =
408 mp.G *
409 (one_gt_pow_m * (s.D + gm * M0) +
410 (mp.alpha * s.I_1 + 4 * boost::math::pow<2>(mp.delta) *
411 boost::math::pow<3>(s.I_1)) *
412 identity2) /
413 (2. * sqrtPhi) +
414 mp.G * (mp.beta + 2 * mp.epsilon_p * s.I_1) * identity2;
415
416 jacobian
417 .template block<1, KelvinVectorSize>(2 * KelvinVectorSize + 2, 0)
418 .noalias() = dF_dsigma.transpose() / mp.G;
419 }
420
421 // G_54
422 jacobian(2 * KelvinVectorSize + 2, 2 * KelvinVectorSize + 1) =
423 -mp.kappa * mp.hardening_coefficient / mp.G;
424
425 // G_52, G_53, G_55 are zero
426 return jacobian;
427}
Eigen::Matrix< double, JacobianResidualSize, JacobianResidualSize, Eigen::RowMajor > JacobianMatrix
Definition Ehlers.h:279
SolidEhlers< DisplacementDim >::KelvinVector plasticFlowDeviatoricPart(PhysicalStressWithInvariants< DisplacementDim > const &s, OnePlusGamma_pTheta const &one_gt, double const sqrtPhi, typename SolidEhlers< DisplacementDim >::KelvinVector const &dtheta_dsigma, double const gamma_p, double const m_p)
Definition Ehlers.cpp:118
double plasticFlowVolumetricPart(PhysicalStressWithInvariants< DisplacementDim > const &s, double const sqrtPhi, double const alpha_p, double const beta_p, double const delta_p, double const epsilon_p)
Definition Ehlers.cpp:105
std::pair< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > thetaSigmaDerivatives(double theta, PhysicalStressWithInvariants< DisplacementDim > const &s)
Definition Ehlers.cpp:149
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
KV::KelvinVectorType< DisplacementDim > KelvinVector
Holds powers of 1 + gamma_p*theta to base 0, m_p, and m_p-1.
Definition Ehlers.cpp:90

References MaterialLib::Solids::Ehlers::MaterialProperties::alpha, MaterialLib::Solids::Ehlers::MaterialProperties::alpha_p, MaterialLib::Solids::Ehlers::MaterialProperties::beta, MaterialLib::Solids::Ehlers::MaterialProperties::beta_p, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::D, MaterialLib::Solids::Ehlers::MaterialProperties::delta, MaterialLib::Solids::Ehlers::MaterialProperties::delta_p, MaterialLib::Solids::Ehlers::MaterialProperties::epsilon_p, MaterialLib::Solids::Ehlers::MaterialProperties::G, MaterialLib::Solids::Ehlers::MaterialProperties::gamma, MaterialLib::Solids::Ehlers::MaterialProperties::gamma_p, MaterialLib::Solids::Ehlers::MaterialProperties::hardening_coefficient, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::I_1, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::J_2, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::J_3, MaterialLib::Solids::Ehlers::MaterialProperties::K, MaterialLib::Solids::Ehlers::MaterialProperties::kappa, MathLib::KelvinVector::kelvin_vector_dimensions(), MaterialLib::Solids::Ehlers::MaterialProperties::m, MaterialLib::Solids::Ehlers::MaterialProperties::m_p, plasticFlowDeviatoricPart(), plasticFlowVolumetricPart(), MaterialLib::Solids::Ehlers::OnePlusGamma_pTheta::pow_m_p, MaterialLib::Solids::Ehlers::OnePlusGamma_pTheta::pow_m_p1, thetaSigmaDerivatives(), and MaterialLib::Solids::Ehlers::OnePlusGamma_pTheta::value.

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ calculatePlasticResidual()

template<int DisplacementDim>
SolidEhlers< DisplacementDim >::ResidualVectorType MaterialLib::Solids::Ehlers::calculatePlasticResidual ( MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const & eps_D,
double const eps_V,
PhysicalStressWithInvariants< DisplacementDim > const & s,
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const & eps_p_D,
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const & eps_p_D_dot,
double const eps_p_V,
double const eps_p_V_dot,
double const eps_p_eff_dot,
double const lambda,
double const k,
MaterialProperties const & mp )

Definition at line 190 of file Ehlers.cpp.

202{
203 static int const KelvinVectorSize =
206 using KelvinVector =
208
209 auto const& identity2 = Invariants::identity2;
210
211 double const theta = s.J_3 / (s.J_2 * std::sqrt(s.J_2));
212
214 // calculate stress residual
215 residual.template segment<KelvinVectorSize>(0).noalias() =
216 s.value / mp.G - 2 * (eps_D - eps_p_D) -
217 mp.K / mp.G * (eps_V - eps_p_V) * identity2;
218
219 auto const [dtheta_dsigma, d2theta_dsigma2] =
220 thetaSigmaDerivatives(theta, s);
221
222 OnePlusGamma_pTheta const one_gt{mp.gamma_p, theta, mp.m_p};
223 double const sqrtPhi = std::sqrt(
224 s.J_2 * one_gt.pow_m_p + mp.alpha_p / 2. * boost::math::pow<2>(s.I_1) +
225 boost::math::pow<2>(mp.delta_p) * boost::math::pow<4>(s.I_1));
226 KelvinVector const flow_D = plasticFlowDeviatoricPart(
227 s, one_gt, sqrtPhi, dtheta_dsigma, mp.gamma_p, mp.m_p);
228 KelvinVector const lambda_flow_D = lambda * flow_D;
229
230 residual.template segment<KelvinVectorSize>(KelvinVectorSize).noalias() =
231 eps_p_D_dot - lambda_flow_D;
232
233 // plastic volume strain
234 {
235 double const flow_V = plasticFlowVolumetricPart<DisplacementDim>(
236 s, sqrtPhi, mp.alpha_p, mp.beta_p, mp.delta_p, mp.epsilon_p);
237 residual(2 * KelvinVectorSize, 0) = eps_p_V_dot - lambda * flow_V;
238 }
239
240 // evolution of plastic equivalent strain
241 residual(2 * KelvinVectorSize + 1) =
242 eps_p_eff_dot -
243 std::sqrt(2. / 3. * lambda_flow_D.transpose() * lambda_flow_D);
244
245 // yield function (for plastic multiplier)
246 residual(2 * KelvinVectorSize + 2) = yieldFunction(mp, s, k) / mp.G;
247 return residual;
248}
Eigen::Matrix< double, JacobianResidualSize, 1 > ResidualVectorType
Definition Ehlers.h:278
double yieldFunction(MaterialProperties const &mp, PhysicalStressWithInvariants< DisplacementDim > const &s, double const k)
Definition Ehlers.cpp:130

References MaterialLib::Solids::Ehlers::MaterialProperties::alpha_p, MaterialLib::Solids::Ehlers::MaterialProperties::beta_p, MaterialLib::Solids::Ehlers::MaterialProperties::delta_p, MaterialLib::Solids::Ehlers::MaterialProperties::epsilon_p, MaterialLib::Solids::Ehlers::MaterialProperties::G, MaterialLib::Solids::Ehlers::MaterialProperties::gamma_p, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::I_1, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::J_2, MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::J_3, MaterialLib::Solids::Ehlers::MaterialProperties::K, MathLib::KelvinVector::kelvin_vector_dimensions(), MaterialLib::Solids::Ehlers::MaterialProperties::m_p, plasticFlowDeviatoricPart(), plasticFlowVolumetricPart(), MaterialLib::Solids::Ehlers::OnePlusGamma_pTheta::pow_m_p, thetaSigmaDerivatives(), MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::value, and yieldFunction().

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ createDamageProperties()

std::unique_ptr< DamagePropertiesParameters > MaterialLib::Solids::Ehlers::createDamageProperties ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
BaseLib::ConfigTree const & config )
inline
Input File Parameter
material__solid__constitutive_relation__Ehlers__damage_properties__alpha_d
Input File Parameter
material__solid__constitutive_relation__Ehlers__damage_properties__beta_d
Input File Parameter
material__solid__constitutive_relation__Ehlers__damage_properties__h_d

Definition at line 16 of file CreateEhlers.h.

19{
21 auto& alpha_d =
22 ParameterLib::findParameter<double>(config, "alpha_d", parameters, 1);
23
24 DBUG("Use '{:s}' as alpha_d.", alpha_d.name);
25
27 auto& beta_d =
28 ParameterLib::findParameter<double>(config, "beta_d", parameters, 1);
29
30 DBUG("Use '{:s}' as beta_d.", beta_d.name);
31
33 auto& h_d =
34 ParameterLib::findParameter<double>(config, "h_d", parameters, 1);
35
36 DBUG("Use '{:s}' as h_d.", h_d.name);
37
38 return std::make_unique<DamagePropertiesParameters>(
39 DamagePropertiesParameters{alpha_d, beta_d, h_d});
40}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:22
OGS_NO_DANGLING Parameter< ParameterDataType > & findParameter(std::string const &parameter_name, std::vector< std::unique_ptr< ParameterBase > > const &parameters, int const num_components, MeshLib::Mesh const *const mesh=nullptr)

References DBUG(), and ParameterLib::findParameter().

Referenced by createEhlers().

◆ createEhlers()

template<int DisplacementDim>
std::unique_ptr< SolidEhlers< DisplacementDim > > MaterialLib::Solids::Ehlers::createEhlers ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
BaseLib::ConfigTree const & config )
Input File Parameter
material__solid__constitutive_relation__type
Input File Parameter
material__solid__constitutive_relation__Ehlers__shear_modulus
Input File Parameter
material__solid__constitutive_relation__Ehlers__bulk_modulus
Input File Parameter
material__solid__constitutive_relation__Ehlers__kappa
Input File Parameter
material__solid__constitutive_relation__Ehlers__beta
Input File Parameter
material__solid__constitutive_relation__Ehlers__gamma
Input File Parameter
material__solid__constitutive_relation__Ehlers__hardening_modulus
Input File Parameter
material__solid__constitutive_relation__Ehlers__alpha
Input File Parameter
material__solid__constitutive_relation__Ehlers__delta
Input File Parameter
material__solid__constitutive_relation__Ehlers__eps
Input File Parameter
material__solid__constitutive_relation__Ehlers__m
Input File Parameter
material__solid__constitutive_relation__Ehlers__alphap
Input File Parameter
material__solid__constitutive_relation__Ehlers__deltap
Input File Parameter
material__solid__constitutive_relation__Ehlers__epsp
Input File Parameter
material__solid__constitutive_relation__Ehlers__mp
Input File Parameter
material__solid__constitutive_relation__Ehlers__betap
Input File Parameter
material__solid__constitutive_relation__Ehlers__gammap
Input File Parameter
material__solid__constitutive_relation__Ehlers__tangent_type
Input File Parameter
material__solid__constitutive_relation__Ehlers__damage_properties
Input File Parameter
material__solid__constitutive_relation__Ehlers__nonlinear_solver

Definition at line 43 of file CreateEhlers.h.

46{
48 config.checkConfigParameter("type", "Ehlers");
49 DBUG("Create Ehlers material");
50
52 auto& shear_modulus = ParameterLib::findParameter<double>(
53 config, "shear_modulus", parameters, 1);
54
55 DBUG("Use '{:s}' as shear modulus parameter.", shear_modulus.name);
56
58 auto& bulk_modulus = ParameterLib::findParameter<double>(
59 config, "bulk_modulus", parameters, 1);
60
61 DBUG("Use '{:s}' as bulk modulus parameter.", bulk_modulus.name);
62
64 auto& kappa =
65 ParameterLib::findParameter<double>(config, "kappa", parameters, 1);
66
67 DBUG("Use '{:s}' as kappa.", kappa.name);
68
70 auto& beta =
71 ParameterLib::findParameter<double>(config, "beta", parameters, 1);
72
73 DBUG("Use '{:s}' as beta.", beta.name);
74
76 auto& gamma =
77 ParameterLib::findParameter<double>(config, "gamma", parameters, 1);
78
79 DBUG("Use '{:s}' as gamma.", gamma.name);
80
82 auto& hardening_modulus = ParameterLib::findParameter<double>(
83 config, "hardening_modulus", parameters, 1);
84
85 DBUG("Use '{:s}' as hardening modulus parameter.", hardening_modulus.name);
86
88 auto& alpha =
89 ParameterLib::findParameter<double>(config, "alpha", parameters, 1);
90
91 DBUG("Use '{:s}' as alpha.", alpha.name);
92
94 auto& delta =
95 ParameterLib::findParameter<double>(config, "delta", parameters, 1);
96
97 DBUG("Use '{:s}' as delta.", delta.name);
98
100 auto& eps =
101 ParameterLib::findParameter<double>(config, "eps", parameters, 1);
102
103 DBUG("Use '{:s}' as eps.", eps.name);
104
106 auto& m = ParameterLib::findParameter<double>(config, "m", parameters, 1);
107
108 DBUG("Use '{:s}' as m.", m.name);
109
111 auto& alphap =
112 ParameterLib::findParameter<double>(config, "alphap", parameters, 1);
113
114 DBUG("Use '{:s}' as alphap.", alphap.name);
115
117 auto& deltap =
118 ParameterLib::findParameter<double>(config, "deltap", parameters, 1);
119
120 DBUG("Use '{:s}' as deltap.", deltap.name);
121
123 auto& epsp =
124 ParameterLib::findParameter<double>(config, "epsp", parameters, 1);
125
126 DBUG("Use '{:s}' as epsp.", epsp.name);
127
129 auto& paremeter_mp =
130 ParameterLib::findParameter<double>(config, "mp", parameters, 1);
131
132 DBUG("Use '{:s}' as mp.", paremeter_mp.name);
133
135 auto& betap =
136 ParameterLib::findParameter<double>(config, "betap", parameters, 1);
137
138 DBUG("Use '{:s}' as betap.", betap.name);
139
141 auto& gammap =
142 ParameterLib::findParameter<double>(config, "gammap", parameters, 1);
143
144 DBUG("Use '{:s}' as gammap.", gammap.name);
145
146 auto tangent_type =
148 makeTangentType(config.getConfigParameter<std::string>("tangent_type"));
149
151 shear_modulus, bulk_modulus, alpha, beta,
152 gamma, delta, eps, m,
153 alphap, betap, gammap, deltap,
154 epsp, paremeter_mp, kappa, hardening_modulus};
155
156 // Damage properties.
157 std::unique_ptr<DamagePropertiesParameters> ehlers_damage_properties;
158
159 auto const& ehlers_damage_config =
161 config.getConfigSubtreeOptional("damage_properties");
162 if (ehlers_damage_config)
163 {
164 ehlers_damage_properties =
165 createDamageProperties(parameters, *ehlers_damage_config);
166 }
167
168 auto const& nonlinear_solver_config =
170 config.getConfigSubtree("nonlinear_solver");
171 auto const nonlinear_solver_parameters =
172 NumLib::createNewtonRaphsonSolverParameters(nonlinear_solver_config);
173
174 return std::make_unique<SolidEhlers<DisplacementDim>>(
175 nonlinear_solver_parameters,
176 mp,
177 std::move(ehlers_damage_properties),
178 tangent_type);
179}
std::unique_ptr< DamagePropertiesParameters > createDamageProperties(std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config)
TangentType makeTangentType(std::string const &s)
Definition Ehlers.h:39
NewtonRaphsonSolverParameters createNewtonRaphsonSolverParameters(BaseLib::ConfigTree const &config)

References BaseLib::ConfigTree::checkConfigParameter(), createDamageProperties(), NumLib::createNewtonRaphsonSolverParameters(), DBUG(), ParameterLib::findParameter(), BaseLib::ConfigTree::getConfigParameter(), BaseLib::ConfigTree::getConfigSubtree(), BaseLib::ConfigTree::getConfigSubtreeOptional(), and makeTangentType().

Referenced by MaterialLib::Solids::createConstitutiveRelation().

◆ makeTangentType()

TangentType MaterialLib::Solids::Ehlers::makeTangentType ( std::string const & s)
inline

Definition at line 39 of file Ehlers.h.

40{
41 if (s == "Elastic")
42 {
44 }
45 if (s == "PlasticDamageSecant")
46 {
48 }
49 if (s == "Plastic")
50 {
52 }
53 OGS_FATAL("Not valid string '{:s}' to create a tangent type from.", s);
54}
#define OGS_FATAL(...)
Definition Error.h:19

References Elastic, OGS_FATAL, Plastic, and PlasticDamageSecant.

Referenced by createEhlers().

◆ plasticFlowDeviatoricPart()

template<int DisplacementDim>
SolidEhlers< DisplacementDim >::KelvinVector MaterialLib::Solids::Ehlers::plasticFlowDeviatoricPart ( PhysicalStressWithInvariants< DisplacementDim > const & s,
OnePlusGamma_pTheta const & one_gt,
double const sqrtPhi,
typename SolidEhlers< DisplacementDim >::KelvinVector const & dtheta_dsigma,
double const gamma_p,
double const m_p )

◆ plasticFlowVolumetricPart()

template<int DisplacementDim>
double MaterialLib::Solids::Ehlers::plasticFlowVolumetricPart ( PhysicalStressWithInvariants< DisplacementDim > const & s,
double const sqrtPhi,
double const alpha_p,
double const beta_p,
double const delta_p,
double const epsilon_p )

Definition at line 105 of file Ehlers.cpp.

109{
110 return 3 *
111 (alpha_p * s.I_1 +
112 4 * boost::math::pow<2>(delta_p) * boost::math::pow<3>(s.I_1)) /
113 (2 * sqrtPhi) +
114 3 * beta_p + 6 * epsilon_p * s.I_1;
115}

References MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::I_1.

Referenced by calculatePlasticJacobian(), and calculatePlasticResidual().

◆ predict_sigma()

template<int DisplacementDim>
SolidEhlers< DisplacementDim >::KelvinVector MaterialLib::Solids::Ehlers::predict_sigma ( double const G,
double const K,
typename SolidEhlers< DisplacementDim >::KelvinVector const & sigma_prev,
typename SolidEhlers< DisplacementDim >::KelvinVector const & eps,
typename SolidEhlers< DisplacementDim >::KelvinVector const & eps_prev,
double const eps_V )

Definition at line 455 of file Ehlers.cpp.

461{
462 static int const KelvinVectorSize =
465 auto const& P_dev = Invariants::deviatoric_projection;
466
467 // dimensionless initial hydrostatic pressure
468 double const pressure_prev = Invariants::trace(sigma_prev) / (-3. * G);
469 // initial strain invariant
470 double const e_prev = Invariants::trace(eps_prev);
471 // dimensioness hydrostatic stress increment
472 double const pressure = pressure_prev - K / G * (eps_V - e_prev);
473 // dimensionless deviatoric initial stress
474 typename SolidEhlers<DisplacementDim>::KelvinVector const sigma_D_prev =
475 P_dev * sigma_prev / G;
476 // dimensionless deviatoric stress
477 typename SolidEhlers<DisplacementDim>::KelvinVector const sigma_D =
478 sigma_D_prev + 2 * P_dev * (eps - eps_prev);
479 return sigma_D - pressure * Invariants::identity2;
480}
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > KelvinVector
Definition Ehlers.h:291

References MathLib::KelvinVector::kelvin_vector_dimensions().

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ sOdotS()

template<int DisplacementDim>
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > MaterialLib::Solids::Ehlers::sOdotS ( MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const & v)

Special product of v with itself: \(v \odot v\). The tensor v is given in Kelvin mapping.

Note
Implementation only for 2 and 3 dimensions.
Attention
Pay attention to the sign of the result, which normally would be negative, but the returned value is not negated. This has to do with \(d(A^{-1})/dA = -A^{-1} \odot A^{-1} \).

Referenced by thetaSigmaDerivatives().

◆ sOdotS< 2 >()

Definition at line 860 of file Ehlers.cpp.

903{
905
906 result(0, 0) = v(0) * v(0);
907 result(0, 1) = result(1, 0) = v(3) * v(3) / 2.;
908 result(0, 2) = result(2, 0) = 0;
909 result(0, 3) = result(3, 0) = v(0) * v(3);
910
911 result(1, 1) = v(1) * v(1);
912 result(1, 2) = result(2, 1) = 0;
913 result(1, 3) = result(3, 1) = v(3) * v(1);
914
915 result(2, 2) = v(2) * v(2);
916 result(2, 3) = result(3, 2) = 0;
917
918 result(3, 3) = v(0) * v(1) + v(3) * v(3) / 2.;
919
920 return result;
921}

◆ sOdotS< 3 >()

Definition at line 860 of file Ehlers.cpp.

865{
867
868 result(0, 0) = v(0) * v(0);
869 result(0, 1) = result(1, 0) = v(3) * v(3) / 2.;
870 result(0, 2) = result(2, 0) = v(5) * v(5) / 2.;
871 result(0, 3) = result(3, 0) = v(0) * v(3);
872 result(0, 4) = result(4, 0) = v(3) * v(5) / std::sqrt(2.);
873 result(0, 5) = result(5, 0) = v(0) * v(5);
874
875 result(1, 1) = v(1) * v(1);
876 result(1, 2) = result(2, 1) = v(4) * v(4) / 2.;
877 result(1, 3) = result(3, 1) = v(3) * v(1);
878 result(1, 4) = result(4, 1) = v(1) * v(4);
879 result(1, 5) = result(5, 1) = v(3) * v(4) / std::sqrt(2.);
880
881 result(2, 2) = v(2) * v(2);
882 result(2, 3) = result(3, 2) = v(5) * v(4) / std::sqrt(2.);
883 result(2, 4) = result(4, 2) = v(4) * v(2);
884 result(2, 5) = result(5, 2) = v(5) * v(2);
885
886 result(3, 3) = v(0) * v(1) + v(3) * v(3) / 2.;
887 result(3, 4) = result(4, 3) =
888 v(3) * v(4) / 2. + v(5) * v(1) / std::sqrt(2.);
889 result(3, 5) = result(5, 3) =
890 v(0) * v(4) / std::sqrt(2.) + v(3) * v(5) / 2.;
891
892 result(4, 4) = v(1) * v(2) + v(4) * v(4) / 2.;
893 result(4, 5) = result(5, 4) =
894 v(3) * v(2) / std::sqrt(2.) + v(5) * v(4) / 2.;
895
896 result(5, 5) = v(0) * v(2) + v(5) * v(5) / 2.;
897 return result;
898}

◆ splitSolutionVector()

template<typename ResidualVector, typename KelvinVector>
std::tuple< KelvinVector, PlasticStrain< KelvinVector >, double > MaterialLib::Solids::Ehlers::splitSolutionVector ( ResidualVector const & solution)

Split the agglomerated solution vector in separate items. The arrangement must be the same as in the newton() function.

Definition at line 486 of file Ehlers.cpp.

487{
488 static auto const size = KelvinVector::SizeAtCompileTime;
489 return std::forward_as_tuple(
490 solution.template segment<size>(size * 0),
491 PlasticStrain<KelvinVector>{solution.template segment<size>(size * 1),
492 solution[size * 2], solution[size * 2 + 1]},
493 solution[size * 2 + 2]);
494}

Referenced by MaterialLib::Solids::Ehlers::SolidEhlers< DisplacementDim >::integrateStress().

◆ thetaSigmaDerivatives()

template<int DisplacementDim>
std::pair< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Ehlers::thetaSigmaDerivatives ( double theta,
PhysicalStressWithInvariants< DisplacementDim > const & s )

Definition at line 149 of file Ehlers.cpp.

151{
152 constexpr int KelvinVectorSize =
155 using KelvinVector =
157 using KelvinMatrix =
159
160 if (theta == 0)
161 {
162 return {KelvinVector::Zero(), KelvinMatrix::Zero()};
163 }
164
165 auto const& P_dev = Invariants::deviatoric_projection;
166
167 // inverse of deviatoric stress tensor
168 if (Invariants::determinant(s.D) == 0)
169 {
170 OGS_FATAL("Determinant is zero. Matrix is non-invertable.");
171 }
172 // inverse of sigma_D
173 KelvinVector const sigma_D_inverse = MathLib::KelvinVector::inverse(s.D);
174 KelvinVector const sigma_D_inverse_D = P_dev * sigma_D_inverse;
175
176 KelvinVector const dtheta_dsigma =
177 theta * sigma_D_inverse_D - 3. / 2. * theta / s.J_2 * s.D;
178 KelvinMatrix const d2theta_dsigma2 =
179 theta * P_dev * sOdotS<DisplacementDim>(sigma_D_inverse) * P_dev +
180 sigma_D_inverse_D * dtheta_dsigma.transpose() -
181 3. / 2. * theta / s.J_2 * P_dev -
182 3. / 2. * dtheta_dsigma / s.J_2 * s.D.transpose() +
183 3. / 2. * theta / boost::math::pow<2>(s.J_2) * s.D * s.D.transpose();
184
185 return {dtheta_dsigma, d2theta_dsigma2};
186}
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > sOdotS(MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &v)
Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > inverse(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)

References MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::D, MathLib::KelvinVector::inverse(), MaterialLib::Solids::Ehlers::PhysicalStressWithInvariants< DisplacementDim >::J_2, MathLib::KelvinVector::kelvin_vector_dimensions(), OGS_FATAL, and sOdotS().

Referenced by calculatePlasticJacobian(), and calculatePlasticResidual().

◆ yieldFunction()

template<int DisplacementDim>
double MaterialLib::Solids::Ehlers::yieldFunction ( MaterialProperties const & mp,
PhysicalStressWithInvariants< DisplacementDim > const & s,
double const k )