Skip to content

Commit f4d183b

Browse files
jduriezjduriez
committed
New MC state variables e_el and w_pl renamed (as such) and included in regression tests
1 parent 2cf71b3 commit f4d183b

File tree

2 files changed

+32
-9
lines changed

2 files changed

+32
-9
lines changed

‎include/materials/mohr_coulomb.tcc‎

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -71,10 +71,10 @@ mpm::dense_map mpm::MohrCoulomb<Tdim>::initialise_state_variables() {
7171
{"theta", 0.},
7272
// Plastic deviatoric strain
7373
{"pdstrain", 0.},
74-
// Elastic energy
75-
{"E_el", 0.},
76-
// Plastic work
77-
{"W_pl", 0.},
74+
// Elastic energy (volume-specific)
75+
{"e_el", 0.},
76+
// Plastic work (volume-specific)
77+
{"w_pl", 0.},
7878
// Flag for current failure because of tensile stress
7979
{"tensile_fail_curr", 0},
8080
// Flag for current failure because of shear stress
@@ -90,7 +90,7 @@ mpm::dense_map mpm::MohrCoulomb<Tdim>::initialise_state_variables() {
9090
template <unsigned Tdim>
9191
std::vector<std::string> mpm::MohrCoulomb<Tdim>::state_variables() const {
9292
const std::vector<std::string> state_vars = {
93-
"phi", "psi", "cohesion", "epsilon", "rho", "theta", "pdstrain", "E_el", "W_pl", "tensile_fail_curr", "shear_fail_curr", "tensile_fail", "shear_fail"};
93+
"phi", "psi", "cohesion", "epsilon", "rho", "theta", "pdstrain", "e_el", "w_pl", "tensile_fail_curr", "shear_fail_curr", "tensile_fail", "shear_fail"};
9494
return state_vars;
9595
}
9696

@@ -492,11 +492,11 @@ Eigen::Matrix<double, 6, 1> mpm::MohrCoulomb<Tdim>::compute_stress(
492492
Vector6d deps_el = ((1+poisson_ratio_)*dstress - poisson_ratio_*tr_dsig*Id) / youngs_modulus_; // Elastic deformation increment (only true strain coefficients)
493493
Vector6d strain_coefs; // coefficients necessary below
494494
strain_coefs << 1, 1, 1, 2, 2, 2;
495-
(*state_vars).at("E_el") += updated_stress.cwiseProduct(deps_el.cwiseProduct(strain_coefs)).sum()*ptr->volume(); // with 2*sigma_xy*eps_xy + .. on the out-of-diagonal terms, this is the correct expression sigma_ij:depsilon^el_ij (times volume)
495+
(*state_vars).at("e_el") += updated_stress.cwiseProduct(deps_el.cwiseProduct(strain_coefs)).sum(); // with 2*sigma_xy*eps_xy + .. on the out-of-diagonal terms, this is the correct expression sigma_ij:depsilon^el_ij
496496

497497
// Update plastic work
498498
Vector6d deps_pl = dstrain - deps_el.cwiseProduct(strain_coefs); // Plastic deformation increment in Voigt notation (with engineering strain coefficents)
499-
(*state_vars).at("W_pl") += updated_stress.cwiseProduct(deps_pl).sum()*ptr->volume();
499+
(*state_vars).at("w_pl") += updated_stress.cwiseProduct(deps_pl).sum();
500500

501501
return updated_stress;
502502
}

‎tests/materials/mohr_coulomb_test.cc‎

Lines changed: 25 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -116,6 +116,12 @@ TEST_CASE("MohrCoulomb is checked in 2D (cohesion only, without softening)",
116116
REQUIRE(state_variables.at("rho") == Approx(0.).epsilon(Tolerance));
117117
REQUIRE(state_variables.at("theta") == Approx(0.).epsilon(Tolerance));
118118
REQUIRE(state_variables.at("pdstrain") == Approx(0.).epsilon(Tolerance));
119+
REQUIRE(state_variables.at("e_el") == 0.);
120+
REQUIRE(state_variables.at("w_pl") == 0.);
121+
REQUIRE(state_variables.at("tensile_fail_curr") == 0);
122+
REQUIRE(state_variables.at("tensile_fail") == 0);
123+
REQUIRE(state_variables.at("shear_fail_curr") == 0);
124+
REQUIRE(state_variables.at("shear_fail") == 0);
119125

120126
const std::vector<std::string> state_vars = {"phi",
121127
"psi",
@@ -124,8 +130,8 @@ TEST_CASE("MohrCoulomb is checked in 2D (cohesion only, without softening)",
124130
"rho",
125131
"theta",
126132
"pdstrain",
127-
"E_el",
128-
"W_pl",
133+
"e_el",
134+
"w_pl",
129135
"tensile_fail_curr",
130136
"shear_fail_curr",
131137
"tensile_fail",
@@ -566,6 +572,11 @@ TEST_CASE("MohrCoulomb is checked in 2D (cohesion only, without softening)",
566572
REQUIRE(updated_stress(3) == Approx(0.).epsilon(Tolerance));
567573
REQUIRE(updated_stress(4) == Approx(0.).epsilon(Tolerance));
568574
REQUIRE(updated_stress(5) == Approx(0.).epsilon(Tolerance));
575+
// Check whether failure flag state variables (after being computed in compute_stress) are correct
576+
REQUIRE(state_variables.at("shear_fail") == 1);
577+
REQUIRE(state_variables.at("shear_fail_curr") == 1);
578+
REQUIRE(state_variables.at("tensile_fail") == 0);
579+
REQUIRE(state_variables.at("tensile_fail_curr") == 0);
569580
}
570581
}
571582
}
@@ -2866,6 +2877,10 @@ TEST_CASE("MohrCoulomb is checked in 3D (c & phi & psi, with softening)",
28662877
// Check plastic strain
28672878
REQUIRE(state_variables.at("pdstrain") ==
28682879
Approx(0.0007751567).epsilon(Tolerance));
2880+
// In the present workflow, all failure flags remain at 0 ?...
2881+
// Check plastic work and elastic energy
2882+
REQUIRE(state_variables.at("e_el") == Approx(15.33614017380).epsilon(Tolerance));
2883+
REQUIRE(state_variables.at("w_pl") == Approx(4.5135130001).epsilon(Tolerance));
28692884
}
28702885

28712886
//! Check for shear failure (pdstrain > pdstrain_residual)
@@ -3373,6 +3388,14 @@ TEST_CASE("MohrCoulomb is checked in 3D (c & phi & psi, with softening)",
33733388
// Check plastic strain
33743389
REQUIRE(state_variables.at("pdstrain") ==
33753390
Approx(0.0002112522).epsilon(Tolerance));
3391+
// Check failure flags
3392+
REQUIRE(state_variables.at("shear_fail") == 1);
3393+
REQUIRE(state_variables.at("shear_fail_curr") == 1);
3394+
REQUIRE(state_variables.at("tensile_fail") == 0);
3395+
REQUIRE(state_variables.at("tensile_fail_curr") == 0);
3396+
// Check plastic work and elastic energy
3397+
REQUIRE(state_variables.at("e_el") == Approx(1.153960429).epsilon(Tolerance));
3398+
REQUIRE(state_variables.at("w_pl") == Approx(1.0964114497).epsilon(Tolerance));
33763399
}
33773400

33783401
//! Check for shear failure (pdstrain > pdstrain_residual)

0 commit comments

Comments
 (0)