Prepring of conference paper "Finite element implementation of finite strain constitutive models formulated in the logarithmic strain space"

Finite element implementation of finite strain constitutive models formulated in the logarithmic strain space A. Moskovkaa,b, M. Frosta, J. Valdmanc, M. Hor´ akc,d aInstitute of Thermomechanics, Czech Academy of Sciences, Dolejˇ skova 5, 18200 Praha, Czechia bDepartment of Mathematics, Faculty of Applied Sciences, University of West Bohemia, Technick´ a 8, 30100 Plzeˇ n, Czechia cInstitute of Theory of Information and Automation, Czech Academy of Sciences, Pod Vod´ arenskou vˇ eˇ z´ ı 4, 18200 Praha, Czechia dDepartment of Mechanics, Faculty of Civil Engineering, Czech Technical University in Prague, Th´ akurova 7, 16629 Praha, Czechia 1. Introduction In many real-world applications, solids undergo severe deformation. In such situations, the continuum mechanics description must extend beyond the infinitesimal strain setting and employ the finite strain theory. It usually requires reformulation of the constitutive model, which poses both theoretical and computational challenges, especially when multiplicative decomposition of strain gradient is applied, e.g., due to inelastic deformation mechanisms. In their seminal papers [2, 3], Miehe and coworkers developed an alternative approach that benefits from constitutive models developed for small deformations and seamlessly integrates them into the generalized logarithmic strain space. The primary concern then shifts from the development of constitutive theories to the adoption of geometric transformations that convert the outputs of the constitutive law from this space into objects appropriate for further computational treatment, e.g., with finite element methods (FEM). This paper exemplifies the procedure and its numerical implementation, including the incorporation into FEM, on a three-dimensional linear elasticity model. The procedure encompasses both the pre-processing phase, in which the logarithmic strain measure is computed from displacements, the resolution of the constitutive model, and the post-processing phase, in which the fundamental quantities (stress tensor, tangent modulus) must be consistently transformed to comply with the particular FEM formulation. Our algorithmic realization of the FEM employing Newton’s iterative method relies on the methodology described in [1]. 2. Constitutive Models in the Logarithmic Strain Space Many FEM software suites enable users to implement specific constitutive behaviours of materials through an independent routine that communicates with the FEM core code through inputs and outputs. A schematic workflow of the constitutive routine takes the form {E,F,P}→Constitutive Model → {T,D, e F},(1) where Edenotes the symmetric strain tensor, Tis the corresponding work-conjugate stress tensor, Dis the tangent (constitutive) modulus, Fand e Fdenote a specific set of internal variables at the input and their updated values at the output, respectively, and Pis a set of auxiliary parameters. Recall that given a displacement U, the corresponding deformation gradient is set as F=∇U+Iand the right Cauchy-Green tensor as C=F⊤F. All the variables are considered as fields over a domain Ω. In classical linear elasticity, the constitutive law, σ=C:ε, returns the Cauchy stress tensor σas the output, T, based on the infinitesimal strain tensor εentering as E;D=Cis the fourth-order stiffness tensor and F, e Fand Pare empty sets. Algorithm 1: Condensed computational workflow 1. FEM assemblies (based on [1]): i) discretization of the computational domain Ωinto elements, ii) computing the values of the local basis functions and their derivatives at the Gauss (integration) points, element volumes, etc. 2. Assembly of the deformation-displacement matrix Band the global elastic stiffness matrix Ke. 3. For-loop over time steps: i) setting the boundary conditions and the initial displacements Ut, ii) while-loop of Newton’s method for solving a nonlinear system: A) Preprocessing: a) determining Ffrom Uand computing C=F⊤F, b) computing the eigenvalues and eigenvectors of Cand then Evia (2), B) Resolution of the constitutive problem (1), C) Postprocessing: a) computing the transformation tensors Pand Lfrom (3), b) computing Pand Musing (4) (and additionally S,σvia (5)), D) assembly of the global tangent matrix Ktutilizing M, E) computing new displacements Uit by solving the system, F) checking the convergence of the Newton’s method. Finite strain reformulation of models according to [2] retains the internal structure of a (rateindependent) “Constitutive Model” in (1) and employs the logarithmic (Hencky) strain as the input strain, i.e. E=1 2ln C.(2) Such a choice is especially beneficial for inelastic models relying on the assumption of additivity of strains, since the additivity is preserved in the logarithmic strain space. However, the stress measure work-conjugate to (2) provided by the model, as well as the tangent modulus, have to be transformed to objects compatible with the formulation of the FEM. Let us consider deformation gradient, F, and the first Piola Kirchhoff (nominal) stress, P, as the conjugate pair for the variational formulation of FEM, and denote the corresponding tangent modulus as M; then the transformations T→Pand D→Mmust be applied. The corresponding transformation tensors Pand Lare obtained with the relations (see [2] for details) P=∂E ∂F,L=∂2E ∂F2,(3) and the transformation formulae take the forms P=T:P,M=PT:D:P+T:L.(4) We also recall the following relation for the (symmetric) second Piola-Kirchhoff stress Sand the (symmetric) Cauchy stress σin terms of P S=F−1P,σ=1 det FP F ⊤.(5) Fig. 1. Reference and computed deformed configurations from the example. The initial state is in green, half tensile stretch in orange, full tensile stretch (five times the initial length) in red, half compression stretch in cyan, and full compression (quarter the initial length) in blue. 3. Example: Uniaxial Stretching of Elastic Cube Let us illustrate the approach and its incorporation into FEM on an example. For simplicity, let us consider a constitutive model without internal variables or parameters, i.e. F=P=∅. In such a case, the complete computational workflow can be summarized as shown in Box “Algorithm 1” above. Let us note that actions A–Cmust be performed at all integration points. The workflow is then implemented and validated on a boundary value problem of uniaxial stretching of an elastic cube. The cube is stretched in one direction up to five times its initial length, then returned to its initial length, and finally compressed to one-quarter of its initial length, see Fig. 1. This is imposed by prescribing displacements (a nonhomogeneous Dirichlet boundary condition) on one face of the cube in the direction of loading/unloading, whereas the opposite face remains fixed in that direction. For the constitutive model, we employ a linear isotropic Hooke’s law i.e. T= Λtr(E)I+ 2µE(6) with Λ, µ being the Lam´ e constants. Motivated by [3] and without loss of generality, we set Λ = µ. Denoting the stretch ratio parallel and perpendicular to the stretching direction as λ1 and λ2, respectively, one gets: F11 =λ1, F22 =F33 =λ2, Fij = 0 for i, j ∈ {1,2,3}∧i=j. Fig. 2. A comparison of analytical and FEM-computed values of component C11 and stretch ratio λ2as functions of λ1for the elastic cube example. Fig. 3. A comparison of the analytical and FEM-computed values of the non-zero component of current (T11), nominal (P11), 2nd Piola-Kirchhoff (S11), and Cauchy (σ11) stresses as functions of λ1for the elastic cube example. It is then a matter of algebraic manipulations employing relations (2)–(5) to successively derive the following analytical formulas: C11 =λ2 1, E11 = ln λ1, λ2=1 4 √λ1 , T11 µ=5 2ln λ1,P11 µ=5 2 ln λ1 λ1 ,S11 µ=5 2 ln λ1 λ2 1 ,σ11 µ=5 2 ln λ1 √λ1 .(7) To verify the numerical implementation, Figs. 2 and 3 compare these analytical formulas with the values computed via in-house FEM. A finite-strain computational analysis of more complex boundary value problems involving large rotations and deflections with a constitutive model featuring a rate-independent inelastic deformation mechanism can be found in [4]. Acknowledgements This work has been financially supported by the Operational Programme Johannes Amos Comenius of the Ministry of Education, Youth and Sport of the Czech Republic, within the frame of project Ferroic Multifunctionalities (FerrMion) [project No. CZ.02.01.01/00/22 008/0004591], co-funded by the European Union. References [1] ˇ Cerm´ ak, M., Sysala, S., Valdman, J., Efficient and flexible MATLAB implementation of 2D and 3D elastoplastic problems, Appl. Math. Comput. 355 (2019), 595. [2] Miehe, C., Apel, N., Lambrecht, M., Anisotropic additive plasticity in the logarithmic strain space: modular kinematic formulation and implementation based on incremental minimization principles for standard materials, Comput. Methods Appl. Mech. Eng. 191 (2002), 5383. [3] Miehe, C., Lambrecht, M., Algorithms for computation of stresses and elasticity moduli in terms of Seth–Hill’s family of generalized strain tensors, Comm. Numer. Methods Eng. 17 (2001), 337. [4] Moskovka, A., Hor´ ak, M., Valdman, J., Knapek, M., Janeˇ cek, M., Sedl´ ak, P. and Frost, M., Finite-Strain Constitutive Model for Shape Memory Alloys Formulated in the Logarithmic Strain Space, Shape Mem. Superelasticity 11 (2025), in print.