MFrontGallery Drucker Prager Triaxial Test

[superflyer11]

Hello,

using MFrontGallery/generic-behaviours/plasticity/DruckerPragerCap.mfront at https://github.com/nagelt/MFrontGallery/blob/master/generic-behaviours/plasticity/DruckerPragerCap.mfront with modifications from @Parameter to @MaterialProperty

Files:
DruckerPragerEdited.zip

Trying to run a triaxial test with the values based on the original file but I am not understanding the ranges of pressure that the model works in

I am referring to the web documentation as well but the +/- of the pa and pb doesnt quite make sense to me

[Asocie13]

Hello,

First, be careful here pa and pb take negative values due to the definition of the hydrostatic pressure. So it is the inverse of the geomechanics notation. The value corresponding to the end of the cap curve is pb, not pa. You can also specify your elastic properties with @MaterialProperty.

The default parameters are quite low and I think you must change them for your application. Furthermore, the DP part could not be stable if you imposed high pressure.

[superflyer11]

if the behaviour is already using @Brick StandardElasticity{young_modulus : 150e3, poisson_ratio : 0.3}; how can I specify E and nu with @MaterialProperty ?

[thelfer]

Don't define them:

@Brick StandardElasticity;

The StandardElasticity automatically defines two material properties with the glossary names YoungModulus and PoissonRatio respectively (see https://thelfer.github.io/tfel/web/glossary.html)

[superflyer11]

@Brick StandardElasticity;
// @Brick StandardElasticity{young_modulus : 150e3, poisson_ratio : 0.3};

@StateVariable strain p[2];
p.setGlossaryName("EquivalentPlasticStrain");

// @Parameter d    = 150;
// @Parameter beta = 0.6;
// @Parameter pa   = -10;
// @Parameter pb   = -75;

// maximum VM stress
@MaterialProperty stress d;
d.setEntryName("d");
// slope of the cap
@MaterialProperty real beta;
beta.setEntryName("beta");
// hydrostatic pressure at the cap
@MaterialProperty stress pa;
pa.setEntryName("pa");
// minimum hydrostatic pressure
@MaterialProperty stress pb;
pb.setEntryName("pb");

@LocalVariable real tg;
@LocalVariable real R;
@LocalVariable real Re;
@LocalVariable bool bpl[2];

@InitLocalVariables {
  tg = tan(beta);
  R = (pa - pb) / (d - pa * tg);
}

@Predictor{
  const auto sigel = computeElasticPrediction();
  // Drucker-Prager
  const auto prel = trace(sigel) / 3;
  const auto qel = sigmaeq(sigel);
  bpl[0] = (qel + prel * tg - d > 0) && (prel > pa);
  // Cap
  const auto seq_c = sqrt(power<2>(prel - pa) +  //
                          power<2>(R * qel));
  bpl[1] = (seq_c > (pa - pb)) && (prel <= pa);
}

@Integrator {
  const auto id = Stensor::Id();
  const auto id4 = Stensor4::Id();
  const auto id_x_id = Stensor4::IxI();
  const auto M = Stensor4::M();
  const auto seps =  real(1.e-12) * young;
  if ((bpl[0]) || (bpl[1])) {
    const auto pr = trace(sig) / 3;
    const auto q = sigmaeq(sig);
    const auto iq = 1 / max(q,seps);
    const auto nq = eval(3 * deviator(sig) * (iq / 2));
    const auto dnq = eval((M - (nq ^ nq)) * iq);
    if (bpl[0]) {
      const auto n   = eval(nq+(tg/3)*id);
      feel += dp[0]*n;
      fp(0) = (q+pr*tg-d)/young;
      // jacobian
      dfeel_ddeel += 2*mu*theta*dp[0]*dnq;

without specifying the quantities shows this error, should I define mu and lambda somewhere?

DruckerPragerEdited.mfront:74:18: error: ‘mu’ was not declared in this scope; did you mean ‘tfel::constants::mu<NumericType, use_quantities>’?
   74 |       dfeel_ddeel += 2*mu*theta*dp[0]*dnq;

[thelfer]

You need to compute it yourself in this case.

@LocalVariable stress mu;

@InitLocalVariables{
  mu = computeMu(young, nu);
}

[superflyer11]

Thank you for the help earlier, I have more questions now

@LocalVariable bool bpl[2];
This is delcared but not used, what is it for?

I am also having trouble understand what variable is automatically declared and evaluated. For a variable like sig, is it computed automatically or I can evaluate it myself by creating a @ComputeStress block?

What does p[0], p[1] stand for and why do I not need to declare them
and is dp[0], dp[1] automatically evaluated from p[0] and p[1]?

Snippet:

@StateVariable strain p[2];
p.setGlossaryName("EquivalentPlasticStrain");

...


@Predictor{
  const auto sigel = computeElasticPrediction();
  // Drucker-Prager
  const auto prel = trace(sigel) / 3;
  const auto qel = sigmaeq(sigel);
  bpl[0] = (qel + prel * tg - d > 0) && (prel > pa);
  // Cap
  const auto seq_c = sqrt(power<2>(prel - pa) +  //
                          power<2>(R * qel));
  bpl[1] = (seq_c > (pa - pb)) && (prel <= pa);
}

@Integrator {
  const auto id = Stensor::Id();
  const auto id4 = Stensor4::Id();
  const auto id_x_id = Stensor4::IxI();
  const auto M = Stensor4::M();
  const auto seps =  real(1.e-12) * young;
  if ((bpl[0]) || (bpl[1])) {
    const auto pr = trace(sig) / 3;
    const auto q = sigmaeq(sig);
    const auto iq = 1 / max(q,seps);
    const auto nq = eval(3 * deviator(sig) * (iq / 2));
    const auto dnq = eval((M - (nq ^ nq)) * iq);
    if (bpl[0]) {
      const auto n   = eval(nq+(tg/3)*id);
      feel += dp[0]*n;
      fp(0) = (q+pr*tg-d)/young;
      // jacobian
      dfeel_ddeel += 2*mu*theta*dp[0]*dnq;
      dfeel_ddp(0) = n;
      dfp_ddp(0,0) = 0;
      dfp_ddeel(0) = (2*mu*theta/young)*n+(tg*theta*lambda/young)*id;
    }

[thelfer]

The following page provides all the relevant information: https://thelfer.github.io/tfel/web/drucker-prager-cap.html

[superflyer11]

I've read through the page and most of my queries above is addressed. For variables like dp, are these automatically declared based on the list of state variables?

Another question irrelevant to above:

When using mtest in python, I am referring to this https://thelfer.github.io/tfel/web/mtest-python.html for the fields inside the MTestCurrentState class. To keep track of the state variables I am using this to view the variables:

s = mtest.MTestCurrentState()
wk = mtest.MTestWorkSpace()
m.completeInitialisation()
m.initializeCurrentState(s)
m.initializeWorkSpace(wk)
...
for i in range(nTime - 1):
    m.execute(s, wk, ltime[i], ltime[i + 1])
    print("===========")
    print(s.iv0)

in which s.iv0 gives out a list of the state variables.
Is there an exposed field to find the order of the state variables, for example the list of entry name or glossary name? I am already referring to the .mfront code but I am finding it difficult to keep track of which one is which

[thelfer]

Hi @superflyer11,

For variables like dp, are these automatically declared based on the list of state variables?

Yes

Is there an exposed field to find the order of the state variables, for example the list of entry name or glossary name?

See the Behaviour class, in particular the expandInternalStateVariablesNames method