new feature request: initialize initial stress in mgis.fenics

[cenguo]

Hi,

I'd like to ask for a new feature in mgis.fenics, which can initialize the initial stress condition before loading is performed (in the context of soil mechanics).

Thanks, Ning

[thelfer]

Hi Ning, This is a discussion that we regularly have with @nagelt for applications in OpenGeoSys. First, you can directly initialize the stress in the mgis.fenics. The real trouble is that most behaviours written with MFront (in particular those based on the StandardElasticity and StandardElastoViscoPlasticity bricks) will discard this information because they rely on the elastic strain. So I expected that what you really need is a way to initialize the elastic strain and maybe other state variables...

Here is a proposal that we could discuss together:

1. I could do is to add a new code block in MFront called @InitStateVariables.
2. By default, this would do nothing. When the StandardElasticity and StandardElastoViscoPlasticity bricks are used, the default behaviour would be to initialize the elastic strain from the stress.
3. This code block is meant to be called explicitly at the beginning of the computation by the calling solver at each integration points. To do this, a new value in the IntegrationType enum could be added.

What do you think ?

[cenguo]

Hi Thomas,

Thanks for your reply. I think your proposal will work perfectly. Previously, I worked with mgis.behaviour and set internal state variables explicitly (see below), and it also worked. I just didn't find a way to do this in mgis.fenics.

epsilon = numpy.dot(compliance, init_sigma.T).T
materData = mgis_bv.MaterialDataManager(behaviour, ngauss)
for s in [materData.s0, materData.s1]:
    s.internal_state_variables[:,:4] = epsilon
    s.internal_state_variables[:,chi_off] = mater['chi0'] # initial value of this state variable is non-zero
    s.thermodynamic_forces[:] = init_sigma

Anyway, I think your solution is more elegant.

Best, Ning

[thelfer]

As a temporary workaroung, I think that you can access the MaterialDataManager as follows:

    material = mf.MFrontNonlinearMaterial("./src/libBehaviour.so",
                                          "IsotropicLinearHardeningPlasticity",
                                          hypothesis=hypothesis,
                                          material_properties=mat_prop)
    problem = mf.MFrontNonlinearProblem(u, material, quadrature_degree=2, bcs=bc)
    for s in [material.data_manager.s0,material.data_manager.s1]:
      s.internal_state_variables[:,:4] = epsilon
      s.internal_state_variables[:,chi_off] = mater['chi0'] # initial value of this state variable is non-zero
      s.thermodynamic_forces[:] = init_sigma

Note that the material data manager is set when the material is assigned to a problem, hence you can't access it before the definition of your problem....

[thelfer]

Hi @cenguo, @nagelt,

I started working of this issue. I came up to this design (not yet fully implemented).

1. In MFront, we should be able to define various @Initialize code blocks as follows:

   
   @Behaviour InitializeTest;
   @Author Thomas Helfer;
   @Date 12/01/2021;
   @Description{
   "First tests of the @Initialize method"
   }

@MaterialProperty stress young; young.setGlossaryName("YoungModulus"); @MaterialProperty real nu; nu.setGlossaryName("PoissonRatio");

@StateVariable StrainStensor eel; @LocalVariable stress K, mu;

@InitLocalVariables { const auto lambda = computeLambda(young, nu); mu = computeMu(young, nu); K = lambda + (2 * mu) / 3; }

@Initialize{ const auto s = deviator(sig); const auto pr = trace(sig) / 3; eel = (pr / (3K)) Stensor::Id() + s / 2mu; }

where `ElasticStrain` is a name used to differentiate this initialisation steps from the others

2. In most interfaces, those code blocks will have no effect
3. In `MGIS`, we could write something like:
````c++
d = mgsi::behaviour MaterialDataManager d(b,100)
// set the stress
..
//
mgis::behaviour::intialize(d,"ElasticStrain");

4. The StandardElasticy and StandardElastoViscoplasticity bricks could provide some predefined initialisation steps, such as the previous example.

What do you think of this proposal ?

Best,

Thomas

[nagelt]

@thelfer great to see this issue is addressed! This way, we can get rid of some of our current work-arounds. It's important also because some soil models have a set of internal variables with non-zero initial values, such as an initial void ratio for example (see Christian's Cam Clay implementation).

So if I understand correctly, all internal variables other than stress can be initialized normally by the solver, and stress is initialized by this inversion eel = inverse(C) * sig or so.

[thelfer]

So if I understand correctly, all internal variables other than stress can be initialized normally by the solver, and stress is initialized by this inversion eel = inverse(C) * sig or so.

Yes, that may be a possible choice. The idea is that the calling solver will have a list of initialize functions that he can call or not :) For example, I plan to automatically declare an initialize function for eel if you use the StandardElasticity brick (or the StandardElastoViscoPlasticity brick. But you can add you own initialisation function, using:

@Initialize<VoidFraction>{
....
}

The user is free to choose an appropriate name.

I also try to figure out if a "global initialisation function" which calls every individual initialisation functions would make sense.

[cbsilver]

Hello Thomas & everybody, I'd like to add something referring to the initialization of state variables (other than the stress): This is very convenient calling a MFront material model directly from mtest or from python like this: @InternalStateVariable 'Porosity' 0.4 m.setInternalStateVariableInitialValue('Porosity', 0.4)

However, calling MFront from OpenGeoSys (OGS) we currently use a workaround, which has the drawback, that the initial value is merely overwritten at t=0, but it is not saved within the state variable for later output. Further, this procedure will fail (re)starting a simulation at some t>0. Considering Thomas Helfer's new (not yet fully implemented) design: So it will then be possible to initialize the state variables manually in (various) @initialize code blocks, prescribing the initial value with some input parameter? That would, indeed, be a nice feature such that we can properly initialize the state variables from OpenGeoSys input files.

[thelfer]

Hi @cbsilver, I would recommend that the initialization of the internal state variables to be handled on the OpenGeoSys side. In some cases, my proposal could be used to provide some complex initialization steps, but I would be cautious about overuse of this feature. Best, Thomas

[thelfer]

Developments started here: https://github.com/thelfer/tfel/issues/25

[thelfer]

@nagelt, @cenguo Would you have a look at Issues https://github.com/thelfer/tfel/issues/24 and https://github.com/thelfer/tfel/issues/25 and make me early feed-backs.

[thelfer]

@nagelt, @cenguo, we have an early working example in the th/work_issue_50. Do not hesitate to make feedbacks on:

• https://github.com/thelfer/MFrontGenericInterfaceSupport/blob/th/work_issue_50/tests/InitializeFunctionTest.mfront
• https://github.com/thelfer/MFrontGenericInterfaceSupport/blob/th/work_issue_50/tests/InitializeFunctionTest.cxx

[thelfer]

@nagelt, @cenguo python bindings are done. See:

• https://github.com/thelfer/MFrontGenericInterfaceSupport/blob/th/work_issue_50/bindings/python/tests/InitializeFunctionTest.py

[thelfer]

The available initialize functions are described by the initialize_functions member of the Behaviour class which associates the name of the initialize function and a small data structure containing:

• the pointer to the initialize function,
• the list of inputs of the initialize function, if any.

The getInitializeFunctionVariablesArraySize function returns the size of an array able to contain the inputs an initialize function for an integration point.

The allocateInitializeFunctionVariables functions return an array able to store the inputs of an initialize function for a given integration point or a set of integrations points.

The executeInitializeFunction functions execute an initialization function on a unique integration point or a set of integration points.

Note about the result of the initialize function

The BehaviourDataView class imposes that the initial state is immutable. Thus, initialize functions can thus only initialize the state of the material at the end of the time step. In most case, the call to the selected initialize functions shall be follow to a call to the update function.

Note about the material properties and the external state variables

The material properties and the external state variables must be set before calling the executeInitializeFunction functions. Only the values at the beginning of the time step are passed to the behaviour.

Example of usage

auto d = BehaviourData{b};
// initialize the material properties and the external state variables
...
// calling an initialize function which requires an input
auto inputs = allocateInitializeFunctionVariables(b, "StressFromInitialPressure");
inputs[0] = pr;
auto v = make_view(d);
executeInitializeFunction(v, b, "StressFromInitialPressure", inputs);