The research project of the Chair of Mechanics - Continuum Mechanics Ruhr-University Bochum in collaboration with Chair of Computational Mathematics - University of Augsburg

The goal of this project is the convexification of some energy function $W$, which becomes non-convex in continuum damage mechanics.


We use the FEM to solve the global structural problem. For this, consider the total potential $\Pi$

\[\Pi = \int_{\mathcal{B}} W(\boldsymbol{F}) \ dV - \Pi^{\text{ext}} \rightarrow \text{stationary}\]

Therefore we set the first variation to zero

\[\delta \Pi = \int_{\mathcal{B}} \delta \boldsymbol{F} \cdot \partial_{\boldsymbol{F}} W(\boldsymbol{F}) \ dV - \Pi^{\text{ext}} = 0\]

Now the question remains open how $W$ is mathematically structured. Since we want to incorporate dissipative processes (history dependent) we need a certain mathematical framework for this, which is called incremental variational formulations.

In this framework, we state the problem, such that we end up with a pseudo elastic potential $\Pi$ that depends only upon the current deformation gradient $\boldsymbol{F}_{k+1}$ Consider the classic continuum damage (1-D) approach, i.e.

\[\psi(\boldsymbol{F}\boldsymbol{F}^T, D(\beta)) = (1-D(\beta)) \psi^0 (\boldsymbol{F}^T\boldsymbol{F})\]

Here, $\psi^0$ denotes a virtually undamaged strain energy, called effective strain energy and \beta denotes an internal variable which is defined by

\[\beta \coloneqq \max_s [\psi^0(s)] \qquad s\in [0,t] \]

where t denotes the current time step. The damage variable $D(\beta)$ is typically some exponential function that never reaches 1. If $D=1$ the (1-D) term becomes 0 and thus reflects a fully damaged state. One example of the function is implemented in src/material.jl at the very top of the file.

Now, we evaluate the Clausius-Duhem inequality and obtain the dissipation potential $\phi$

\[\phi \coloneqq \psi^0 \dot{D} = \beta \frac{\partial D}{\partial \beta} \dot{\beta} > 0\]

Now, we define the incremental stress potential $W$ as

\[W(\boldsymbol{F}_{k+1}) = \inf_{\beta}[W(\boldsymbol{F}_{k+1}, \beta_{k+1})] \qquad \text{with} \qquad \mathcal{W}(\boldsymbol{F}_{k+1}, \beta_{k+1}) \coloneqq \int_{t_k}^{t_{k+1}} \dot{\psi} + \phi \ dt\]

Luckily, we can resolve the minimizing w.r.t. the internal variabel $\beta$ analytically, which is done in Balzani and Ortiz, 2012. Analytic minimization yields the following formula for the incremental stress potential

\[W(\boldsymbol{F}_{k+1}) = \psi(\boldsymbol{F}_{k+1}, D_{k+1}) - \psi(\boldsymbol{F}_k,D_k)+\beta D_{k+1} - \beta_k D_k - \tilde{D}_{k+1} + \tilde{D}_k\]

note that $\tilde{D}$ denotes the anti derivative of the damage function, which is implemented one line below the explicit example in src/materials.jl

Since we have now the incremental stress potential, we are good to go and can take the first and second derivative to obtain the first Piola-Kirchhoff stress tensor and the tangent modul, respectively.

\[\boldsymbol{P}_{k+1} = \frac{\partial W(\boldsymbol{F}_{k+1})}{\partial \boldsymbol{F}_{k+1}} \qquad \text{and} \qquad \mathbb{A}_{k+1} = \frac{\partial^2 W(\boldsymbol{F}_{k+1})}{\partial \boldsymbol{F}_{k+1}^2}\]

In this Julia package the heavy lifting w.r.t. the derivatives of $W$ is done by automatic differentiation, however, in Balzani and Ortiz, 2012 are explicit formulae as well.

The last undiscussed fact is the convexification of $W$. There exists different approaches but right now we are focusing on the one in Balzani and Ortiz, 2012. Hopefully, soon there will be an approach from the book of Sören Bartels.

In general we substitute $W$ by its rank-one convex hull $W_c$ in originally nonconvex regimes For an one dimensional case the convexified problem is given as

\[F\coloneqq \xi F^+ + (1-\xi)F^- \]

based thereon, we can compute $F^+$ and $F^-$

\[F^+ \coloneqq F(1+(1-\xi)d) \qquad \text{and} \qquad F^- \coloneqq F(1-\xi d)\]

Here, $F^+$ and $F^-$ can be interpreted as two virtually homogenised phases where the prior corresponds to a strongly damaged phase and the latter to a weakly damage phase, respectively. Now, in the nonconvex regime we keep those two states fixed, i.e. interpolate between the two states, and increase the volume fraction $\xi$ of those states. This mixture of states seems to be energetically favourable.

\[W \leftarrow W_C(F) = \inf_{\xi,d} [\overline{W}(F,d,\xi)] \qquad \text{with} \qquad \overline{W} = \xi W(F^+) + (1-\xi) W(F^-)\]

This yields a nonconvex, global two dimensional optimization problem on the material point (integration point) level, which is implemented in constitutive_driver(Material Routine section). The major problem of this approach is, that we cannot observe strain softening, which is the main motivation behind the project. One possible trick is to reconvexify the problem in each time step as depicted below

convexification_new relaxed