3. The IMPLEX method

The general elements of the IMPLEX method are shown here. For more details, we can refer in particular to [1].

The proposed method is based on two successive steps carried out to determine all the unknowns at the « pseudo-instant » (load) \({t}_{n\text{+}1}\).

The first step consists of an explicit extrapolation of the internal variables, then of the constraints, based on the quantities calculated previously (at load \({t}_{n}\)) and the load step \(\Delta {t}_{n\text{+}1}\). Thanks to this extrapolation, the local tangent matrix is evaluated and fixed; the resolution of the equilibrium equation (3) makes it possible to determine the displacements, which are considered fair and which are therefore fixed in turn.

The second step consists in achieving the implicit integration of behavior, based on the degrees of freedom evaluated in the previous step.

The first stage is of explicit type, while the second is of the implicit type, hence the name IMPLEX. At the end of these two steps, the equilibrium equation is not verified (it is true in terms of the fields extrapolated to step (1) only); however, the smaller the load increment, the smaller the error committed must be; it is therefore advisable to carry out the calculations with several different load steps in order to verify that the difference between the responses obtained is small (converged response in terms of step of load). In the following, the two steps are described. A summary of the key points is drawn up.

3.1. A 2-step method

Here we will present the 2 fundamental steps of method IMPLEX in succession, summarized in Table 1.

3.1.1. Explicit extrapolation and determination of degrees of freedom

This first step concerns in particular the internal variable \(\alpha\), whose evolution is governed by a type law (2).

At the start of load step \({t}_{n\text{+}1}\), we have all the information from load step \({t}_{n}\) and previous load steps. It is then possible to write the following Taylor developments (the evolution of the internal variable is considered sufficiently regular to be able to do so):

\(\begin{array}{}\{\begin{array}{cccc}{\alpha }_{n\text{+}1}& \text{=}& {\alpha }_{n}\text{+}\Delta {t}_{n\text{+}1}{\dot{\alpha }}_{n}+O(\Delta {t}_{n\text{+}1}^{2})& \\ {\alpha }_{n}& \text{=}& {\alpha }_{n\text{-}1}\text{+}\Delta {t}_{n}{\dot{\alpha }}_{n\text{-}1}+O(\Delta {t}_{n}^{2})& \\ \Delta {t}_{n}{\dot{\alpha }}_{n}& \text{=}& \Delta {t}_{n}{\dot{\alpha }}_{n\text{-}1}\text{+}\Delta {t}_{n}(\Delta {t}_{n}{\ddot{\alpha }}_{n\text{-}1})+\Delta {t}_{n}O(\Delta {t}_{n}^{2})& \Rightarrow \Delta {t}_{n}{\dot{\alpha }}_{n}\text{=}\Delta {t}_{n}{\dot{\alpha }}_{n\text{-}1}\text{+}O(\Delta {t}_{n}^{2})\end{array}\\ \Rightarrow \{\begin{array}{ccc}{\dot{\alpha }}_{n}& \text{=}& \frac{\Delta {\alpha }_{n}}{\Delta {t}_{n}}\text{-}O(\Delta {t}_{n}^{2})\\ {\alpha }_{n\text{+}1}& \text{=}& {\alpha }_{n}\text{+}\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\Delta {\alpha }_{n}-\Delta {t}_{n\text{+}1}O(\Delta {t}_{n}^{2})+O(\Delta {t}_{n\text{+}1}^{2})\end{array}\end{array}\) (5)

with \(\Delta {X}_{i}={X}_{i}-{X}_{i\text{-}1}\). By truncation, neglecting the terms of order two, we obtain the following prediction \(\tilde{{\alpha }_{n\text{+}1}}\) for the variable \({\alpha }_{n\text{+}1}\):

\({\tilde{\alpha }}_{n\text{+}1}={\alpha }_{n}+\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\Delta {\alpha }_{n}\) (6)

Through (6), we note that \({\tilde{\alpha }}_{n\text{+}1}\) is indeed obtained explicitly, at load \({t}_{n\text{+}1}\), as a function of the quantities obtained implicitly at load \({t}_{n}\). FIG. 1 schematizes this extrapolation method.

The extrapolation error can be defined as the difference between the value of the extrapolated internal variable and its real final value; subject to a sufficiently regular evolution of the internal variable to carry out the Taylor developments in sufficient order, this error can be evaluated at:

\({e}_{{\alpha }_{n\text{+}1}}=∣{\tilde{\alpha }}_{n\text{+}1}-{\alpha }_{n\text{+}1}∣\approx ∣{\ddot{\alpha }}_{n}∣\Delta {t}_{n\text{+}1}^{2}\) (7)

This shows that in the case of a sufficiently regular evolution of the internal variable, the error decreases quadratically as a function of the load step.

Once the internal variable has been extrapolated, we can determine the constraint \({\tilde{\sigma }}_{n\text{+}1}({\varepsilon }_{n\text{+}1},{\tilde{\alpha }}_{n\text{+}1})\) and the local tangent operator \({\tilde{C}}_{{e}_{n\text{+}1}}\).

In the particular case of isotropic damage laws, it is the local secant operator that is used; the constraints and this extrapolated operator are then written:

\(\mathrm{\{}\begin{array}{ccc}{\tilde{\sigma }}_{n\text{+}1}({\varepsilon }_{n\text{+}1},{\tilde{\alpha }}_{n\text{+}1})& \text{=}& (1\mathrm{-}{\alpha }_{n\text{+}1}){C}^{\text{elas}}\mathrm{:}{\varepsilon }_{n\text{+}1}\\ {\tilde{C}}_{{e}_{n\text{+}1}}\text{=}\frac{\mathrm{\partial }{\tilde{\sigma }}_{n\text{+}1}}{\mathrm{\partial }{\tilde{\varepsilon }}_{n\text{+}1}}& \text{=}& (1\mathrm{-}{\tilde{\alpha }}_{n\text{+}1}){C}^{\text{elas}}\end{array}\) (8)

These extrapolated local quantities are used to assemble the global tangent matrix \({\tilde{K}}_{n+1}\) via (4), and to determine the internal force.

The equilibrium equation in terms of extrapolated field \({F}_{\text{ext}}({t}_{n\text{+}1})-{F}_{\text{int}}({U}_{n\text{+}1},{\tilde{\sigma }}_{n\text{+}1},{t}_{n\text{+}1})=0\) is then solved with this extrapolated matrix and the displacement field \({\tilde{U}}_{n\text{+}1}\) obtained. It is solved by a classical Newton-Raphson method; it is therefore linearized by \({\tilde{K}}_{n\text{+}1}\delta {\tilde{U}}_{n\text{+}1}={F}_{\text{ext}}({t}_{n\text{+}1})-{F}_{\text{int}}({\tilde{\sigma }}_{n\text{+}1},{t}_{n\text{+}1})\). In the case of the laws of behavior treated here (isotropic damage with secant operator and linear isotropic plasticity), the behavior operator is constant for a time step (i.e. independent of the current deformation state \({\varepsilon }_{n\text{+}1}\)); the Newton-Raphson resolution process becomes linear by time step and therefore converges in an iteration. However, this step is iterative in some cases (non-linear plasticity for example), due to the non-constant nature of the tangent operator. At the first iteration, we would then take \({\varepsilon }_{n\text{+}1}^{0}={\varepsilon }_{n}\).

At the end of this first step, the displacement field \({U}_{n\text{+}1}\) is considered to be equal to the field obtained by solving the equilibrium equation written in terms of the extrapolated fields: \({U}_{n\text{+}1}={\tilde{U}}_{n\text{+}1}\). This travel field will not be modified later.

++——————————————————————————————————————-+ || | ++——————————————————————————————————————-+ || | ++ .. image:: images/10041D74000069F000002EEC1E1CD17971BC585E.svg + || :width: 256 | ++ :height: 174 + || | ++ + || | ++——————————————————————————————————————-+

Figure 1: Diagram of the extrapolation method

3.1.2. Implicit determination of the elements of the law of behavior

After the first phase, the field of movement is known. It is fixed for the load step and will therefore not be modified.

The deformation \({\varepsilon }_{n\text{+}1}({U}_{n\text{+}1})\) is then determined, then the equations (2) of the law of behavior are solved implicitly in order to obtain the stress fields \({\sigma }_{n\text{+}1}\) and the internal variable \({\alpha }_{n\text{+}1}\).

This step is identical to the standard solution of equations for laws of behavior. The only difference is the need to store \(\frac{\Delta {\alpha }_{n\text{+}1}}{\Delta {t}_{n\text{+}1}}\) to extrapolate the internal variable to the next step.

At the end of this step, the displacement field \({U}_{n\text{+}1}\) and the constraint fields \({\sigma }_{n\text{+}1}\) and the internal variable \({\alpha }_{n\text{+}1}\) fields are therefore known. A major difference compared to an implicit classical calculation is the fact that the real equilibrium equation is not verified; it is only true in terms of the extrapolated fields.

Step 1: Explicit Extrapolation	Step 2: Implicit Integration
Extrapolation of internal variables: \({\tilde{\alpha }}_{n\text{+}1}={\alpha }_{n}+\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\Delta {\alpha }_{n}\)	Freeze moves: \({\dot{U}}_{n\text{+}1}=0\)
Calculation of the extrapolated stress: \({\tilde{\sigma }}_{n\text{+}1}({\varepsilon }_{n\text{+}1},{\tilde{\alpha }}_{n\text{+}1})\)	Implicit resolution of laws of behavior: State law: \({\sigma }_{n\text{+}1}\mathrm{=}\mathrm{\sum }(\varepsilon ({U}_{n\text{+}1}),{\alpha }_{n\text{+}1})\) Law of evolution: \(\mathrm{\{}\begin{array}{c}f({\alpha }_{n\text{+}1},{\sigma }_{n\text{+}1})\le 0\\ {\dot{\alpha }}_{n\text{+}1}\mathrm{=}{\dot{\lambda }}_{n\text{+}1}\ge 0\\ {\dot{\lambda }}_{n\text{+}1}f({\alpha }_{n\text{+}1},{\sigma }_{n\text{+}1})\mathrm{=}0\end{array}\)
Solving the equilibrium equation in extrapolated fields: \({F}_{\text{ext}}({t}_{n\text{+}1})-{F}_{\text{int}}({U}_{n\text{+}1},{\tilde{\sigma }}_{n\text{+}1},{t}_{n\text{+}1})=0\)
Stage 1 exit: \({U}_{n\text{+}1}\)	Stage 2 exit: \({U}_{n\text{+}1}\), \({\sigma }_{n\text{+}1}\), \({\alpha }_{n\text{+}1}\)

Table 1: Summary of method IMPLEX

3.2. Automatic time step management

Like all explicit methods, method IMPLEX introduces an intrinsic error that must decrease almost quadratically with the load step. The solution may therefore depend slightly on the load step selected by the user.

If he wishes, he can use automatic management of the load step via the command DEFI_LIST_INST, by specifying METHODE = “AUTO” and MODE_CALCUL_TPLUS =” IMPLEX “(cf. [9]). This method makes it possible to control the error while optimizing the calculation time provided that the user has chosen a well-calibrated first time step.

The aim is therefore to minimize the error defined by equation (7). To do this, we will maximize the variable increment extrapolated by a quantity marked \(\mathit{Tol}\):

\(\Delta {\tilde{\alpha }}_{n\text{+}1}=\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\Delta {\alpha }_{n}\le \mathrm{Tol}\) (9)

From where:

\(\Delta {t}_{n\text{+}1}\le \frac{\mathrm{Tol}\mathrm{.}\Delta {t}_{n}}{\Delta {\alpha }_{n}}\) (10)

As the internal variable increments depend on the point in question, the time step will be chosen as the minimum value of the expression (10) over the whole structure, i.e. finally:

\(\Delta {t}_{n\text{+}1}=\mathrm{Tol}\mathrm{.}\Delta {t}_{n}{\mathrm{MIN}}_{x\in \Omega }(\frac{1}{\Delta {\alpha }_{n}(x)})\) (11)

Conditions limiting acceleration and deceleration are then added, as well as minimum and maximum limits for the time increment. The acceleration conditions \((\mu )\) and deceleration are not modifiable and have been calibrated to concrete cases, while the minimum and maximum limits can be modified by the user (the default values are given here).

\(\begin{array}{c}\Delta {t}_{n\text{+}1}\le \mu \Delta {t}_{n}=\mathrm{1,2}\Delta {t}_{n}\\ \Delta {t}_{n\text{+}1}\ge 0.5\Delta {t}_{n}\\ \Delta {t}_{n\text{+}1}\le 10\Delta {t}_{0}\\ \Delta {t}_{n\text{+}1}\ge 0.001\Delta {t}_{0}\end{array}\) (12)

Thus, the first time step \(\Delta {t}_{0}\), provided by the user, defines the time increment limits. To choose it, the user is advised to have carried out preliminary calculations with the Newton method and to know the elastic limit of the structure; it would seem that the choice of a first time step equal to half the elastic limit allows a good efficiency of the method.

Let us then show that the error committed is well controlled. To do this, the assumption is made that equation 9 is true. So we have \(\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\Delta {\alpha }_{n}\le \mathrm{Tol}\), and if we consider that the acceleration and deceleration limit factors allow us to write \(\frac{\Delta {t}_{n\text{+}1}}{\Delta {t}_{n}}\approx 1\), we have \(\Delta {\alpha }_{n}\le \mathrm{Tol}\); and progressively:

\(\forall n\in \mathbb{N},\Delta {\alpha }_{n}\le \mathrm{Tol}\) (13)

This finally makes it possible to write:

\({e}_{{\alpha }_{n\text{+}1}}\approx ∣\ddot{{\alpha }_{n}}∣\Delta {t}_{n\text{+}1}^{2}\le \mathrm{Tol}{\mu }^{2}(1\text{+}\mu )\) (14)

The error is therefore under control.

3.3. Key points of the method

After this summary presentation, some remarks are made here to understand the value of the method and also its limitations.

First of all, it is based on an extrapolation of internal variables, determined from Taylor development. Since Taylor’s developments are only valid for sufficiently regular functions, so will the method. Thus, the crossing of the elastic limit, or the transition from a state of charge to discharge, are points for which the method is, from a point of view, not suitable: the derivative of the damage is zero on the discharge side or the elastic limit side, and not zero on the load side. However, if the load steps are sufficiently small, the approximation can be made, insofar as the implicit correction takes place and therefore the extrapolation errors are partly erased. In unstable situations, however, for example when the damaged area grows abruptly (which is generally characterized by a significant snap-back of the overall force/displacement response), the method, although robust, cannot guarantee a reliable response, regardless of the load increment used: it is not adapted to this type of situation.

According to equation (8), in the case of isotropic damage laws, due to the limitation of \(\alpha\) to 1 and the symmetric and positive definite nature of the local elastic matrix \({C}^{\text{elas}}\), the local secant matrix \({\tilde{C}}_{{e}_{n\text{+}1}}\) is always symmetric definite positive. Therefore, by assembly, the global tangent matrix \({\tilde{K}}_{n\text{+}1}\) remains well conditioned: the robustness problems mentioned in the introduction, linked to the increasingly singular nature of \(K\), must therefore be eliminated.

Moreover, for a given load step and for the laws of behavior developed here (ENDO_FRAGILE confer [3], ENDO_ISOT_BETON confer [], and VMIS_ISOT_LINE confer [5]), the local tangent matrix is known by explicit extrapolation and remains constant throughout the load step. 4 The linearization of the equilibrium equation (1) by a Newton-Raphson method results in a global tangentconstant matrix \({\tilde{K}}_{n\text{+}1}\); in other words, the equilibrium equation becomes linear at each load step, and its resolution requires only one iteration.

Then, and at the risk of being redundant, the method leads to unbalanced external and internal forces at the end of each load step; they are only balanced at the end of the first phase, i.e. only in terms of extrapolated fields:

\(\{\begin{array}{ccc}{F}_{\text{ext}}({t}_{n\text{+}1})-{F}_{\text{int}}({\sigma }_{n\text{+}1}({U}_{n\text{+}1}),{t}_{n\text{+}1})& \ne & 0\\ {F}_{\text{ext}}({t}_{n\text{+}1})-{F}_{\text{int}}({U}_{n\text{+}1},{\tilde{\sigma }}_{n\text{+}1},{t}_{n\text{+}1})& \text{=}& 0\end{array}\) (15)

Therefore, during the calculation, the algorithm should not be asked to check the residue to within a tolerance, as is usual to do (implicitly). Again, from this point of view, robustness is guaranteed, since the classical stopping criterion is irrelevant.

Finally, this method is only intended to increase the robustness of the calculation, and not the quality of the response obtained. Thus, it introduces an intrinsic error through extrapolation. At best, it is only possible to obtain results that are almost as good as those obtained by a conventional implicit method of local resolution. This error must decrease in a near quadratic manner as a function of the imposed load step, subject to a sufficiently regular evolution of the internal variable (which in fact excludes the unstable propagation of damaged areas).

The use of this method therefore requires a certain amount of critical attention. In order to avoid major errors, it is recommended to perform calculations with different and small load increments, in order to verify that the solution does not differ too much from one increment to another (in other words, that the solution is close to convergence in terms of load increments). Moreover, in the case of overcoming unstable situations, the results obtained should only be considered for their qualitative aspect.