Machine Learning flow rule¶
Machine Learning (ML) algorithms provide a great flexibility to describe aribitrary mathematical functions. At the same time they offer the possibility to handle large data sets and multi-dimensional features as input. Hence, using ML algorithms as constitutive rules for plastic material behavior offers the possibility to explicitly take into account microstructural information of the material in the constitutive modeling. Furthermore, data resulting from experiment and micromechanical simulations can be hybridized to generate training data sets. The present example is using Support Vector Classification (SVC) as yield function. The SVC algorithm is trained by using deviatoric stresses as input data and the information whether a given stress state leads to purely elastic or rather to elastic-plastic deformation of the material as result data. In this way, a ML yield function is obtained, which can determine whether a given stress state lies inside or outside of the elastic regime of the material. Furthermore, the yield locus, i.e., the hyperplane in stress space on which plastic deformation occurs, can be reconstructed from the SVC, and the gradient on this yield locus can be conveniently calculated. Therefore, the standard formulations of continuum plasticity, as the return mapping algorithm, can be applied in Finite Element Analysis (FEA) in the usual way. Thus, it is demonstrated that the new ML yield function can replace conventional FEA yield functions. In a next step, microstructural information will be considered directly as feature. Machine learning algorithms have been adopted from the scikit-learn platform (https://scikit-learn.org/stable/).
This package demonstrates the training and application of ML flow rules in FEA in form of a simple example with data synthetically produced from standard flow rules, like isotropic J2 (von Mises) and Hill-type anisotropic plasticity. It uses the pyLabFEM module. This notebook uses the matplotlib (https://matplotlib.org/) library for the visualization of results, and NumPy (http://www.numpy.org) for mathematical operations with arrays.
The scientific background of this work is desribed in mode detail in the article by A. Hartmaier "Data-Oriented Constitutive Modeling of Plasticity in Metals" Materials 2020, 13(7), 1600. This open access article is available here.
Author: Alexander Hartmaier, ICAMS, Ruhr-Universtität Bochum, Germany, March 2020
This work is licensed under a Creative
Commons Attribution-NonCommercial-ShareAlike 4.0 International License
(CC-BY-NC-SA)
The pyLabFEA package comes with ABSOLUTELY NO WARRANTY. This is free software, and you are welcome to redistribute it under the conditions of the GNU General Public License (GPLv3)
Theoretical background¶
The yield function of a material is defined as \begin{equation} f = \sigma_{eq} - \sigma_y\,, \end{equation} where plastic deformation sets in at $f=0$, i.e. when the the equivalent stress $\sigma_{eq}$ equals the yield strength $\sigma_y$ of the material. The equivalent stress used here is based on the pricipal stresses $\sigma_i$ with $(i=1, 2, 3)$, as
\begin{equation} \sigma^\mathrm{J2}_{eq} = \sqrt{ \frac{1}{2}\left[ H_1 \left(\sigma_1-\sigma_2\right)^2 + H_2 \left(\sigma_2-\sigma_3\right)^2 + H_3 \left(\sigma_3-\sigma_1\right)^2 \right] } . \end{equation}In this yield function, the anisotropy of the material's flow behavior is described in a Hill-like approach for orthotropic materials. If the axes of the pricipal stresses do not match the symmetry axes of the material, the material axes and with it the parameters $H_1, H_2$ and $H_3$ must be rotated into the coordinate system of the eigenvectors of the stress tensor.
The gradient of the yield function with respect to the pricipal stresses is needed for calculating the plastic strain increments in the return mapping algorithm of continuum plasticity, and can be evaluated analytically as \begin{equation} \frac{\partial f}{\partial \sigma_1} = \frac{\partial \sigma_{eq}}{\partial \sigma_1} = \frac{\left( H_1+H_3 \right) \sigma_1 - H_1 \sigma_2 - H_3 \sigma_3}{\sigma_{eq}} \\ \frac{\partial f}{\partial \sigma_2} = \frac{\partial \sigma_{eq}}{\partial \sigma_2} = \frac{\left( H_2+H_1 \right) \sigma_2 - H_1 \sigma_1 - H_2 \sigma_3}{\sigma_{eq}} \\ \frac{\partial f}{\partial \sigma_3} = \frac{\partial \sigma_{eq}}{\partial \sigma_3} = \frac{\left( H_3+H_2 \right) \sigma_3 - H_3 \sigma_1 - H_2 \sigma_2}{\sigma_{eq}} \\ \end{equation}
Note that in the case of isotropic plasticity, i.e. $H_1=H_2=H_3=1$, the gradient takes the simple form \begin{equation} \frac{\partial f}{\partial \sigma_i} = 3 \frac{\sigma_i - p}{\sigma_{eq}} \hspace{2em} (i=1,2,3) , \end{equation} where $p = 1/3 \mbox{Tr}(\sigma)$.
Deviatoric stress space¶
Since plastic deformation in most metals does not depend on hydrostatic stress components, it is useful to transform the principal stresses from the representation as Cartesian (3D) vector $\sigma=(\sigma_1, \sigma_2, \sigma_3)$ in the principal stress space into a vector $s=(\sigma_{eq}, \theta, p)$ in the cylindrical coordinate system, with the equivalent stress $\sigma_{eq}$ representing the norm of the deviator of $\sigma$ and the polar angle $\theta$ lying in the deviatoric stress plane normal to the hydrostatic axis $p$. This coordinate transformation improves the efficiency of the training, because only 2D data for the equivalent stress and the polar angle $\theta$ need to be used as training features, whereas the hydrostatic component is disregarded. The coordinate transformation is performed by introducing a complex-valued deviatoric stress \begin{equation} \sigma'_c = \pmb{\sigma}\cdot \pmb{a} + i\,\pmb{\sigma}\cdot \pmb{b} = \sqrt{2/3}\sigma_{eq}\; e^{i\theta}, \end{equation} where $i$ is the imaginary unit, such that the polar angle is obtained as \begin{equation} \theta = \mathrm{arg}\, \sigma'_c = -i \, \ln \frac{\pmb{\sigma}\cdot \pmb{a} + i\,\pmb{\sigma}\cdot \pmb{b}} {\sqrt{2/3}\sigma_{eq}} \,, \end{equation}
where $a=(2,-1,-1)/\sqrt{6}$ (real axis) and $b=(0,1,-1)/\sqrt{2}$ (imaginary axis) are the unit vectors that span the deviatoric stress plane normal to the hydrostatic axis $c=(1,1,1)/\sqrt{3}$.
An advantage of this coordinate transformation is that the gradient of the yield function w.r.t. the cylindrical coordinates has only one non-constant component. The complete gradient w.r.t. the cylindrical coordinates reads \begin{eqnarray*} \frac{\partial f}{\partial \sigma_{eq}} &=& 1 \\ \frac{\partial f}{\partial \theta} &=& \frac{\partial \sigma_y}{\partial \theta} \\ \frac{\partial f}{\partial p} &=& 0 . \end{eqnarray*}
To transform the gradient of the flow rule from this cylindrical coordinate system back to the principle stress space, in which form it is used later on to calculate the direction of the plastic strain increments in the return mapping algorithm of the plasticity model, we introduce the Jacobian matrix for this coordinate transformation as \begin{equation} J=\frac{\partial s}{\partial \sigma} = \left( \begin{array}{ccc} \frac{\partial \sigma_{eq}}{\partial \sigma_1} & \frac{\partial \theta}{\partial \sigma_1} & \frac{\partial p}{\partial \sigma_1} \\ \frac{\partial \sigma_{eq}}{\partial \sigma_2} & \frac{\partial \theta}{\partial \sigma_2} & \frac{\partial p}{\partial \sigma_2} \\ \frac{\partial \sigma_{eq}}{\partial \sigma_3} & \frac{\partial \theta}{\partial \sigma_3} & \frac{\partial p}{\partial \sigma_3} \\ \end{array} \right) \end{equation} with ${\partial \sigma_{eq}}/{\partial \sigma_i}$ as given above, ${\partial p}/{\partial \sigma_i}=1/3$, and \begin{equation} \frac{\partial \theta}{\partial \sigma_j} = -i\, \left( \frac{\pmb{a} + i\,\pmb{b}}{\pmb{\sigma}\cdot \pmb{a} + i\,\pmb{\sigma}\cdot \pmb{b}} - \frac{3\, \pmb{\sigma}'}{\sigma_{eq}^2} \right) \hspace{2em} (j=1,2,3) . \end{equation} Finally, the gradient in the 3D principle stress space is obtained as \begin{equation} \frac{\partial f}{\partial \sigma_j} = \sum\limits_{k=1}^{3} J_{jk} \frac{\partial f}{\partial s_k} \hspace{2em} (j=1,2,3) . \end{equation}
Material definition¶
In the following, the Python class Material is invoked to be used as a material card in FEA, demonstrating the application of standard flow rules and machine learning (ML) flow rules. This class Material contains the attributes and methods to defining the elastic and plastic material behavior and to evaluate the materials constitutive behavior. Furthermore, all necessary subroutines for plotting the results are defined.
Two identical materials are defined,mat_h will be used to apply the standard Hill-like flow rules in FEA and to generate synthetical data that is used to train a ML flow rule in mat_ml. Consequently,mat_h and mat_ml should reveal identical properties, which is verified in simple FEA demonstrations for 2D plane stress cases. As reference, a further material with the same yield strength, but isotropic J2 plasticity is defined as mat_iso. Here, we consider ideal plasticity, i.e. no work hardening, and no dependence on the hydrostatic stress.
