Karlsruher Institut für Technologie

**Schriftenreihe Kontinuumsmechanik im Maschinenbau**

Johannes Ruck

**18**

J. Ruck

**Modeling martensite transformation based on a sharp interface**

Modeling martensitic phase transformation in dual phase steels based on a sharp interface theory

Johannes Ruck

**Modeling martensitic phase transformation in dual phase steels based on a sharp interface theory**

#### **Schriftenreihe Kontinuumsmechanik im Maschinenbau Band 18**

Karlsruher Institut für Technologie (KIT) Institut für Technische Mechanik Bereich Kontinuumsmechanik

Hrsg. Prof. Dr.-Ing. habil. Thomas Böhlke

Eine Übersicht aller bisher in dieser Schriftenreihe erschienenen Bände finden Sie am Ende des Buchs.

## **Modeling martensitic phase transformation in dual phase steels based on a sharp interface theory**

by Johannes Ruck

Karlsruher Institut für Technologie Institut für Technische Mechanik Bereich Kontinuumsmechanik

Modeling martensitic phase transformation in dual phase steels based on a sharp interface theory

Zur Erlangung des akademischen Grades eines Doktor-Ingenieurs von der KIT-Fakultät für Maschinenbau des Karlsruher Instituts für Technologie (KIT) genehmigte Dissertation

von Dipl.-Ing. Johannes Ruck

Tag der mündlichen Prüfung: 14. Mai 2020 Hauptreferent: Prof. Dr.-Ing. Thomas Böhlke Korreferent: Prof. Dr.-Ing. Andreas Menzel

#### **Impressum**

Karlsruher Institut für Technologie (KIT) KIT Scientific Publishing Straße am Forum 2 D-76131 Karlsruhe

KIT Scientific Publishing is a registered trademark of Karlsruhe Institute of Technology. Reprint using the book cover is not allowed.

www.ksp.kit.edu

*This document – excluding the cover, pictures and graphs – is licensed under a Creative Commons Attribution-Share Alike 4.0 International License (CC BY-SA 4.0): https://creativecommons.org/licenses/by-sa/4.0/deed.en*

*The cover page is licensed under a Creative Commons Attribution-No Derivatives 4.0 International License (CC BY-ND 4.0): https://creativecommons.org/licenses/by-nd/4.0/deed.en*

Print on Demand 2021 – Gedruckt auf FSC-zertifiziertem Papier

ISSN 2192-693X ISBN 978-3-7315-1072-7 DOI 10.5445/KSP/1000128076

# **Modeling martensitic phase transformation in dual phase steels based on a sharp interface theory**

**Zur Erlangung des akademischen Grades**

**Doktor der Ingenieurwissenschaften**

**der KIT-Fakultät für Maschinenbau**

**Karlsruher Institut für Technologie (KIT)**

**vorgelegte**

**Dissertation**

**von**

**Dipl.-Ing. Johannes Ruck**

**Tag der mündlichen Prüfung: 14.05.2020 Hauptreferent: Prof. Dr.-Ing. Thomas Böhlke Korreferent: Prof. Dr.-Ing. Andreas Menzel**

## **Zusammenfassung**

Die vorliegende Arbeit behandelt die thermomechanische martensitische Phasentransformation in austenitischen Stählen mit niedrigem Kohlenstoffgehalt, die z.B. für Strukturbauteile in der Automobilindustrie von großer Bedeutung sind. Diese Stahlsorten bilden bei der Umwandlung von Austenit in Martensit vorwiegend Lattenmartensit. Die martensitische Umwandlung wird durch eine kontinuumsmechanische scharfe Grenzflächentheorie im Kontext großer Deformationen modeliert. Die treibende Kraft der scharfen Grenzfläche wird auf eine Laminatsubstruktur aufgebracht. Plastische Effekte werden durch ein Kristallplastizitätsmodel in allen Laminatphasen abgebildet. Die Vererbung von Versetzungen im Austenit wird durch einen einfachen Ansatz berücksichtigt. Dieser erlaubt die Vererbung eines konstanten Versetzungsdichteanteils vom Austenit in den sich neu bildenden Martensit. Zusätzlich wird zur Modellierung der Dynamik der scharfen Grenzfläche eine thermodynamisch konsistente kinetische Beziehung eingeführt. Die Beziehung von scharfen Grenzflächenmodellen, die auf eine Laminatstruktur angewendet werden und Modellen basierend auf einer Rang-1 Energierelaxation wird diskutiert. Um Nichtgleichgewichtszustände im Rahmen einer Energieminimierung zu beschreiben, wird das Prinzip der minimalen dissipativen Leistung verwendet. Dies ermöglicht die Plastizität und die kinetische Beziehung für die Austenit Martensit Grenzfläche über eine Minimierung auszudrücken. Die numerische Implementierung des entwickelten Models als User-Subroutine in die kommerzielle FEM Software ABAQUS wird vorgestellt. Numerische Anwendungen des Models zeigen die Möglichkeiten das experimentell beobachtbare thermische Transformationsverhalten in Dualphasenstählen qualitativ zu beschreiben. Darüberhinaus lassen sich basierend auf dem Model drei verschiedene Transformationsstadien, die von verschiedenen physikalisch Effekten dominiert werden, ableiten.

## **Summary**

The present works deals with the thermomechanical martensitic phase transformation of low carbon steels, which is of high significance in e.g., automotive structural applications. In this category of steels predominantly lath martensite forms. The martensite transformation is modeled based on a continuum mechanical sharp interface theory in the context of finite deformations. The sharp interface driving force is applied to a laminate substructure. Plastic effects are captured by a crystal plasticity model in all laminate phases. Inheritance of austenite dislocations is accounted for by a very simple approach, allowing a constant dislocation fraction to be inherited by the newly forming martensite. Additionally, a kinetic relation using a simple thermodynamical consistent approach is introduced. The relation of sharp interface models applied to a laminate and models based on rank 1 energy relaxation is discussed. To also formulate non-equilibrium states in a minimization framework the principle of minimum dissipative power is used. This allows to include plasticity and the kinetic relation for the austenite − martensite interface. The numerical implementation of the model as a user subroutine in the commercial FEM software ABAQUS is presented. Numerical applications of the model show its ability to qualitatively capture the experimentally observed thermal transformation behavior of low carbon austenite into martensite in dual phase steels. Furthermore, the model predicts a martensite transformation, where three different stages with different dominating mechanisms can be distinguished.

## **Acknowledgments**

This thesis originated during five challenging and interesting years as an accademic employee at the Chair of Continuum Mechanics. Throughout this time many people contributed in various ways to this thesis and I am deeply grateful to them.

First and foremost I am grateful to Prof. Thomas Böhlke for attracting my interest to the field of continuum mechanics during my early undergraduate studies, for offering me the possibility to work in his research group at the Institute of Engineering Mechanics later on, for his advise, the fruitful discussions, the various interesting topics I could work on and his enduring support during and after my time at the institute.

Furthermore, I would like to greatly thank Prof. Andreas Menzel for the co-review of my thesis and his valuable comments.

I would also like to thank the members of the Graduate School 1483 for the extensive and interesting discussions. Moreover, I greatly appreciate the funding of the German Research Foundation (DFG) for my project within the Research Training Group GRK 1483.

A big thank goes to my colleagues at the institute for sharing countless discussions as well as enjoyable non-scientific conversations and unforgettable moments at conferences and summer schools. Especially, I would like to thank my office mate Hannes for many scientific and private conversations and for supporting the same - recently not so successfull - football team.

Next, I thank Ute, Helga and Tom for taking care of organizational

matters be it conferences, travel, teaching or the IT infrastructure. Also, you always had an open ear and contributed a lot to the warm and pleasant atmosphere at the institute.

Last but not least I would like to thank my family and friends for their support and the necessary distraction they provided me outside the scientific field.

## **Contents**




# **Chapter 1 Introduction**

### **1.1 Scope of the work and outline**

The scope of this work is to model the thermomechanical phase transformation behavior of low carbon austenite into martensite as it is observed in dual phase steels. The martensite transformation model is based on a continuum mechanical sharp interface theory accounting for finite deformations. Additionally, to model non-equilibrium processes the kinetic of the interface needs to be resolved in time. Since plasticity plays an important role in the formation of lath martensite, which is the preferred morphology of martensite in low carbon steels, the model has to account for plastic effects inside the austenite and martensite.

A short introduction into the martensitic phase transformation in carbonated steels is given in chapter 1.3, followed by an overview of existing modeling approaches in chapter 2. In a continuum mechanical sense the phase transformation can be understood as the motion of a singular surface, separating two different solid materials, through a body. Therefore, the fundamental principles of continuum mechanics and the treatment of bodies with singular surfaces is discussed in chapter 3. The martensitic microstructure can often be approximated by laminates, which is a simple method to homogenize a heterogeneous material. Hereby, the different phases are separated by planar interfaces. Therefore, chapter 4 deals with the properties of laminates

with and without mobile interfaces. The formation of martensite from a physical point of view can be understood as the product of an energy minimizing process inside the material. Mathematically, this fact is related to relaxation which deals with the minimization of non-convex potentials by minimizing sequences. The minimizing sequences can be interpreted as the microstructure of the material. A short insight into this topic is given in chapter 5. Interestingly, the sharp interface theory can be related to the rank 1 relaxation. This and a phase transformation model based on a laminate approach is discussed in chapter 6. In chapter 7 the basic phase transformation model is specified to the characteristics of the martensitic phase transformation in low carbon steels. To distinguish between martensite which has already been transformed and transforming martensite a three phase rank 1 laminate is introduced. The numerical implementation of the governing equations on the material point level and into the FEM package Abaqus as a UMAT is discussed in chapter 8. In chapter 9 at first the convex and different rank 1 convex envelopes are compared with respect to the rate independent elastic phase transformation. This gives an insight into the difference between the convex and approximations of rank 1 convex envelopes. Afterwards, numerical examples of the transformation model discussed in chapter 7, which adds dissipation to the phase transformation and plasticity to all phases, are presented in the context of the laminate theory. The sensibility of the transformation model with respect to the material parameters is discussed. This is followed by finite element application of the phase transformation model to a laminate structure comprising ferrite and austenite, representing a dual phase steel. In this context the effect of the sorrounding ferrite matrix on the austenite martenstie transformation is investigated. At last the finite element results of the transformation behavior of a synthetic ferrite austenite microstructure subjected to rapid cooling are presented.

**Originality of the work** To the author's knowledge the following points represent aspects which have not been published before.


## **1.2 Notation**

A direct symbolic tensor notation is preferred throughout the text. Tensor components are expressed by latin indices and Einstein's summation convention is applied. Vectors and second-order tensors are denoted by lowercase and uppercase bold letters, *a* and *A*, respectively. Fourth-order tensors are written as double uppercase letters A. In general, any second and higher-order tensor can be written as *A*h*r*<sup>i</sup> , where *r* indicates the tensor rank. The second-order identity tensor is written

as *I*. The fourth-order identity on second-order tensors is denoted by I = *I*✷*I* and its symmetric part by I S . Completely symmetric and traceless, i.e. irreducible tensors are denoted with a prime, e.g., *A* ′ . Spherical counterparts are denoted by *A* ◦ .

Arbitrary vectors *a* and *b*, second-order tensors *A*, *B* and *C* and the fourth-order tensor C are used in the following definitions. The composition of two second-order or two fourth-order tensors is formulated by *C* = *AB* and C = AB. A linear mapping of second-order tensors by a fourth-order tensor is written as *A* = C[*B*]. The scalar product is denoted by *A* · *B* and A · B, respectively. We define the composition ✷ via (*A*✷*B*)[*C*] = *ACB*, the dyadic product ⊗ as (*A* ⊗ *B*) [*C*] = (*B* · *C*) *A*, and the contraction · with (*a* ⊗ *b*) · (C[*a* ⊗ *b*]) = (*a* ⊗ *a*) · (C[*b* ⊗ *b*]). The gradient operator with respect to the reference position is denoted as Grad (·) while the operator with respect to the current position is grad (·). Analogously the divergence operators are defined as Div (·) and div (·). The trace of a second order tensor is denoted by tr(·) and the determinant by det(·). The cofactor of a tensor cof(·) is given by cof(*A*) = det(*A*)*A* -T. The symmetric and skew symmetric parts of a tensor are written as sym(·) and skw(·). The Rayleigh product *⋆* for a second order tensor *F* and fourth order tensor C is defined as *F ⋆*C = *FimFjnFkoFlpCmnopei*⊗*e<sup>j</sup>* ⊗*e<sup>k</sup>* ⊗*e<sup>l</sup>* . The product *⋆* 2 (·) between a second order tensor *F* and fourth order tensor C is defined as *F ⋆* <sup>2</sup>C = *FjmFlnCimkne<sup>i</sup>* ⊗ *e<sup>j</sup>* ⊗ *e<sup>k</sup>* ⊗ *e<sup>l</sup>* .

Throughout the work the right gradient of a tensor valued function is used which increases the order of the respective function by one

$$\frac{\partial f(x)}{\partial x} = \frac{\partial f(x\_i)}{\partial x\_j} \mathbf{e}\_i \otimes \mathbf{e}\_j. \tag{1.1}$$

In a similar way the divergence is defined by

$$\operatorname{div}(f(x)) = \frac{\partial f(x)}{\partial x} \cdot I = \frac{\partial f(x\_i)}{\partial x\_i} \mathbf{e}\_i \otimes \mathbf{e}\_i. \tag{1.2}$$

The jump of a quantity *A* is defined as [*A*] = *A*<sup>+</sup> − *A*<sup>−</sup>. The average value of a quantity *A* is given by h*A*i = 1*/*2(*A*<sup>+</sup> + *A*<sup>−</sup>).

### **1.3 Martensite in ferritic alloys**

**General characteristics of martensitic transformations** The ability of carbonated steels to transform upon rapid cooling into a material with a comparable higher strength and hardness has been known and used for over a thousand years. This phenomenon is due to the transformation of a high temperature phase austenite into a low temperature phase martensite, accompanied by a local change of the microstructure in terms of the crystal symmetry. The attractive mechanical property of such a material makes steel still today an important material, especially for structural and highly loaded components, e.g., as TRIP and dual phase steels in cars.

Martensite is a thermodynamical non-equilibrium phase (Gottstein, 2013), which forms in carbonated, austenitic steels through rapid cooling after falling below the martensite start temperature *MS*. The necessary cooling rates exceed around 100 K/s (see, e.g., Olasolo et al., 2011). The rapid change in temperature prevents diffusion of the interstitial carbon inside the austenite and thus the diffusion controlled transformation into ferrite, bainite, perlite or cementite. This, in consequence leads to a distorted lattice of the parent austenite which then transforms into martensite. While the parent phase austenite is in most carbonated steels face centered cubic (fcc) the crystal symmetry of the transformed martensite depends on the chemical composition of the steel or the chosen material for non-ferrous materials. These are in Fe-C alloys for example fcc to body centered cubic (bcc) or body centered tetragonal distorted (bct), depending on the carbon content, and fcc to bcc in Fe-Ni allyos (see, e.g., Nishiyama, 1978; Umemoto et al., 1983). Non ferritic martensitic transformation are observed in e.g. Co as fcc to hexagonal closed packed (hdp), in Cu-Ni-Al shape memory alloys as cubic to orthorombic, cubic to monoclinic in Ni-Ti and cubic to tetragonal in In-Tl Bhattacharya (2004).

Despite their different crystallographic properties, all martensitic phase transformations share some common features. Martensitic phase transformations are characterized by a cooperative movement of atoms over a range less than an interatomic distance, (see, e.g., Clapp, 1995; Pereloma and Edmonds, 2012). Such a phase transformation is termed displacive (Christian et al., 1995). The phase transformation process occurs in a very short period of time, thus suppressing diffusion. The chemical composition of austenite and martensite are therefore practically identical and the transformation process itself can be considered as adiabatic (Nishiyama, 1978; Fischer et al., 1994). To allow such a cooperative movement of atoms with high speeds the interface between austenite and martensite, also called habit plane, must be highly mobile and coherent. Hence, it should by an unrotated and undistorted plane (Wechsler et al., 1953; Mackenzie and Bowles, 1954; Bhadeshia et al., 2009).

The change in crystal symmetry from austenite to martensite is characterized by a homogeneous deformation, converting the austenite into the martensite lattice (Christian et al., 1995). The transformation strain, which is mainly deviatoric with only small changes in volume, leads to a macroscopic shape change of the transformed prior austenite regions (Bhadeshia et al., 2009; Pereloma and Edmonds, 2012). The transformation strain is unique for each material with a particular chemical composition. A first attempt to describe this characteristic transformation strain was made by Bain and Dunkirk (1924) in the case of an fcc to bct transformation. Bain identified the transformation strain by a stretch tensor *U*, which transforms the austenite lattice vectors into the martensite lattice vectors with minimal atomic movement. For the considered fcc to bct transformation three different stretch tensors, related to each other by appropriate rotations, can be distinguished (see

also Fig. 1.1)

$$\begin{aligned} \mathbf{U}\_1 = \begin{pmatrix} \gamma\_1 & 0 & 0 \\ 0 & \gamma\_1 & 0 \\ 0 & 0 & \gamma\_2 \end{pmatrix}, \mathbf{U}\_2 = \begin{pmatrix} \gamma\_1 & 0 & 0 \\ 0 & \gamma\_2 & 0 \\ 0 & 0 & \gamma\_1 \end{pmatrix}, \mathbf{U}\_3 = \begin{pmatrix} \gamma\_2 & 0 & 0 \\ 0 & \gamma\_1 & 0 \\ 0 & 0 & \gamma\_1 \end{pmatrix}. \end{aligned} \tag{1.3}$$

The elongation is represented by *γ*<sup>1</sup> while *γ*<sup>2</sup> is identified with the compression of a lattice vector. The transformation of the lattice vectors is then identified with

$$
\mathbf{e}\_M = \mathbf{U}\_i \mathbf{e}\_A. \tag{1.4}
$$

#### **1.4 Crystallography of martensite**

**Lattice orientation relationships** As martensite forms within austenite specific crystallographic relations between the respective lattices can be established. They relate two parallel planes and directions of the austenite and martensite lattice to each other (Bain and Dunkirk, 1924; Nishiyama, 1978; Kurdjumow and Sachs, 1930). The most prominent orientation relationships for face centered cubic austenite (*A*) to body centered martensite (*M*) are summarized in Table 1.1 and visualized in Fig. 1.1, where the corresponding planes and directions of the austenite and martensite are colored in red and blue, respectively. As opposed to


**Table 1.1:** Orientation relationships for austenite/martensite

**Figure 1.1:** Transformation mechanism according to Bain and Dunkirk (1924) a), Nishiyama (1978) b) and Kurdjumow and Sachs (1930) c). On the left two austenite cells and a corresponding martensite cell are depicted. The parallel planes and directions are highlighted in blue, respectively red. In the center the martensite cell and the corresponding transformation mechanism is visualized. The transformed martensite is shown on the right.

the Bain mechanism which accomplishes the transformation through particular stretches, Kurdjumow and Sachs (1930) and Nishiyama (1978) proposed a shear controlled transformation along the {111}*<sup>A</sup>* plane in directions *<* 01¯1 *><sup>A</sup>* (KS) and *<* 11¯2 *><sup>A</sup>* (NW), see also Fig. 1.1. Another distinct feature between the Bain, Nishiyama-Wassermann and Kurdjomov-Sachs orientation is the number of martensite variants. The Bain construction predicts three crystallographic distinguishable variants. KS predicts 24 of them due to four parallel closest packed planes and six parallel directions, see also Table 1.2. NW distinguishes

between twelve different variants due to the same four closest packed planes as KS with three closest packed parallel directions. Other orientation relationships have been introduced by Greninger and Troiano (1949) and Pitsch (1959). While the Bain orientation is experimentally not observed, most steels exhibit orientation relationships between KS and NW<sup>1</sup> . Therefore, a model resulting in a continuous description of orientation relationships between NW and KS has been proposed by Cayron (2015). The interface or habit plane between austenite and martensite is another distinct feature which takes on specific orientations for different kinds of martensite, see also Section 1.5. Under the assumption that the interface is coherent and represents an invariant plane, i.e. an undistorted and unrotated plane, the Phenomenological theory of martensite crystallography (PMTC) (Wechsler et al., 1953; Mackenzie and Bowles, 1954) derives transformation deformations and interface orientations for particular martensite transformations. The transformation deformation is considered as an invariant plane strain, i.e. a shear deformation along the habit plane. The transformation deformation is based on the Bain stretch tensor, which is not an invariant plane strain. This is corrected by a rotation *R*. A lattice invariant shear *S* is introduced, produced by either slip or twinning, to account for the correct macroscopic shape change

$$F\_T = RUS.\tag{1.5}$$

The transformation deformation can be equivalently expressed as

$$F\_T = I + \gamma d \otimes N = I + a \otimes N,\tag{1.6}$$

which allows to determine the habit plane orientation *N*. This representation is nothing more than the Hadamard compatibility condition

<sup>1</sup> The KS and NW orientation relationships assume the transformation path along closest packed planes (KS) and directions (NW) which are a much better choice in terms of a coherent interface between austenite and martensite compared to the Bain orientation.


**Table 1.2:** 24 KS martensite variants (Mishiro 2013)

between martensite and austenite in terms of the inelastic deformation contributions<sup>2</sup> , excluding plastic deformation inside the austenite. Modifications of the PMTC have been introduced by Kelly (1992) in form of a

<sup>2</sup> Since coherency of a planar interface, or compatibility of two adjacent deformation gradients is defined by the jump of the total deformation gradient, Eq. (1.6) can be interpreted as a special case in the absence of elastic deformations. This also renders the jump of the elastic stress tensor zero. In the case of a geometrical linear theory and no inelastic austenite strains it can be shown that Eq. (1.6) satisfies force equilibrium and compatibility at the interface for identical elastic moduli of austenite and martensite.

second lattice invariant shear. This allows to determine the habit plane orientation of lath martensite, which could not be reproduced before.

### **1.5 Morphology of martensite**

Martensite can appear in various different microstructural shapes and arrangements in ferrous alloys. The two most common ones are namely lath and plate martensite.

**Lath martensite** forms in steels with carbon contents up to 0.6-0.8% (Krauss and Marder, 1971; Stormvinter et al., 2011). Lath martensite shows a characteristic multilevel substructure, see Fig. 1.2, which on the lowest scale is identified with martensite laths. Martensite laths consist of martensite single crystals with a typical width of less than 1 *µ*m and show high dislocation densities (Kitahara et al., 2006). Several laths with the same crystallographic orientation, i.e. the same Bain variant, aggregate to blocks. In turn, several blocks with the same {111}*<sup>A</sup>* plane with respect to the parent austenite constitute a packet. Since there are only four crystallographic equivalent planes of martensite with respect to the {111}*<sup>A</sup>* plane in austenite, a maximum of four different orientated packets inside a prior austenite grain is possible. The orientation of lath martensite is approximately KS, Maki (1990). The packets and blocks show a tendency to become finer with increased carbon content (Maki et al., 1980). The large transformation deformations are accommodated by plastic slip, which results in very high dislocation densities of the order 10<sup>14</sup> − 10<sup>15</sup> 1/m<sup>2</sup> (Morito et al., 2003). The observed habit planes for lath martensite are most commonly of type {111}*<sup>A</sup>* and {557}*<sup>A</sup>* (Pereloma and Edmonds, 2012).

**Plate martensite** forms in high-carbon as well as high-nickel ferrous alloys (Krauss and Marder, 1971). The first plates to form partition the prior austenite grain, limiting the size of the subsequently forming plates. The partitioning of the prior austenite grain is visible in the so

**Figure 1.2:** Substructure of martensite on the left, KS variants of lath martensite inside a former austenite grain (Kitahara et al., 2006)

called midrib, which is evident for most plate martensite (Krauss and Marder, 1971). In contrast to lath martensite, where the martensite laths inside a packet are arranged in parallel layers, adjacent plates form at angles to each other. The transformation deformations are thus accommodated by twinning, also referred to as self accommodation in contrast to plastic accommodation for lath martensite. Hence, plate martensite exhibits a fine microstructure of twins (Kelly and Nutting, 1961). For moderate martensite start temperatures lenticular shaped plates arise. In that case only the area close to the midrib is twinned. Further away from the midrib partially twinned regions with dislocations and untwinned regions with high dislocation densities are observed. In the case of low martensite start temperatures and high carbon concentration, like Fe-Ni-C shape memory alloys , thin plate martensite forms. Thin plate martensite is characterized by a highly mobile and planar interface (Maki et al., 1973). This allows also the reverse martensiteaustenite transformation under heating or unloading as observed in shape memory alloys. The dislocation density inside the martensite and parent austenite is low. The orientation of the habit plane for plate martensite is in most cases of type {259}*A*, {225}*<sup>A</sup>* or {3 10 15}*A*.


**Table 1.3:** Particular martensite types (Maki, 1990)

**Chapter 2**

# **State of the art modeling martensitic phase transformation**

## **2.1 Classification of martensite modeling approaches**

**Modeling martensite transformation** Martensite transformation is mainly driven by the mechanical, e.g. stress as well as strain, and thermal state of a material. It appears with roughly two different characteristics. Plate martensite, which exhibits a fine scale microstructure consisting of twin related martensite variants. Due to the in most cases insignificant plastic flow the change in material properties is mainly attributed to the transformation strains of the respective martensite variants (Cherkaoui and Berveiller, 2000). This allows to model the overall material behaviour as thermoelastic. On the contrary, lath martensite in high strength steels is accompanied by a significant amount of plastic strain, caused by the locally large magnitude of the transformation strains. This establishes on the local, micromechanical level significantly different stress and strain states compared to the macroscopic fields. Furthermore, the transformational and plastic strain fields interact which leads to the Magee (Magee, 1970) (variant selection) and Greenwood-Johnson effect (Greenwood and Johnson, 1965) (plastic flow despite macroscopic stresses which are

lower than the yield stress). To adequately capture the characteristics of the martensitic phase transformation a large variety of models have been developed. A first attempt of a classification is given by the choice of the scale on which these models are established. In the framework of continuum mechanics these are the microscopic and macroscopic scale. Models based on the microscale take into account material properties and physical processes on the crystallographic level. Micromechanical phase transformation models can be further subdivided into phenomenological or thermodynamical approaches related to internal variables and models based on the notion of a driving force, governing the motion of an interface separating two solid phases with different elastic properties and eigendeformations. In contrast, macroscopic models are concerned with the description of the material properties and physical processes on the continuum level. Therefore, they are more suitable to model large structural applications where a detailed understanding of the micromechanical fields and processes is redundant. Models based on atomistic (e.g., Maresca and Curtin, 2017) and molecular dynamic approaches deliver a detailed understanding of the local processes or actions, which induce or accompany martensitic phase transformation. For the latter a detailed review in the context of shape memory alloys can be found in Cisse et al. (2016). However, the present work deals with the continuum mechanical description of the austenite to martensite transformation, which puts these models out of scope. They will not be considered in the following discussion.

## **2.2 Macroscopic or phenomenological phase transformation models**

**Modeling stress or strain induced martensite transformation** Macroscopic or phenomenological models often distinguish between thermal

and stress or strain induced martensitic phase transformation. Actually, as will be discussed later on, those contributions are not separable in a clear way and the phase transformation is rather driven by a combination or interaction of them. One of the earliest models in the isothermal, mechanical context was proposed by Olson and Cohen (1975). It accounts for strain induced martensite transformation. The model is based on the physically motivated assumption that the probability of martensite nucleation is related to plastic shear bands which act as nucleation sites. A later modification of Stringfellow et al. (1992) allowed to account for the present stress state by a self-consistent approach using spherical martensite inclusions. These models often serve as the starting point for various other macroscopic phase transformation models (see, e.g., Olson and Cohen, 1982; Tomita and Iwamoto, 1995). Tomita and Iwamoto (1995) included the strain rate influence on toughness and ductility as well as deformation mode dependence and effects due to thermomechanical coupling, grain size and stress triaxiality (Iwamoto et al., 1998; Iwamoto and Tsuta, 2000; 2002). Application of such models in the context of Finite Element simulations can be found in (Papatriantafillou et al., 2006; Sierra and Nemes, 2008). Both publications investigate the effect of retained austenite on the work hardening behaviour of TRIP steels, where the mechanical induced martensite transformation is modeled according to Stringfellow et al. (1992). Martensite transformation in high strength steels is accompanied by plastic deformation, even in cases where the macroscopic stress is well beyond the yield limit. This is due to the locally large transformation deformations and is known as the Greenwood-Johnson effect (Greenwood and Johnson, 1965). Martensitic phase transformation models accounting for the Greenwood-Johnson effect have been proposed by e.g. (Greenwood and Johnson, 1965; Leblond et al., 1985; 1989; Hallberg et al., 2007). The original one dimensional model by Greenwood and Johnson (1965) was substantially improved by Leblond et al. (1989) taking into account three dimensional stress, respectively strain contributions assuming ideal plastic phases and a small strain theory.

**Modeling temperature induced martensite transformation** A model for the temperature induced martensite transformation in iron based alloys has been proposed by Koistinen (1959), who modeled the martensite evolution solely as a function of the undercooling below the martensite start temperature. It is often applied in the simulation of the hot stamping or press hardening process as e.g. (Lee et al., 2010; Hippchen et al., 2016; Neumann, 2017). However, experiments show that the martensite start temperature depends on the applied stress and moreover on the chemical composition of the prior austenite (Lee and Van Tyne, 2012). Therefore, several extensions of the Koistinen-Marburger model have been suggested to include additionally the chemical composition (see, e.g., Van Bohemen and Sietsma, 2009) and a stress dependent martensite start temperature (see, e.g., Fischer et al., 1998). Additionally, to capture the S-shaped martensite evolution non linear extensions of the Koistinen-Marburger model have been suggested (see, e.g., Neumann and Böhlke, 2016). Further applications of the Koistinen-Marburger model can be found in a finite strain multi phase transformation model by Mahnken et al. (2012) and a finite strain martensitic phase transformation model by Hallberg et al. (2007). Both models are set in the thermodynamic, energetic framework of Generalized Standard Materials.

### **2.3 Micromechanical transformation models**

**Sharp interface approaches** Models considering an interface with a vanishing thickness are denoted as sharp interface approaches. An expression of a sharp interface driving force based on considerations of conservation laws across the singular surface and the second law

of thermodynamics can be found in (Abeyaratne and Knowles, 1990; Silhavy, 1997; Levitas, 1998; Maugin, 1998), where the sharp interface driving force *fSI* is identified with *fSI* = [*ρ*ˇΨ] − h*σ*ˇi · [*F*]. From this representation it is already visible that the sharp interface driving force predicts the motion of the interface by a combination of stress and thermal contributions, which cannot be separated. The driving force is an energetic quantity which also yields a local transformation criterion. Early approaches to model phase transformation based on an energetic quantity or driving force were introduced by Patel and Cohen (1953); Roitburd and Temkin (1986); Fischer et al. (1994) and Levitas (1998). Patel and Cohen (1953) derived in the small strain context an energetic phase transformation criterion given by *f* = [*ρ*Ψ] − *σ*¯ · [*ε*]. This expression also identifies the driving force with the jump of the Helmholtz free energy and a second contribution by the product of a stress and the jump of the strain tensor. In contrast to Abeyaratne and Knowles (1990) the stress is hereby assumed as the effective stress. Hence, this representation coincides with the sharp interface driving force for identical stresses in both phases. Roitburd and Temkin (1986) defined a transformation condition for an elastic two phase material without any external stresses as the difference of the Helmholtz free energy between the two phases. The driving force is thus identified with *f* = [*ρ*ˇΨ] and coincides with the driving force of Patel and Cohen (1953) in the special case of a vanishing external stress and with that of Abeyaratne and Knowles (1990) in the case of a vanishing average stress inside the phases. Fischer et al. (1994) presented a derivation of a driving force as the result of the variation of the volume integrated Gibbs free energy, accounting for the jump of the integrands across the interface separating two different phases. Under the assumption of identical elastic moduli of the phases and small strains the driving force reads *f* = [*ρ*Ψ*chem*] − h*σ*i · [*εin*], where *εin* characterizes inelastic strain contributions. Due to the identical elastic moduli the jump of the mechanical, elastic Helmholtz free energy equals the scalar product

of the jump of the elastic strain and the average stress. Thus, the driving force by Fischer et al. (1994) coincides with the sharp interface driving force, in the case of small strains. Further decomposing the jump of the inelastic strains [*εin*] into transformational and plastic parts Fischer et al. (1994) associated the term driving force only with the transformational contributions. A phase transformation criterion was accordingly identified with *f* = [*ρ*Ψ*chem*] − h*σ*i · [*εt*] = *f<sup>c</sup>* + h*σ*i · [*εt*], equating the driving force to a reaction force comprising a threshold energy *f<sup>c</sup>* and a contribution of the average stress and the jump of the plastic strain. In a later publication this concept was dropped and the driving force was identified with the full expression.

Besides a thermodynamical consistent description of the motion of the phase boundary the microstructure of the material has to be accounted for. While low alloyed, low carbon high strength steels in most cases exhibit lath martensite, high carbon steels as well as shape memory alloys produce (multiple) twinned martensite microstructures. These microstructures can approximately be described by laminate models. Especially in the field of shape memory alloys a large class of models uses laminate based approaches to represent the microstructure (see, e.g., Kruzík, 1998; Stupkiewicz and Petryk, 2002; Kružík and Luskin, 2003; Aubry et al., 2003; Bartel and Hackl, 2009; Bartel et al., 2011). In these models the laminate substructure acts as a tool to regularize the, due to the different eigendeformations of the individual phases, non-convex effective free energy of the laminate through an approximation of the rank 1 convex envelope. This relates laminate models to the mathematical theory of relaxation, which deals with the convexification of non-convex potentials by finding suitable convex envelopes, establishing a lower bound of the non-convex potential, (see, e.g., Dacorogna, 2008; Schröder and Hackl, 2014). The development of a microstructure, respectively a transformation of phases, is hereby connected with the loss of quasiconvexity of the corresponding effective Helmholtz free energy. The evolution of phases with different eigen-

deformations represents minimizing sequences of the corresponding boundary value problem as e.g., shown in Ball and James (1987) in the case of stress free transformations. The main goal of the relaxation is to find an approximation as close as possible to the quasiconvex envelope, which is known to be the weakest notion of convexity in finite hyperelasticity (Morrey et al., 1952). The analytical computation of the generalized convex envelopes is only possible in rare cases. Examples are Le Dret et al. (1985) for the quasiconvex envelope of the St. Venant Kirchhoff energy, Conti and Theil (2005); Conti et al. (2013) for the application to single crystal plasticity and DeSimone and Dolzmann (2002) in case of the phase transition of nematic polymers. The usual strategy hereby is to find with the polyconvex envelope a lower and the rank 1 convex envelope an upper bound for the quasiconvex envelope. These can then be used to verify the quasiconvex envelope. In most other cases this task is accomplished numerically. This procedure has been applied to a large class of continuum mechanical problems. Application to damage mechanics can be found in (Gürses and Miehe, 2011; Balzani and Ortiz, 2012; Junker et al., 2017), where most relaxations result in the convex envelope. A one dimensional model using the quasiconvex envelope was proposed by Schmidt-Baldassari and Hackl (2003). In crystal plasticity the rank 1 relaxation is often chosen to account for a deformation compatible microstructure of the material, which provides an approximation of the observed microstructure. Investigations of rank 1 relaxation in single slip crystal plasticity were done by (Ortiz and Repetto, 1999; Carstensen et al., 2002; Miehe and Lambrecht, 2003) for fixed laminate orientations. Bartels et al. (2004) investigated relaxation with respect to first and second order laminates and evolving laminate orientations for finite strains. Remarkably, a comparison between the second order laminate and the polyconvex envelope yielded, apart from discretization errors, the same energy. A non-local extension to include interface energy was introduced by Conti and Theil (2005). To model the evolution of

the considered microstructure in time Kochmann and Hackl (2011) proposed a relaxation where the variables are split into an elastic (equilibrium) and a non-elastic (time-dependent) part. To model the martensitic phase transformation laminates constitute a suitable choice to construct a deformation compatible upper bound of the quasiconvex envelope. Furthermore, the relaxation conditions lead to the equilibrium conditions of a sharp interface on each laminate level. Thus models based on laminate, respectively rank 1 relaxation satisfy the sharp interface driving force, as will also be discussed later.

First investigations of elastic relaxation approaches with respect to phase transformation problems in the finite strain context can be found in Kruzík (1998); Kružík and Luskin (2003) where a double well energy is relaxed by first and second order laminates. Aubry et al. (2003) applied an elastic relaxation to the cubic to orthorombic phase transformation of a Cu-Ni-Al shape memory alloy using sequential laminates. By considering the interface energy between austenitemartensite interfaces and martensite twins, an infinitely fine microstructure is precluded. Furthermore, taking into account interfacial energy contributions different microstructures in loading and unloading are realized. This reproduces the typical stress hysteresis of shape memory alloys. Stupkiewicz and Petryk (2002) presented a laminate approach to model shape memory alloys using small strains and a rank 2 laminate, representing martensite twins. The orientation of the twin laminate as well as the volume fraction of each martensite variant was prescribed based on the crystallographic theory of martensite. The evolution of the averaged martensite twin and the austenite phase on the upper laminate level is modeled according to the driving force by Abeyaratne and Knowles (1990). Hence this approach is equivalent to a partial relaxation of the effective Helmholtz free energy with respect to the interface jump vector and martensite volume fraction on the upper laminate scale. Bartel and Hackl (2009) proposed a rank 1 relaxation approach using small strains. To model the phase transformation of

shape memory alloys with a theoretically arbitrary number of variants, a special structured rank 2 laminate was introduced. By splitting the local variables into elastic and inelastic contributions the effective Helmholtz free energy is partially relaxed with respect to the elastic variables, while keeping the laminate orientation fixed after initiation of the phase transformation. The evolution of the remaining inelastic variables is determined by dissipation potentials. This allows to model the stress hystereses of shape memory alloys in a different way to Aubry et al. (2003). An extension to also include plasticity was proposed in Bartel et al. (2011) using a rank 1 laminate and the same strategy in terms of decomposing the internal variables in elastic and inelastic contributions. Additionally, by treating the laminate orientation as an inelastic variable a history dependent spatial evolution of the laminate is made possible. Using the convex envelope Hackl and Heinen (2008a) proposed a polycrystal, small strain model for multivariant shape memory alloys. The convex relaxation still preserves the general form of the driving force, but leads to a non deformation compatible microstructure. The stress hystresis of shape memory alloys was reproduced by a similar elastic-inelastic split of the internal variables as in Bartel and Hackl (2009) and Bartel et al. (2011).

A more realistic approach in terms of an accurate microstructure is given by the small strain martensitic phase transformation model of (Cherkaoui et al., 1998; Cherkaoui and Berveiller, 2000). The evolution of martensite plates in an austenite matrix is modeled by homogenization of a spherical martensite inclusion subjected to the sharp interface driving force. Plasticity inside the austenite phase as well as martensite phase is taken into account, which allows to reproduce the Greenwood-Johnson effect. The plastic strain of the austenite phase is upon transformation fully inherited by the corresponding martensite variant.

**Diffuse interface models** Despite being the most prominent approach to model micromechanical martensitic phase transformation the contribution of phase field models to martensitic phase transformation is kept rather short at this place, since the present work focuses on sharp interface models. Phase field models assign the interface a finite thickness, which makes the corresponding fields and internal variables continuous across the interface. This is usually achieved by one or multiple order parameter(s) which vary continuously between 0 and 1, where the lower and upper bound identify different, pure phases. This makes an explicit tracking of the position of the interface redundant and allows to use a regular grid. Phase field models are commonly based on a Ginzburg-Landau type polynomial which chracterizes a multiwell potential for the order parameter (Landau, 1965). The order parameter is then coupled to the field variable of the global problem, which usually is the strain or the displacement, respectively. Additionally the gradient of the order parameter is introduced which represents an interface energy contribution and regularizes the due to the multiwell potential otherwise ill-posed problem. First approaches to model phase transformation with phase field models in the small strain elastic case have been proposed by Chen et al. (1992) and Wang and Khachaturyan (1997). A first approach to regularize a sharp interface theory by using an order parameter was presented in a large strain setting by Grach and Fried (1997) based on the theory developed in Fried and Gurtin (1996). Later on, plasticity was included and large strain frameworks were introduced (see, e.g., Yamanaka et al., 2008; Levitas et al., 2009). While in most models only the three Bain variants are considered, Du (2017) modeled martensitic phase transformation using the 24 Kurdjomow-Sachs variants accounting for plasticity, which allowed to reproduce the hierarchical substructure of lath martensite. Depending on the formulation of the order parameter phase field models can recover the sharp interface driving force in the limit of a vanishing interface (see, e.g., Hildebrand and Miehe, 2012; Schneider et al., 2017). **Internal variable and phenomenological approaches** Micromechanical approaches modeling the martensite evolution based on energetic,

thermodynamical principles or physically motivated stochastic equations can be found in (Turteltaub and Suiker, 2006a; Tjahjanto et al., 2008; Yadegari et al., 2012; Ostwald et al., 2010; 2012; 2015; Brocca et al., 2002; Mehrabi and Kadkhodaei, 2013). Turteltaub and Suiker (2006a;b) proposed a finite strain model for thin plate martensite as it is found in high carbon steels or shape memory alloys undergoing cubic to tetragonal martensitic transformation. The martensite variants are treated as transformation systems analogously to the slip systems in crystal plasticity. The transformation deformation of each martensite variant is determined based on the crystallographic theory of martensite. For each transformation system the Taylor assumption of a uniform deformation gradient is made. In Tjahjanto et al. (2008) the model was later extended to include plastic effects of the austenite, as well as interface energy contributions to the driving force. A fully thermomechanical coupled version was presented in Yadegari et al. (2012). One drawback of the fully resolved micromechanical models is the computational cost. An effort to circumvent this issue is given by microplane (Brocca et al., 2002; Mehrabi and Kadkhodaei, 2013) and microsphere (Ostwald et al., 2010; 2012; 2015) approaches, where the dimensionality of the problem for solving the local equations is reduced. Microplane approaches describe the multiaxial macroscopic behavior as a superposition of uniaxial responses within several planes of different orientation, assuming suitable static and kinematical constraints. Microsphere approaches are based on the projection of the macroscopic stress and strain onto normal directions of a sphere. In each direction the local equations are evaluated in one dimension. In Ostwald et al. (2010; 2012; 2015) a microsphere approache to model martensitic phase transformation in shape memory alloys and TRIP steels was proposed, where the phase transformation is modeled based on statistical physics and plasticity is accounted for.

#### **Chapter 3**

# **Basics of continuum mechanics and thermodynamics of singular surfaces**

### **3.1 Kinematics and strains**

**Motion and deformation of a body** To model the behavior of materials and the underlying processes altering its microstructure, a suitable scale of modeling has to be established. In this work a continuum mechanical formulation is used. This assumes that the considered length scale allows to model the material as a continuous material without taking into account the discrete nature of the real material (atoms, electrons,...). In the following the continuum mechanical description of a body is introduced, where the focus is laid on materials with a singular surface separating different material phases.

A material body B, embedded in the Euclidean space, consists of a set of material points *P*. For a given time *t* = 0 the body is said to occuppy its reference placement - not necessarily a state of zero deformation. For any actual time *t >*= 0 the body is in its current placement. The motion of a material body is described by a function *χ* with

$$x = \chi(X, t),\tag{3.1}$$

which identifies the current position of a material point *x* with the corresponding position vector *X* in the reference placement<sup>1</sup> . The mapping function *χ* is assumed to be continuous, invertible, twice differentiable and orientation preserving, Bertram (2008). The displacement of a material point is determined by

$$
\boldsymbol{w}\_{\rm L} = \boldsymbol{\chi}(\boldsymbol{X}, t) - \boldsymbol{X}.\tag{3.2}
$$

The material time derivative (·)˙ of a Lagrangian quantity ΞL(*X, t*) is defined by

$$
\dot{\Xi}\_{\rm L} = \frac{\partial \Xi\_{\rm L}}{\partial t}. \tag{3.3}
$$

Moreover, the velocity of a material point is given by

$$v\_{\mathbb{L}} = \dot{\chi} = \frac{\partial \chi(\mathbf{X}, t)}{\partial t},\tag{3.4}$$

With the displacement *u<sup>L</sup>* and the mapping function *χ* at hand deformation measures are defined. The deformation gradient *F* and the displacement gradient *H* are defined by

$$F = \frac{\partial \chi(X, t)}{\partial X}, \quad H = \frac{\partial u\_{\mathcal{L}}}{\partial X} = F - I,\tag{3.5}$$

respectively.

The polar decomposition allows to decompose the deformation gradient uniquely into stretch tensors *U* or *V* and an average rotation tensor *R* by

$$F = RU = VR, \qquad \mathcal{R} \in \mathcal{Or} \\ \text{f} \\ \text{h}^+ \qquad \mathcal{U}, \mathcal{V} \in \text{Psym.} \tag{3.6}$$

<sup>1</sup> If a field quantity Ξ is parametrized with respect to the current configuration Ξ<sup>E</sup> = ΞE(*x, t*) it is called Eulerian. If it is parametrized with respect to the reference configuration Ξ<sup>L</sup> = ΞL(*X, t*) it is termed Lagrangian.

Hence, the deformation can be interpreted as a composition of a stretch and a rotation or vice versa. Furthermore, the deformation gradient transforms infinitesimal line elements, infinitesimal area elements and infinitesimal volume elements from the reference into the current configuration, (see, e.g., Silhavy, 1997)

$$\mathbf{dx} = F \mathbf{d}X,\tag{3.7}$$

$$
\mathbf{da} = \det(\mathbf{F}) \mathbf{F}^{\mathbf{T}} \mathbf{d}A,\tag{3.8}
$$

$$\mathrm{d}v = \det(F)\mathrm{d}V.\tag{3.9}$$

By making use of the deformation gradient, the right and left Cauchy-Green tensors are defined as

$$C = F^T F = U^T U, \quad B = F F^T = V V^T,\tag{3.10}$$

which are symmetric and only dependent on the left or right stretch tensor *V* , respectively *U*. Taking the right Cauchy-Green tensor *C* a special class of the generalized Hill strains, the Seth strains (Seth, 1961), are defined by

$$k(x) = \begin{cases} \frac{1}{2m} \left( x^m - 1 \right), & m \neq 0, \\\frac{1}{2} \ln \left( x \right), & m = 0 \end{cases} \tag{3.11}$$

The Seth strain measures provide a zero-strain for an undeformed body *k*(1) = 0. They are monotonuous growing and twice differentiable. The property *k* ′ (1) = 1*/*2 guarantees compatibility with the infinitesimal strain tensor *ε*, see Eq. (3.13). The well known Green-Lagrange strain Tensor *E* is hereby recovered for *m* = 1, while the Biot strain *E* B and Hencky strain *E* H are identified by *m* = 1*/*2 and *m* = 0, respectively

$$E = \frac{1}{2} \left( C - I \right), \quad E^{\rm B} = U - I, \quad E^{\rm H} = \ln \left( U \right). \tag{3.12}$$

The previously defined strain measures are valid for large stretches and rotations. If one considers only small stretches and rotations, which satisfy ||*H*|| ≪ 1, a linearization of the finite strain measures yields the infinitesimal strain tensor

$$\varepsilon = \frac{1}{2} \left( H + H^T \right). \tag{3.13}$$

The velocity gradient is given as the partial derivative of the velocity in the current configuration

$$L = \frac{\partial v\_E}{\partial x} = \dot{F}F^{-1}, \quad D = \text{sym}(L), \quad W = \text{skw}(L). \tag{3.14}$$

The symmetric tensor *D* is the strain rate, whereas the unsymmetric tensor *W* is the spin tensor.

### **3.2 Balance laws for bodies with singular surfaces**

In the following, the index *<sup>L</sup>* of a Lagrangian quantity is omitted, since all considerations take place in the reference configuration.

**A body with a singular surface** Consider a body B in the reference configuration, containing the same set of material points P(*X*) at any instant in time. The body moves with a velocity *v*, has density *ρ*ˇ(*X*), and absolute temperature *θ*(*X, t*).The body B is composed of two domains B <sup>+</sup> and B <sup>−</sup>. The aforementioned fields are assumed to be continuous in each domain. The domains are separated by a generally non material singular surface S(*X*<sup>∗</sup> *, t*) = S*t*. The singular surface as introduced here

**Figure 3.1:** Body B with singular surface

is considered to have zero thickness and no thermodynamic properties like surface energy or surface tension. The orientation of the singular surface S*<sup>t</sup>* in the reference placement at a material point *X*<sup>∗</sup> is given by a unit vector *N*<sup>∗</sup> = *N*<sup>∗</sup> (*X*<sup>∗</sup> *, t*). The singular surface S*<sup>t</sup>* propagates through B with a velocity *v* ∗ <sup>⊥</sup> in normal direction *N*<sup>∗</sup> of the singular surface.

Certain fields, like for example the stress tensor *σ*ˇ, the velocity *v*, the deformation gradient *F* or the temperature *θ*(*X, t*) can suffer jumps across the singular surface, while remaining continuous everywhere else. Using the left (-) and right (+) limit of a general field quantity Ξ, see Fig. 3.1, with respect to the singular surface the jump [Ξ] and the average value hΞi are defined as

$$\{\Xi\} = \Xi^{+} - \Xi^{-} \qquad , \qquad \langle \Xi \rangle = \frac{1}{2} (\Xi^{+} + \Xi^{-}) . \tag{3.15}$$

The jump of a product of two quantities Ξ<sup>1</sup> and Ξ<sup>2</sup> can be decomposed into

$$\{\Xi\_1 \Xi\_2\} = \{\Xi\_1 \mathbf{\bar{\xi}} \langle \Xi\_2 \rangle + \{\Xi\_2 \mathbf{\bar{\xi}} \langle \Xi\_1 \rangle \,. \tag{3.16}$$

A singular surface that is associated with the same set of material points P(*X*<sup>∗</sup> ) at any time is a material singular surface S*<sup>t</sup>* = S*mt*. A singular surface that propagates relative to a material point is called a non-material singular surface S*t*, e.g. a shock wave or a moving phase front during a solid to solid phase transformation<sup>2</sup> .

**Compatibility of the deformation gradient at a singular surface** The deformation gradient at the left and right side of a singular surface cannot be arbitrary, if continuity of the displacement field *u* across a singular surface S*<sup>t</sup>* is required (for a proof see Appendix 11)

$$\{F\} = a \otimes \mathcal{N}^\*.\tag{3.17}$$

Eq. (3.17) is the first Hadamard jump condition, (see, e.g., Silhavy, 1997), which ensures geometric compatibility at the singular surface. This can also be seen by considering the normal projection of the cofactor of the jump of the deformation gradient (see, e.g., Silhavy, 1997)

$$\begin{split} \text{cof}(\{F\})N^\* &= \left(\text{cof}(F\_1 + a \otimes N^\*) - \text{cof}(F\_1)\right)N^\* \\ &= \left(\text{cof}(F\_1)(I + N^\* \cdot F\_1^{-1} aI - N^\* \otimes F\_1^{-1} a) - \text{cof}(F\_1)\right)N^\* \\ &= 0. \end{split} \tag{3.18}$$

Thus, the area of the singular surface S*<sup>t</sup>* approaching it from one side or the other is identical, which characterizes a geometrical compatible, coherent interface. Likewise, the jump of the velocity is constrained to

$$\{v\} = -v\_{\perp}^\* \{F\} N^\*,\tag{3.19}$$

<sup>2</sup> While both, phase fronts and shock waves, are concerned with discontinuities of first gradients or time derivatives of the displacement, they are subject to different jump conditions at the singular surface. A shock wave is associated with an adiabatic, dissipative transition. A phase front is characterized by a in general non-adiabtic, dissipative motion with a spatially continuous temperature field at the singular surface (Maugin, 1998).

where the normal velocity of the singular surface is given by *v* ∗ ⊥. Equation (3.19) is the second Hadamard jump condition, (see, e.g., Silhavy, 1997), at the interface relating the jump of the deformation gradient via the normal of the singular surface to the jump of the velocity field. By inserting the first Hadamard jump condition Eq. (3.17) for a geometrically compatible interface into the second, one finds that

$$\{v\} = -v\_{\perp}^\* a.\tag{3.20}$$

Thus, the direction of the jump of the velocity is related to the negative interface jump vector *a*, whereas the magnitude is given by product of the normed interface jump vector *a* and the normal velocity of the singular surface *v* ∗ ⊥.

**Balance laws** If a body B is subjected to a thermomechanical process the change over time of displacement, temperature and density at any material point is defined by balance laws. In the following, balance relations for mass, linear momentum, angular momentum, energy and entropy at regular material points inside the volume of the body B, the bulk, are stated as well as balance equations on the singular surface.

All balance relations and corresponding equations are established in the reference configuration. Therefore, the geometry of the volume of the body V*m*<sup>0</sup> and boundary *∂*V*m*<sup>0</sup> is time independent. With the transformation formulas, Eq. (3.8) and Eq. (3.9), the push-back of volume and area related Eulerian quantities, into the reference configuration can be performed. The integration can alternatively be performed with respect to the reference configuration. To distinguish such a quantity from a Lagrangian parametrized quantity, it is identified by Ξˇ. For volume

densities and Lagrangian surface densities the following equations hold

$$\int\_{\mathcal{V}\_{\rm mt}} \Xi\_{\rm E} \mathrm{d}v = \int\_{\mathcal{V}\_{\rm m0}} \det(\mathbf{F}) \Xi\_{\rm E} \mathrm{d}V = \int\_{\mathcal{V}\_{\rm m0}} \breve{\Xi}(\chi(\mathbf{X}, t) \mathrm{d}V, \tag{3.21}$$

$$\int\_{\partial \mathcal{V}\_{\rm m}} q\_{\rm E} \mathrm{d}a = \int\_{\partial \mathcal{V}\_{\rm m0}} \det(\boldsymbol{F}) \boldsymbol{F}^{-\rm T} q\_{\rm E} \mathrm{d}A = \int\_{\partial \mathcal{V}\_{\rm m0}} \check{\boldsymbol{q}}(\boldsymbol{\chi}(\boldsymbol{X}, t) \mathrm{d}A. \tag{3.22}$$

Before starting with the balance equations two fundamental theorems are introduced.

**Divergence Theorem** The divergence theorem (Liu, 2013) for bodies with singular surfaces relates the divergence of a vector field *q*ˇ<sup>Ξ</sup> with the flux of the field over a closed surface and the jump over the singular surface in normal direction

$$\int\_{\mathcal{V}\_{m0}} \mathrm{Div}\left(\check{\mathbf{q}}\_{\Xi}\right) \mathrm{d}V = \oint\_{\partial \mathcal{V}\_{m0}} \check{\mathbf{q}}\_{\Xi} \cdot \mathbf{N} \mathrm{d}A - \int\_{S\_t} \{\check{\mathbf{q}}\_{\Xi}\} \cdot \mathbf{N}^\* \mathrm{d}A. \tag{3.23}$$

**Reynold's Transport Theorem** Reynold's transport theorem (Liu, 2013) gives the time derivative of a quantity Ξ = ˇ Ξ( ˇ *X, t*) integrated over a (e.g., material) time dependent volume. While in the reference configuration the placement and hence volume V*m*<sup>0</sup> of the body B is constant, the position of the singular surface is in general not and therefore determines time dependent volumes V + *mt* and V − *mt*

$$\frac{\mathrm{d}}{\mathrm{d}t} \int\_{\mathcal{V}\_{m0}} \check{\Xi} \mathrm{d}V = \int\_{\mathcal{V}\_{m0}} \frac{\partial \check{\Xi}}{\partial t} \mathrm{d}V - \int\_{S\_{mt}} \{\check{\Xi}\} v^\* \cdot \mathbf{N}^\* \mathrm{d}A. \tag{3.24}$$

**General balance equation** The general form of a global, integral balance equation specifies the change in time of a quantity Ξˇ integrated over the volume V*m*<sup>0</sup> of the body B. The change in time of Ξˇ is balanced by the sum of a production term, *p*ˇ<sup>Ξ</sup> = ˇ*p*Ξ(*χ*(*X, t*)*, t*), and a supply *s*ˇ<sup>Ξ</sup> = ˇ*s*Ξ(*χ*(*X, t*)*, t*) inside the volume as well as the flux over the boundaries, see Eq. (3.22), *q*ˇ<sup>Ξ</sup> = *q*ˇΞ(*χ*(*X, t*)*, t*) in direction of the outward normal (Müller, 2007)

$$\frac{\mathrm{d}}{\mathrm{d}t} \int\_{\mathcal{V}\_{m0}} \check{\Xi} \mathrm{d}V = \int\_{\mathcal{V}\_{m0}} (\check{p}\Xi + \check{s}\Xi) \,\mathrm{d}V + \oint\_{\partial\mathcal{V}\_{m0}} \check{q}\_{\Xi} \cdot \mathbf{N} \mathrm{d}A. \tag{3.25}$$

Applying the transport theorem, the divergence theorem and considering contributions inside the bulk and on the singular surface separately, local forms of the general balance equation are obtained (Müller, 2007)

$$\frac{\partial \check{\Xi}}{\partial t} = \check{p}\_{\Xi} + \check{s}\_{\Xi} + \text{Div}\left(\check{q}\_{\Xi}\right),\tag{3.26}$$

$$\{\check{\Xi}\}v\_{\perp}^{\*} = -\{\check{q}\_{\Xi}\} \cdot \mathbf{N}^{\*}.\tag{3.27}$$

With the general balance equation at hand, balance equations for mass, linear momentum, angular momentum, energy and entropy can be formulated.

**Balance of mass** The balanced quantity is the density Ξ = ˇ ˇ *ρ*, where the production and supply of mass vanish *p*ˇ*ρ*<sup>ˇ</sup> = ˇ*sρ*<sup>ˇ</sup> = 0 as well as the flux *q*ˇ*ρ*<sup>ˇ</sup> = **0**. The balance of mass states that the total mass of a body B is constant in time *ρ*ˇ = ˇ*ρ*(*X*), while the jump of the density over the singular surface multiplied with the normal speed of the singular surface vanishes

$$
\dot{\tilde{\rho}}(\mathcal{X}, t) = 0,\tag{3.28}
$$

$$\{\not p(\mathcal{X})\} v\_{\perp}^{\*} = 0.\tag{3.29}$$

Consequently, for a material singular surface S*mt* (*v* ∗ <sup>⊥</sup> = 0) the jump of *ρ*ˇ is arbitrary, while for a moving singular surface the density has to be spatially continuous, since the mass flux is assumed to be zero.

**Balance of linear momentum and angular momentum** As stated by Newton and Euler the change of mass specific momentum in time *ρ*ˇ*v*˙ of a body B is balanced by the resulting externally applied forces. The supply is hereby given as external body forces *s*ˇ*ρv*<sup>ˇ</sup> = *b* while the flux is identified with the First Piola-Kirchhoff stress tensor or nominal stress, respectively *q*ˇ*ρv*<sup>ˇ</sup> = *σ*ˇ, representing tractions applied at the boundary. The production of momentum is zero, *p*ˇ*ρv*<sup>ˇ</sup> = 0. The balance of linear momentum, taking into account the balance of mass and identities (3.15) and (3.16), reads

$$
\Diamond \dot{\rho} \mathbf{b} + \text{Div}\left(\check{\sigma}\right) = \check{\rho} \dot{v}, \tag{3.30}
$$

$$-\langle \check{\rho} \rangle \{ \underline{v} \} \underline{v}\_{\perp}^{\*} = \{ \check{\sigma} \} \underline{N}^{\*},\tag{3.31}$$

where *ρ*ˇ*v* is the momentum vector. Considering a stationary state the balance of linear momentum inside the bulk reduces to the balance of internal and external forces given by

$$
\check{\rho}\mathbf{b} + \text{Div}\left(\check{\sigma}\right) = \mathbf{0}.\tag{3.32}
$$

For a material singular surface Eq. (3.31) yields the traction continuity or force equilibrium at the singular surface

$$\{\check{\sigma}\}N^\* = 0.\tag{3.33}$$

For a non-zero normal velocity of the singular surface, under consideration of the Hadamard conditions Eq. (3.17) and Eq. (3.19), Eq. (3.31) states that the jump of the traction vector is collinear to the jump vector *a* of the deformation gradient at the singular surface

$$
\langle \check{\rho} \rangle v\_{\perp}^{\*} \,^2 \mathbf{a} = \{ \check{\sigma} \} \mathbf{N}^\*. \tag{3.34}
$$

The balance of angular momentum states that the Kirchoff stress tensor is symmetric

$$
\check{\sigma}\boldsymbol{F}^T = \boldsymbol{F}\check{\sigma}^T.\tag{3.35}
$$

The jump of angular momentum over the singular surface imposes the same condition as the jump of momentum. Therefore, it does not yield an additional, independent equation.

**Balance of energy** The balance of total energy, also called the first law of thermodynamics, is given by

$$\not{\rho}\dot{e} + \frac{\not{\rho}}{2}(\boldsymbol{v}\cdot\boldsymbol{v})\dot{\boldsymbol{\epsilon}} = \not{\rho}\boldsymbol{w} - \text{Div}\left(\not{\dot{q}}\right) + \not{\rho}\boldsymbol{b}\cdot\boldsymbol{v} + \text{Div}\left(\not{\dot{\sigma}}^T\boldsymbol{v}\right), \tag{3.36}$$

$$-\{\check{\rho}e + \frac{\check{\rho}}{2}v \cdot v\} v\_{\perp}^{\*} = -\{\check{q}\} \cdot \mathcal{N}^{\*} + \{\check{\sigma}^{\mathrm{T}}v\} \cdot \mathcal{N}^{\*}.\tag{3.37}$$

The change in time of the total energy is given by the sum of the massspecific internal and kinetic energy *ρe*ˇ + ˇ*ρv* · *v*. It is balanced by contributions related to the internal and kinetic energy. The supply of kinetic energy is given by the power of external forces and the supply of the internal energy by the mass specific heat generation of a material point *s*ˇ*<sup>e</sup>* = ˇ*ρ*(*b* · *v* + *w*). The flux of the kinetic energy is represented by the power of tractions applied at the boundary, whereas the non-convective flux of internal energy is identified by the negative heat flux vector, <sup>ˇ</sup>*f<sup>e</sup>* <sup>=</sup> Div *σ*ˇ T *v* − *q* . The rate of mechanical work added to the body by the power of external forces is identified with the contributions of the kinetic supply and flux. By using the divergence theorem one can express the power of tractions by the sum of the power of the internal stresses and the flux of kinetic energy over the boundary

$$\operatorname{Div}\left(\dot{\boldsymbol{\sigma}}^{\mathsf{T}}\boldsymbol{v}\right) = \dot{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} + \operatorname{Div}\left(\dot{\boldsymbol{\sigma}}\right) \cdot \boldsymbol{v}.\tag{3.38}$$

For a body with vanishing velocity, *v* = **0**, the balance of total energy inside the bulk reduces to the balance of internal energy

$$
\check{\rho}\dot{e} = \check{\rho}w - \operatorname{Div}\left(\check{q}\right) + \check{\sigma} \cdot \dot{\bar{F}},\tag{3.39}
$$

where the balance of mass Eq. (3.28) and Eq. (3.38) have been accounted for. The power of internal stresses *σ*ˇ · *F*˙ couples the balance of momentum with the balance of internal energy. Therefore, a stress accompanied deformation of a body leads to a change of temperature inside the body. The jump of the total energy can be further simplified by using the jump of momentum Eq. (3.31). By expanding first the power due to the jump of the stress vector using identity Eq. (3.16) and then applying the jump of momentum Eq. (3.31) one finds

$$\{\dot{\sigma}^T v\} \cdot \mathbf{N}^\* = \{\dot{\sigma}\} \mathbf{N}^\* \cdot \langle v \rangle + \langle \dot{\sigma} \rangle \mathbf{N}^\* \cdot \{v\} \tag{3.40}$$

$$=-\frac{\check{\rho}}{2}v\_{\perp}^{\*}\left\{\boldsymbol{v}\cdot\boldsymbol{v}\right\}+\left\{\boldsymbol{v}\right\}\cdot\left\langle\check{\boldsymbol{\sigma}}\right\rangle\boldsymbol{N}^{\*}.\tag{3.41}$$

Therefore, the power of the jump of the stress vector decomposes into the jump of the kinetic energy and the mechanical power due to the average stress vector. Plugging Eq. (3.41) back into the jump of the total energy and applying the Hadamard jump conditions, Eq. (3.19) and Eq. (3.20), yields

$$-\{\check{\rho}e\}v\_{\perp}^{\*} = -\{\check{q}\} \cdot \mathbf{N}^{\*} - \mathbf{a} \cdot \langle \check{\sigma} \rangle \mathbf{N}^{\*} v\_{\perp}^{\*}.\tag{3.42}$$

Hence, by accounting for the jump of linear momentum and a continuous displacement field across the singular surface, the jump of the total energy reduces to the jump of the internal energy. The jump of the internal energy is then equal to the jump of the heat flux vector projected onto the normal of the interface and a mechanical power comprising the average stress vector. Considering a non-material singular surface, in the special case of an either spatially constant temperature field or an adiabatic singular surface, the jump of the internal energy is equal to the power of the average stress. If additionally, the stress vector on either side of the interface is zero, the jump of the internal energy vanishes, implying that the internal energy on either side is identical. In the purely thermal case, i.e. zeros stress vectors on both sides of the interface, the jump of the internal energy is proportional to the jump of the heat flux vector.

For a vanishing normal velocity of the singular surface the jump condition implies an adiabatic singular surface and therefore either a spatially constant temperature field in B or temperature gradients, which are identical approaching *S* in the limit from B <sup>+</sup> and B −.

**Balance of Entropy** The above mentioned balance equations are complemented by the balance of entropy. The balance of entropy in its local form states that the change in time of the mass specific entropy *ρ*ˇ*η*˙ is equal to the sum of mass specific entropy supply *s*ˇ*<sup>η</sup>* = *w/θ*, the flux of entropy over the boundaries *q*ˇ*<sup>η</sup>* = *q*ˇ*/θ* and the mass specific entropy production *p*ˇ*η*. The form of the entropy supply and flux as given above are of constitutive nature, and therefore not universally valid (Müller, 2007). Further, introducing the assumption that the entropy production for each material and thermokinematic process is non-negative (Müller, 2007), *p*ˇ*<sup>η</sup>* ≥ 0, the second law of thermodynamic is recovered

$$
\check{\rho}\dot{\eta} - \frac{\check{\rho}w}{\theta} + \text{Div}\left(\frac{\check{q}}{\theta}\right) \ge 0,\tag{3.43}
$$

$$-\{\check{\rho}\eta\}v\_{\perp}^{\*} + \{\frac{\check{q}}{\theta}\} \cdot \mathbf{N}^{\*} \geq 0.\tag{3.44}$$

For a reversible process, e.g. a non-dissipative process, the second law of thermodynamics is equal to zero. For an irreversible process or a dissipative process the entropy is increasing.

## **3.3 Implications of the second law of thermodynamics**

**Clausius Duhem inequality at regular points** The second law of thermodynamics is of essential importance in continuum mechanics, imposing restrictions on the considered processes and constitutive relations. To this end, a relation between the internal energy *e*, and the

Helmholtz free-energy Ψˇ <sup>3</sup> per unit mass and the entropy is established (see, e.g., Abeyaratne and Knowles, 2006)

$$
\Psi = e - \eta \theta,\tag{3.45}
$$

Substituting the Helmholtz free-energy Eq. (3.45) back into the internal energy balance Eq. (3.39) and plugging it into the entropy inequality Eq. (3.43) yields the so called Clausius-Duhem inequality at all regular points (Truesdell, 1964)

$$-\not{\rho}\dot{\Psi} - \not{\rho}\eta \dot{\theta} + \not{\sigma} \cdot \dot{\mathbf{F}} - \frac{\not{q}}{\theta} \cdot \text{Grad}\left(\theta\right) \ge 0. \tag{3.46}$$

On substituting *e*˙ −*ηθ*˙ = Ψ +˙ *η* ˙*θ* it can be seen that the Clausius-Duhem inequality implies for the bulk, that the rate of internal energy minus the rate of entropy must not exceed the sum of the internal stress power and the heat flux over the boundaries. In the special case of an isothermal, adiabatic process the rate of mechanically stored energy does not exceed the stress power

$$
\check{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} - \check{\rho} \dot{\Psi} \ge 0. \tag{3.47}
$$

**Clausius-Duhem inequality on a singular surface** In the following the implications of the jump of the entropy balance for a general process are discussed. The only restrictions which are imposed are the continuity of the displacement and mean density across the singular surface. Multiplying (3.44) by h*θ*i the entropy jump balance for a non-continuous

<sup>3</sup> The corresponding Gibbs free energy Ψ∗(*σ*ˇ*, θ*) or Ψ∗(*p, θ*) can be found by performing a Legendre transformation of the Helmholtz free energy Ψ(*F, θ*) with respect to the deformation gradient *F* or the volumetric part of the deformation gradient det(*F*). The drawback of the Legendre transformation with respect to *F* is that Ψ∗(*σ*ˇ*, θ*) cannot generally be defined as function of *σ*ˇ, since *σ*ˇ is not necessarily invertible in *F* (Silhavy, 1997).

temperature across the interface reads

$$-\{\check{\rho}\eta\}\langle\theta\rangle v\_{\perp}^{\*} + \langle\theta\rangle\{\frac{\check{q}}{\theta}\} \cdot \mathbf{N}^{\*} \geq 0. \tag{3.48}$$

With identity (3.16)

$$
\langle \theta \rangle \{ \frac{\check{q}}{\theta} \} \cdot \mathbf{N}^\* = \{ \check{q} \} \cdot \mathbf{N}^\* - \{ \theta \} \langle \frac{\check{q}}{\theta} \rangle \cdot \mathbf{N}^\* \tag{3.49}
$$

holds. Solving Eq. (3.42) for [*q*ˇ] · *N*<sup>∗</sup> and plugging the result into Eq. (3.49) gives

$$
\langle \theta \rangle \{ \frac{\check{\mathbf{q}}}{\theta} \} \cdot \mathbf{N}^\* = \{ \check{\rho} e \} v\_\perp^\* - \langle \check{\sigma} \rangle \cdot \{ F \} v\_\perp^\* - \{ \theta \} \langle \frac{\check{\mathbf{q}}}{\theta} \rangle \cdot \mathbf{N}^\*. \tag{3.50}
$$

Substituting the above result into Eq. (3.48) yields an expression for the jump of the entropy, where the jump of momentum and continuity of the displacement field across the singular surface have been accounted for

$$\check{\rho}\left(\{e\}-\{\eta\}\langle\theta\rangle\right)v\_{\perp}^{\*}-\langle\check{\sigma}\rangle\cdot\{F\}v\_{\perp}^{\*}-\{\theta\}\langle\frac{\check{q}}{\theta}\rangle\cdot\mathbf{N}^{\*}\geq 0.\tag{3.51}$$

Decomposing the internal energy into the Helmholtz free energy and the entropy using Eq. (3.45), the jump of the internal energy reads

$$\{e\} = \{\Psi + \eta\theta\} = \{\Psi\} + \langle \eta \rangle \{\theta\} + \langle \theta \rangle \{\eta\}.\tag{3.52}$$

Finally, substituting Eq. (3.52) into Eq. (3.51) yields a Clausius-Duhemlike inequality on the singular surface (see, e.g., Silhavy, 1997)

$$\begin{split} \left( \langle \dot{\rho} (\{\Psi\} + \langle \eta \rangle \{\theta\} \rangle - \langle \check{\sigma} \rangle \cdot \{F\} \right) v\_{\perp}^{\*} - \{\theta\} \langle \frac{\check{\underline{\theta}}}{\theta} \rangle \cdot \mathbf{N}^{\*} &= \\ f^{\*} v\_{\perp}^{\*} - \{\theta\} \langle \frac{\check{\underline{\theta}}}{\theta} \rangle \cdot \mathbf{N}^{\*} &\geq 0. \end{split} \tag{3.53}$$

The term

$$f^\* = \check{\rho}(\{\Psi\} + \langle \eta \rangle \{\theta\}) - \langle \check{\sigma} \rangle \cdot \{F\} \tag{3.54}$$

is often denoted as the driving force *f* <sup>∗</sup> and represents an energetic quantity conjugated to the motion of the interface identified by *v* ∗ ⊥. Besides the continuity of displacement and density across the interface no restrictions on the fields and constitutive relations have been imposed. Thermal as well as inertial effects are taken into account, where the inertial effects actually cancel in the reduced balance of energy, compare Eq. (3.42), and do not contribute to the motion of the interface. The dissipation at the singular surface is then given by two contributions. The first is related to the motion of the interface, given by the product of the driving force and the normal velocity of the singular surface. The second is given by a purely thermal contribution, expressed by the jump of the temperature across the singular surface and the average of the product of the heat flux and the temperature.

It is now assumed that the singular surface represents a phase boundary of a body with two different, solid phases occupying B <sup>+</sup>, B <sup>−</sup> and the transformation of phase <sup>+</sup> into <sup>−</sup> or vice versa is a diffusionless process. The phase transformation is then characterized by a local motion of the singular surface accompanied by dissipation. Such a phase transformation is termed displacive (Christian et al., 1995). To be thermodynamically consistent the motion of the phase boundary must be in accordance with Eq. (3.55). This implies that the jump balances Eqs. (3.29), (3.31), (3.37) and (3.44) are satisfied for the phase transformation process.

Phase transformations with large velocities of the phase front can often be considered to be adiabatic, i.e. no heat is exchanged at the interface *q* <sup>+</sup>*/*<sup>−</sup> = **0**. The Clausius Duhem inequality on the interface then reads

$$\left(\dot{\rho}\left(\{\Psi\} + \langle \eta \rangle \{\theta\} - \langle \check{\sigma} \rangle \cdot \{F\}\right)v\_{\perp}^{\*} = f\_{\text{ad}}^{\*}v\_{\perp}^{\*} \ge 0. \tag{3.55}$$

The driving force *f* ∗ ad is still the same as in the non-adiabtic case, but the dissipation is now only related to the motion of the singular surface. In the special case of a continuous temperature across the interface, *θ* <sup>+</sup> = *θ* <sup>−</sup>, Eq. (3.55) reduces to

$$\left(\left\langle\dot{\rho}\{\Psi\}-\langle\dot{\sigma}\rangle\cdot\{F\}\right)v\_{\perp}^{\*} = f\_{\theta\_{con}}^{\*}v\_{\perp}^{\*} \ge 0.\tag{3.56}$$

In the case of an adiabatic singular surface or continuous temperature across the singular surface the second law of thermodynamics trivially implies that the driving force and normal velocity are pointing into the same direction. Hence, the direction in which the phase boundary propagates is established.

Reformulating the second term of Eq. (3.56) by use of identity (3.16) and the jump of momentum, Eq. (3.29), yields

$$\langle \check{\boldsymbol{\sigma}} \rangle \cdot \{ \boldsymbol{F} \} = \boldsymbol{N}^\* \cdot \{ \boldsymbol{F}^\mathsf{T} \} \langle \check{\boldsymbol{\sigma}} \rangle \boldsymbol{N}^\* = \boldsymbol{N}^\* \cdot \left( \{ \boldsymbol{F}^\mathsf{T} \check{\boldsymbol{\sigma}} \} - \langle \boldsymbol{F}^\mathsf{T} \rangle \{ \check{\boldsymbol{\sigma}} \} \right) \boldsymbol{N}^\*$$

$$= \boldsymbol{N}^\* \cdot \left( \{ \boldsymbol{F}^\mathsf{T} \check{\boldsymbol{\sigma}} \} - \rho \langle \boldsymbol{F}^\mathsf{T} \rangle \{ \boldsymbol{F} \} v\_\bot^{\*2} \right) \boldsymbol{N}^\* \tag{3.57}$$

$$=\boldsymbol{N}^\* \cdot \left( \{ \boldsymbol{F}^\mathsf{T} \check{\boldsymbol{\sigma}} \} - \frac{\rho}{2} \{ \boldsymbol{F}^\mathsf{T} \boldsymbol{F} \} \boldsymbol{v}\_\perp^\mathsf{T} \right) \boldsymbol{N}^\*. \tag{3.58}$$

Thus, Eq. (3.56) can alternatively be expressed by (Abeyaratne and Knowles, 1990)

$$\mathbf{N}^\* \cdot \{\boldsymbol{\beta}\Psi I - \dot{\boldsymbol{\sigma}}^T \boldsymbol{F} + \frac{\langle \boldsymbol{\beta} \rangle v\_{\perp}^\*}{2} \mathbf{F}^T \mathbf{F} \} \mathbf{N}^\* v\_{\perp}^\* = f\_{\theta\_{con}}^\* v\_{\perp}^\* \ge 0. \tag{3.59}$$

This result allows the interpretation of the driving force *f* <sup>∗</sup> as the projection of the jump of the Eshelby momentum tensor *ρ*ˇΨ*I* − *F* <sup>T</sup>*σ*ˇ and a kinetic term onto the normal of the singular surface. However, this representation is misleading in the sense that it pretends the existence of a kinetic term. This term is actually not present, see Eq. (3.56), and is just evoked by expressing the result using the Eshelby-energy momentum tensor.

An equilibrium state of the phase boundary or interface is characterized by a vanishing normal velocity and a continuous temperature across the singular surface, resulting in zero dissipation. Likewise, the jump of momentum reduces to the traction continuity Eq. (3.31). Consequently, the set of equations describing a phase boundary in equilibrium are

$$f\_{eq}^{\*} = \{\check{\rho}\Psi\} - \mathbf{a} \cdot \langle \check{\sigma} \rangle \mathbf{N}^{\*} = 0,\tag{3.60}$$

$$\{\check{\sigma}\}N^\* = 0.\tag{3.61}$$

The equilibrium driving force is given as the jump of the mass specific Helmholtz free-energy minus the power of the average stress vector. Eqs. (3.60) and (3.61) are also known as the Weierstrass-Erdmann relations which amongst other conditions are necessary for a strongly locally stable interface of gradient discontinuity (Grabovsky and Truskinovsky, 2014). Note, that the orientation of the phase boundary is not specified by the above equations and requires an additional equation, if it is not determined by the microstructure of the material in advance. For a stress free phase transformation the driving force reduces to the jump of the Helmholtz free-energy, which is identical on both sides of the interface in an equilibrium state. Thus, the phase transformation is only driven by the difference of the Helmholtz free energies. For a non stress-free phase transformation the equilibrium position of the interface also depends on the average stress and therefore may also favor a phase with a larger Helmholtz free energy.

If one considers a non-equilibrium state in a dynamic phase transformation, the speed of the phase boundary or phase front can not be determined as a consequence of the above mentioned balance equations. For this reason a kinetic relation between the driving force and speed of the phase boundary needs to be established. Phenomenological approaches are given in Abeyaratne and Knowles (1990) which result in a linear relation of the normal velocity and the driving force with a temperature

dependent function *h*(*θ*)

$$v\_{\perp}^{\*} = h(\theta)f^{\*}, \qquad h(\theta) \ge 0 \tag{3.62}$$

where the non-negativity of *h*(*θ*) is a consequence of Eq. (3.55). A physically more sound approach is proposed by Berezovski and Maugin (2007) who use non-equilibrium balance equations to derive a kinetic relation in the isothermal case, supplemented by a nucleation criterion in the following general form

$$\left\|\boldsymbol{\rho}\boldsymbol{v}\_{\perp}^{\*}\right\|^{2} = \frac{D(f^{\*}-f\_{\rm cr})}{ADf^{\*}-B}.\tag{3.63}$$

*A, B* and *D* are parameters which depend on the thermoelastic properties of the material, the temperature, deformation and possible transformation strains of a phase. *f*cr is a threshold value for the driving force, which has to be overcome to initiate interface motion.

#### **3.4 Stress measures and constitutive relations**

**Stress measures** To characterize the load state of a body B subjected to external forces stress measures are defined. Considering a body in its current configuration, the symmetric Cauchy stress tensor *σ* is defined as

$$\int\_{a} \boldsymbol{\sigma} \cdot \boldsymbol{n} \mathrm{d}a = \int\_{a} \boldsymbol{t} \mathrm{d}a = \boldsymbol{f} \tag{3.64}$$

relating a resultant force *f* in the current configuration to a force density per unit area. The force density is given by the Cauchy stress tensor, projected onto the outward normal *n* of the respective sectional plane inside the considered body B. The resulting vector *t* is called the Cauchy stress vector. If the area in the current configuation is pushed back into

the reference configuration, see Eq. (3.8), the in general unsymmetric first Piola-Kirchhoff stress tensor *σ*ˇ is obtained as

$$\int\_{A} \check{\sigma} \mathrm{d}A = \int\_{A} \sigma \, \det(F) F^{\mathrm{T}} \mathrm{d}A = f. \tag{3.65}$$

Additionally, the symmetric second Piola-Kirchhoff stress tensor *S* <sup>2</sup>*PK*, acting in the reference configuration, is defined by

$$S^{2PK} = F^{-1} \sigma \det(F) F^{\mathrm{T}}.\tag{3.66}$$

**Constitutive relations** The stress measures defined above can also be related to the internal deformation state by constitutive relations. The constitutive relations are hereby established in a sense that a generalized stress measure *S <sup>G</sup>* is a function of a power conjugated generalized strain measure *E G*

$$\mathbf{S}^G = h(\mathbf{E}^G) \tag{3.67}$$

where for an elastic, simple material the current stress only depends on the current, local generalized strain. A generalized stress and strain measure are called power conjugated if

$$S^G \cdot E^G = \sigma \cdot D. \tag{3.68}$$

holds. In the special case of a hyperelastic material the considered generalized stress measure is related to a potential, given by the mass specific Helmholtz free energy Ψ

$$\mathbf{S}^G = \check{\rho} \frac{\partial \Psi}{\partial \mathbf{E}^G}.\tag{3.69}$$

In the case of a thermoelastic material the Helmholtz free energy is only dependent on the external observable quantities *E <sup>G</sup>* and *θ*

$$
\Psi = \Psi(\boldsymbol{E}^G, \boldsymbol{\theta}).\tag{3.70}
$$

To be a reasonable choice the Helmholtz free energy has to satisfy physically consistent, positive energy states in the case of zero or infinite deformation. Moreover, the principle of material symmetry and in the context of geometrical non-linear theories material frame-indifference have to be adequately reflected. Furthermore, the Helmholtz free energy has to treatable by an existence theory in the case of finite hyperelasticity, see section 5. These requirements are satisfied by the following conditions, (see, e.g., Silhavy, 1997, p. 192) and (Ball, 2002, p. 2),


$$W(F, \theta) \to \infty \quad \forall \, |\det(F)| \to 0 \tag{3.71}$$

$$W(F, \theta) \to \infty \quad \forall \, |\det(F)| \to \infty \tag{3.72}$$


A possible choice for a Helmholtz free energy satisfying the above conditions is a Neo-Hookean energy formulated with invariants of the deformation gradient *F* in a general additive form (see e.g, Hesch et al., 2017)

$$W(F, \theta) = W\_1(||F||^2) + W\_3(\det(F), \theta),\tag{3.73}$$

$$W = \frac{\mu}{2} (\text{tr}(\mathbf{C}) - 3 - 2 \ln(\det(\mathbf{F}))) + (\frac{\lambda}{2} - \beta \triangle \theta)(\det(\mathbf{F}) - 1)^2,\tag{3.74}$$

relating the first Piola-Kirchhoff stress to the deformation gradient via the hyperelastic relation

$$
\check{\sigma} = \frac{\partial W(F, \theta)}{\partial F}. \tag{3.75}
$$

To account for the microstructure evolution of an inelastic material accompanied by dissipation, internal variables are introduced (Coleman and Gurtin, 1967). Internal variables are classified as on the macroscopic scale non observable quantities (see, e.g., Maugin, 2015). The collection of all of them is in the following denoted by the power conjugated stress-like vector quantity *β* and strain-like quantity *α*. To describe the state of an inelastic material the Helmholtz free energy is extended by the vector of inelastic variables

$$W = W(\mathbf{F}, \theta, \underline{\alpha}).\tag{3.76}$$

In the sense of a Generalized Standard Material the vector of the stress-like internal variables *β* is related to the strain-like quantites *α* by the negative derivative of the Helmholtz free energy with respect to *α* (Halphen and Son Nguyen, 1975)

$$
\underline{\beta} = -\frac{\partial W}{\partial \underline{\alpha}}.\tag{3.77}
$$

Evaluating the Clausius-Duhem inequality Eq. (3.46) for an inelastic material in the framework of a Generalized Standard Material yields

$$-\check{\rho}\left(\frac{\partial\Psi}{\partial F}\cdot\dot{F} + \frac{\partial\Psi}{\partial\theta}\cdot\dot{\theta} + \frac{\partial\Psi}{\partial\underline{\alpha}}\cdot\dot{\underline{\alpha}}\right) - \check{\rho}\eta\dot{\theta} + \check{\sigma}\cdot\dot{F} - \frac{\check{\mathbf{q}}}{\theta}\cdot\text{Grad}\,(\theta) \ge 0,\tag{3.78}$$

where the entropy *η* and the first Piola-Kirchhoff stress tensor *σ*ˇ are due to the linearity and arbitrariness of the rates *F*˙ and ˙*θ* identified by

$$
\check{\sigma} = \check{\rho} \frac{\partial \Psi}{\partial \mathbf{F}}, \quad \eta = -\frac{\partial \Psi}{\partial \theta}. \tag{3.79}
$$

Such relations confine the stress and entropy to energetic quantities. The remaining inequality, under the assumption of a spatially homogeneous temperature, implies for the non thermal related entropy production rate or dissipation D

$$\mathcal{D} = -\not{\rho} \frac{\partial \Psi}{\partial \underline{\alpha}} \cdot \dot{\underline{\alpha}} = \underline{\beta} \cdot \dot{\underline{\alpha}} \ge 0. \tag{3.80}$$

An inelastic material that obeys the above inequality for all thermokinematical processes is called thermodynamic consistent. The evolution of the internal variables is usually established as constitutive equations. In most cases these constitutive relations are not available and have to be derived based on physically motivated principles. In mechanics these are potential based approaches as e.g. the Hamilton principle, the principle of minimum potential energy and in the case of equilibrium thermodynamics the minimum of the Gibbs free energy. This renders the task to describes the evolution of the considered internal variables a variational problem. Amongst the first to use a variational formulation was Onsager (1931) introducing the principle of maximum dissipation or the principle of maximum entropy production. The principle of maximum dissipation states that in the isothermal case the difference of the dissipation D and a non-negative dissipation potential Φ ∗ ( ˙*α*), depending on the rates of the strain like variables, is maximized

$$\max\_{\underline{\dot{\alpha}}} \left\{ \mathcal{D} - \Phi = \underline{\beta} \cdot \underline{\dot{\alpha}} - \Phi \right\}. \tag{3.81}$$

This yields evolution equations of the type

$$
\underline{\beta} = \frac{\partial \Phi}{\partial \underline{\dot{\alpha}}}.\tag{3.82}
$$

Onsager restricted his work to quadratic dissipation potentials leading to linear evolution equations. This concept was later generalized by Ziegler (1963) to non linear evolution equations derived from dissipation potentials with an arbitrary exponent. The principle of maximum dissipation has been widely accepted as a tool to derive evolution equations for internal variables, which is also due to the fact that in plasticity the experimentally observed normality of the yield surface and the plastic strain can be captured (see, e.g., Lubliner, 1984).

Another possibility to derive evolution equations following a variational formulation is the minimum principle for the dissipation potential (see, e.g., Gill et al., 2001). The principle of minimum dissipation power states that the total power *P*, composed of the rate of the Helmholtz free energy and the dissipation potential Φ ∗ , is minimized by the rates of the internal strain like variables

$$\min\_{\underline{\dot{\alpha}}} \left\{ P \right\} = \min\_{\underline{\dot{\alpha}}} \left\{ \dot{W} + \Phi^\* \right\} \tag{3.83}$$

resulting in a stationarity condition which defines the evolution of the stress-like internal variables *β*

$$\frac{\partial W}{\partial \dot{\underline{\alpha}}} + \frac{\partial \Phi^\*}{\partial \dot{\underline{\alpha}}} = \mathbf{0} \quad \underline{\beta} = \frac{\partial \Phi^\*}{\partial \dot{\underline{\alpha}}} \tag{3.84}$$

This allows to include the evolution of inelastic variables in a minimization framework, which is especially useful to study the evolution of materials with microstructure (see e.g., Ortiz and Repetto, 1999; Miehe and Lambrecht, 2003; Hackl and Heinen, 2008b; Kochmann and Hackl, 2011; Bartel et al., 2011). In cases where the dissipation potential based on the rates is unknown, as e.g. in the case of crystal plasticity, a LegendreFenchel transformation of Φ ∗ ( ˙*α*) is applied

$$\Phi(\underline{\beta}) = \sup\_{\underline{\dot{\alpha}}} \left\{ \underline{\beta} \cdot \underline{\dot{\alpha}} - \Phi^\*(\underline{\alpha}) \right\}. \tag{3.85}$$

In the case of the principle of minimum dissipation power the Legendre transformation of Φ ∗ ( ˙*α*) yields the same stationarity conditions as the principle of minimum dissipation power. Hence, the rates of the inelastic variables can be determined by the dissipation potential Φ(*β*)

$$
\dot{\underline{\alpha}} = \frac{\partial \Phi(\underline{\beta})}{\partial \beta} \tag{3.86}
$$

which holds per definition of the Legendre-Fenchel transformation. For dissipation potentials homogeneous of degree one with respect to *α*˙ the principle of minimum dissipation power is equivalent to the principle of maximum dissipation, see Hackl and Fischer (2008). The differences and similarities between the principle of maximum dissipation and the principle of minimum dissipation power are also discussed in this publication. Further details concerning different approaches to determine evolution equations for non-equilibrium thermodynamics are covered by an extensive review of Fischer et al. (2014).

#### **Chapter 4**

## **Laminate homogenization**

## **4.1 Laminate approach in the context of homogenization**

**Effective behavior of heterogeneous materials** When dealing with materials with different microstructural constituents, especially in engineering applications, often only the effective or macroscopic material behavior is of primary interest. The effective material behavior can be derived by homogenization which denotes the calculation of effective material properties based on the micromechanical material properties. This task is accomplished by a large variety of models, which can be roughly classified into mean field and full field approaches. Mean field approaches comprise on the one hand bounds derived from energetic principles (Nadeau and Ferrari, 2001; Voigt, 1889; Reuss, 1929; Hashin and Shtrikman, 1963; Kröner, 1977; Willis, 1977) representing a range of feasible effective material properties. On the other hand, there exist a large number of models based on the analytical solution of an inclusion problem (Eshelby, 1957). This allows to identify the local stress fields of the inclusion as well as matrix and consequently the effective material behavior. Such approaches deliver estimates for the effective material behaviour (see, e.g., Mori and Tanaka, 1973; Christensen and Lo, 1979; Zheng and Du, 2001). Full field models are usually applied in the context of numerical approaches operating on a grid of a fully resolved

microstructure. Examples in this class are Fast Fourier Transform (FFT) based and reduced order model approaches (Moulinec and Suquet, 1995; Fritzen and Böhlke, 2013).

A simple, analytical way to homogenize a two phase material is given by laminate approaches. A laminate comprises different phases with constant deformation gradients, separated by a planar interface (see, e.g., Francfort and Murat, 1986; Glüge and Kalisch, 2014). The laminate is used as a specific microstructure at each material point, where a scale separation between the material point and the corresponding phases of the laminate substructure is assumed. Thus, the microscale represented by the laminate has to be much smaller than the observed structure at the material point. Examples of such microstructures are the martensite laths of a block or the twinned martensite variants of plate martensite. A laminate is characterized by its rank, where in the following the focus is first laid on a rank 1 laminate in the context of a geometrical nonlinear setting.

### **4.2 Rank 1 laminates**

**Defintion of a rank 1 laminate** A rank 1 laminate comprises a region Ω with constant volume *V* , a singular surface with constant cross section *A* and periodic, pairwise homogeneous (one dimensional) deformations (Silhavy, 1997). A deformation is referred to as one dimensional if the deformation gradient outside of the interface can be parametrized by a scalar projection along the interface normal *ξ* = *X* · *N*

$$F(X) = F(\xi). \tag{4.1}$$

If the deformation is additionally *p*-periodic, the average deformation

**Figure 4.1:** Rank 1 laminate, *p* = 1

gradient is identical to the effective deformation gradient

$$\bar{F} = \langle F \rangle = \frac{1}{p} \int\_0^p F(\xi) d\xi \tag{4.2}$$

with the periodicity *p*. The corresponding deformation and deformation gradient can be shown to be of the following form (Silhavy, 1997, p. 41)

$$\chi(X) = \bar{F}X + u(\xi) \qquad F = \bar{F} + u'(\xi) \otimes N \tag{4.3}$$

where *u*(*ξ*) is a *p*-periodic, continuous and differentiable function such that

$$u'(\xi) = \mathbf{a} \begin{cases} -c\_2 & 0 \le \xi < c\_1 p \\ c\_1 & c\_1 p \le \xi < p \end{cases} \tag{4.4}$$

A phase of the laminate is now identified by its volume fraction *c<sup>i</sup>* , possible corresponding internal variables *α<sup>i</sup>* and a homogeneous deformation gradient *F<sup>i</sup>*

$$F\_1 = \dot{F} - c\_2 a \otimes \mathcal{N} \qquad F\_2 = \dot{F} + c\_1 a \otimes \mathcal{N}.\tag{4.5}$$

The resulting phase deformations can be interpreted as the sum of an average macroscopic deformation gradient and a periodic perturbation. The perturbation is given by the rank 1 connection of the interface normal *N* and a vector *a*, scaled by the respective volume fraction. It represents the fluctuation of the individual phase deformation gradients due to the microstructure. The phase deformation gradients satisfy the Hadamard jump condition for a coherent, stress free interface, Eq. (3.17),

$$\{F\} = F\_2 - F\_1 = \mathbf{a} \otimes \mathbf{N},\tag{4.6}$$

where *a* represents the magnitude of the jump in normal direction at the interface. The deformation of the rank 1 laminate is thus characterized by periodic, piecewise homogeneous deformation gradients, separated by a planar interface with normal *N*. The volume fractions of the laminate phases can be related to the position of the interface *ξ<sup>S</sup>* through

$$c\_1 = \frac{V\_1}{V} = \frac{\xi\_S}{L}, \qquad c\_2 = \frac{V\_2}{V} = 1 - \frac{\xi\_S}{L}, \quad c\_1 + c\_2 = 1,\tag{4.7}$$

introducing the characteristic length *L* of the laminate, see Fig. 4.1. The normal velocity of the interface is identified with *v* ∗ <sup>⊥</sup> <sup>=</sup> ˙*ξS*. Therefore, the evolution of the phase fractions can be expressed by the normal velocity of the interface

$$
\dot{c}\_1 = -\dot{c}\_2 = \frac{v\_\perp^\*}{L}.\tag{4.8}
$$

#### **4.3 Macroscopic laminate behaviour**

Due to the periodic deformation gradients *F*(*ξ*) the macroscopic deformation gradient is equal to the volume averaged deformation gradient

of the laminate. Since *F* is piecewise constant, the effective deformation gradient *F*¯ = h*F*(*X*)i is given by

$$\bar{F} = \frac{1}{|V|} \int\_{V} F(X)dV = \sum\_{i=1}^{2} c\_{i}F\_{i}.\tag{4.9}$$

Due to the fact of piecewise constant values of *F* the macroscopic free energy reads

$$\bar{W} = \sum\_{i=1}^{2} c\_i W\_i(F\_i, \underline{\alpha}\_i, c\_i). \tag{4.10}$$

Analogously to the effective deformation gradient, the effective first Piola-Kirchhoff stress tensor, accounting for a stress free interface as well as material singular surface, is given by

$$\bar{\tilde{\sigma}} = \langle \check{\sigma} \rangle = \frac{1}{|V|} \int\_{V} \check{\sigma}(\mathbf{X}) \mathrm{d}V = \sum\_{i=1}^{2} c\_{i} \check{\sigma}\_{i}. \tag{4.11}$$

Assuming a hyperelastic material behaviour (see Eq. (3.75)), it can be shown that the effective first Piola-Kirchhoff stress tensor is also recovered as the partial derivative of the effective Helmholtz free energy by the effective deformation gradient in the same form as above. Accounting for Eq. (4.5)

$$\bar{\bar{\sigma}} = \frac{\partial \bar{W}}{\partial \bar{F}} = \sum\_{i=1}^{2} c\_{i} \left( \frac{\partial F\_{i}}{\partial \bar{F}} + \frac{\partial F\_{i}}{\partial a} \frac{\partial a}{\partial \bar{F}} + \frac{\partial F\_{i}}{\partial N} \frac{\partial N}{\partial \bar{F}} \right)^{\text{T}\_{\text{H}}} [\frac{\partial W\_{i}}{\partial F\_{i}}] \tag{4.12}$$

$$\mathcal{L} = \sum\_{i=1}^{2} c\_i \check{\sigma}\_i + c\_1 c\_2 (\{\check{\sigma}\} \mathcal{N}) \frac{\partial a}{\partial \bar{F}} + c\_1 c\_2 (\{\check{\sigma}\} \mathcal{)^\top a \frac{\partial \mathcal{N}}{\partial \bar{F}}.\tag{4.13}$$

The second term vanishes due to traction continuity across the interface. The third term vanishes for a constant interface orientation *N*. In the case of a non-constant interface orientation, the normal vector has to

be chosen in a way that [*σ*ˇ] T *a* = **0** is satisfied. This condition, as will be shown later, arises when minimizing the effective Helmholtz free energy of a laminate material with respect to the interface orientation. In the above cases relation Eq. (4.11) holds.

**Hill Mandel condition in the context of a rank 1 laminate** With ¯*σ*ˇ and *F*¯ at hand it can be shown that the Hill-Mandel condition (Hill, 1963; Mandel, 1971) is satisfied by the laminate homogenization. The Hill-Mandel condition states that for either periodic or homogeneous displacement, respectively stress boundary conditions the macroscopic energy of the average quantities is equal to the energy, found by integrating over the micromechanical phases

$$\frac{1}{|V|}\int\_{V}\check{\sigma}(x)\cdot\dot{F}(x)dV=\langle\check{\sigma}\rangle\cdot\langle\dot{F}\rangle.\tag{4.14}$$

As a result the energy rate is conserved when bridging the microscopic scale (in the following represented by the laminate) towards the macroscopic scale through the process of homogenization. A similar form is given by a work like expression

$$\frac{1}{|V|}\int\_{V}\check{\sigma}(x)\cdot F(x)dV = \bar{\check{\sigma}}\cdot\bar{F}.\tag{4.15}$$

Considering a rank 1 laminate in mechanical equilibrium, i.e. constant volume fractions and fixed interface orientation, with piecewise constant deformation gradients (Eq. (4.5)) and stress tensors the volume average over the micromechanical phases yields

$$\frac{1}{|V|}\int\_{V}\check{\sigma}(x)\cdot F(x)dV = c\_{1}\check{\sigma}\_{1}\cdot F\_{1} + c\_{2}\check{\sigma}\_{2}\cdot F\_{2}$$

$$= c\_{1}\check{\sigma}\_{1}\cdot(\bar{F}-c\_{2}\boldsymbol{a}\otimes\boldsymbol{N}) + c\_{2}\check{\sigma}\_{2}\cdot(\bar{F}+c\_{1}\boldsymbol{a}\otimes\boldsymbol{N}).\tag{4.16}$$

Reordering the terms and accounting for traction continuity (Eq. (3.31)), Eq. (4.16) can be simplified to

$$(c\_1 \check{\sigma}\_1 + c\_2 \check{\sigma}\_2) \cdot \bar{F} + c\_1 c\_2 \mathbf{a} \cdot (\check{\sigma}\_2 - \check{\sigma}\_1) \mathbf{N} = \langle \check{\sigma} \rangle \cdot \langle F \rangle = \bar{\check{\sigma}} \cdot \bar{F}. \tag{4.17}$$

Therefore, the Hill-Mandel condition Eq. (4.15) is satisfied for a rank 1 laminate in equilibrium.

It can be further shown that in the case of a rank 1 laminate with evolving phases and a variable interface orientation, a minimum state of the corresponding effective laminate Helmholtz free energy implies the Hill Mandel condition. To that end the effective stress power, in the case of a hyperelastic material, is identified with

$$
\langle \check{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} \rangle = \frac{1}{|V|} \int\_{V} \check{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} \mathrm{d}V = \frac{1}{|V|} \frac{\mathrm{d}}{\mathrm{d}t} \int\_{V} \check{\rho} \Psi \mathrm{d}V \tag{4.18}
$$

Assuming a coherent and hence stress free interface, applying the transport theorem and decomposing the volume integral into the respective phases yields

$$\frac{1}{|V|}\frac{d}{dt}\int\_{V}\dot{\rho}\Psi \mathbf{d}V = \frac{1}{|V|}\left(\int\_{V\_1}\check{\sigma}\cdot\dot{\mathbf{F}} \mathrm{d}V\_1 + \int\_{V\_2}\check{\sigma}\cdot\dot{\mathbf{F}} \mathrm{d}V\_2 - \{\check{\rho}\Psi\}v\_\perp^\*\,\mathrm{d}A\right). \tag{4.19}$$

Due to the phasewise homogeneous stresses and deformation gradients (Eq. (4.5)) as well as a constant laminate geometry, Eq. (4.18) can be expressed as

$$
\langle \check{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} \rangle = \frac{1}{|V|} (V\_1 \check{\boldsymbol{\sigma}}\_1 \cdot (\dot{\bar{\boldsymbol{F}}} - \dot{c}\_2 \boldsymbol{a} \otimes \boldsymbol{N} - c\_2 \dot{\boldsymbol{a}} \otimes \boldsymbol{N} - c\_2 \boldsymbol{a} \otimes \dot{\boldsymbol{N}})
$$

$$
+ V\_2 \check{\boldsymbol{\sigma}}\_2 \cdot (\dot{\bar{\boldsymbol{F}}} + \dot{c}\_1 \boldsymbol{a} \otimes \boldsymbol{N} + c\_1 \dot{\boldsymbol{a}} \otimes \boldsymbol{N} + c\_1 \boldsymbol{a} \otimes \dot{\boldsymbol{N}}))
$$

$$
$$

59

Using the relations Eq. (4.7) and Eq. (4.8), accounting in the special case of a laminate for the identities *c*1*σ*ˇ <sup>1</sup> + *c*2*σ*ˇ <sup>2</sup> = ¯*σ*ˇ as well as ¯*σN*ˇ = h*σ*ˇi*N* and rearranging the terms yields

$$\frac{1}{|V|}\int\_{V}\check{\boldsymbol{\sigma}}\cdot\dot{\boldsymbol{F}}dV = \overset{\ddot{\boldsymbol{\sigma}}}{\boldsymbol{\sigma}}\cdot\dot{\check{\boldsymbol{F}}} + c\_{2}c\_{1}\{\check{\boldsymbol{\sigma}}\}N\cdot\dot{\boldsymbol{a}} + c\_{2}c\_{1}\{\check{\boldsymbol{\sigma}}^{\mathsf{T}}\}\boldsymbol{a}\cdot\dot{N}$$

$$+ \left(\langle\check{\boldsymbol{\sigma}}\rangle\cdot\left(\boldsymbol{a}\otimes N\right) - \{\hat{\boldsymbol{\rho}}\Psi\}\right)v\_{\perp}^{\ast}\frac{A}{|V|}.\tag{4.21}$$

In the case of a rank 1 laminate, the effective stress tensor and deformation gradient is identical to the volume averaged quantities assuming homogeneous or periodic displacement, respectively stress boundary conditions. Hence the Hill-Mandel condition holds if the time derivatives of the second and third term, as well as the normal velocity of the interface vanish. In any other case traction continuity, a zero driving force or vanishing velocity of the interface, characterizing an infinitely slow process close to equilibrium, and an additional equation

$$\{\check{\sigma}^{\mathsf{T}}\}a = \mathbf{0} \tag{4.22}$$

have to hold to satisfy the Hill- Mandel condition for a rank 1 laminate

$$
\langle \check{\boldsymbol{\sigma}} \cdot \dot{\boldsymbol{F}} \rangle = \langle \check{\boldsymbol{\sigma}} \rangle \cdot \langle \dot{\boldsymbol{F}} \rangle = \earrow \dot{\bar{\boldsymbol{F}}}.\tag{4.23}
$$

It will be shown later that these three conditions also arise when considering minimum states of the effective Helmholtz free energy - as necessary conditions of the rank 1 convex envelope - for a two phase material. Thus, an equilibrium state of a two phase hyperelastic material, where the phases are identified by laminates, equally satisfies the Hill-Mandel condition. The same can also be shown to be true for a multiphase effective Helmholtz free energy using a laminate substructure.

**Effective balance of angular momentum** Furthermore, the laminate homogenization preserves the symmetry of the macroscopic Kirchhoff stress tensor *τ*¯ and in consequence does not violate the macroscopic balance of angular momentum. The macroscopic angular momentum, in terms of the effective quantities, is identified with

$$
\bar{\tau} = \stackrel{\bar{\circ}}{\bar{\sigma}} \bar{\bar{F}}^T = \bar{F} \,\, \bar{\bar{\sigma}}^T,\tag{4.24}
$$

where for a Rank-1 laminate ¯*σ*ˇ and *F*¯ are given by Eqs. (4.13) and (4.9). Eq. (4.24) then yields

$$
\check{\sigma}\_2 \check{\mathbf{F}}\_1^T + \check{\sigma}\_1 \check{\mathbf{F}}\_2^T = \boldsymbol{F}\_1 \check{\sigma}\_2^T + \boldsymbol{F}\_2 \check{\sigma}\_1^T \tag{4.25}
$$

Adding *σ*ˇ <sup>1</sup>*F* T <sup>1</sup> − *σ*ˇ <sup>1</sup>*F* T 1 to the right and *F*2*σ*ˇ T <sup>2</sup> − *F*2*σ*ˇ T 2 to the left term as well as replacing *F*<sup>2</sup> − *F*<sup>1</sup> by [*F*] one obtains

$$
\check{\sigma}\_1 \{ \boldsymbol{F}^\mathsf{T} \} + (\check{\sigma}\_2 + \check{\sigma}\_1) \boldsymbol{F}\_1^\mathsf{T} = -\{ \boldsymbol{F} \} \check{\sigma}\_2^\mathsf{T} + \boldsymbol{F}\_2 (\check{\sigma}\_1^\mathsf{T} + \check{\sigma}\_2^\mathsf{T}) \Leftrightarrow
$$

$$
\check{\sigma}\_1 \{ \boldsymbol{F}^\mathsf{T} \} + \{ \boldsymbol{F} \} \check{\sigma}\_2^\mathsf{T} = -(\check{\sigma}\_2 + \check{\sigma}\_1) \boldsymbol{F}\_1^\mathsf{T} + \boldsymbol{F}\_2 (\check{\sigma}\_1^\mathsf{T} + \check{\sigma}\_2^\mathsf{T}) \tag{4.26}
$$

Performing a similar operation as above, by adding *σ*ˇ <sup>1</sup>*F* T <sup>1</sup> − *σ*ˇ <sup>1</sup>*F* T <sup>1</sup> + *σ*ˇ <sup>2</sup>*F* T <sup>2</sup> − *F*2*σ*ˇ T 2 to the right equation yields

$$
\check{\sigma}\_1 \{ F^\Gamma \} + \{ F \} \check{\sigma}\_2^\Gamma = \{ F \} \check{\sigma}\_1^\Gamma + \check{\sigma}\_2 \{ F^\Gamma \}. \tag{4.27}
$$

Due to the Hadamard condition, Eq. (3.17), and traction continuity *σ*ˇ <sup>1</sup>*N* = *σ*ˇ <sup>2</sup>*N* = ˇ*t* Eq. (4.27) reads

$$
\check{t}\otimes\mathfrak{a}+\mathfrak{a}\otimes\check{t} = \mathfrak{a}\otimes\check{t}+\check{t}\otimes\mathfrak{a},\tag{4.28}
$$

which in conclusion shows that the symmetry of the effective Kirchoff stress tensor and equally the local macroscopic angular momentum balance is preserved for the laminate.

**Geometrical properties of laminates** Equation (4.9) fulfills the correct mapping of line elements from the macroscale to the microscale of the laminate phases. Additionally, it can be shown that the corresponding

area and volume elements are mapped correctly, (see, e.g., Hildebrand, 2013) and references therein

$$\text{cof}(\bar{F}) = c\_1 \text{cof}(F\_1) + c\_2 \text{cof}(F\_2), \tag{4.29}$$

$$\det(\bar{F}) = c\_1 \det(F\_1) + c\_2 \det(F\_2). \tag{4.30}$$

A detailed derivation of these results can also be found in the appendix 11.

**Equilibrium state** Static, mechanical equilibrium of a laminate is characterized by traction continuity and in the case of a material capable of undergoing phase transformation, a vanishing driving force, see Eqs. (3.61) and (3.60), which is equal to force- and phase equilibrium at the interface. By satisfying the force or momentum balance at the interface the angular momentum balance holds likewise and consequently the laminate is in static, mechanical equilibrium.

**Laminates of higher order** Laminates of higher order can be derived by introducing a rank 1 laminate in each phase (*j*) of a rank (*k* − 1) laminate, subsequently leading to a rank *k* laminate. A rank *k* laminate can be identified with a tree structure, see Fig. 4.2, where the nodes of the tree represent the individual *j* phases of the *k*-th order laminate. Starting from the top which represents a homogeneous material or rank 0 laminate, the laminate of rank *k* is found through branching of the respective laminates, respectively nodes of the lower order laminate. Hereby, the number of sublaminates grows exponentially which limits the applicability of higher order laminates. Equivalently to a rank 1 laminate, the average deformation of a sublaminate on the *k*-th level is given by the deformation of the corresponding *k* − 1 laminate

$$F\_{k-1,j} = c\_{k,j} F\_{k,j,1} + (1 - c\_{k,j}) F\_{k,j,2} \tag{4.31}$$

where *j* denotes the total number of laminates on the *k* − 1th laminate level. Due to the coherent laminate interfaces the Hadamard compati-

**Figure 4.2:** Tree-like structure for a rank 3 laminate

bility condition, Eq. (3.17), holds for each sublaminate

$$\{F\}\_{k,j} = a\_{k,j} \otimes \mathbf{N}\_{k,j}.\tag{4.32}$$

This imposes a fluctuation on the deformation of the respective phases of the *k*-th order sublaminate and in turn implies a separation of scales between the laminates of Rank *k* − 1 and *k*.

#### **Chapter 5**

# **Variational principle in elasticity and energy relaxation**

#### **5.1 Existence of solutions in hyperelasticity**

**Variational principle for non-linear hyperelasticity** Performing a static, mechanical structural analysis of a body B subjected to displacement and traction boundary conditions one is usually interested in finding the displacement field of that body. With the displacement field at hand, deformation and consequently stress measures inside the body are obtained. Considering a hyperelastic material the Helmholtz freeenergy *W*(*F*) serves as a potential for the stress which, with respect to the deformation gradient *F*, is given as the first Piola-Kirchhoff stress tensor *σ*ˇ = *∂W*(*F*)*/∂F*. In the mechanical, geometrically nonlinear case computing the displacement field amounts to solving the momentum balance supplemented by traction *t*<sup>0</sup> and displacement boundary conditions *u*<sup>0</sup> on the respective boundaries *∂*Ω*t* and *∂*Ω*u* of B

$$\begin{aligned} \text{Div}\left(\check{\sigma}\right) + \check{\rho}b &= 0 & \text{in } \Omega, \\ u(X) &= u\_0(X) & X \in \partial\Omega\_L, \\ \frac{\partial W(F)}{\partial F}N &= t(X) = t\_0(X) & X \in \partial\Omega\_t. \end{aligned}$$

The above equations (5.1) can be equivalently expressed by the weak solutions of the formal variation of a total potential Π(*u*) comprising the Helmholtz free-energy *W*(*F*(*u*)) and a potential *U*(*u*) representing the volume forces *ρ*ˇ*b* and tractions *t*<sup>0</sup> on the boundary (see e.g., Ball, 2002). In this case the volume forces and tractions must be conservative forces, which means that the total work done by moving the corresponding forces between two points is independent of the taken path. The corresponding respective generalized force *f* is then related to the negative derivative of that potential *f* = − *∂U*(*u*) *∂u* . Π(*u*) is then given as

$$
\Pi(u) = \int\_{\Omega} W(F) + U(u) \,\mathrm{d}V. \tag{5.2}
$$

Looking for an equilibrium position of the body implies that the critical points of the functional Π(*u*) are replaced by minimizers. Existence of a minimizing displacement field of a particular space requires lower semicontinuity of the Helmholtz free-energy *W*(*F*) with respect to the weak topology over the Banach space together with a coercivity assumption (Morrey et al., 1952; Ball, 1976; Ciarlet, 1988; Ball, 2002; Dacorogna, 2006; 2008). The considered space is usually taken as the Sobolev space<sup>1</sup> *W*<sup>1</sup>*,p*. However, the above requirements have to satisfy

$$L^p(\Omega) := \{ f : \Omega \subset \mathbb{R}^n \to \mathbb{R}, \int\_{\Omega} |f^p| \mathrm{d}V < \infty \} \quad 1 < p < \infty \tag{5.3}$$

for an open, bounded set Ω ⊂ R*n*, with a norm ||*f*||*<sup>p</sup> <sup>L</sup>p*(Ω) := (R Ω |*f <sup>p</sup>*|d*V* ) (1*/p*) . The Sobolev space *Wk,p* is defined as

$$W^{k,p} := \left\{ f \in L^p(\Omega), \ f^{'\alpha} \in L^p(\Omega), \ \forall \alpha < |k| \right\},\tag{5.4}$$

where for the generalized *α*-th weak derivative *f* ′*α*

$$(\Phi, \left. f \right.^{\prime} \alpha \right)\_{L^{p}} = (-1)^{\alpha} (\partial^{\alpha} \Phi, f)\_{L^{p}} \, \forall \Phi \in C\_{0}^{\infty} \tag{5.5}$$

holds.

<sup>1</sup> The space *Lp*(Ω) of Lebegues-integrable functions is defined as

the physically motivated assumptions on the elastic free-energy, see section 3.4.

### **5.2 Generalized notions of convexity**

If one considers a non-linear hyperelastic material the notion of weak lower semicontinuity is linked with convexity of the corresponding free-energy. In the case of a scalar functional depending on vectorial or tensorial quantities different generalized notions of convexity are defined (Dacorogna, 2008). These are convexity, polyconvexity, quasiconvexity and rank 1 convexity, given in the order of decreasing generality. The chosen type of convexity poses constraints on the constitutive relations as will be discussed in the following.

**Convexity** In the one dimensional case convexity translates to the graphical illustration of a function, where the connecting line between two function values is always above or equal to the lowest function value in that segment. With regard to a hyperelastic material convexity states that a homogeneous deformation is energetically favorable (strictly convex) or equal (convex) to an inhomogeneous deformation state. Furthermore, it is ensured that the stress *σ*ˇ increases monotonously with increasing deformation and that the tangent in the intermediate configuration has positive eigenvalues (see e.g., Silhavy, 1997, p. 257-258).

Convexity of the free-energy *W*(*F*) implies weak lower semicontinuity of the energy functional Π(*u*) in *W*<sup>1</sup>*,p* if supplemented by a suitable coercivity condition and upper bound (Dacorogna, 2008, p. 105-108). The drawback of a convex free-energy is that due to a unique minimizer, material instabilities like buckling are excluded. Furthermore, convexity in combination with the principle of material frame indifference implies unphysically restriction on the principal stresses of the Kirchhoff stress tensor (Truesdell, 1964; Ciarlet, 1988), see also

appendix Section 11. Therefore, a convex free-energy is not desirable from a physical point of view. Hence, generalized (weaker) convexity notions are necessary, which still allow for weak lower semicontinuity.

**Quasiconvexity** The weakest known generalized convexity which still leads to weak lower semicontinuity of the energy functional Π(*u*) is quasiconvexity. Quasiconvexity has been introduced by Morrey et al. (1952), see also Meyers (1965). Quasiconvexity states that the energy of a homogeneously deformed body is less or equal than that of a body with a perturbed displacement field and consequently non homogeneous deformation. The concept of quasiconvexity can also be linked with the principle of material stability, see (Krawietz, 1986, p. 315-316). A homogeneous body constrained to homogeneous boundary conditions is in a state of material stability, if for an isothermal process starting from this state no work or energy can be released. From a practical point of view, quasiconvexity is hard to verify due to its non-local definition (see e.g., Morrey et al., 1952).

**Polyconvexity** Polyconvexity (Ball, 1976) is stronger than quasiconvexity, yet allowing for multiple possible minimizers. In contrast to the integral inequality of quasiconvexity it is formulated as a pointwise, local condition. This allows for easier investigation of free-energy functionals *W*(*F*). The free energy functional is required to be representable as

$$W(F) = g(F, \text{cof}(F), \det(F)),\tag{5.6}$$

(Ball, 1976; Dacorogna, 2008). Polyconvexity implies quasiconvexity and therefore weak lower semicontinuity of the energy functional Π(*u*) (Ball, 1976, p. 377). Corresponding growth conditions can be found in (Ball, 2002; Müller et al., 1994; Fonseca et al., 2005) and are less restrictive than those associated with quasiconvexity. Polyconvexity relates the free-energy to the deformations of line-, area- and volume elements of the body B. It also satisfies ellipticity and consequently only allows

for real wave speeds of the considered material.

**Rank 1 convexity** The weakest notion of convexity is rank 1 convexity. For dimensions greater than *n* = *m* = 2 Šverák (1992) showed that in general rank 1 convexity does not imply quasiconvexity, which hints that the minimizing deformation fluctuations due to the microstructure are more complex, than those imposed by homogeneous rank 1 connected deformation gradients. Consequently, weak lower semicontinuity of the functional Π(*u*) cannot be established in general. As a result the existence of minimizers cannot be guaranteed. Nevertheless, rank 1 convexity plays an important role in modeling phase transitions and stable coexisting phases. A rank 1 convex function is convex along rank 1 connected lines, e.g. the difference of the initial and the end point of the line segment are given by a rank 1 tensor. Analogously to convexity, the stress *σ*ˇ is monotonously increasing, but now with respect to a rank 1 connected tensor, as opposed to convexity. Similarly as in the convex case the material favors a homogeneous deformation state over a heterogeneous deformation state where the deformation gradients differ by rank 1. Conversely, if rank 1 convexity is violated the material breaks up into heterogeneous deformation states. This heterogeneous deformation states can be identified by different (deformation) compatible phases.

The mentioned generalized notions of convexity can be related to each other in a sense that the stronger implies the weaker. The reverse is in general not true (Dacorogna, 2008).

$$\text{Convex} \to \text{Polyconv} \times \text{Quasisinvover} \to \text{Rank 1 convex} \qquad (5.7)$$

The discussion above solely focused on (geometrical) non-linear elasticity. However, a wide class of problems in continuum mechanics also deals with physical non-linearities like plasticity, damage and phase transition. The corresponding free-energies or generalized potentials Π(*F*) often lack any convexity at all. The consequence are mesh dependent results, localization phenomena, fine scale oscillating minimizing sequences of the solution field with infinitely local solutions or no solution at all. Remedies for these problems are:


## **5.3 Incremental energy minimization for inelastic solids**

The investigation of non-convex inelastic problems with respect to the generalized convexity notions is made possible by introducing an incremental stress potential Π*n*+1(*F <sup>n</sup>*+1). This allows to analyze inelastic problems using a variational formulation. The corresponding boundary value-problem can then be formulated as a sequence of incremental minimization problems. In their seminal work concerning this subject Ortiz and Repetto (1999) derived an expression for the incremental stress potential Π*n*+1(*F <sup>n</sup>*+1) depending on the rate of the free energy and the power of the internal variables *α* and *β* respectively

$$\Pi^{n+1}(\mathcal{F}^{n+1}) = \inf\_{\underline{\alpha}(t)} \int\_{t\_0}^{t\_1} \dot{W}(\mathcal{F}, \underline{\alpha}) + \underline{\beta} \cdot \dot{\underline{\alpha}} \,\mathrm{d}t \tag{5.8}$$

and initial condition

$$
\underline{\alpha}(t^n) = \underline{\alpha}^n \tag{5.9}
$$

This idea was extended by Carstensen et al. (2002) and Miehe and Lambrecht (2003) by including a dissipation potential Φ, which is often chosen as the dissipation distance, in the sense of a generalized standard material for rate independent materials

$$\Pi^{n+1}(F^{n+1}) = \inf\_{\underline{\alpha}(t)} \int\_{t^n}^{t\_{n+1}} \dot{W}(F, \underline{\alpha}) + \Phi(\underline{\alpha}, \dot{\underline{\alpha}}) \, \mathrm{d}t \tag{5.10}$$

where *∂*Φ/*∂α*˙ = *β* holds. It can be shown that Eq. (5.10) satisfies the Euler and Biot equation for small time steps and hence provides a pointwise approximation (see e.g., Miehe and Lambrecht, 2003). With the incremental stress potential at hand related convexity and growth conditions of Π*<sup>n</sup>*+1 as well as weak lower semicontinuity of the corresponding integrand can be investigated in an incremental sense. The above introduced conditions for solutions in elasticity and relaxation strategies, which will be discussed in the following section, can then be incrementally applied to inelastic materials.

#### **5.4 Relaxation**

**Relaxation for materials with microstructure** Relaxation is often applied to treat non-convex problems, for materials with microstructure. In the relaxed problem the development of a microstructure regularizes the originally non-convex problem. The microstructure can be understood as a physical realization of minimizing sequences of the nonconvex integrand (see e.g., Ball and James, 1987). Fine microstructures arise in materials which undergo phase transition like shape-memory alloys, dual -, and trip-steels as wells as in plasticity where twinning

is often observed in metallic materials. A direct solution of the nonconvex free-energy is possible, see e.g. Tadmor et al. (1999), but practically unfeasible. This is due to high computational costs related to the (re)meshing required to capture the fine scale oscillations of the solution representing the microstructure. Therefore, the original non-convex problem is often replaced by its relaxed convex-, polyconvex-, quasi-, or rank 1 convex envelope. The respective envelope is defined as their corresponding pointwise supremum which is bounded from above by the respective non-convex potential. The corresponding envelope then describes the macroscopic response of the material. The fractions of the different microstructural constituents can be interpreted in terms of a probability measure. The exact arrangement and number of alternating deformation gradients is only reflected in a statistical sense. Thus the microstructure is condensed via a statistical averaging process. The optimal microstructure is in general not unique, see e.g. Kohn (1991). In the following, the definitions for the respective envelopes are given with respect to the free-energy function *W*(*F*) depending on the deformation gradient *F*. These definitions are also valid in the general context of the incremental stress potential Π*<sup>n</sup>*+1(*F <sup>n</sup>*+1) or any other arbitrary function *f*(*A*) depending on a general tensor valued argument *A*.

**Convex envelope** The convex envelope is defined as (Dacorogna, 1985, p. 418) or (Schröder and Hackl, 2014, p. 87-88)

$$CW(F) = \sup\{f(F) \le W(F) \mid f(F) \text{ convex}\},\tag{5.11}$$

and can be computed by

$$CW(F) = \inf\_{\lambda\_{\alpha}, F\_{\alpha}} \left\{ \sum\_{\alpha} \lambda\_{\alpha} W(F\_{\alpha}), \sum\_{\alpha} \lambda\_{\alpha} = 1, \, \lambda\_{\alpha} > 0, \, \sum\_{\alpha} \lambda\_{\alpha} F\_{\alpha} = F \right\} \tag{5.12}$$

Regarding a mechanical problem the convex envelope yields the lowest possible energy of the original non-convex free-energy or potential. Apart from the drawbacks of a convex free-energy mentioned above, another one becomes obvious if one considers *λ<sup>i</sup>* as the volume fractions of different material phases with corresponding phase deformations *F<sup>i</sup>* and *F* as an average macroscopic deformation. Although the macroscopic deformation is recovered as the volume averaged phase deformations *F<sup>i</sup>* , the geometric compatibility condition Eq. (3.17) between the individual phases is in general not met. However, the jump of momentum between the phases, represented by the continuity of the traction vector, is satisfied since the above minimization results in equal stresses for both phases. The convex envelope is therefore identical to the Reuss-bound (see also chapter 9), which in terms of homogenization represents a lower bound with respect to the free energy of the material (see e.g. Castaneda, 1992) and thus generalized tangent moduli.<sup>2</sup> **Polyconvex envelope** The polyconvex envelope is defined as (Dacorogna, 1985, p. 418) or (Schröder and Hackl, 2014, p. 87-88)

$$PW(F) = \sup\{f(F) \le W(F) \mid f(F) \text{ polynover}\},\tag{5.13}$$

and can be computed by

$$PW(F) = \inf\_{\lambda\_{\alpha}, F\_{\alpha}} \left\{ \sum\_{\alpha} \lambda\_{\alpha} W(F\_{\alpha}), \ \sum\_{\alpha} \lambda\_{\alpha} = 1, \ \lambda\_{\alpha} > 0, \ \sum\_{\alpha} \lambda\_{\alpha} F\_{\alpha} = F, \ \sum\_{\alpha} \lambda\_{\alpha} > 0 \ \text{ret}(F) \right\}. \tag{5.14}$$

$$\sum\_{\alpha} \lambda\_{\alpha} \text{cof}(F\_{\alpha}) = \text{cof}(F), \quad \sum\_{\alpha} \lambda\_{\alpha} \text{det}(F\_{\alpha}) = \text{det}(F), F \in \text{Inv}^{+} \begin{cases} \text{cof} \\ \end{cases}. \tag{5.14}$$

<sup>2</sup> The corresponding upper bound (Voigt bound) is found by minimization of the Gibbs energy with respect to the stresses, which results in constant phase deformations. In this case compatibility of the individual phases is satisfied at the expense of traction continuity, which in general does not hold.

The energy of the polyconvex envelope is larger than the convex. Yet, it preserves a correct mapping of the line-, area- and volume elements from the different phases on the micro- to the macroscopic level.

**Quasiconvex envelope** The quasiconvex envelope is given by (Dacorogna, 1982, p. 104)

$$QW(F) = \sup \{ f(F) \le W(F) \: \mid \: f(F) \text{ quasiconvex} \}, \tag{5.15}$$

and can be computed by

$$QW(\mathbf{F}) = \inf\_{\varphi} \left\{ \int\_{\Omega} W(\mathbf{F} + \text{Grad}\,(\varphi)) \, \text{d}V, \, \varphi \in C\_0^{\infty}(\Omega, \mathbb{R}^n), \, \varphi = \mathbf{0} \in \partial\_{\Omega\_u^\*} \right\}.\tag{5.16}$$

Again, the quasiconvex envelope is hard to determine since it requires an infinite dimensional minimization with respect to the displacement perturbation *ϕ*. The displacement perturbation *ϕ* reflects the fluctuation of the deformation field due to the different phases on the microscale. Therefore, the quasiconvex envelope determines an energetically optimal microstructure by displacement fluctuations in *C*<sup>∞</sup> 0 . The resulting energy is larger than or equal to the polyconvex envelope. **Rank 1 convex envelope** The rank 1 convex envelope is given by (Dacorogna, 1985, p. 418) or (Schröder and Hackl, 2014, p. 87-88)

$$RW(F) = \sup \{ f(F) \le W(F) \: \mid \: f(F) \text{ rank } 1 \text{ convex} \},\tag{5.17}$$

where an approximation by first order laminates can be computed by

$$R\_1 W(\mathcal{F}) = \inf\_{\lambda\_\alpha, \mathcal{F}\_\alpha} \left\{ \sum\_{\alpha=1}^2 \lambda\_\alpha W(\mathcal{F}\_\alpha), \sum\_{\alpha=1}^2 \lambda\_\alpha = 1, \lambda\_\alpha > 0, \\ \sum\_{\alpha=1}^2 \lambda\_\alpha \mathcal{F}\_\alpha = \mathcal{F}, \text{Rank}(\mathcal{F}\_2 - \mathcal{F}\_1) \le 1 \right\}. \tag{5.18}$$

In general one arrives at the rank 1 convex envelope by repeating the laminate averaging process for each phase infinitely often (see e.g., Schröder and Hackl, 2014)

$$R\_k(F) \xrightarrow[k \to \infty]{} RW(F). \tag{5.19}$$

How often one has to repeat the lamination process for a close approximation of the rank 1 convex envelope crucially depends on the applied macroscopic deformation and material intrinsic inelastic deformations. The rank 1 convex envelope can be related to a relaxation with respect to a specific fixed microstructure. While the convex envelope creates an isotropic microstructure, the rank 1 convex envelope exhibits an anisotropic, oriented microstructure, characterized by the rank 1 convex connection of homogeneous phase deformation gradients inside each laminate. The rank 1 relaxation process leads to a geometrically compatible microstructure. In addition, the force equilibrium between each phase is satisfied. From a physical point of view the rank 1 convex envelope is very well suited to model materials with nonconvex potentials, where laminate microstructures arise.

Such laminate microstructures are observed in materials undergoing martensitic phase transformation like e.g. shape memory alloys, dual phase and trip steels (Ball and James, 1987; Bhattacharya, 2004; Pereloma and Edmonds, 2012). They are also observed in crystal plasticity due to motion of dislocations (see e.g., Christian and Mahajan, 1995). In crystal plasticity the phenomenom of latent hardening is attributed to the evolution of laminate microstructures (see e.g., Ortiz and Repetto, 1999).

#### **Chapter 6**

## **Modeling phase transformation with a sharp interface theory**

## **6.1 Basic thermomechanical, finite deformation phase transformation model**

**A laminate approach to model phase transformation** This section deals with a finite deformation approach to model thermomechanical, solid to solid phase transformation. Solid to solid phase transformation is e.g., observed in steels and shape memory alloys, where a high temperature phase austenite, transforms under suitable thermal and mechanical boundary conditions into the low temperature phase martensite. At first, the general concept of the model is outlined, constituting a basic model for the description of thermomechanical phase transformation. In the following the basic model is modified to account for the characteristics of the martensitic phase transformation in low carbon steels. The basic model is built on a rank 1 laminate, see chapter 4, and the sharp interface theory of moving strain discontinuities (Abeyaratne and Knowles, 1990; Silhavy, 1997) in a finite deformation context. All equations are formulated with respect to the reference placement. This renders the density , the volume and surface of the body as well as the dimensions of the singular surface, separating the participating phases, independent of the deformation. The stress in the reference configura-

**Figure 6.1:** Microstructure of DP-steel (Böhlke et al., 2014) and corresponding rank 1 laminate model

tion is given by the first Piola-Kirchhoff stress tensor *σ*ˇ. Now, a body B is considered consisting of austenite and possible other phases, which in the following are assumed not to be subjected to phase transformation, like e.g. ferrite. Each austenite material point which is undergoing phase transformation is represented by a laminate substructure, see Fig. 6.1. The body B is subjected to external loads, temperature fields and kinematical constraints. Each material point of that body is exposed to an effective macroscopic deformation gradient *F*¯ and an effective, homogeneous temperature ¯*θ*. The effective deformation gradient *F*¯ results from the mechanical boundary value problem and consequently varies in space. In contrast, the effective temperature ¯*θ* is assumed to be homogeneous inside the body B or material point as well as in the laminate substructure. Thus, a change of the effective temperature is considered as a slow equilibrium process, which implies stationarity of the balance of energy, respectively balance of internal energy in the case of a vanishing velocity of the body B. The most basic laminate substructure is given in terms of a rank 1 laminate. The rank 1 laminate comprises an austenite phase (*A*) and martensite phase (*M*), separated by an interface with normal *N*. In the elastic case the austenite and respectively martensite phase is identified by volume fractions *cA*, *cM*, total homogeneous deformation gradients *F <sup>A</sup>*, *F <sup>M</sup>* as well as possible transformation deformation gradients, constant in space and time, *F At*, *F Mt*. With the above assumption the temperature of both

phases reads

$$
\theta\_A = \theta\_M = \bar{\theta} = \theta. \tag{6.1}
$$

Additionally, the following relations, expressing a constant laminate volume, deformation compatibility of the two phases and compatibility of the phase deformation gradients with the effective deformation gradient, see section 4 Eq. (4.7) and Eq. (4.9), hold

$$c\_A + c\_M = 1,\tag{6.2}$$

$$c\_A F\_A + c\_M F\_M = \bar{F},\tag{6.3}$$

$$\{F\} = \mathfrak{a} \otimes \mathcal{N}.\tag{6.4}$$

With these identites the phase deformations follow as

$$F\_A = \bar{F} - c\_M \mathbf{a} \otimes \mathbf{N}, \quad F\_M = \bar{F} + c\_A \mathbf{a} \otimes \mathbf{N}. \tag{6.5}$$

The normal *N* characterizing the orientation of the interface is interpreted as the habit plane between austenite and martensite (see, e.g., Ball and James, 1987). It can either be prescribed based on experimental observations or determined by consideration of energy minimization states, see section 6 or Ball and James (1987).

The effective free-energy *W*¯ and effective first Piola-Kirchhoff stress tensor ¯*σ*ˇ of the rank 1 laminate are given by, compare section 4 Eq. (4.10) and Eq. (4.11),

$$
\bar{W} = c\_A W\_A(\mathcal{F}\_A, \mathcal{F}\_{TA}, \theta) + c\_M W\_M(\mathcal{F}\_M, \mathcal{F}\_{MT}, \theta), \tag{6.6}
$$

$$
\bar{\breve{\sigma}} = c\_A \breve{\sigma}\_A + c\_M \breve{\sigma}\_M. \tag{6.7}
$$

The mechanical equilibrium of the rank 1 laminate in the isothermal case is characterized by traction continuity and a vanishing driving force, see section 3.1 Eq. (3.31) and Eq. (3.60),

$$\{\check{\sigma}\}N = 0,\tag{6.8}$$

$$f = W\_M - W\_A - \mathbf{a} \cdot \langle \check{\sigma} \rangle \mathbf{N} = 0. \tag{6.9}$$

A non-equilibrium state of the laminate promotes phase transformation, which on the laminate level is related to motion of the interface, separating the two phases, see chapter 4, Eq. (4.7). The motion of the interface is governed by the driving force *f*, see section 3.1. A new equilibrium state is subsequently characterized by a new position of the interface and hence different volume fractions of austenite and martensite, see Fig. 4.1. Accordingly, the deformation and stress state inside the laminate phases changes.

**Kinetic relation interface** A vanishing driving force is identified with an equilibrium state of the laminate. The corresponding phase transformation process leading to that equilibrium state is dissipationless and occurs instantaneous in time.

Usually, in materials undergoing solid to solid phase transformation the phase transformation process is accompanied by dissipation as e.g. visible in the stress hysteresis of shape memory alloys. Moreover, it is inherently time dependent. The velocity, with which the phase interface propagates, is a specific material property. It can reach speeds, e.g. in steels, approaching the speed of sound (Bhadeshia and Honeycombe, 2017). To determine the propagation speed of the phase interface a relation of the driving force and the normal velocity of the interface has to be established, also denoted as kinetic relation. As discussed in section 3.1 the second law of thermodynamics requires the product of the normal velocity of the interface *v* ∗ <sup>⊥</sup> and the driving force to be non-negative

$$v\_{\perp}^{\*}f \ge 0,\tag{6.10}$$

where the normal velocity of the interface can be related to the evolution of the martensite volume fraction by Eq. (4.8). In the framework of Generalized Standard Materials one can now choose a phenomenological dissipation potential, depending on the driving force, of the form

$$
\Phi = \frac{\dot{c}\_0 f\_0}{m+1} \left\langle \frac{|f| - f\_c}{f\_0} \right\rangle^{m+1}.\tag{6.11}
$$

Treating the driving force *f* and the martensite volume fraction *c<sup>M</sup>* as power-conjugated quantities, the evolution of the martensite volume fraction reads

$$\dot{c}\_{M} = \frac{\partial \Phi}{\partial f} = \dot{c}\_{0} \left\langle \frac{|f| - f\_{c}}{f\_{0}} \right\rangle^{m} \text{sgn}(f). \tag{6.12}$$

The critical driving force *f<sup>c</sup>* can be interpreted as an energy barrier which has to be overcome to propagate the interface. *f*<sup>0</sup> represents a normalizing constant and is coupled with the intrinsic interface velocity *c*˙0, which is a material constant. The exponent *m* reflects the rate sensitivity. Decomposing *c*˙*<sup>M</sup>* = *v* ∗ <sup>⊥</sup>*/L* and *c*˙<sup>0</sup> = *v* ∗ <sup>⊥</sup><sup>0</sup>*/L*, introducing *L* as the length of the laminate which is equal to a characteristic length of the considered microstructure, see Fig. 4.1, Eq. (6.12) can be formulated equivalently by

$$v\_{\perp}^{\*} = v\_{\perp 0}^{\*} \left\langle \frac{|f| - f\_c}{f\_0} \right\rangle^m \text{sgn}(f). \tag{6.13}$$

Eq. (6.13) implies a kinetic relation for the motion of the interface. For a non-negative material intrinsic velocity *v* ∗ <sup>⊥</sup><sup>0</sup> and positive normalizing constant *f*<sup>0</sup> Eq. (6.13) satisfies the second law of thermodynamics on the interface

$$f v\_{\perp}^{\*} = v\_{\perp 0}^{\*} \left\langle \frac{|f| - f\_c}{f\_0} \right\rangle^m |f| \ge 0. \tag{6.14}$$

Hence, applying a dissipation potential reflects a possibility to introduce a simple, thermodynamically consistent kinetic relation. The above kinetic relation constitutes a non-linear relationship between the driving force *f* and the normal velocity of the interface *v* ∗ <sup>⊥</sup>. In the special case of *m* = 1 a linear relation of the driving force and the interface velocity is recovered. This compares to the likewise linear relation between the normal velocity of the interface and the driving force *f* with a temperature dependent parameter *ν*(*θ*) as a proportionality factor proposed byAbeyaratne and Knowles (1990)

$$
v\_{\perp}^{\*} = \nu(\theta) f.\tag{6.15}$$

The influence of temperature in the above non-linear kinetic relation can be included by assuming a temperature dependency of *f<sup>c</sup>* = *f*(*θ*) and *v* ∗ <sup>⊥</sup> = *v* ∗ ⊥(*θ*). This reflects changes due to temperature in the material specific normal velocity of the interface *v* ∗ <sup>⊥</sup><sup>0</sup> and the energetic barrier *fc*.

Applying the above kinetic relation can also be interpreted as a Perzyna type viscous evolution of the martensite volume fraction. In this case the viscous damping parameter is identified with *η* = 1*/c*˙<sup>0</sup>

$$
\dot{c} = \frac{1}{\eta} \left\langle \frac{|f| - f\_c}{f\_0} \right\rangle^m \text{sgn}(f). \tag{6.16}
$$

This, from a numerical point of view, can be interpreted as a viscous regularization for the martensite evolution equation. As discussed, rank 1 convexity and in particular the chosen approach using rank 1 laminates, is in general not equal to quasiconvexity. A viscous regularization, where *η* is now interpreted as a numerical damping parameter, is able to recover for small time steps and suitable chosen damping parameters quasiconvexity (Needleman, 1988; Friesecke and Dolzmann, 1997).

## **6.2 Equivalence of energy minimization and sharp interface theory**

**General aspects of energy minimization** The transformation from austenite to martensite allows a material point to locally reduce its energy due to given load and temperature boundary conditions. The reduction of energy is thus the main aspects which drives phase transformation. Therefore, martensitic phase transformation is often considered in terms of energy minimizing principles (Ball and James, 1987; Wang and Khachaturyan, 1997; Levitas, 1998). An energetic minimum describes an equilibrium state of a body B subjected to given boundary conditions. Energy minimization also satisfies thermodynamic consistency (see, e.g., Müller, 2007), if for a given set of boundary conditions the minimization is performed with respect to a suitable energy, as discussed in the following. Dealing with an adiabatic, closed system an equilibrium state of that system for given boundary conditions is identified with a maximum of the entropy (Müller, 2007). When analyzing a body in continuum mechanics, the notion of an adiabatic system is in general not satisfied. However, accounting for the energy balance and entropy inequality (Eqs. (3.36) and (3.43)) it can be shown that for given constant temperature and displacement boundary conditions the Helmholtz free energy is minimized. Therefore, a nonadiabatic material body B with a spatially homogeneous temperature field *θ* 6= *f*(*x*) is considered. It is assumed that the constitutive relation for the entropy flux *η<sup>S</sup>* = *q*ˇ*/θ* is valid for the considered material. Multiplying Eq. (3.43) with *θ* 6= *θ*(*x*) and substituting the result into Eq. (3.36) yields

$$\rho \not\models \left(\dot{\eta}\theta - \dot{e} - \frac{1}{2}(v \cdot v) \cdot \right) + \check{\rho}\mathbf{b} \cdot v + \text{Div}\left(\check{\sigma}\right) \cdot v + \check{\sigma} \cdot \dot{F} \ge 0,\tag{6.17}$$

where Div *σ*ˇ T *v* = Div (*σ*ˇ) · *v* + *σ*ˇ · *F*˙ has been used. Additionally, considering a static problem, *v* = **0**, and the Legendre transformation of the Helmholtz free energy Ψ = *e* − *θη*, Eq. (6.17) results in

$$
\dot{\rho}\dot{\Psi} \le \dot{\boldsymbol{\sigma}} \cdot \dot{\mathbf{F}} - \dot{\rho}\dot{\theta}\eta.\tag{6.18}
$$

Hence, static equilibrium of a body B, subjected to given constant temperature and deformation boundary conditions, *F*˙ = 0 and ˙*θ* = 0, is identified by a minimum of the Helmholtz free energy<sup>1</sup> (Müller, 2007). Considering a two phase material comprising e.g. austenite and martensite subjected to displacement and temperature boundary conditions, an equilibrium state is characterized by a minimum of the effective Helmholtz free energy *W*¯ . However, due to the (homogeneous) eigendeformation of the martensite phase the effective Helmholtz free energy *W*¯ is non-convex, respectively non quasi-convex. Consequently, for given temperature and displacement boundary conditions, a homogeneous deformation state of the material might become energetically unfavorable. In that case, the material can further reduce its energy by breaking up into a mixture of phases, exhibiting a microstructure in a way that it satisfies compatibility with the imposed effective deformation. This implies, for prescribed temperature and displacement boundary conditions, a minimization of the effective Helmholtz free energy with respect to the non-constant phase fractions, respectively phase deformations, and yields a point on the convex envelope. For that point an energetically optimal arrangement of the considered microstructure, characterized by the volume fractions and phase deformation gradients subjected to the constraints discussed in section 5,

<sup>1</sup> In the context of the solution of a non-linear incremental problem by use of a displacement controlled Finite-Element analysis, each converged increment implies an equilibrium state and hence a minimum of the Helmholtz free energy for a fixed (converged) deformation and temperature field. Thus, the set of solutions at the minimum represents a sequence of equilibrium states, where each increment corresponds to an equilibrium state along the deformation and temperature path, prescribed by the respective boundary condition.

is found. To account for a deformation compatible microstructure of the two phases, minimization of the effective Helmholtz free energy is performed with an additional constraint, imposing the rank 1 connection of the phase deformation gradients (see section 5). This results in an approximation of the rank 1 convex envelope by a first order laminate, comprising two phases.

**Energy minimization and sharp interface theory for the elastic case** By constructing an approximation of the rank 1 convex envelope an upper bound for the quasiconvex envelope can be found, see section 5. The rank 1 convex envelope with respect to first order laminates is given by (compare Eq. (5.18))

$$W\_{R1} = \inf\_{\mathbb{Z}\_{R1}} \{ \bar{W}, \ c\_A + c\_M = 1, \ c\_A \mathcal{F}\_A + c\_M \mathcal{F}\_M = \bar{\mathcal{F}},$$

$$\mathcal{F}\_M - \mathcal{F}\_A = \mathfrak{a} \otimes \mathcal{N} \}, \tag{6.19}$$

$$\underline{x}\_{R1} = \left(c\_A, c\_M, \mathbf{F}\_A, \mathbf{F}\_M\right)^T. \tag{6.20}$$

The relaxed energy, given by *WR*<sup>1</sup>, is equal or lower than the Helmholtz free energy of a homogeneous deformation state with a single austenite or martensite phase. Therefore, it allows the material to reduce its energy at a material point by a locally inhomogeneous deformation state. To compute the relaxed energy *WR*<sup>1</sup> three Lagrange parameters *λ*1*,λ*2*,λ*<sup>3</sup> are introduced taking into account the constraints with respect to a constant laminate volume, deformation compatibility between the austenite and martensite as well as compatibility with the imposed effective deformation *F*¯ . The Lagrange functional which is to be minimized now reads

$$\mathcal{L} = \bar{W} + \lambda\_1 (c\_A + c\_M - 1) + \lambda\_2 \cdot (\mathcal{F}\_M - \mathcal{F}\_A - \mathbf{a} \otimes \mathcal{N})$$

$$+ \lambda\_3 \cdot (c\_A \mathcal{F}\_A + c\_M \mathcal{F}\_M - \bar{\mathcal{F}}).\tag{6.21}$$

The relaxed energy *WR*<sup>1</sup> is given by the infimum of the Lagrange functional

$$\bar{W}\_{R1} = \inf\_{\underline{\mathcal{L}}\_{\mathcal{L}}} \left\{ \mathcal{L} \right\}, \quad \underline{\mathcal{L}}\_{\mathcal{L}} = \left( c\_A, c\_M, \mathcal{F}\_A, \mathcal{F}\_M, a, \mathcal{N}, \lambda\_1, \lambda\_2, \lambda\_3 \right)^{\mathrm{T}}.\tag{6.22}$$

with minimizing conditions

$$\frac{\partial \mathcal{L}}{\partial c\_A} = W\_A + \lambda\_1 + \lambda\_3 \cdot F\_A = 0,\tag{6.23}$$

$$\frac{\partial \mathcal{L}}{\partial c\_M} = W\_M + \lambda\_1 + \lambda\_3 \cdot F\_M = 0,\tag{6.24}$$

$$\frac{\partial \mathcal{L}}{\partial F\_A} = c\_A \check{\sigma}\_A - \lambda\_2 + c\_A \lambda\_3 = \mathbf{0},\tag{6.25}$$

$$\frac{\partial \mathcal{L}}{\partial F\_M} = c\_M \check{\sigma}\_M - \lambda\_2 + c\_M \lambda\_3 = \mathbf{0},\tag{6.26}$$

$$\frac{\partial \mathcal{L}}{\partial a} = -\lambda\_2 N = 0,\tag{6.27}$$

$$\frac{\partial \mathcal{L}}{\partial \mathbf{N}} = -\lambda\_2^\mathrm{T} \mathbf{a} = \mathbf{0},\tag{6.28}$$

subjected to the constraints

$$\frac{\partial \mathcal{L}}{\partial \lambda\_1} = c\_A + c\_M - 1 = 0,\tag{6.29}$$

$$\frac{\partial \mathcal{L}}{\partial \lambda\_2} = F\_M - F\_A - a \otimes \mathcal{N} = 0,\tag{6.30}$$

$$\frac{\partial \mathcal{L}}{\partial \lambda\_3} = c\_A F\_A + c\_M F\_M - \bar{F} = \mathbf{0}.\tag{6.31}$$

The tensorial Lagrange multipliers *λ*<sup>2</sup> and *λ*<sup>3</sup> can be analytically solved for. *λ*<sup>2</sup> is obtained by multiplying Eq. (6.26) with *c<sup>A</sup>* and subtracting Eq. (6.25) times *cM*. *λ*<sup>3</sup> follows by adding Eq. (6.25) to Eq. (6.26).

$$
\lambda\_2 = -c\_M c\_A (\check{\sigma}\_M - \check{\sigma}\_A),
\tag{6.32}
$$

$$
\lambda\_3 = -c\_A \check{\sigma}\_A - c\_M \check{\sigma}\_M = -\bar{\check{\sigma}}.\tag{6.33}
$$

Hence, *λ*<sup>3</sup> is identical to the negative effective stress and *λ*<sup>2</sup> is the jump of the phase stress tensors weighted by the product of the volume fractions *cAcM*. Substituting the Lagrange parameters *λ*<sup>2</sup> and *λ*<sup>3</sup> back into Eqs. (6.23), (6.24), (6.27) and (6.28) and subtracting Eq. (6.23) from Eq. (6.24), eliminating the Lagrange multiplier *λ*1, yields the necessary equations to compute the relaxed energy *WR*<sup>1</sup>

$$W\_M - W\_A - \mathbf{a} \cdot \bar{\breve{\sigma}} \,\mathrm{N} = W\_M - W\_A - \mathbf{a} \cdot \langle \breve{\sigma} \rangle \,\mathrm{N} = 0,\tag{6.34}$$

$$\{\ddot{\sigma}\}N = 0,\tag{6.35}$$

$$\{\check{\sigma}^{\mathrm{T}}\}a = 0,\tag{6.36}$$

with

$$F\_A = \dot{F} - c\_M a \otimes \mathcal{N}, \quad F\_M = \bar{F} + (1 - c\_M) a \otimes \mathcal{N}.\tag{6.37}$$

Due to Eq. (6.35) the effective stress vector ¯*σN*ˇ is equal to the average stress vector h ¯*σ*ˇ i*N* = 1*/*2(*σ*ˇ*<sup>M</sup>* + *σ*ˇ *<sup>A</sup>*)*N* of the austenite and martensite phase. Thus, Eq. (6.34) is identical to the driving force in the case of an isothermal phase transformation in the sharp interface setting, see Eq. (3.60). Eq. (6.35) satisfies traction continuity at a material singular surface, given by the austenite-martensite interface. Hence, the approximation of the rank 1 convex envelope by a two phase laminate satisfies the interface equilibrium conditions of a two phase material, Eqs. (3.60) and (3.61), under the assumption of a spatially homogeneous density. The latter is due to the small volume change between 2−3% of most martensitic phase transformations approximately valid, despite large directional transformation deformations (compare, e.g., Moyer and Ansell (1975)). Furthermore, Eq. (6.36) represents an additional equation to determine the interface orientation in terms of an energetically optimal arrangement of the microstructure. The conditions (6.34)- (6.36) are often referred to as configurational phase, force and torque

equilibrium (see, e.g., Aubry et al., 2003)<sup>2</sup> . As discussed in section 4 a laminate substructure fulfills the correct mapping of differential line, area and volume elements, Eqs. (4.9), (11.31), (11.19), from the microto the macroscale, and therefore also the relaxed energy *WR*<sup>1</sup>. This in turn implies that the rank 1 relaxation with respect to two phases implicitly satisfies the constraints for the polyconvex envelope of a two phase material, see Eq. (5.14).

**Minimization accounting for dissipative effects of inelastic variables** In the above minimization all quantities describing the microstructure through the two phase rank 1 laminate are treated elastically, i.e. a change of these variables is not accompanied by dissipation. Since the evolution of inelastic variables is accompanied by energy dissipation an extension of the elastic relaxation with respect to the inelastic quantities is not reasonable. To include the evolution of inelastic variables in a minimization framework the principle of minimum dissipation power is applied, see section 3.4

$$\min\_{\underline{\dot{\alpha}}} \left\{ P \right\} = \min\_{\underline{\dot{\alpha}}} \left\{ \dot{W} + \Phi^\* \right\}. \tag{6.38}$$

Therefore, all internal variables of austenite and martensite - without further specification yet - are collected in the vectors *α<sup>A</sup>* and *αM*. *α<sup>A</sup>* and *α<sup>M</sup>* describe separate inelastic processes in austenite and martensite, as for example plasticity, where the dissipation is confined to the respective phase. Since a change of the microstructure in terms of the martensite volume fraction *c<sup>M</sup>* and the interface orientation *N* is also accompanied by dissipation, they are considered as additional inelastic variables. In contrast to *α<sup>A</sup>* and *α<sup>M</sup>* changes of these quantities cause

<sup>2</sup> The notion of torque equilibrium seems to be a little bit artificial. The balance of moments which is given by the symmetry of the Kirchhoff stress tensor *τ* , is satisfied at the interface due to the balance of linear momentum. The symmetry of the effective Kirchhoff stress holds likewise as a consequence of the laminate homogenization, see Eq. (4.28). Therefore, Eq. (6.36) does not represent a 'classical' balance of angular momentum.

dissipation inside both phases, as a change of the martensite volume fraction causes a change of the austenite volume fraction. Reorientation of the interface affects both phases, too. Hence the dissipation of these two quantities is not locally limited to one phase, as can be also seen by the dependency of the Helmholtz free energy for each phase on them. The vector of strain-like internal variables, describing inelastic quantities inside the laminate, is now introduced as

$$\underline{\alpha} = \left( \underline{\alpha}\_{A}, \underline{\alpha}\_{M}, c\_{M}, \mathbf{N} \right)^{\mathrm{T}} \tag{6.39}$$

with corresponding stress-like variables defined as

$$\underline{\beta} = -\frac{\partial W}{\partial \underline{\alpha}} = \left(\underline{\beta}\_A, \underline{\beta}\_M, f, \underline{\zeta}\right)^T. \tag{6.40}$$

The effective dissipation potentials comprising contributions of both phases then read

$$\bar{\Phi}(\underline{\beta}) = (1 - c\_M)\Phi\_A(\underline{\beta}\_A) + c\_M\Phi\_M(\underline{\beta}\_M) + \Phi\_f(f) + \Phi\_N(\zeta), \tag{6.41}$$

$$\bar{\Phi}^\*(\underline{\dot{\alpha}}) = (1 - c\_M)\Phi\_A^\*(\underline{\dot{\alpha}}\_A) + c\_M\Phi\_M^\*(\underline{\dot{\alpha}}\_M) + \Phi\_f^\*(\dot{c}\_M) + \Phi\_N^\*(\dot{\mathcal{N}}).\tag{6.42}$$

For the current work, the total power is defined as the sum of the rate of the partially relaxed effective Helmholtz free energy *<sup>W</sup>*¯˙ *<sup>R</sup>*1*,P T* and the conjugated dissipation potential Φ¯ <sup>∗</sup> ( ˙*α*). The effective Helmholtz free energy is now relaxed with respect to the interface jump vector *a* and the phase deformation gradients *F <sup>A</sup>* and *F <sup>M</sup>*, as for example also applied in Bartel and Hackl (2009) and Bartel et al. (2011)

$$\bar{W}\_{R1, PT} = \inf\_{\mathbf{a}, \lambda\_1, \lambda\_2, \mathbf{F}\_A, \mathbf{F}\_M} \left\{ \bar{W} + \lambda\_1 \cdot (\mathbf{F}\_M - \mathbf{F}\_A - \mathbf{a} \otimes \mathbf{N}) \right\} $$
 
$$+ \lambda\_2 \cdot (c\_A \mathbf{F}\_A + c\_M \mathbf{F}\_M - \bar{\mathbf{F}}) \Big| . \tag{6.43}$$

The total power of the two phase laminate is then defined as

$$P = \dot{\bar{W}}\_{R1, PT}(\bar{F}, \theta, \mathbf{a}, \underline{\alpha}) + \bar{\Phi}^\*(\dot{\underline{\alpha}}).\tag{6.44}$$

The rates of the dissipative variables *α*˙ are assumed to minimize the total power and therefore yield the following stationarity conditions

$$\min\_{\underline{\dot{\alpha}}} \{ P \} \to -\frac{\partial \bar{W}\_{R1, PT}}{\partial \underline{\alpha}} + \frac{\partial \bar{\Phi}^\*}{\partial \underline{\dot{\alpha}}} = -\underline{\beta} + \frac{\partial \bar{\Phi}^\*}{\partial \underline{\dot{\alpha}}} = \mathbf{0}.\tag{6.45}$$

Since dissipation potentials based on the rate of the strain-like internal variables are e.g. in the case of crystal plasticity in general not known, a Legendre transformation of Φ¯ <sup>∗</sup> ( ˙*α*) is applied (see, e.g., Bartel et al., 2011)

$$\bar{\Phi}(\underline{\beta}) = \sup\_{\underline{\alpha}} \left\{ \underline{\beta} \cdot \underline{\alpha} - \bar{\Phi}^\*(\underline{\alpha}) \right\}. \tag{6.46}$$

The Legendre transformation of Φ¯ <sup>∗</sup> ( ˙*α*) yields the same stationarity conditions as the principle of minimum dissipation power (Eq. (6.45)). Hence, the rates of the inelastic variables can be determined by the dissipation potential Φ( ¯ *β*)

$$
\dot{\underline{\alpha}} = \frac{\partial \bar{\Phi}(\underline{\beta})}{\partial \beta} \tag{6.47}
$$

which holds per definition of the Legendre transformation. The complete set of equations to describe the inelastic minimization, using the principle of minimum dissipation power, are then given by the rates of the strain-like inelastic variables for each phases complemented by the traction continuity, as the minimizing condition for the partial relaxed

energy

$$
\dot{\underline{\alpha}}\_A = \frac{\partial \bar{\Phi}(\underline{\beta})}{\partial \underline{\beta}\_A}, \quad \dot{\underline{\alpha}}\_M = \frac{\partial \bar{\Phi}(\underline{\beta})}{\partial \underline{\beta}\_M}, \tag{6.48}
$$

$$\dot{c}\_{M} = \frac{\partial \bar{\Phi}(\underline{\beta})}{\partial f}, \quad \dot{\mathbf{N}} = \frac{\partial \bar{\Phi}(\underline{\beta})}{\partial \underline{\zeta}},\tag{6.49}$$

$$\{\dot{\sigma}\}\_{A/M} \mathbf{N} = 0.\tag{6.50}$$

By choosing a suitable dissipation potential Eq. (6.49)1, as discussed above, yields again a kinetic relation for the propagation of the interface.

#### **Chapter 7**

# **Modeling martensitic phase transformation in low carbon dual phase steel**

## **7.1 Modeling of the single phases austenite and martensite**

**General aspects of martensite in low carbon steel** Due to its low carbon content and high martensite start temperature low carbon steels mainly form lath martensite. Lath martensite reduces the elastic energy by slip or plastic accommodation. Thus, lath martensite as well as the surrounding austenite exhibit high dislocation densities. Additionally, it forms within a multilevel substructure consisting of packets, blocks and laths, see also chapter 1.3.

**Decomposition of the deformation gradient** Before the details of the phase transformation model are outlined the constitutive relations describing the single phases are established. This is, as already mentioned, done in a finite strain framework. To describe elastic and possible inelastic contributions to the total deformation the multiplicative decomposition of the latter according to (Bilby et al., 1955; Kröner, 1959; Lee, 1969) is introduced. This defines for each inelastic contribution to the total deformation an intermediate configuration between the reference B<sup>0</sup> and the current configuration B of a body B. In the case of the

*F <sup>A</sup>* **Figure 7.1:** Austenite configurations

austenite phase which only experiences elastic and plastic deformations the multiplicative decomposition reads

$$F\_A = F\_{eA} F\_{pA} \tag{7.1}$$

defining a plastic intermediate configuration, see Fig. 7.1. The velocity gradient with corresponding elastic and plastic contribution is given by

$$L\_A = \dot{F}\_{eA} F\_{eA}^{-1} + F\_{eA} \dot{F}\_{pA} F\_{pA}^{-1} F\_{eA}^{-1} \tag{7.2}$$

$$L\_{eA} = \dot{F}\_{eA} F\_{eA}^{-1} \qquad L\_{pA} = \dot{F}\_{pA} F\_{pA}^{-1}.\tag{7.3}$$

If a plastically deformed austenite single crystal transforms into martensite, the martensite develops in the plastically distorted austenite lattice. In the following, it is now assumed that the plastic slips and therefore the plastic distortion of the austenite lattice is fully inherited by the newly formed martensite. Hence, the austenite plastic deformation is the first contribution to the martensite deformation gradient, which in the following is denoted by *Fip M* (inherited plasticity). The transformation deformation given by the corresponding Bain deformation represents the second contribution, followed by plastic deformation of the martensite lattice and elastic deformation. Therefore, the total

deformation gradient of the martensite phase is decomposed into

$$F\_M = F\_{e\,M} F\_{p\,M} F\_{t\,M} F\_{ip\,M} \tag{7.4}$$

with intermediate configurations as shown in Fig. 7.2. The martensite transformation deformation is constant in time as well as the inherited plastic deformation, representing a lattice predeformation, which does not evolve in time. Now the inelastic contributions to the martensite deformation gradient as the transformation deformation gradient *F*˜ *t M* are defined as,

$$
\tilde{F}\_{t\,M} = F\_{p\,M} F\_{t\,M} F\_{ip\,M} = R\_{p\,M} U\_{p\,M} F\_{t\,M} F\_{ip\,M}.\tag{7.5}
$$

This shows a similar structure as the transformation deformation gradient defined by the crystallographic theory of martensite Eq. (1.5). The rotation *R* is now given by the rotation of the plastic deformation gradient inside the martensite lattice. The stretch tensor *U* is in this case composed of the stretch of the martensite plastic deformation gradient, *Up M*, and the Bain deformation *Ft M*. The lattice invariant shear *S* is identified with the inherited austenite plastic deformation gradient *Fip M*. The martensite velocity gradient follows as

$$L\_M = \dot{\mathbf{F}}\_{eM} \mathbf{F}\_{eM}^{-1} + \mathbf{F}\_{eM} \dot{\mathbf{F}}\_{pM} \mathbf{F}\_{pM}^{-1} \mathbf{F}\_{eM}^{-1} \tag{7.6}$$

$$L\_e M = \dot{F}\_{eM} F\_{eM}^{-1} \qquad L\_p M = \dot{F}\_{pM} F\_{pM}^{-1} . \tag{7.7}$$

**Crystal plasticity** The evolution of the plastic deformation gradient of the austenite and martensite is modeled by crystal plasticity models. Plasticity in metals is related to motion of dislocations on crystallographic defined closest packed planes in closest packed directions. The closest packed planes are identified by their normal *n<sup>β</sup>* while the closest packed directions on those planes are denoted by a unit vector *dβ*. *n<sup>β</sup>* and *d<sup>β</sup>* are defined in the intermediate configuration B*p*. The set of a

**Figure 7.2:** Martensite configurations

normal and corresponding slip direction (on that plane) constitutes a slip system *d<sup>β</sup>* ⊗ *nβ*. Using the Orowan equation the propagation of dislocations on that slip system can be related to the rate of the plastic slip by

$$
\dot{\gamma}\_{\beta} = \rho\_{\beta} b v\_{\beta}. \tag{7.8}
$$

*ρ<sup>β</sup>* is the dislocation density, *b* the magnitude of the Burgers vector and *v<sup>β</sup>* the averaged dislocation velocity on a slip system. The velocity gradient of the plastic deformation is assumed to be related to the slip rate and the orientation of the corresponding slip systems by (see, e.g., Peirce et al., 1983)

$$\mathbf{L}\_p = \sum\_{\beta=1}^N \dot{\gamma}\_\beta \mathbf{d}\_\beta \otimes \mathbf{n}\_\beta,\tag{7.9}$$

where *N* denotes the number of slip systems, which depends on the crystal symmetry. As can be seen the plastic velocity gradient identifies plasticity with shear in crystallographic defined directions, relating dislocations to a macroscopic shear like deformation. To reduce the

numerical effort an accumulated plastic slip is often introduced

$$\gamma\_{ac} = \int\_{t\_0}^{t\_1} \sum\_{\beta=1}^{N} |\dot{\gamma}\_{\beta}| \mathbf{d}t. \tag{7.10}$$

**Constitutive properties of the single phases** The Helmholtz free energy for each phase *α* consists of a mechanical thermoelastic part, a plastic contribution, related to stored dislocations, and a chemical part

$$W\_{\alpha} = W\_{e\,\alpha} + W\_{p\,\alpha} + W\_{chem\,\alpha}.\tag{7.11}$$

The thermoelastic part of the Helmholtz free energy is modeled by an elastic and thermal isotropic Neo-Hookean energy

$$W\_{e\alpha} = \frac{\mu\_{\alpha}}{2} \left( \text{tr}(\mathbf{C}\_{e\alpha}) - 2\text{ln}(J\_{e\alpha}) - 3 \right) + \frac{\lambda\_{\alpha}}{2} \text{ln}(J\_{e\alpha})^2 - \beta\_{\alpha} (\text{ln}(J\_{e\alpha})) \Delta\theta \tag{7.12}$$

with the Lame constants *µ<sup>α</sup>* and *λα*, the thermal expansion coefficient *α<sup>α</sup>* and related quantity *β<sup>α</sup>* = (3*λ<sup>α</sup>* + 2*µα*)*αα*. Additionally, △*θ* = *θ* −*θ*<sup>0</sup> is given with respect to to a reference temperature *θ*0. The plastic part is now defined as

$$W\_{p\,\alpha} = \int\_0^{\gamma\_{ac\,\alpha}} \tau\_{C\,\alpha}(\gamma\_{ac\,\alpha}) \mathbf{d}\gamma\_{ac\,\alpha} \tag{7.13}$$

with an isotropic hardening, represented by the critical shear stress *τc α*. The critical shear stress represents a threshold value which has to be overcome to allow further propagation of dislocations. Due to the accumulated plastic slip, *τC α* is in the above case identical on each slip system. The chemical part reads (see, e.g., Fischer et al., 1994)

$$W\_{chem\,\alpha} = \check{\rho}\_{\alpha} \left( e\_{0\,\alpha} - \eta\_{0\,\alpha} \triangle \theta - c\_{\alpha,\epsilon} (\theta \ln(\frac{\theta}{\theta\_0}) - \triangle \theta) \right). \tag{7.14}$$

It can be interpreted as the energy to arrange atoms in a specific lattice for a given temperature (see, e.g., Fischer et al., 1994). The chosen Helmholtz free energy is polyconvex for modest volume changes, see definition 5, and satisfies the principle of material frame indifference and material symmetry, see section 3.4. Additionally, the Helmholtz free energy is normalized such that,

$$W\_{\alpha}(F\_{\alpha} = I, \theta = \theta\_0, \underline{\alpha}\_{\alpha} = \mathbf{0}) = 0 \tag{7.15}$$

holds. Collecting all inelastic deformation contributions in *Fin α*, the first Piola-Kirchhoff stress *σ*ˇ *<sup>α</sup>* and the Cauchy stress tensor *σ<sup>α</sup>* are then given by

$$\check{\sigma}\_{\alpha} = \frac{\partial W\_{\alpha}}{\partial F\_{\alpha}} = \mu\_{\alpha} (\boldsymbol{F}\_{e\,\alpha} \boldsymbol{F}\_{\alpha\,\alpha}^{\mathrm{T}}) + \left(\lambda\_{\alpha} \ln \left(J\_{e\,\alpha}\right) - \mu\_{\alpha} - \beta\_{\alpha} \triangle \theta\right) \boldsymbol{F}\_{\alpha}^{\mathrm{T}},\tag{7.16}$$

$$\boldsymbol{\sigma}\_{\alpha} = \frac{1}{J\_{e\,\alpha}} \check{\sigma}\_{\alpha} \boldsymbol{F}\_{\alpha}^{\mathrm{T}} = \frac{1}{J\_{e\,\alpha}} \left(\mu\_{\alpha} \boldsymbol{B}\_{e\,\alpha} + \left(\lambda\_{\alpha} \ln \left(J\_{e\,\alpha}\right) - \mu\_{\alpha} - \beta\_{\alpha} \triangle \theta\right) \boldsymbol{I}\right). \tag{7.17}$$

**Evolution of inelastic variables of the single phases** The set of strainlike internal variables for austenite and martensite comprises the (accumulated) plastic slip and the plastic deformation gradient contribution

$$\underline{\alpha}\_{A} = \left(\gamma\_{ac}\boldsymbol{A}, \boldsymbol{F}\_{p\boldsymbol{A}}\right)^{\mathrm{T}}, \qquad \underline{\alpha}\_{M} = \left(\gamma\_{ac\,M}, \boldsymbol{F}\_{p\,M}\right)^{\mathrm{T}}.\tag{7.18}$$

For the martensite phase the transformation deformation and inherited plasticity is not treated as an internal variables since they are constant in time, see above. The same holds for the inherited slip from the prior austenite phase. Conjugated stress-like internal variables are then derived by

$$-\frac{\partial W\_{\alpha}}{\partial \gamma\_{ac\,\alpha}} = -\tau\_{c\,\alpha}, \qquad -\frac{\partial W\_{\alpha}}{\partial F\_{p\,\alpha}} = M\_{e\,\alpha} F\_{p\,\alpha}^{\mathrm{T}},\tag{7.19}$$

where *M<sup>e</sup>* = *F* T *<sup>e</sup>σ*ˇ *<sup>e</sup>* denotes the elastic Mandel stress tensor. With this quantities at hand dissipation potentials can be derived to model the evolution of the plastic internal variables for austenite and martensite. In the case of austenite, which possesses an fcc lattice, plasticity is driven by edge dislocations. Due to the planar structure of an edge dislocation core, the motion of these type of disclocations is driven by the stress inside the corresponding slip system. In the case of bcc lattices plasticity is characterized by screw dislocations, which show a nonplanar core structure. Therefore, additionally to the stress inside the slip system of the dislocation, stresses from other slip systems contribute to the motion of screw dislocations. In a first step those non-planar stress contributions are not taken into account. The motion of dislocations in a slip system is then related to the projected stress onto that slip system *τβ α* with

$$
\tau\_{\beta\,\alpha} = \sigma\_{\alpha} \cdot \mathbf{F}\_{e\,\alpha} \mathbf{d}\_{\beta\,\alpha} \otimes \mathbf{F}\_{e\,\alpha}^{\mathrm{T}} \mathbf{n}\_{\beta\,\alpha} = \mathbf{F}\_{e\,\alpha}^{\mathrm{T}} \breve{\sigma}\_{e\,\alpha} \cdot \mathbf{d}\_{\alpha} \otimes \mathbf{n}\_{\alpha},\tag{7.20}
$$

where small elastic volume changes, det(*Fe α*) ≈ 1, have been assumed. This is also known as Schmid behavior (Schmid and Boas, 1937). The evolution of the austenite and martensite internal variables can then be modeled by the following phenomenological dissipation potentials

$$\Phi\_{\alpha} = \frac{\dot{\gamma}\_{0\,\alpha} \tau\_{0\,\alpha}}{m+1} \sum\_{\beta=1}^{N} \left\langle \frac{|\tau\_{\beta\,\alpha}| - \tau\_{c\,\alpha}(\gamma\_{ac\,\alpha})}{\tau\_{0\,\alpha}} \right\rangle^{m+1},\tag{7.21}$$

$$\dot{\gamma}\_{ac\,\alpha} = \frac{\partial \Phi\_{\alpha}}{\partial - \tau\_{c\,\alpha}} = \dot{\gamma}\_{0\,\alpha} \sum\_{\beta=1}^{N} \left\langle \frac{|\tau\_{\beta\,\alpha}| - \tau\_{c\,\alpha} \langle \gamma\_{ac\,\alpha} \rangle}{\tau\_{0\,\alpha}} \right\rangle^{m} = \sum\_{\beta=1}^{N} \dot{\gamma}\_{\beta\,\alpha}, \tag{7.22}$$

$$\dot{\boldsymbol{F}}\_{p\alpha} = \frac{\partial \Phi\_{\alpha}}{\partial \mathbf{M}\_{e\alpha} \mathbf{F}\_{p\alpha}^{\rm T}} = \sum\_{\beta=1}^{N} \dot{\gamma}\_{\beta\alpha} \mathbf{d}\_{\beta\alpha} \otimes \mathbf{n}\_{\beta\alpha} \mathbf{F}\_{p\alpha} \text{sgn}\left(\tau\_{\beta\alpha}\right),\tag{7.23}$$

$$L\_{p\,\alpha} = \sum\_{\beta=1}^{N} \dot{\gamma}\_{\beta\,\alpha} \mathbf{d}\_{\beta\,\alpha} \otimes \mathfrak{n}\_{\beta\,\alpha} \operatorname{sgn} \left( \tau\_{\beta\,\alpha} \right). \tag{7.24}$$

99

## **7.2 Laminate approach for low carbon steel**

**Phase transformation model for low carbon steel** The basic phase transformation model, see section 6, is now applied to model the martensitic phase transformation in low carbon steels at the level of a martensite block, see Fig. 7.3. To reduce the computational effort only the three Bain-variants are taken into account, the misorientation of the corresponding Kurdjomov-Sachs variants is therefore not considered. With this approximation a martensite block consists of one Bain variant. The transformation of an austenite single crystal into a martensite crystal of such a block is modeled by a rank 1 laminate. Plasticity inside the martensite and austenite phase plays an important role in the transformation process, due to the plastic accommodation. Martensite shows a heterogeneous distribution of dislocations Takebayashi et al. (2010), originating from plastic inheritance of the prior austenite as well as stress induced dislocations in the regions, which have already transformed into martensite. Now an austenite single crystal or material point at the phase front which transforms into a martensite single crystal is considered. It is assumed that during the process of transformation this material point inherits dislocations of the prior austenite and exhibits plasticity inside the martensite lattice but has no plastic predeformation. The martensite phase transformation is in that case defined between an originally elastic martensite variant, inheriting dislocations of the parent austenite crystal, and a (possibly) plastically predeformed austenite phase.

A laminate only accounts for homogeneous deformation states inside the laminate phases. Therefore, when considering the evolution in a time incremental sense, a two phase rank 1 laminate, as discussed in section 6, is not able to reflect the different levels of intrinsic martensite plasticity, i.e. non-inherited dislocations at the phase front and in regions which have already been transformed. Thus, using a two phase rank 1 laminate would result in a phase transformation of a plastically

**Figure 7.3:** Relation of the laminate model to the martensite substructure, left image according to Kitahara et al. (2006)

predeformed martensite, due to the homogeneous plastic deformation state inside the martensite laminate. Hence, the two phase rank 1 laminate is modified in an incremental sense by a third phase. The third phase represents a transition zone or the region, where during the considered time interval austenite is transformed into martensite. It is denoted by *M*△. Since the transition zone is only an extension of the previously transformed martensite *M*0, it comprises the same transformation deformation and does not represent a microstructure on a lower scale compared to the austenite phase. Therefore, no separation of scales is present and the transition zone is located on the same laminate level as the previously transformed martensite and the untransformed austenite. The previously transformed martensite *M*<sup>0</sup> represent the martensite which has been transformed before the considered time interval. Its volume fraction is constant and only the elastic and plastic deformation and correspondingly the stress state evolve. Considering now a time interval *t <sup>n</sup>* → *t <sup>n</sup>*+1, the martensitic phase transformation is modeled by the evolution of the transition zone. At the beginning of the time interval *t <sup>n</sup>* the laminate only contains the previously transformed martensite phase *M*<sup>0</sup> with volume fraction *c n M*<sup>0</sup> and the austenite phase *A* with volume fraction *c n <sup>A</sup>*. The volume

fraction *cM*△ of the transition zone is zero, see Fig. 7.4 a). If the new

**Figure 7.4:** Laminate at the beginning a), middle b) and end c) of a time interval *t <sup>n</sup>* → *t n*+1

deformation and temperature state promotes phase transformation inside the laminate, the phase transition zone is introduced, where the current transformation from austenite to martensite takes place, see Fig. 7.4 b). At the end of the considered time interval the transition zone and the previously transformed martensite are homogenized into one martensite phase, see Fig. 7.4 c). This martensite phase represents at the beginning of the next time interval the previously transformed martensite *M*0. Due to the inherent non-linearity in the finite strain measures, a simple volume fraction based homogenization of the plastic deformations and plastic slips is not reasonable. The chosen approach is outlined below.

**Laminate relations and effective laminate quantities** In a similar way as outlined in section 4, the following relations hold for the three phase rank 1 laminate. A constant laminate volume implies

$$c\_{M\_0} + c\_{M\_\Delta} + c\_A = 1.\tag{7.25}$$

Compatibility between the neighboring phases and with respect to the effective deformation *F*¯ requires

$$\{F\}\_{M\_0/M\_\Delta} = F\_{M\_0} - F\_{M\_\Delta} = \mathbf{b} \otimes \mathbf{N},\tag{7.26}$$

$$\{F\}\_{M\_{\Delta}/M\_{A}} = F\_{M\_{\Delta}} - F\_{A} = \mathbf{a} \otimes \mathbf{N},\tag{7.27}$$

$$c\_A \mathbf{F}\_A + c\_{M\_0} \mathbf{F}\_{M\_0} + c\_{M\_\Delta} \mathbf{F}\_{M\_\Delta} = \bar{\mathbf{F}}.\tag{7.28}$$

With these conditions at hand localization of the effective deformation into the austenite and martensite phases yields

$$F\_A = \ddot{F} - c\_M a \otimes \mathbf{N} - c\_{M\_0} \mathbf{b} \otimes \mathbf{N},\tag{7.29}$$

$$F\_{M\_0} = \bar{\mathbf{F}} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} + (1 - c\_{M\_0})\mathbf{b} \otimes \mathbf{N},\tag{7.30}$$

$$F\_{M\_{M\_{\Delta}}} = \bar{\mathbf{F}} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} - c\_{M\_0}\mathbf{b} \otimes \mathbf{N},\tag{7.31}$$

where *c<sup>M</sup>* = *cM*<sup>0</sup> + *cM*△. The effective or macroscopic Helmholtz free energy is given by the volume average of the Helmholtz free energy of the phases

$$\bar{W} = c\_A W\_A(\mathbf{F}\_A, \theta, \underline{\alpha}\_A) + c\_{M\_0} W\_{M\_0}(\mathbf{F}\_{M\_0}, \theta, \mathbf{F}\_T, \mathbf{F}\_{M\_0, ip}, \underline{\alpha}\_{M\_0})$$

$$+ c\_{M\_\triangle} W\_{M\_\triangle}(\mathbf{F}\_{M\_\triangle}, \theta, \mathbf{F}\_T, \mathbf{F}\_{M\_\triangle, ip}, \underline{\alpha}\_{M\_\triangle}).\tag{7.32}$$

and the effective first Piola-Kirchhoff stress as

$$
\bar{\tilde{\sigma}} = \frac{\partial \bar{W}}{\partial \bar{F}} = c\_A \check{\sigma}\_A + c\_{M0} \check{\sigma}\_{M\_0} + c\_{M\_\triangle} \check{\sigma}\_{M\_\triangle}.\tag{7.33}
$$

Again it can be shown that the Hill-Mandel condition for a laminate in equilibrium is satisfied.

$$\begin{split} &\frac{1}{|V|}\int\_{V}\dot{\sigma}(\mathbf{z})\cdot\mathbf{F}(\mathbf{z})\mathrm{d}V = c\_{A}\dot{\sigma}\_{A}\cdot\mathbf{F}\_{A} + c\_{M\_{0}}\dot{\sigma}\_{M\_{0}}\cdot\mathbf{F}\_{M\_{0}} + c\_{M\_{\Delta}}\dot{\sigma}\_{M\_{\Delta}}\cdot\mathbf{F}\_{M\_{\Delta}} \\ &= c\_{A}\dot{\sigma}\_{A}\cdot\left(\bar{\mathbf{F}} - c\_{M}\mathbf{a}\otimes N - c\_{M\_{0}}\mathbf{b}\otimes N\right) \\ &+ c\_{M\_{0}}\dot{\sigma}\_{M\_{0}}\cdot\left(\bar{\mathbf{F}} + (1 - c\_{M})\mathbf{a}\otimes N + (1 - c\_{M\_{0}})\mathbf{b}\otimes N\right) \\ &+ c\_{M\_{\Delta}}\dot{\sigma}\_{M\_{\Delta}}\cdot\left(\bar{\mathbf{F}} + (1 - c\_{M})\mathbf{a}\otimes N - c\_{M\_{0}}\mathbf{b}\otimes N\right) \\ &= \bar{\mathbf{\sigma}}\cdot\bar{\mathbf{F}} + (1 - c\_{M})\left(c\_{M\_{0}}\dot{\sigma}\_{M\_{0}} + c\_{M\_{\Delta}}\dot{\sigma}\_{M\_{\Delta}} - c\_{M}\dot{\sigma}\_{A}\right)\cdot\mathbf{a}\otimes N \\ &+ c\_{M\_{0}}\left(\left(1 - c\_{M\_{0}}\dot{\sigma}\_{M\_{0}} - c\_{M\_{\Delta}}\dot{\sigma}\_{M\_{\Delta}} - (1 - c\_{M})\dot{\sigma}\_{A}\right)\cdot\mathbf{b}\otimes N \\ &= \bar{\mathbf{\sigma}}\cdot\bar{\mathbf{F}}. \end{split} \tag{7.34}$$

Moreover, by recursive application of the procedure outlined in section 11 it follows that for the three phase rank 1 laminate the line, see Eq. (7.28), area and volume elements are mapped correctly from the micro- to the macroscale

$$\text{cof}(\bar{F}) = c\_A \text{cof}(F\_A) + c\_{M\_0} \text{cof}(F\_{M\_0}) + c\_{M\_\triangle} \text{cof}(F\_{M\_\triangle}),\tag{7.35}$$

$$\det(\bar{F}) = c\_A \det(F\_A) + c\_{M\_0} \det(F\_{M\_0}) + c\_{M\_\triangle} \det(F\_{M\_\triangle}), \tag{7.36}$$

see also appendix, chapter 11.

#### **7.3 Energy minimization for lath martensite**

**Elastic energy relaxation** As shown before in section 6, energy minimization is equivalent to the equilibrium conditions of a two phase material with a deformation compatible microstructure, represented by a

laminate. In this section it is first shown that the elastic relaxation again coincides with the equilibrium conditions. Afterwards, the evolution of the inelastic variables, i.e. plastic slips and deformations, is included by applying the principle of minimum dissipation power. To compute an approximation of the rank 1 convex envelope of the effective energy for the three phase laminate the following Lagrange functional is introduced

$$\mathcal{L} = \bar{W} + \lambda (c\_A + c\_M - 1) + \lambda\_1 \cdot (\mathcal{F}\_{M\_0} - \mathcal{F}\_{M\_\Delta} - \mathbf{b} \otimes \mathcal{N})$$

$$+ \lambda\_2 \cdot (\mathcal{F}\_{M\_\Delta} - \mathcal{F}\_A - \mathbf{a} \otimes \mathcal{N})$$

$$+ \lambda\_3 \cdot (c\_A \mathcal{F}\_A + c\_M \mathcal{F}\_M + c\_{M\_\Delta} \mathcal{F}\_{M\_\Delta} - \bar{\mathcal{F}}).\tag{7.37}$$

The Lagrange parameters *λ,λ*1*,λ*<sup>2</sup> and *λ*<sup>3</sup> account for a constant volume, deformation compatibility at the interfaces of the three phases and compatibility with respect to the effective deformation gradient as the volume averaged phase deformation gradients. Since *c<sup>M</sup>*<sup>0</sup> is constant, the relaxation is only performed with respect to *c<sup>A</sup>* and *c<sup>M</sup>*△, which can evolve due to phase transformation for given deformation and temperature boundary conditions. The rank 1 convex envelope is then identified by

$$\bar{W}\_{R1} = \inf\_{\underline{x}\_{\mathcal{L}}} \left\{ \mathcal{L} \right\}, \ \underline{x}\_{\mathcal{L}} = \left( c\_A, c\_{M\_{\Delta}}, \mathcal{F}\_A, \mathcal{F}\_{M\_0}, \mathcal{F}\_{M\_{\Delta}}, a, b, N, \lambda, \lambda\_1, \lambda\_2, \lambda\_3 \right)^T \tag{7.38}$$

with minimizing conditions

$$\frac{\partial \mathcal{L}}{\partial c\_A} = W\_A + \lambda + \lambda\_3 \cdot F\_A = 0,\tag{7.39}$$

$$\frac{\partial \mathcal{L}}{\partial c\_M} = W\_{M\_\triangle} + \lambda + \lambda\_3 \cdot \mathbf{F}\_{M\_\triangle} = 0,\tag{7.40}$$

$$\frac{\partial \mathcal{L}}{\partial F\_A} = c\_A \check{\sigma}\_A - \lambda\_2 + c\_A \lambda\_3 = \mathbf{0},\tag{7.41}$$

105

$$\frac{\partial \mathcal{L}}{\partial F\_{M\_0}} = c\_{M\_0} \check{\sigma}\_{M\_0} - \lambda\_1 + c\_{M\_0} \lambda\_3 = \mathbf{0},\tag{7.42}$$

$$\frac{\partial \mathcal{L}}{\partial F\_{M\_{\triangle}}} = c\_{M\_{\triangle}} \check{\sigma}\_{M\_{\triangle}} + \lambda\_2 - \lambda\_1 + c\_{M\_{\triangle}} \lambda\_3 = \mathbf{0},\tag{7.43}$$

$$\frac{\partial \mathcal{L}}{\partial a} = -\lambda\_2 N = \mathbf{0},\tag{7.44}$$

$$\frac{\partial \mathcal{L}}{\partial \mathbf{b}} = -\lambda\_1 \mathbf{N} = \mathbf{0},\tag{7.45}$$

$$\frac{\partial \mathcal{L}}{\partial N} = -\lambda\_2^\mathrm{T} a - \lambda\_1^\mathrm{T} b = 0,\tag{7.46}$$

subjected to the constraints

$$\frac{\partial \mathcal{L}}{\partial \lambda} = c\_A + c\_{M\_0} + c\_{M\_\Delta} - 1 = 0,\tag{7.47}$$

$$\frac{\partial \mathcal{L}}{\partial \lambda\_1} = F\_{M\_0} - F\_{M\_\Delta} - \mathbf{b} \otimes \mathbf{N} = \mathbf{0},\tag{7.48}$$

$$\frac{\partial \mathcal{L}}{\partial \lambda\_2} = F\_{M\_\Delta} - F\_A - a \otimes \mathcal{N} = \mathbf{0},\tag{7.49}$$

$$\frac{\partial \mathcal{L}}{\partial \lambda\_2} = c\_A \mathbf{F}\_A + c\_{M\_0} \mathbf{F}\_{M\_0} + c\_{M\_\triangle} \mathbf{F}\_{M\_\triangle} - \bar{\mathbf{F}} = \mathbf{0}.\tag{7.50}$$

The Lagrange multipliers *λ*1*,λ*<sup>2</sup> and *λ*<sup>3</sup> can again be analytically solved for. By adding Eqs. (7.41) and (7.42) to Eq. (7.43) one finds

$$
\lambda\_3 = -(c\_A \check{\sigma}\_A + c\_{M\_0} \check{\sigma}\_{M\_0} + c\_{M\_\triangle} \check{\sigma}\_{M\_\triangle}) = -\bar{\check{\sigma}}.\tag{7.51}
$$

Substituting *λ*<sup>3</sup> back into Eq. (7.41) and Eq. (7.42) the two remaining Lagrange parameters are identified by

$$
\lambda\_1 = c\_A (c\_M \check{\sigma}\_A - c\_{M\_0} \check{\sigma}\_{M\_0} - c\_{M\_\Delta} \check{\sigma}\_{M\_\Delta}), \tag{7.52}
$$

$$
\lambda\_2 = -c\_{M\_0}(-c\_A \check{\sigma}\_A - (c\_{M\_0} - 1)\check{\sigma}\_{M\_0} - c\_{M\_\Delta} \check{\sigma}\_{M\_\Delta}).\tag{7.53}
$$

Adding Eq. (7.41) and Eq. (7.42), accounting for the solutions of *λ*<sup>1</sup> and *λ*2, yields

$$\{\check{\sigma}\}\_{M\_0/A}N = \mathbf{0}.\tag{7.54}$$

Subtracting Eq. (7.39) from Eq. (7.40), considering the solution for *λ*3, and using Eq. (7.54) to modify Eq. (7.44) and Eq. (7.45) the set of equations corresponding to the relaxed energy for the three phase laminate are given by

$$W\_{M\triangle} - W\_A - \mathbf{a} \cdot \bar{\hat{\sigma}} \mathbf{N} = W\_{M\triangle} - W\_A - \mathbf{a} \cdot \langle \check{\sigma} \rangle \mathbf{N} = 0,\tag{7.55}$$

$$(\check{\sigma}\_{M\_{\triangle}} - \check{\sigma}\_{A})\mathbf{N} = \mathbf{0},\tag{7.56}$$

$$(\check{\sigma}\_{M\_0} - \check{\sigma}\_{M\_\Delta})\mathbf{N} = \mathbf{0},\tag{7.57}$$

$$\left(c\_{M\_0}(-c\_A\check{\sigma}\_A - (c\_{M\_0} - 1)\check{\sigma}\_{M\_0} - c\_{M\_\triangle}\check{\sigma}\_{M\_\triangle})\right)^\mathrm{T}\mathfrak{b} -$$

$$\left(c\_A(c\_M\check{\sigma}\_A - c\_{M\_0}\check{\sigma}\_{M\_0} - c\_{M\_\triangle}\check{\sigma}\_{M\_\triangle})\right)^\mathrm{T}\mathfrak{a} = 0,\tag{7.58}$$

$$F\_A = \dot{F} - c\_M a \otimes \mathbf{N} - c\_{M\_0} \mathbf{b} \otimes \mathbf{N},\tag{7.59}$$

$$F\_{M\_0} = \bar{\mathbf{F}} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} + (1 - c\_{M\_0})\mathbf{b} \otimes \mathbf{N},\tag{7.60}$$

$$F\_{M\_{\triangle}} = \bar{\mathbf{F}} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} - c\_{M\_0}\mathbf{b} \otimes \mathbf{N} \,. \tag{7.61}$$

**Inelastic minimization** Since lath martensite reduces its energy through plastic accommodation or slip, the plasticity related internal variables have to be taken into account. Therefore, as introduced above single crystal plasticity models featuring an accumulated plastic strain are used to model the plastic behavior of each phase. The vector of the strain-like internal variables, describing the three phase rank 1 laminate is partitioned as follows

$$\underline{\alpha} = (\underline{\alpha}\_A, \underline{\alpha}\_{M\_0}, \underline{\alpha}\_{M\_\Delta}, c\_{M\_\Delta}, \mathbf{N})^\mathrm{T}. \tag{7.62}$$

The internal variables, related to the local dissipation of each phase, are the corresponding accumulated plastic slip and deformation gradient

$$
\underline{\alpha}\_A = (\gamma\_{ac \, A}, \mathbf{F}\_{p \, A})^\mathrm{T}, \tag{7.63}
$$

$$\underline{\alpha}\_{M\_0} = (\gamma\_{ac \, M\_0}, \mathbf{F}\_{p \, M\_0})^{\mathrm{T}},\tag{7.64}$$

$$\underline{\alpha}\_{M\_{\triangle}} = \left(\gamma\_{ac}\underline{M}\_{\triangle}, \boldsymbol{\mathcal{F}}\_{p\,M\_{\triangle}}\right)^{\rm T},\tag{7.65}$$

with corresponding stress like internal variables

$$\underline{\beta}\_{A} = (1 - c\_{M})(-\tau\_{cA}, \mathbf{M}\_{eA} \mathbf{F}\_{pA}^{\cdot \mathbf{T}})^{\mathbf{T}},\tag{7.66}$$

$$\underline{\beta}\_{M\_0} = c\_{M\_0} (-\tau\_{c\,M\_0}, \mathbf{M}\_{e\,M\_0} \mathbf{F}\_{p\,M\_0}^{\cdot \mathbf{T}})^{\mathbf{T}},\tag{7.67}$$

$$
\underline{\beta}\_{M\_{\triangle}} = c\_{M\_{\triangle}} (-\tau\_{cM\_{\triangle}}, \mathbf{M}\_{e\,M\_{\triangle}} \mathbf{F}\_{p\,M\_{\triangle}}^{\rm T})^{\rm T}.\tag{7.68}
$$

The dissipation of the volume fraction *c<sup>M</sup>*△ as well as *N* are again related to the whole laminate. However due to the specific orientation of the austenite-martensite interface or habit plane, *N* is fixed after initiation of the transformation process, using the values given in Table 1.3. The effective conjugated dissipation potential describing the three phase rank 1 laminate reads

$$\bar{\Phi} = c\_A \Phi\_A(\underline{\beta}\_A) + c\_{M\_0} \Phi\_{M\_0}(\underline{\beta}\_{M\_0}) + c\_{M\_\triangle} \Phi\_{M\_\triangle}(\underline{\beta}\_{M\_\triangle}) + \Phi(f). \tag{7.69}$$

A corresponding primal dissipation potential can be found by a Legendre-Fenchel transformation

$$\bar{\Phi}^\*(\underline{\dot{\alpha}}) = \sup\_{\underline{\beta}} \left( \underline{\beta} \cdot \dot{\underline{\alpha}} - \Phi(\underline{\beta}) \right). \tag{7.70}$$

The individual contributions to the effective conjugated dissipation potential are given by, compare Eq. (7.21) for plasticity of each phase and Eq. (6.11) to model the kinetic relation,

$$\Phi\_A = \frac{\dot{\gamma}\_{0\,A} \tau\_{0\,A}}{m+1} \sum\_{\beta=1}^N \left\langle \frac{|\tau\_{\beta\,A}| - \tau\_{c\,A}(\gamma\_{ac\,A})}{\tau\_{0\,A}} \right\rangle^{m+1},\tag{7.71}$$

$$\Phi\_{M\_0} = \frac{\dot{\gamma}\_{0\,M\_0} \tau\_{0\,M\_0}}{m+1} \sum\_{\beta=1}^{N} \left\langle \frac{|\tau\_{\beta\,M\_0}| - \tau\_{c\,M\_0}(\gamma\_{ac\,M\_0})}{\tau\_{0\,M\_0}} \right\rangle^{m+1},\tag{7.72}$$

$$\Phi\_{M\triangle} = \frac{\dot{\gamma}\_{0\,M\triangle} \tau\_{0\,M\triangle}}{m+1} \sum\_{\beta=1}^{N} \left\langle \frac{|\tau\_{\beta\,M\triangle}| - \tau\_{c\,M\triangle} \left(\gamma\_{ac\,M\triangle}\right)}{\tau\_{0\,M\triangle}} \right\rangle^{m+1},\tag{7.73}$$

$$
\Phi\_f = \frac{\dot{c}\_0 f\_0}{m+1} \left\langle \frac{|f| - f\_c}{f\_0} \right\rangle^{m+1} \,. \tag{7.74}
$$

With the above definition of the internal variables a partially relaxed energy *WR*1*,P T* is defined by

$$W\_{R1,PT} = \inf\_{\underline{\mathcal{Z}}\_{R1,PT}} \{ \bar{W} + \lambda\_1 \cdot (\boldsymbol{F}\_{M\_{\triangle}} - \boldsymbol{F}\_A - \boldsymbol{a} \otimes \boldsymbol{N}) \}$$

$$\quad + \lambda\_2 \cdot (\boldsymbol{F}\_{M\_{\triangle}} - \boldsymbol{F}\_{M\_0} - \boldsymbol{b} \otimes \boldsymbol{N})$$

$$\quad + \lambda\_3 \cdot (c\_A \boldsymbol{F}\_A + c\_{M0} \boldsymbol{F}\_{M0} + c\_{M\triangle} \boldsymbol{F}\_{M\triangle} - \bar{\boldsymbol{F}}) \},\tag{7.75}$$

$$\quad - \left(\boldsymbol{F}\_A \cdot \boldsymbol{F}\_A - \bar{\boldsymbol{F}}\_A \cdot \boldsymbol{b} \right) \quad \forall \quad \lambda\_1, \lambda\_2, \lambda\_3 \in \mathbb{Z}\_{\triangle}.\tag{7.76}$$

$$\underline{x}\_{R1,PT} = \left(\boldsymbol{F}\_A, \boldsymbol{F}\_{M\_0}, \boldsymbol{F}\_{M\_\triangle}, \boldsymbol{a}, \boldsymbol{b}, \lambda\_1, \lambda\_2, \lambda\_3\right)^T. \tag{7.76}$$

The total power *P* and the stationarity conditions using the principle of minimum dissipation power are then given by

$$\min\_{\underline{\dot{\alpha}}} \left\{ P \right\} = \min\_{\underline{\dot{\alpha}}} \left\{ \dot{W}\_{R1, PT} + \bar{\Phi} \right\}. \tag{7.77}$$

Applying as outlined in section 6 a Legendre-transformation, evolution equations for the strain-like internal variables can be derived. The set of equations describing the three phase rank 1 laminate are now given

by

$$\{\ddot{\sigma}\}\_{M\_{\triangle}/M\_0} \mathbf{N} = \mathbf{0}, \quad \{\ddot{\sigma}\}\_{M\_{\triangle}/A} \mathbf{N} = \mathbf{0},\tag{7.78}$$

$$\dot{\gamma}\_{ac\,A} = \frac{\partial \bar{\Phi}}{\partial(-c\_A \tau\_{c\,A})} = \sum\_{\beta=1}^N \dot{\gamma}\_{\beta\,A},\tag{7.79}$$

$$\dot{\boldsymbol{F}}\_{p\boldsymbol{A}} = \frac{\partial \bar{\boldsymbol{\Phi}}}{\partial (\boldsymbol{c}\_{A} \boldsymbol{M}\_{e\boldsymbol{A}} \boldsymbol{F}\_{p\boldsymbol{A}}^{\rm T})} = \sum\_{\beta=1}^{N} \gamma\_{\beta\,A} \mathbf{d}\_{\beta\,A} \otimes \boldsymbol{n}\_{\beta\,A} \boldsymbol{F}\_{p\boldsymbol{A}},\tag{7.80}$$

$$\dot{\gamma}\_{ac \, M\_0} = \frac{\partial \bar{\Phi}}{\partial - (c\_{M\_0} \tau\_{c \, M\_0})} = \sum\_{\beta=1}^{N} \dot{\gamma}\_{\beta \, M\_0},\tag{7.81}$$

$$\dot{\boldsymbol{F}}\_{p\,M\_0} = \frac{\partial \bar{\Phi}}{\partial (\boldsymbol{c}\_{M\_0} \boldsymbol{M}\_{\boldsymbol{e}\,M\_0} \boldsymbol{\mathcal{F}}\_{p\,M\_0}^T)} = \sum\_{\beta=1}^N \gamma\_{\beta\,M\_0} \mathbf{d}\_{\beta\,M} \otimes \boldsymbol{n}\_{\beta\,M} \boldsymbol{F}\_{p\,M\_0},\tag{7.82}$$

$$\dot{\gamma}\_{ac \, M\_{\Delta}} = \frac{\partial \bar{\Phi}}{\partial(-c\_{M\_0} \tau\_{c \, M\_{\Delta}})} = \sum\_{\beta=1}^{N} \dot{\gamma}\_{\beta \, M\_{\Delta}},\tag{7.83}$$

$$\dot{\mathbf{F}}\_{p\,M\_{\triangle}} = \frac{\partial \bar{\Phi}}{\partial (c\_{M\_{\triangle}} \mathbf{M}\_{e\,M\_{\triangle}} \mathbf{F}\_{p\,M\_{\triangle}}^{\mathrm{T}})} = \sum\_{\beta=1}^{N} \gamma\_{\beta\,M\_{\triangle}} \mathbf{d}\_{\beta\,M} \otimes \mathbf{n}\_{\beta\,M} \mathbf{F}\_{p\,M\_{\triangle}},\tag{7.84}$$

$$\dot{c}\_{M} = \frac{\partial \bar{\Phi}}{\partial f} = \dot{c}\_{0} \left\langle \frac{|f| - f\_{c}}{f\_{0}} \right\rangle^{m} \text{sgn}(f),\tag{7.85}$$

and

$$F\_A = \dot{\mathbf{F}} - c\_M \mathbf{a} \otimes \mathbf{N} - c\_{M\_0} \mathbf{b} \otimes \mathbf{N},\tag{7.86}$$

$$F\_{M\_0} = \bar{\mathbf{F}} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} + (1 - c\_{M\_0})\mathbf{b} \otimes \mathbf{N},\tag{7.87}$$

$$F\_{M\_{\triangle}} = \bar{F} + (1 - c\_M)\mathbf{a} \otimes \mathbf{N} - c\_{M\_0}\mathbf{b} \otimes \mathbf{N}.\tag{7.88}$$

**Homogenization of the martensite phases at the end of a time interval** As mentioned above the previously transformed martensite *M*<sup>0</sup> and the currently transformed martensite inside the transition zone *M*△ are homogenized at the end of a time interval, yielding the previously

transformed martensite *M*<sup>0</sup> at the beginning of the next time interval. Therefore, at the end of each time interval a homogeneous phase *M* is constructed which yields the same effective behavior as the volume fraction weighted sum of the individual phases *M*<sup>0</sup> and *M*△. The volume fraction of the homogeneous phase *M* and total deformation gradient of the homogeneous phase *M* determined by the laminate kinematics read

$$c\_M = c\_{M\_0} + c\_{M\_\Delta},\tag{7.89}$$

$$F\_M = \frac{1}{c\_M} \left( c\_{M\_0} F\_{M\_0} + c\_{M\_\triangle} F\_{M\_\triangle} \right). \tag{7.90}$$

Apart from the total deformation the homogenization has to account for the plastic deformation and the inherited plastic deformation, which are in general different for the previously transformed martensite and the currently transformed martensite. Additionally, the accumulated plastic slip inside both martensite laminate phases has to be homogenized. A simple volume averaging of these quantities is not possible because of the multiplicative decomposition of the total deformation gradient and the non-linearity of the plasticity related part of the Helmholtz free energy. Since the elastic and plastic properties of *M*<sup>0</sup> and *M*△ are identical, they can be described by a Helmholtz free energy of the same type as the individual phases, see Eq. (7.12). A Helmholtz free energy for the homogeneous phase *M* is then constructed by

$$W\_M = W\_M(\mathbf{F}\_M, \theta, \mathbf{F}\_T, \mathbf{F}\_{ip\ M}, \mathbf{F}\_{p\ M}, \gamma\_{ac\ M}). \tag{7.91}$$

*γac M* can be determined by requiring equivalence of the hardening contributions to the Helmholtz free energy

$$c\_M W\_{p\,M} = c\_{M\_0} W\_{p\,M\_0} + c\_{M\_\triangle} W\_{p\,M\_\triangle}.\tag{7.92}$$

The quantities *Fip M,Fp M* are determined by the requirement that the Helmholtz free energy for the homogeneous phase *M* is minimized subjected to the constraint, that it is identical to the volume fraction weighted sum of the phase energies *W<sup>M</sup>*<sup>0</sup> and *W<sup>M</sup>*△

$$\min\_{\mathbf{F}\_{p\ M}, \mathbf{F}\_{ip\ M}} \left\{ W\_M + \lambda (W\_M - \frac{c\_{M\_0}}{c\_M} W\_{M\_0} - \frac{c\_{M\_\Delta}}{c\_M} W\_{M\_\Delta}) \right\}. \tag{7.93}$$

The minimization conditions and constraint read

$$(1+\lambda)\frac{\partial W\_M}{\partial \mathbf{F}\_{p\,M}} = \mathbf{0},\tag{7.94}$$

$$(1+\lambda)\frac{\partial W\_M}{\partial \mathbf{F}\_{ip\,M}} = \mathbf{0},\tag{7.95}$$

$$
\delta W\_M - \frac{c\_{M\_0}}{c\_M} W\_{M\_0} - \frac{c\_{M\_\triangle}}{c\_M} W\_{M\_\triangle} = 0. \tag{7.96}
$$

# **Chapter 8 Numerical implementation**

#### **8.1 Material point level**

**Integration of the evolution equations** To implement the governing equations for the three phase rank 1 laminate, Eqs. (7.78)- (7.85), the evolution equations (7.79)- (7.85) are discretized in time. The evolution of the martensite volume fraction of the transition zone Eq. (7.85) and the accumulated plastic slip of the three phases austenite *A*, Eq. (7.79), transformed martensite *M*0, Eq. (7.81), and transforming martensite *M*△, Eq. (7.83) are integrated by a backward Euler scheme yielding

$$c^{n+1} - c^n - \triangle t \left(\frac{|f^{n+1}| - f\_c}{f\_0}\right)^m \text{sgn}(f^{n+1}) \dot{c}\_0 = 0,\tag{8.1}$$

$$
\gamma\_{ac}^{n+1} - \gamma\_{ac}^{n}{}\_{A} - \triangle t \dot{\gamma}\_{ac}{}\_{A} = 0,\tag{8.2}
$$

$$
\gamma\_{ac}^{n+1}{}\_{M\_0} - \gamma\_{ac}^{n}{}\_{M\_0} - \triangle t \dot{\gamma}\_{ac}{}\_{M\_0} = 0,\tag{8.3}
$$

$$
\gamma\_{ac}^{n+1}{}\_{M\triangle} - \gamma\_{ac}^{n+1}{}\_{M\triangle} - \triangle t \dot{\gamma}\_{ac}{}\_{M\triangle} = 0,\tag{8.4}
$$

where the accumulated plastic slip of the transforming martensite *MM*△ at time *t <sup>n</sup>* is assumed to be the current accumulated slip of the austenite *γ n ac M*△ = *γ n*+1 *ac A*. To preserve the incompressible nature of the plastic deformation gradients the corresponding evolution equations are integrated using the exponential map (see, e.g., Miehe, 1996). Since the plastic deformation gradient is unsymmetric, a direct calculation of the

exponential function is numerically expensive. An approximation of the tensor exponential can be given in form of a midpoint rule by (Stein and Steinmann, 1995)

$$F\_p^{n+1} = \exp(\triangle tL\_p) F\_p^n \tag{8.5}$$

$$\exp(\triangle tL\_p) \approx (I - \frac{\triangle t}{2}L\_p)^{-1}(I + \frac{\triangle t}{2}L\_p). \tag{8.6}$$

This is equivalent to the Pade approximation of the tensor exponential, which has also been applied and discussed in Baaser (2004). As a result, the evolution equations for the plastic deformation gradients follow as

$$F\_{p\,A}^{n+1} - (I - \frac{\triangle t}{2} L\_{p\,A})^{-1} (I + \frac{\triangle t}{2} L\_{p\,A}) F\_{p\,A}^n = \mathbf{0},\tag{8.7}$$

$$\mathbf{F}\_{p\,M\_0}^{n+1} - (\mathbf{I} - \frac{\triangle t}{2} \mathbf{L}\_{p\,M\_0})^{-1} (\mathbf{I} + \frac{\triangle t}{2} \mathbf{L}\_{p\,M\_0}) \mathbf{F}\_{p\,M\_0}^n = \mathbf{0},\tag{8.8}$$

$$(F\_{p\ M\_{\triangle}}^{n+1} - (I - \frac{\triangle t}{2} \mathbf{L}\_{p\ M\_{\triangle}})^{-1} (I + \frac{\triangle t}{2} \mathbf{L}\_{p\ M\_{\triangle}}) F\_{p\ M\_{\triangle}}^n = \mathbf{0}.\tag{8.9}$$

The remaining equations corresponding to traction continuity across the phase interfaces are minimizing equations and therefore algebraic. They are now listed for the sake of completeness

$$\{\check{\sigma}^{n+1}\}\_{M\_{\triangle}/M\_0} \mathbf{N} = \mathbf{0},\tag{8.10}$$

$$\{\check{\sigma}^{n+1}\}\_{M\triangle/A}N = \mathbf{0}.\tag{8.11}$$

Eqs. (8.1)-(8.4) and (8.7)-(8.11) constitute a nonlinear system of equations in 37 unknowns, which in the following are collected in the vector of local unknowns

$$s = \left(a^{n+1}, b^{n+1}, c^{n+1}, \gamma\_{ac\,A}^{n+1}, \gamma\_{ac\,M\_0}^{n+1}, \gamma\_{ac\,M\_\Delta}^{n+1}, \mathbf{F}\_{p\,A}^{n+1}, \mathbf{F}\_{p\,M\_0}^{n+1}, \mathbf{F}\_{p\,M\_\Delta}^{n+1}\right)^{\mathrm{T}}.\tag{8.12}$$

The unknowns are solved for by a Newton scheme. To guarantee a numerical stable and robust scheme a staggered Newton scheme is used. Thus the Newton scheme is partitioned into two steps 1) and 2).

1) For a given total deformation gradient *F <sup>n</sup>*+1 and temperature *θ n*+1 the crystal plasticity and traction continuity equations are solved for the unknowns *a n*+1 *, b n*+1 *, γ<sup>n</sup>*+1 *ac α ,F n*+1 *p α* , corresponding to a laminate with stationary interfaces. To avoid numerical instabilities related to the power law of the plastic evolution equations, the exponent *p* is first set to 1. The equations (8.2)-(8.11) are then solved to a predefined, coarse tolerance. Afterwards, the exponent *p* is increased until a prescribed tolerance for the residuals of the corresponding equations is exceeded. The system is then solved again until the prescribed tolerance is reached. This is repeated until the desired exponent *p* is reached.

2) Upon convergence of Eqs. (8.2)-(8.11) a second Newton scheme is initiated to solve for all 37 unknowns, in case that the driving force is larger than the critical driving force *fc*. Since changes of the volume fraction can lead to large changes of the phase deformation gradients, a homotopy-scheme is applied for Eq. (8.1). Therefore, the residual of Eq. (8.1) at the end of the first Newton scheme is saved as *g* ∗ . Denoting Eq. (8.1) by *g*(1)

$$g^\* = g(1)(c^n, \mathbf{a}^{n+1,1}, \mathbf{b}^{n+1,1}, \gamma\_{ac\cdot\alpha}^{n+1,1}, F\_{P\alpha}^{n+1,1})\tag{8.13}$$

with the intermediate results for the jump vectors and plasticity related internal variables after the first Newton scheme (·) *n*+1*,*1 . The first equation *g*(1) is then modified to

$$c^{n+1} - c^n - \triangle t \left(\frac{f^{n+1} - f\_c}{f\_0}\right)^m \dot{c}\_0 - \lambda\_H g^\* = 0. \tag{8.14}$$

At the beginning of the second Newton scheme the parameter *λ<sup>H</sup>* is chosen close to 1 and the full system of equations is solved to a prescribed tolerance. In the following *λ<sup>H</sup>* is decreased by a prescribed factor and

the system is solved again until *λ<sup>H</sup>* reaches zero and the prescribed tolerance for the residual of all 37 equations is met.

#### **8.2 Consistent algorithmic tangent**

**Stress integration for finite strains in Abaqus** To solve the global Finite-Element equations, related to the weak form of the linear momentum balance, the Jacobian of the material law at each material or integration point is required. In the case of a finite deformation anaylsis Abaqus uses a Jaumann rate form to integrate the stress equations (see, e.g., Nguyen and Waas, 2016). Therefore, the consistent algorithmic tangent provided by the user has to be derived from the Jaumann stress rate. Using the relation between the Kirchhoff stress tensor *τ* and the second Piola-Kirchhoff stress tensor *S* the rate of the Kirchhoff stress tensor can be expressed as

$$
\dot{\tau} = \dot{F}SF^T + F\dot{S}F^T + F\dot{S}\dot{F}^T.\tag{8.15}
$$

Rearranging the terms and accounting for *L* = *F F*˙ <sup>−</sup><sup>1</sup> yields

$$
\dot{\tau} - L\tau - \dot{\tau}L^T = F \frac{\partial S}{\partial E^{GL}} [\dot{E}^{GL}] F^\Gamma,\tag{8.16}
$$

where the left hand side is the Oldroyd rate of the Kirchhoff stress. The Jaumann rate of the Kirchhoff stress tensor, decomposing the velocity gradient into *L* = *D* + *W*, reads

$$
\dot{\tau} - W\tau + \tau W = F \frac{\partial S}{\partial E^{GL}} [\dot{E}^{GL}] F^{\mathrm{T}} + D\tau + \tau D. \tag{8.17}
$$

Using the definition of the Green-Lagrangian strain tensor one finally obtains

$$
\dot{\tau} - W\tau + \tau W = F \star \frac{\partial S}{\partial E}[D] + D\tau + \tau D. \tag{8.18}
$$

The algorithmic tangent C*algo* which has to be provided to Abaqus is additionally normalized by the volume and thus given by

$$\mathbb{C}\_{algo} = \frac{1}{\det(F)} \left( F \star \frac{\partial S}{\partial E^{GL}} + \frac{1}{2} \left( I \Box \tau + (I \Box \tau)^{\text{T}\_{\text{R}}} + \tau \Box I + (\tau \Box I)^{\text{T}\_{\text{R}}} \right) \right) \tag{8.19}$$

Due to the laminate interface equations, Eqs. (8.1)-(8.3), which are set up in the intermediate configuration and are related to the Piola-Kirchhoff stress tensor, a derivative of the local equations with respect to the symmetric Green-Lagrange *E GL* strain tensor is not possible. Therefore, an approximation of the Jaumann rate is used. The rate of the Kirchhoff stress tensor is thus related to the first Piola-Kirchhoff stress tensor by

$$
\dot{\boldsymbol{\tau}} = \dot{\boldsymbol{\sigma}} \dot{\boldsymbol{F}}^{\mathrm{T}} + \frac{\partial \check{\boldsymbol{\sigma}}}{\partial \boldsymbol{F}} [\dot{\boldsymbol{F}}] \boldsymbol{F}^{\mathrm{T}}, \tag{8.20}
$$

$$=\boldsymbol{\tau}\boldsymbol{L}^{\mathrm{T}} + \frac{\partial \check{\boldsymbol{\sigma}}}{\partial \boldsymbol{F}}[\dot{\boldsymbol{F}}]\boldsymbol{F}^{\mathrm{T}}.\tag{8.21}$$

Rearranging the terms and expanding the left and right side by −*Lτ* yields the Lie derivative which is then used to formulate the Jaumann stress rate using the first Piola-Kirchhoff stress tensor

$$
\dot{\tau} - L\tau - \tau L^{\mathrm{T}} = \frac{\partial \check{\sigma}}{\partial F} [LF] F^{\mathrm{T}} - \tau W + D\tau. \tag{8.22}
$$

Decomposing *L* = *D*+*W* and neglecting on the right side the terms related to the spin tensor *W* an approximation of the algorithmic tangent

is found by

$$\mathbb{C}\_{algo} \approx \frac{1}{\det(F)} \left( \frac{1}{2} F \star^2 \left( \frac{\partial \check{\sigma}}{\partial F} + \frac{\partial \check{\sigma}}{\partial F}^{T\_R} \right) + I \Box \tau \right) \tag{8.23}$$

The partial derivative on the local level can now be expressed by the total derivative of the effective Piola-Kirchhoff stress tensor ¯*σ*ˇ *<sup>n</sup>*+1 with respect to to the deformation gradient *F* = *F*¯ *<sup>n</sup>*+1, provided as input by Abaqus,

$$\mathbf{d}\bar{\boldsymbol{\sigma}}^{n+1} = \sum\_{\alpha} \left( \frac{\partial \check{\sigma}^{n+1}\_{\alpha}}{\partial \bar{\boldsymbol{F}}^{n+1}} \mathbf{d}\bar{\boldsymbol{F}}^{n+1} + \frac{\partial \check{\sigma}\_{\alpha}}{\partial \mathbf{s}^{n+1}} \mathbf{d}\mathbf{s}^{n+1} \right), \tag{8.24}$$

with the vector of unknowns *s* as introduced in Eq. (8.12). Note that the phase stresses do not depend on all of the unknowns. The total derivative of the first Piola-Kirchhoff stress with respect to the effective deformation gradient is then given by

$$\frac{\mathbf{d}\bar{\bar{\sigma}}^{n+1}}{\mathbf{d}\bar{F}^{n+1}} = \sum\_{\alpha} \frac{\partial \check{\sigma}^{n+1}\_{\alpha}}{\partial \bar{F}^{n+1}} + \frac{\partial \check{\sigma}\_{\alpha}}{\partial s^{n+1}} \frac{\mathbf{d}s^{n+1}}{\mathbf{d}\bar{F}^{n+1}}.\tag{8.25}$$

The total derivative d*s <sup>n</sup>*+1*/*d*F*¯ *<sup>n</sup>*+1 can be found by static condensation of the local equations (8.1)-(8.4) and (8.7)-(8.11), which are now collected in the vector *g*, by

$$\frac{\mathbf{d}\mathbf{s}^{n+1}}{\mathbf{d}\bar{F}^{n+1}} = -\left(\frac{\partial \mathbf{g}^{n+1}}{\partial \mathbf{s}^{n+1}}\right)^{-1} \frac{\partial \mathbf{g}^{n+1}}{\partial \bar{F}^{n+1}}.\tag{8.26}$$

# **Chapter 9 Applications**

### **9.1 Elastic laminate theory**

**Convex vs rank 1 convex envelope with respect to first order laminates** A material point subjected to an effective uniaxial deformation gradient at constant temperature is considered. The microstructure of the material comprises an austenite phase and one martensite variant. The material is assumed to be isotropic elastic. The Helmholtz freeenergy is chosen as a Neo-Hookean free energy of the form given in Eq. (7.12) without hardening and thermal contributions. Thus, the chemical energy reduces to a constant. The material parameters for austenite and martensite can be found in Table 9.1, and have been taken from Bartel and Hackl (2009). They represent a Cu-Ni shape-memory alloy. The eigendeformation of the austenite is equal to the identity


**Table 9.1:** Material parameters for shape memory alloy

while that of the martensite variants is identified with

$$\begin{aligned} \boldsymbol{F}\_{M1,T} &= \begin{pmatrix} \gamma\_{M,1} & 0 & 0\\ 0 & \gamma\_{M,2} & 0\\ 0 & 0 & \gamma\_{M,2} \end{pmatrix}, & \boldsymbol{F}\_{M2,T} &= \begin{pmatrix} \gamma\_{M,2} & 0 & 0\\ 0 & \gamma\_{M,1} & 0\\ 0 & 0 & \gamma\_{M,2} \end{pmatrix},\\ \boldsymbol{F}\_{M3,T} &= \begin{pmatrix} \gamma\_{M,2} & 0 & 0\\ 0 & \gamma\_{M,2} & 0\\ 0 & 0 & \gamma\_{M,1} \end{pmatrix}. \end{aligned} \tag{9.1}$$

The Helmholtz free energy for both phases and the effective laminate energy is given by

$$W\_{\alpha} = \frac{\mu\_{\alpha}}{2} \left( \text{tr}(\mathbf{C}\_{\alpha,e}) - 2\text{ln}(J\_{\alpha,e}) - 3 \right) + \frac{\lambda\_{\alpha}}{2} \text{ln}(J\_{\alpha,e})^2 + W\_{\alpha,ch}, \tag{9.2}$$

$$
\bar{W} = \sum\_{\alpha=1}^{2} c\_{\alpha} W\_{\alpha} \tag{9.3}
$$

The convex envelope for a two phase material, see also Eq. (5.12), is defined by

$$W\_{CON} = \inf\_{c\_{\alpha}, \mathbf{F}\_{\alpha}} \left\{ \sum\_{\alpha=1}^{2} c\_{\alpha} W(\mathbf{F}\_{\alpha}) + \lambda (\sum\_{\alpha=1}^{2} c\_{\alpha} - 1) + \lambda\_1 \cdot (\sum\_{\alpha=1}^{2} c\_{\alpha} \mathbf{F}\_{\alpha} - \mathbf{F}) \right\} \tag{9.4}$$

and yields the following minimization conditions and constraint

$$W\_M - W\_A - \langle \check{\sigma} \rangle \cdot \{F\} = 0,\tag{9.5}$$

$$
\check{\sigma}\_M - \check{\sigma}\_A = \mathbf{0},\tag{9.6}
$$

$$\sum\_{\alpha=1}^{2} c\_{\alpha} F\_{\alpha} - \bar{F} = \mathbf{0}.\tag{9.7}$$

As can be seen, it provides a driving force, which is basically of the same form as the sharp interface driving force. However, while compatibility of the phase deformation gradients with the effective deformation gradient is satisfied, the jump of the phase deformation gradients is in general not a rank 1 tensor. Thus, compatibility of the phase deformation gradients is in general not satisfied. The second equation corresponds to equal stresses inside both phases. Thus, it renders the convex envelope equivalent to the Reuss bound of the effective energy and represents an energetically lower bound for the energy of the effective material for a given deformation gradient.

The rank 1 convex envelope with respect to first order laminates is given by Eq. (6.34), Eq. (6.35) and Eq. (6.36). The phase deformation gradients are defined by Eq. (6.5). To investigate the influence of the laminate orientation the rank 1 convex envelope is investigated both for an evolving laminate orientation (R1) and a fixed laminate orientation (R1P). The latter corresponds to a partial relaxation of the effective Helmholtz free energy , since the minimization is not performed with respect to the laminate normal *N*. In the case of an evolving laminate orientation an initial guess for the interface normal *N* is found in each increment by scanning a unit sphere accounting for a constant volume fraction of the phases and traction continuity across the interface, Eq. (6.35). With this initial guess Eq. (6.34)-Eq. (6.36) are solved for the jump vector *a*, the normal *N* and the martensite volume fraction *cM*. In the case of a fixed laminate orientation the normal orientation is chosen such, that it is colinear to the imposed effective uniaxial deformation gradient, which is given by *F*¯ = *γe*<sup>1</sup> ⊗ *e*<sup>1</sup> + *I*, *γ* ∈ [0*,* 0*.*025]. The chosen interface normal in the case for a fixed laminate orientation is thus *N* = *e*<sup>1</sup> ⊗ *e*1. For the following comparison only the first martensite variant is considered, which in terms of a homogeneous material yields the lowest overall Helmholtz free energy.

The relaxed energies, together with the energies for a homogeneous austenite, *WA hom*, and martensite, *WM hom*, material subjected to the same effective deformation gradient, are depicted in Fig. 9.1 a). It can be seen that all relaxed energies constitute an envelope to the homogeneous phase energies. The convex envelope represents the lowest possible energy, while the rank 1 convex envelope for an evolving interface results in a slightly higher effective Helmholtz free energy. The

**Figure 9.1:** Comparison of the convex, rank 1 convex and partially rank 1 relaxed envelope for an effective uniaxial deformation gradient. a) Relaxed and homogeneous phase energies for uniaxial deformation, b) evolution of martensite volume fraction for different relaxations, c) effective stress component coaxial to effective deformation gradient, d) rank 1 laminate

rank 1 laminate with a constant interface orientation shows the highest effective energy. The difference between the rank 1 laminate with fixed interface orientation and evolving interface orientation is in comparison

larger than that of the latter compared with the convex envelope. The connecting points of the corresponding relaxed energies with the energies of the homogeneous phases are such, that for lower relaxed energies it is shifted to a smaller effective deformation gradient for the start of the transformation and a larger effective deformation gradient for the end of the transformation process. This results in a phase transformation over a longer deformation range for more relaxed energies. Or in other words, the less relaxed the effective energy is, the faster the phase transformation proceeds, see also Fig. 9.1 b). The extremal example is a non relaxed energy, for which the phase transformation occurs instantaneously at the point where the homogeneous phase energies intersect. Interestingly the relaxed effective energy for the rank 1 laminate with fixed interface orientation shows the characteristic Maxwell line, which is observed in phase transformation of fluids or gases. The corresponding stress in direction of the effective deformation gradient is consequently constant, Fig. 9.1 c). The reason is, that the chosen interface normal leads to a colinear jump vector *a* = *kN*, with a constant *k*. Thus, the driving force and the traction continuity basically degenerate to scalar equations, since the normal and the jump of the phase deformation gradients has only one non-zero component. This can be compared to the phase transformation of two gases, where the pressure and volume are the only, scalar values. By choosing the volume and volume fractions in a proceeding phase transformation the pressure is fixed.

The stresses in direction of the effective deformation gradient corresponding to the convex and rank 1 convex envelope show a similar characteristic, with lower values for the stress of the convex envelope. The initiation of the transformation process is hereby visible in the transition of the linear stress regime to the nonlinear stress regime, respectively vice versa for the end of the transformation process.

**Rank 1 convex envelope with respect to higher order laminates** Next the rank 1 convex envelope of an austenite-martensite material is compared with respect to first, second and third order laminates for an effective uniaxial deformation gradient *F*¯ = *te*<sup>1</sup> ⊗ *e*<sup>1</sup> + *I*, *t* ∈ [0*,* 0*.*07]. The material parameters are the same as above, except for the transformation stretch tensors, which have been increased to better distinguish the results for the rank 2 and rank 3 laminate.

$$\begin{aligned} \mathbf{F}\_{M1,T} &= \begin{pmatrix} 1.03 & 0 & 0 \\ 0 & 1.03 & 0 \\ 0 & 0 & 0.95 \end{pmatrix}, & \mathbf{F}\_{M2,T} &= \begin{pmatrix} 1.03 & 0 & 0 \\ 0 & 0.95 & 0 \\ 0 & 0 & 1.03 \end{pmatrix}, \\ \mathbf{F}\_{M3,T} &= \begin{pmatrix} 0.95 & 0 & 0 \\ 0 & 1.03 & 0 \\ 0 & 0 & 1.03 \end{pmatrix}. \end{aligned} \tag{9.8}$$

For the first order laminate the first martensite variant is chosen. The second order laminate comprises the first two martensite variants, since they are from an energetic point of view the most favorable with respect to the imposed macroscopic deformation gradient. The third order laminate comprises two second order laminates with variants one, two and one, three respectively, see Fig. 9.2 d). For all considered laminates the laminate orientation is evolving. To identify good candidates for the corresponding initial values a similar procedure as above is applied. In this respect, it turns out that the energetically optimal interface orientation of the martensite sublaminates is only slightly dependent on the applied effective deformation gradient and mainly on the respective transformation deformation gradients<sup>1</sup> .

The results show, that the relaxed energy of the rank 1 laminate, *WR*<sup>1</sup> constitutes an envelope for the homogeneous phases *WA, WM*<sup>1</sup>*, WM*<sup>3</sup>.

<sup>1</sup> To that end, a rank 2 or rank 3 laminate is considered in the case of small strains and identical elastic moduli for all martensite variants. It can now be shown, that the minimizing condition with respect to the corresponding martensite normals is independent of the imposed effective strain and the martensite volume fractions. It only depends on the transformation deformation stretch tensors of each variant. Thus, suitable initial values for finite strains can be found, without the numerically expensive task of scanning a unit sphere for each normal

**Figure 9.2:** Comparison of the rank 1 convex envelope with respect to first (*R*1), second (*R*2) and third order laminates (*R*3) for a uniaxial effective deformation gradient. a) Phase energies, b) effective stress coaxial to macroscopic deformation gradient, c) volume fractions, d) structure of rank 2 and rank 3 laminates

However, the rank 2, *WR*<sup>2</sup>, and rank 3 laminate, *WR*<sup>3</sup>, yield lower energies but do not represent an envelope with respect to the homogeneous phase of one martensite variant. They represent envelopes of homogeneous phases, for which the transformation deformation gradient is a mixture of martensite variants. The homogeneous energy, which the envelope of the rank 2 laminate aims at, can be approximated by a homogeneous phase with an eigendeformation equal to a transformation deformation gradient consisting of those of the martensite variants one

and two with equal volume fractions, compare Fig. 9.3 c). This is due to the fact, that those two variants constitute two energetic equivalent phases with respect to to the imposed effective deformation gradient, see Fig. 9.2a). In the case of the rank 3 laminate all three variants are accounted for, with an unequal distribution of the volume fractions. Due to the different hypothetical homogeneous target energies for each laminate, the transformation rates of austenite into martensite are significantly faster for the higher order laminates, see Fig. 9.2c). Therefore, it can be concluded that higher order laminates do not only lower the effective Helmholtz free energy , but also shorten the 'transformation distance'.

The rank 2 laminate yields over the considered range of the effective deformation gradient an energy, which is always lower than that of the rank 1 laminate. The energy of the rank 3 laminate is during the transformation process lower than that of the rank 2 laminate, see Fig. 9.2a). This is made possible through a finer mixture of phases using the third variant, even though it is energetically unfavorable as a homogeneous phase. After the point, at which all austenite of the rank 3 laminate has transformed into martensite, its effective Helmholtz free energy increases over that of the rank 2 laminate. This marks the macroscopic deformation state, for which a rank 2 laminate becomes energetically preferable to the rank 3 laminate. Hence, the rank 3 laminate would actually degenerate into a rank 2 laminate. This case has not been considered here. After all austenite has transformed the volume fractions are kept constant.

Comparing the stresses in direction of the uniaxial effective deformation gradient, the rank 1 laminate yields the highest followed by the rank 2 laminate and the lowest for the rank 3 laminate. While the rank 1 and 2 laminate show significantly curved stress distributions during the transformation process, the rank 3 laminate nearly reaches a stress plateau shortly after the start of the phase transformation.

**Figure 9.3:** Evolution of the individual volume fraction of each martensite phase for the rank 2 a) and rank 3 laminate b). c) Effective free energy of rank 2 laminate, compared with homogeneous austenite and martensite phases as well as a homogeneous phase containing a mixture of martensite variant 1 and 2. Effective stress state for the rank 1 d), rank 2 e) and rank 3 laminate f).

The stress increases just slightly in the given deformation range. Thus, higher order laminates reduce the slope of the stress curve during the transformation, Fig. 9.2b).

A more detailed investigation of the effective stress state shows that for the given effective deformation gradient the rank 3 laminate reaches a stress state, for which all normal stresses are approximately equal with comparatively small shear stresses, Fig. 9.3 f). The rank 2, Fig. 9.3 e), and rank 1 laminate, Fig. 9.3 d), show a tendency towards a larger deviation of the normal stress components with even smaller shear stresses. This can be explained by the appearance of more martensite variants in higher order laminates, which in combination lead to a volumetric inelastic deformation and a spherically dominated stress state.

Comparing the evolution of the individual volume fractions for the rank 2 laminate, Fig. 9.3 a), shows, that the martensite variants one and two are coexisting throughout the considered deformation range. At the beginning variant one has a slightly larger volume fraction, while at the end both phases have identical volume fractions. The difference at the beginning might be due to the reason that compatibility of the phase deformation gradients prevents equal volume fractions of both variants. The rank 3 laminate shows over the whole range of the transformation process identical volume fractions for the first and second variant, Fig. 9.3 b). Compared to the rank 2 laminate the third martensite variants enables the laminate to construct a mixture of identical volume fractions for the first and second variant, while still remaining compatible. The third variant has a lower volume fraction, which is still comparatively high, bearing in mind that it is macroscopically undesirable with respect to a homogeneous strain energy.

### **9.2 Inelastic laminate theory - lath martensite**

**Material parameters and load scenario** In the following the transformation behavior of the lath martensite model, as introduced in Chapter 6 is investigated at the level of a material point. Because no complete set of temperature dependent material parameters was available for austenite and martensite, the considered steel is more of hypothetical nature. However, they are closely related to those presented in Neumann and Böhlke (2016). Both phases, austenite and martensite are assumed to be isotropic with a temperature dependent Young' Modulus *E*(*θ*) and constant Poisson's ratio *ν*. The Voce hardening law is adopted for both


**Table 9.2:** Material parameters for a hypothetical, low carbon steel

phases, with a temperature dependent initial hardening modulus *h*0(*θ*), final hardening modulus *h*∞(*θ*), initial yield stress *τ*0(*θ*) and final yield stress *τ*∞(*θ*). The rate sensitivity parameter *m* for the corresponding plastic dissipation potentials Φ*<sup>α</sup>* is assumed to be temperature inde-

**Figure 9.4:** Temperature dependent material parameters for austenite and martensite

pendent. The thermal expansion coefficient *α* and the specific heats *c<sup>v</sup>* are assumed to be constant. The difference of the internal energy and entropy constant in the Helmholtz free energy has also been assumed temperature independent Du (2017), with [*e*0] *A/M* = 340 N/mm<sup>2</sup> and [*η*0] *A/M* = 1*.*333 N/(K mm<sup>2</sup> ). The above material parameters constitute the set of reference parameters and are illustrated in Fig. 9.5 a) -c). The chemical energies of martensite and austenite are depicted in Fig. 9.5 d), and show an equilibrium temperature of around 550°C for both phases in the absence of mechanical energy contributions. The

**Figure 9.5:** Evolution of laminate quantities for a constant effective deformation gradient and prescribed temperature path. a) martensite volume fraction, b) accumulated plastic slip, c) elastic, d) plastic and e) chemical free energy, f) absolute difference of the respective energy contributions)

first load scenario for the material point is a fixed effective deformation gradient *F*¯ = *I* and a prescribed cooling rate *θ*(*t*) = 505°C − *t* · 100K/s for *t* ∈ [0*,* 4*.*8]s. The exponent for the dissipation potential, describing the kinetic relation Eq. (6.12) is chosen as *p* = 1 resulting in a linear relation between the driving force and the normal velocity of the interface. The growth rate of lath martensite is given between 10<sup>−</sup><sup>6</sup> and 10<sup>−</sup><sup>1</sup> m/s (Nishiyama, 1978; Villa et al., 2015). Choosing the size of a martensite block as 10 *µ*m the parameter *c*˙<sup>0</sup> = *v* ∗ ⊥0 */L* follows as *c*˙<sup>0</sup> = 1000/s. In Fig. 9.5 a) the martensite volume fraction is depicted versus the imposed temperature and a fixed effective deformation gradient. The martensite evolution can be roughly split into three different stages. In the first stage a slight increase in martensite volume fraction is visible,

followed by a second stage, where the martensite evolution increases rapidly and a third stage, in which a saturation of the martensite evolution can be observed. To investigate this behavior in detail, the sharp interface driving force in the adiabatic case Eq. (3.55) is decomposed in the respective energy contributions. These are the elastic part, the hardening contribution and the chemical part constituting the Helmholtz free energy of each phase

$$f = W\_{e\,M} - W\_{e\,A} + W\_{p\,M} - W\_{p\,A} + W\_{ch\,M} - W\_{ch\,A} - \mathbf{a} \cdot \langle \dot{\sigma} \rangle \mathbf{N} = 0\tag{9.9}$$

In the following the impact of each of them on the phase transformation is discussed. In the first stage martensite forms within a prior untransformed austenite grain. Due to its large eigendeformation the martensite deforms plastically, while the austenite is still able to elastically compensate the deformation due to the forming martensite, see Fig. 9.5 b). This makes the martensite from an energetic point of view unpreferable. The transformation process is quite slow and is only driven by the difference in the chemical free energies Fig. 9.5 e). The plastic deformation of the martensite at the phase front stays

approximately constant after initiation of the transformation process Fig. 9.5 b). The beginning of the second stage is characterized by plastic yielding of the austenite. This penalizes the austenite Helmholtz free energy and leads to an increased martensite transformation rate, which in turn leads to further plastic yielding inside the austenite. In this stage the martensite at the phase front keeps a constant level of martensite intrinsic plasticity and only experiences an increase in plasticity due to inherited dislocations of the prior austenite regions, see Fig. 9.5 b).

The increased accumulated plastic slip inside the martensite at the phase front, due to inherited plasticity, leads to a saturation of the martensite transformation in the third stage. The reason for this can be found in the larger interaction energy of dislocations inside the martensite lattice compared to the austenite lattice, reflected by the higher hardening energy of the martensite. Thus, a considerable amount of energy is needed to accommodate the austenite dislocations inside the martensite lattice. This energy is consequently not available to drive the martensite phase front. The significant increase in the difference of the hardening energy of martensite at the phase front and the austenite at the onset of plastic inheritance can be seen in Fig. 9.5 d). In this respect, the advantage of the three phase rank 1 laminate, compared to a two phase rank 1 laminate becomes obvious. While the three phase rank 1 laminate is able to resolve the different levels of plasticity at the phase front and already transformed martensite regions, a two phase laminate is only able to provide the average plastic quantities inside the laminate. Since the predicted averaged plastic slips and accumulated plastic slip at the phase front and transformed regions show a significant difference, Fig. 9.5 b), the two phase rank 1 laminate would predict a much faster martensite transformation in the third stage, lacking a saturation mechanism.

The influence of the elastic energies, Fig. 9.5 c), f), is small compared to the chemical and plastic energy. Thus, in conclusion the current model predicts a martensite evolution dominated by the difference in plastic

**Figure 9.6:** a) Contributions to the driving force, b) relaxed effective free energy vs homogeneous austenite free energy, c) effective Cauchy stress components for the considered constant effective deformation gradient and prescribed temperature path

and chemical energies of both phases, Fig. 9.5 f). For the considered load scenario, the elastic and plastic energy as well as the energy due to the jump of the deformation gradient and the average stress all act in favor of the austenite Fig. 9.6 a). The phase transformation is in this particular case only driven by the difference in chemical energy. Fig. 9.6 a) shows that the mixture of austenite and martensite relaxes the energy and leads to a significantly reduced effective energy, compared to a homogeneous austenite phase. The effective stress state is mainly hydrostatic, which is due to the imposed effective deformation gradient and the particular transformation deformation

gradient of the martensite. At the end of the transformation process, where the transformation rate slows significantly, a decreasing stress magnitude is observed. This is due to the thermal strains which slightly compensate for the deformations imposed by the martensite transformation gradient inside the laminate.

#### **9.3 Investigation of material parameters**

In the following, the influence of the material parameters on the martensite evolution is investigated for the above load scenario of a constant effective deformation gradient *F*¯ = *I* and a prescribed cooling rate *θ*(*t*) = 505°C − *t* · 100K/s for *t* ∈ [0*,* 4*.*8]s.

**Kinetic relation** The dynamic of the phase front is characterized by the evolution for the martensite volume fraction Eq. (6.12). Throughout the following, the exponent is assigned a constant value *m* = 1. The normalizing constant *f*<sup>0</sup> is chosen as 1. The rate *c*˙<sup>0</sup> is now varied in 4 steps, corresponding to different material intrinsic velocities of the martensite phase front. First the minimizing condition *f* = 0 corresponding to *c*˙<sup>0</sup> = ∞ and in the following *c*˙<sup>0</sup> = 1 , 0*.*1 and 0*.*01 1/s are investigated. The results are depicted in Fig. 9.7, where the inverse relation *η* = 1*/c*˙<sup>0</sup> has been used. Varying *η* between 0 and 1 shows no significant difference on the martensite evolution, Fig. 9.7 a). A slightly delayed martensite evolution is observed for *η* = 10. For *η* = 100 the martensite evolution is visibly slowed down at the beginning, but catches up with decreasing temperatures. The reason for that can be found in the higher levels of plasticity inside the austenite with larger martensite volume fractions. Through inheritance of dislocations the martensite evolution is slowed down and the influence of the phase front velocity reduces significantly. In the following, the plastic evolution of both phases, the phase energies as well as the 11 stress component are only compared for *c*˙ = ∞ and *c*˙ = 0*.*01, since the other

**Figure 9.7:** Comparison of different intrinsic velocities of the martensite phase front *c*˙<sup>0</sup> = 1*/η* and influence on the martensite volume fraction a), accumulated plastic slip b), elastic energy c), plastic energy d) and effective stress for a constant effective deformation gradient and prescribed cooling rate.

values show only minor differences. The accumulated plastic slip inside the austenite is lower for smaller velocities of the phase front in the early stages of the martensite transformation due to a smaller transformation rate. Analogously to the martensite volume fraction the accumulated plastic slip inside the austenite catches up with decreasing temperatures. The accumulated plastic slip inside the martensite at the phase front differs for the two considered velocities of the phase front only through the inherited slip from the austenite, Fig. 9.7 b). With a slower interface velocity the elastic energy of the martensite at the phase front and the austenite is decreased at high temperatures, due to a smaller martensite volume fractions and consequently less deformation inside the laminate. The same holds for the plastic energies, where less accumulated slip at high temperatures equals smaller plastic energies for lower interface velocities. The relaxed effective energy is still smaller for smaller velocities of the phase front, since the chemical energy of the austenite is larger than that of the martensite and thus a lower martensite transformation rate increases the effective laminate energy.

**Hardening parameters austenite** First the initial hardening modulus *hA,*<sup>0</sup> is varied in three steps, see Table 9.3. The impact of the initial hardening modulus is relatively minor for the considered range, see Fig. 9.8 a), b). The difference of the martensite volume fraction compared to the reference parameters is at all temperatures less than 0.1 %. As can be expected a higher hardening modulus increases the austenite hardening energy and therefore slightly accelerates the transformation. The opposite holds for a smaller initial hardening modulus. Second, the initial and final yield stress are varied separately in the same manner as the hardening modulus. The variation of the initial and final yield stress of the Voce hardening law shows a noticeable impact on the martensite evolution. As discussed in the beginning, the onset of plastic deformation inside the austenite marks the beginning of the second transformation stage, where a rapid increase of the martensite volume fraction takes place. Therefore, by increasing the initial yield stress the


**Table 9.3:** Variation of austenite hardening parameters

beginning of the second stage is shifted to lower temperatures, Fig. 9.8 c). However, an increasing yield stress results in a larger hardening energy, Fig. 9.8 d). This causes the transformation to proceed faster in the following, due to an energetically less favorable austenite phase. Again, the accumulated plastic slip inside the austenite evolves analogously to the martensite, Fig. 9.8 d). For lower values of the initial and final yield stress the martensite start temperature is shifted to higher temperatures and the martensite transformation proceeds slower, for the opposite reasons as pointed out above.

**Hardening parameters martensite** In the same way as for the austenite the hardening parameters for the martensite phase are varied, see Table 9.4. The variation of the martensite hardening parameters shows,


**Table 9.4:** Variation of martensite hardening parameters

that a higher hardening modulus yields a lower martensite start temperature and slower evolution rate. This is due to the higher hardening energy, rendering the martensite phase energetically less favorable and thus slowing its evolution. In contrast, a smaller value leads to a higher martensite start temperature and faster overall transformation, see Fig. 9.9 a). Again, the accumulated slip inside the austenite evolves simultaneously to the martensite and is therefore larger for faster martensite transformation rates, Fig. 9.9 b). This increases the

**Figure 9.8:** Variation of the austenite initial hardening parameter *h*0*,M*. Influence on martensite volume fraction a) and austenite accumulated plastic slip b). Variation of the austenite initial and final critical shear stress *τ*0*,A* and *τ*∞*,A*. Influence on martensite volume fraction c), austenite hardening energy and austenite accumulated plastic slip e).

**Figure 9.9:** Variation of the martensite initial hardening parameter *h*0*,M*. Influence on martensite volume fraction a), martensite accumulated plastic slip b), martensite hardening energy c) and austenite hardening energy d)

hardening energy of the martensite, due to more inherited plasticity at lower temperatures, Fig. 9.9 c). However, the austenite plastic energy is also increased, which cancels the saturation effect of the higher martensite plastic energy. Hence, a lower hardening modulus still leads to an overall faster martensite transformation. The variation of the initial and final yield stress of the martensite shows that for smaller values a faster martensite transformation rate is achieved, Fig. 9.10 a). A smaller value in this regard means less martensite hardening energy. Additionally, the accumulated slip inside the austenite increases due to the faster martensite evolution which penalizes the austenite phase, Fig. 9.10 b). Again, this implies that the martensite phase inherits more

**Figure 9.10:** Variation of the martensite initial and final critical shear stress *σ*0*,M* and *σ*∞*,M*. Influence on martensite volume fraction a), martensite accumulated plastic slip b) and martensite hardening energy c)

austenite plasticity and the effect on the martensite evolution is less pronounced, but still significant. Higher values of the initial and final yield stresses lead to a delayed start of the martensite transformation, since the martensite hardening energy increases. This implies that the austenite exhibits less accumulated slip at higher temperatures and the transformation rate is slower.

**Inheritance of plasticity from austenite to martensite** As already discussed the proposed model predicts a crucial impact of plastic inheritance on the martensite evolution. In the following, it is distinguished between the plastic slips or plastic deformation gradient of a phase *Fp α*, which represents the lattice distortion caused by dislocations and the

accumulated plastic slip which characterizes the self hardening contribution due to mobile dislocations. It is now assumed that martensite forms in a plastically distorted austenite lattice and thus austenite plastic slips are fully incorporated into the martensite. This implies *Fihp M* = *Fp A*. Since not all dislocations of the prior austenite region might be situated on martensite slip planes and are thus not mobile, the inheritance of the accumulated plastic slip is gouverend by a constant parameter *ξih*.

Different levels of plastic inheritance with respect to the accumulated plastic slip are now considered. Therefore, a very simplistic inheritance approach, in a way that the inherited accumulated plastic slip is multiplied by a constant factor *ξih*, is investigated. *ξih* is now varied between 0, 0*.*5 and 1. As can be seen, the inheritance factor has a significant impact on the phase transformation, 9.11a) . A small inheritance factor accelerates drastically the transformation process, leading even to an unstable, respectively instantaneous transformation behavior shortly before all austenite has been transformed. Since in this case martensite at the phase front does not inherit plasticity from austenite, the accumulated plastic slip at the martensite phase front stays approximately constant during the transformation, see Fig. 9.11 c). However, the austenite level of plasticity steadily increases and penalizes the austenite Helmholtz free energy, Fig. 9.11 b), d), while the plastic contribution to the martensite Helmholtz free energy does not increase by much, Fig. 9.11 e). For larger values of the inheritance factor the accumulated plastic slip at the martensite phase front increases, which penalizes the martensite Helmholtz free energy and slows the transformation considerably. The relaxed effective energies are increased for a higher level of plastic inheritance.

**Influence of the deformation state on the transformation** Apart from the temperature, the applied effective deformation gradient plays an important role regarding the martensite transformation. The influence of different effective deformation gradients on the rank 1 laminate are investigated. Hereby, the applied effective deformation gradient is kept constant. The considered effective deformations gradients are

$$
\bar{F}\_0 = I, \quad \bar{F}\_1 = I + \gamma \mathbf{e}\_1 \otimes \mathbf{e}\_1, \quad \bar{F}\_2 = I + \gamma \mathbf{e}\_1 \otimes \mathbf{e}\_2
$$

$$
\bar{F}\_3 = I(1 + \gamma), \quad \bar{F}\_4 = -\bar{F}\_3, \quad \gamma = 0.02. \tag{9.10}
$$

As can be seen in Fig. 9.12 a) the effective deformation gradient significantly changes the transformation initiation and behavior, compared to the macroscopically undeformed case. Applying a uniaxial and spherical effective deformation gradient, case (2) and (3), shifts the martensite start temperatures to higher temperatures. The reason is that these effective deformation gradients partially relax the martensite transformation gradient given by the Bain strains. This facilitates the martensite transformation and thus leads to higher martensite start temperatures, because a smaller difference in chemical energies is necessary to drive the phase front. A weak point of the current model is seen, when applying an effective shear deformation gradient. This, in contrast to experimental observations, shifts the martensite initiation to cooler temperatures. The reason herefore can also be found in the Bain transformation strains, which do not include any shear component. However, with lower temperatures the transformation curve for the effective shear deformation gradient approaches the one for the effectively undeformed laminate. This is due to the additional plastic deformation contributions to the martensite deformation which comprise shear deformations. The application of a compressive spherical effective deformation gradient leads to a significantly delayed start of the martensite transformation, for the opposite reasons as above. This compares quite well to experimental results, where an applied pressure partially suppresses martensite transformation. The effective compressive deformation gradient leads to a similar pressure dominated stress state.

**Figure 9.11:** Comparison of different inheritance factors *ξIH* with respect to the martensite volume fraction a), austenite and martensite accumulated plastic slip b), c), austenite and martensite hardening energies d), e) and effective free energy f) for a constant effective deformation gradient and prescribed cooling rate

**Figure 9.12:** Influence of a constant effective deformation gradient on martensite volume fraction a) and effective stress for deformation scenarios *F*¯ <sup>0</sup> a), *F*¯ <sup>1</sup> b), *F*¯ <sup>2</sup> c), *F*¯ <sup>3</sup> d), *F*¯ <sup>4</sup> e) for a prescribed cooling rate.

#### 9 Applications

**Figure 9.13:** a) Ferrite (blue) austenite/martensite (red) laminate, b) macroscopic shear deformation in 23 and 12 plane at the end of transformation

## **9.4 Application of the inelastic laminate theory in FEM**

**A ferrite austenite laminate - elastoplastic ferrite** As a first step to simulate the martensite evolution in a dual phase steel, a microstructure comprising an austenite and a ferrite phase is modeled as a simple laminate, see Fig. 9.13 a). The laminate is discretized with C3D8 elements and periodic boundary conditions are applied to the faces and edges of the structure. To simulate a macroscopic stress free microstructure one node is fixed with all translational degree of freedoms inside the volume. This, in contrast to the material point computations should lead to much more realistic stress fields inside the austenite and martensite phase as well as for the effective stress. A homogeneous temperature field *θ*(*t*) = 505 °C−*t* · 100 K/s, *t* ∈ [0*,* 4*.*8]s is imposed on the structure, which leads to the formation of martensite inside the austenite regions. The martensite variant, which forms inside the austenite regions is the same for all elements, such that the laminate has two effective deformation gradients each for the ferrite and austenite, respectively martensite region. The material parameters correspond to the reference parameters, see Table 9.2. The ferrite is modeled elastoplastic, with a bcc crystal plasticity model. The material parameters for ferrite are listed in Table 9.5.



**Table 9.5:** Material parameters for ferrite

**Results** Comparing the transformation curve, Fig. 9.14 a), for the austenite ferrite laminate with the material point computations of the rank 1 laminate shows a significantly faster martensite evolution. This is due to the less constrained structure with a lower overall stress state. Consequently the, in the considered cases, unfavorable contribution to the driving force *a* · h*σ*ˇi*N* is reduced. Furthermore, an earlier martensite start temperature at around 501◦C can be observed. The lower stress level is also reflected by lower elastic energies, which still show the austenite as the energetically preferable phase, Fig. 9.14 c). The level of accumulated plastic slips in both phases is also smaller compared to the investigations of the material point rank 1 laminate, Fig. 9.14 b). Due to similar initial yield stresses, plastic deformation inside the ferrite matrix starts at the same time as the austenite starts to deform plastically. However, the level of plasticity inside the ferrite is much smaller compared to the austenite, due to the higher final yield stress. Again, the relaxed effective energy of the laminate at a material point shows the significant reduction of energy compared to a homogeneous austenite elastic phase, Fig. 9.14 e). The effective stress state, Fig. 9.15 a), shows much smaller stress components than the rank 1 laminate with a fixed effective deformation gradient. Since the Bain variant 1 was used for this simulation, stresses in 33-direction develop.

**Figure 9.14:** Comparison of elastic and plastic-ferrite austenite laminate. Evolution of martensite volume fraction a), accumulated plastic slips b), *a* · h*σ*ˇi*N* c). Plastic ferrite laminate, elastic energy d), hardening energy e) and effective free energies f).

**Figure 9.15:** Effective and laminate stresses for an elastoplastic ferrite and a prescribed cooling rate. a) Effective stress *σ*¯, b) austenite laminate stress, c) stress in the transformed martensite region, d) stress inside the martensite transformation zone.

This is due to the large compressive transformation deformation in this direction, which is constrained by the ferrite matrix. All three phases, the austenite, the already transformed martensite and the transforming martensite exhibit shear stresses in 12 direction, see Fig. 9.15 b), c), d), which on a macroscopic level cancel each other due to opposite signs. These shear stresses are caused by the martensite transformation deformation in combination with the plastic deformation inside the laminate phases, which are constrained by the periodic boundary

conditions of the macroscopic laminate. The effective deformation of the laminate shows the typical shear deformation of martensite inside a matrix Fig. 9.13 b).

**A ferrite austenite laminate - elastic ferrite** To identify the effect of the ferrite matrix on the transformation behavior, the ferrite phase is now modeled elastically. The cooling profile, boundary conditions and material parameters are the same as above. An elastic ferrite phase slows the evolution of the martensite phase significantly, especially at lower temperatures, as can be seen in Fig. 9.14 a). This is due to the higher stresses and consequently larger contribution *a* · h*σ*ˇi*N* to the driving force, see Fig. 9.14 f). Furthermore, the stress level inside the laminate is slightly increased which results in slightly larger accumulated plastic slips. This slows the martensite evolution even further in the stages, where plastic inheritance plays an important role. **A synthetic ferrite austenite microstructure** As a last example the martensite transformation of austenite grains in a ferritic matrix, similar to that of a dual phase steel is considered. To that end a 15<sup>3</sup> (representative) volume element with a periodic microstructure and random grain orientation is created with the free available software DREAM3D (Groeber and Jackson, 2014). The small numbers of elements is on the one hand due to the high computational cost of the proposed model, including three crystal plasticity models inside the rank 1 laminate and another crystal plasticity model to simulate plasticity inside the ferrite matrix. On the other hand choosing reasonable values for the kinetic relation of the martensite phase front yield especially for a larger number of elements severe convergence issues. Thus, the size of the volume element is not really representative. However, the current investigation just aims at presenting some basic features of the proposed phase transformation model and the impact of the martensitic phase transformation on the surrounding ferrite of a dual phase steel. The synthetic microstructure and corresponding grain orientation is depicted in Fig. 9.16 b). Two configurations have been investigated,

the first containing 35 % austenite prior to the transformation and the second containing 10 % austenite. Both RVEs are subjected to the same temperature profile and the material parameters are in both cases those listed in Table 9.2, except for the difference in the energy constant, which has been chosen to approximately match the martensite start temperature of the experiments presented in da Silva et al. (2014).

**Variant selection** To determine which martensite variants evolves at a material point a simple energy based approach is chosen. Until initiation of the phase transformation at a material point, it is checked which Bain variant yields the lowest elastic energy with respect to the effective deformation gradient at that point. This is also equivalent to a maximization of the driving force with respect to the Bain variants. The chosen approach yields for the considered microstructure a distribution of Bain variants as depicted in Fig. 9.16 a). The Bain variants are distributed with approximately equal amounts, which is also often observed experimentally.

**Martensite transformation** The effective martensite volume fraction for both configurations is shown in Fig. 9.16 e). Again the three stages of the martensite evolution are present, where for the first stage a martensite start temperature of around 458 °C can be observed. A more detailed investigation of the martensite volume fraction in each grain, Fig. 9.17 a) and Fig. 9.18 a) , shows that at the beginning and in the early stages of the second transformation stage the individual martensite volume fraction is quite inhomogeneously distributed. A comparatively low martensite volume fraction is observed in regions where different Bain variants are adjacent. The reason for this are different eigendeformations of the Bain variants. As a consequence different effective deformation fields evolve in the region of each Bain variant, which are unfavorable for the respective neighboring Bain variants. With ongoing transformation the relative difference between this regions and regions, where a homogeneous distribution of Bain variants is present becomes smaller. This can be explained by the domi-

**Figure 9.16:** Synthetic microstructure with 35 % a) and % 10 b) austenite, distribution of Bain variants for 35 % c), and 10 % d) prior austenite. Comparison of martensite evolution for both levels auf initial austenite e) and comparison with experiments of da Silva et al. (2014) f).

nance of plastic effects with advanced martensite volume fractions, which outweigh the influence of the elastic deformation fields.

Comparing the martensite transformation of the RVE with 35 % initial austenite and 10 % initial austenite shows, that the martensite transformation proceeds faster with less initial austenite. This is due to the lower overall stress state for less martensite, see also Fig. 9.17 b), e) and Fig. 9.18 b), e).

The overall characteristic of the martensite evolution with respect to the imposed temperature - again - shows qualitatively the experimentally observed course. Comparing it with experiments presented in da Silva et al. (2014), see Fig. 9.16 e), f), which cover the martensite transformation in a a low carbon steel, shows that the first transformation stage is much more pronounced in the observed experiment. The simulations predict a much smaller temperature range for the first stage. This leads to a higher martensite volume fraction for higher temperatures. The initial increase of the martensite transformation in the second stage is quite similar in the experiment and the simulation. The third stage or saturation of the martensite evolution sets in earlier for the simulated martensite transformation and is also more pronounced. This leads to a slightly increased overall temperature range, necessary to transform all austenite into martensite, compared to the experiment. The effective martensite evolution requires a temperature difference of about 300 K, which is larger than that of the experiment.

All in all the second and the third stage of the transformation (rapid growth due to austenite plasticity and saturation due to inheritance of austenite dislocations) match qualitatively comparable experiments relatively well. However, the saturation effect of the martensite evolution, due to plastic inheritance is too strong, which slows the transformation. This might be due to the quite simplistic inheritance model. A more advanced model which models the plastic inheritance as a function of the martensite volume fraction has been proposed by Kundin et al. (2015).

**Figure 9.17:** a) Martensite volume fraction, b) effective v. Mises stress laminate, c) accumulated plastic slip austenite , d) accumulated plastic slip martensite, e) v. Mises stress ferrite, f) accumulated plastic slip ferrite

9.4 Application of the inelastic laminate theory in FEM

**Figure 9.18:** a) Martensite volume fraction, b) effective v. Mises stress laminate, c) accumulated plastic slip austenite , d) accumulated plastic slip martensite, e) v. Mises stress ferrite, f) accumulated plastic slip ferrite

The first stage of the martensite transformation is captured rather poorly. The reasons for this can for example be found in the fact that no interface energy is considered in the proposed model. Especially in the early stages, where the martensite blocks are small the effect of the interface energy on the transformation is higher compared to stages, where more martensite has already formed. This would lead to slower transformation rates, as observed in the experiment by da Silva et al. (2014).

The v. Mises stress inside the martensite reaches values up to 1700 MPa, Fig. 9.17 b) for both configurations , until the end of the transformation. The accumulated plastic slip shows again in both cases a maximum of about 4*.*7 · 10<sup>−</sup><sup>1</sup> , Fig. 9.17 d). This values seem to be a little bit too high, but qualitatively represent the high dislocation density, which is observed in lath martensite upon cooling. The remaining austenite exhibits significant higher plastic deformation with an accumulated plastic slip of about 2, Fig. 9.17 d).

The values of the respective austenite and martensite quantities do not differ significantly with respect to different levels of initial austenite inside the RVE. This is due to the comparatively similar size of the austenite grains and the fact that the laminate behavior is mostly dominated by the properties of the austenite and martensite phase. In comparison, the stresses inside the ferrite matrix are much lower as well as the accumulated plastic slip for a smaller overall volume fraction of initial austenite, Fig. 9.17 e) and Fig. 9.17 f). The accumulated plastic slip inside the ferrite matrix still shows relatively high values around the martensite regions, where the large transformation deformation of the martensite causes large stress fields. For both, the stress and the accumulated plastic slip inside the ferrite matrix smaller levels of accumulated plastic slip are observed for less initial austenite. This yields lower overall stresses in the ferrite matrix due to less martensite. For a martensite fraction of around 35 % plastic deformation is caused in all ferrite crystals. For the considered RVE with 10 % martensite not

**Figure 9.19:** Uniaxial tension in 33 direction of the RVE. a) Effective Cauchy stress in 33-direction, b) accumulated plastic slip inside the ferrite matrix for 10 % initial austenite, c) accumulated plastic slip inside the ferrite matrix for 35 % initial austenite

all ferrite crystals are deformed plastically. Upon cooling to 25 ◦C a uniaxial tensile test of the RVE is performed, Fig. 9.19. While the 10 % austenite RVE exhibits an elastic deformation for small deformations, the RVE with 35 % austenite deforms immediately plastically since all ferrite and martensite, respectively martensite grains have already experienced plasticity during cooling. However, due to its higher overall martensite volume fraction the 35 % austenite RVE exhibits higher stresses during the subsequent uniaxial tensile test. For both RVEs the ferrite matrix shows a significant increase in plastic deformation due to its lower yield stress compared to the martensite grains.

# **Chapter 10 Conclusions and outlook**

A finite strain thermomechanical martensitic phase transformation model based on a sharp interface theory has been proposed. The model is based on the application of the sharp interface driving force to a laminate substructure. Plastic effects are captured by a crystal plasticity model in all laminate phases. Additionally, a kinetic relation using a simple thermodynamical consistent approach is introduced. The relation of sharp interface models and models based on rank 1 energy relaxation has been discussed. To extend the relaxation to nonequilibrium states of the material the principle of minimum dissipative power has been used to formulate all equations in a minimization framework. The numerical implementation on the material point level as well as in the FE program Abaqus as a UMAT has been presented. An analytical expression for the algorithmic tangent has been derived. Numerical examples for different load and temperature scenarios have been discussed.

In comparison to the phenomenological Koistinen Marburger model the s-shaped trend of the martensite evolution is captured. Furthermore, the local stress and plastic fields inside the laminate are resolved. The sharp interface driving force also accounts for the actual stress state with respect to the phase transformation. The model is able to capture the characteristic of the martensite evolution qualitatively well. In comparison to the experiments of da Silva et al. (2014) the first

transformation stage and the saturation of the martensite evolution are under- respectively overpredicted. The deviation in the first transformation stage might be due to non-comparable material parameters or the lack of an austenite martensite interface energy. The martensite evolution in the third stage can be improved by a physically based model for the inheritance of austenite dislocation into martensite. An exact comparison to the experiments was not possible, since the applied material parameters only represent a rough estimate. However, the results are expected to improve quantitatively for a more realistic set of material parameters.

Since this model represents the first step towards a micromechanical sharp interface model for the thermomechanical phase transformation of martensite a lot of interesting points are still left for future works. The model can be extended to include a full thermomechanical coupling. This allows to investigate the influence of e.g. heat generated by plastic dissipation on the phase transformation. Additionally, the interaction of the martensite phase front with austenite dislocations could be investigated based on a more physical approach. Moreover, the non-Schmidt behaviour of bcc materials as well as temperature dependent activity of bcc slip systems might be incorporated for the martensite and ferrite to reach more realistic results concerning the plastic effects. Additionally, to improve global convergence on the one hand and on the other hand include interface energy a gradient extension of the free energy with respect to the gradient of the martensite volume fraction could be introduced.

# **Chapter 11**

## **Appendix**

#### **11.1 Geometrical compatibility at interfaces**

Parametrization of a singular surface with respect to the reference configuration, using particular points *ξ*1*, ξ*<sup>2</sup> on the singular surface gives

$$\mathbf{X}^\* = \mathbf{y}^\xi \left(\xi\_1, \xi\_2, t\right). \tag{11.1}$$

Choosing the parametrization in a way that the tangential vectors *t*<sup>1</sup> and *t*<sup>2</sup>

$$t\_1 = \frac{\partial X^\*}{\partial \xi\_1}, \quad t\_2 = \frac{\partial X^\*}{\partial \xi\_2}, \tag{11.2}$$

are linear independent one can show (Abeyaratne and Knowles, 2006) that the normal velocity of the singular surface is independent of the chosen parametrization, *ξ*<sup>1</sup> and *ξ*2. The parametrization of the singular surface in the current configuration is defined by the mapping function *χ* ∗

$$x^\* = \chi^\*\left(\chi^\xi\left(\xi\_1, \xi\_2, t\right), t\right). \tag{11.3}$$

Accounting for the continuity of the displacement field, the gradient of *x* <sup>∗</sup> with respect to *ξ*<sup>1</sup> and *ξ*<sup>2</sup> must be identical approaching the interface

in tangent direction from the <sup>+</sup> and <sup>−</sup> side

$$\frac{\partial x^{\*,+}}{\partial \xi\_{\alpha}} = \frac{\partial x^{\*,-}}{\partial \xi\_{\alpha}} = F^{+}t\_{\alpha} = F^{-}t\_{\alpha}, \quad \alpha = 1,2,\tag{11.4}$$

which yields, independent of the chosen parametrization,

$$\{F\}t\_{\alpha} = 0.\tag{11.5}$$

This result implies that the set *t<sup>α</sup>* is a 2-dimensional null space of [*F*], while the deformation gradient itself is defined in a 3-dimensional vector space. Consequently, the jump of the deformation gradient is of rank one (Abeyaratne and Knowles, 2006)

$$\{F\} = a \otimes \mathbf{N}^\* \tag{11.6}$$

where *N*<sup>∗</sup> is perpendicular to *t<sup>α</sup>* and represents the normal of the singular surface. Repeating the same steps for the time derivative gives

$$\frac{\partial \mathbf{x}^{\*,+}}{\partial t} - \frac{\partial \mathbf{x}^{\*,-}}{\partial t} = \mathbf{F}^{+} \frac{\partial \chi^{\xi}}{\partial t} + \mathbf{v}^{+} - \mathbf{F}^{-} \frac{\partial \chi^{\xi}}{\partial t} - \mathbf{v}^{-} = \mathbf{0}.\tag{11.7}$$

Accounting for the continuity of the displacement and the fact that the null space of the jump of the deformation gradient is a set of tangential vectors only the normal contribution to the velocity of the singular surface is relevant and Eq. (11.7) gives

$$\{v\} = -v\_{\perp}^\* \{F\} N^\* \tag{11.8}$$

where the normal velocity of the singular surface *v* ∗ <sup>⊥</sup> is given by

$$
v\_{\perp}^{\*} = \boldsymbol{v}^{\*} \cdot \mathbf{N}^{\*} = \frac{\partial \chi^{\xi}}{\partial t} \cdot \mathbf{N}^{\*}.\tag{11.9}$$

#### **11.2 Convexity and material frame indifference**

Assume a convex function *W*(*F*). The principle of material frame indifference states that

$$W(QF) = W(F). \tag{11.10}$$

Using the following convexity condition for *F*<sup>1</sup> = *F,F*<sup>2</sup> = *QF*, see e.g., Ciarlet (1988), it holds that

$$W(F) + \frac{\partial W(F)}{\partial F} \cdot (QF - F) \le W(QF). \tag{11.11}$$

Applying the principle of material frame indifference and performing some algebraic manipulations one finds

$$\frac{\partial W(F)}{\partial F} F^\top \cdot (Q - I) = \tau \cdot (Q - I) \le 0 \tag{11.12}$$

with the Kirchhoff stress tensor *τ* . Since *τ* ∈ *Sym* can be decomposed into a diagonal part *D* with eigenvalues *τ<sup>i</sup>* and an orthogonal tensor *Q<sup>τ</sup>* ∈ *Orth*<sup>+</sup>. Eq. (11.12) then yields

$$Q\_\tau^\mathrm{T} DQ\_\tau \cdot (Q - I) = D \cdot (R - I) \qquad R = Q\_\tau Q Q\_\tau^\mathrm{T} \in Orb^+. \tag{11.13}$$

Specific choices of *R* as *R*<sup>1</sup> = diag(1*,* −1*,* −1), *R*<sup>2</sup> = diag(−1*,* 1*,* −1) and *R*<sup>3</sup> = diag(−1*,* −1*,* 1) show that

$$
\tau\_1 + \tau\_2 \ge 0 \quad \tau\_1 + \tau\_3 \ge 0 \quad \tau\_2 + \tau\_3 \ge 0 \tag{11.14}
$$

for arbitrary *F* ∈ *Inv*<sup>+</sup>. This is in general clearly not true considering for example the case of a spherical, compressive deformation state.

## **11.3 Mapping of differential area and volume elements for the three phase rank 1 laminate**

**Mapping of volume elements for a two phase rank 1 laminate** The determinant of the sum of a second order tensor *A* and two Rank-1 connected tensors *a* and *b* is given by (Silhavy, 1997, p. 11)

$$\det(\mathbf{A} + \mathbf{a} \otimes \mathbf{b}) = \det(\mathbf{A})(1 + \mathbf{b} \cdot \mathbf{A}^{-1}\mathbf{a}).\tag{11.15}$$

The change of volume with respect to the effective macroscopic deformation of the laminate, accounting for the geometrical compatibility condition (3.17) and identity (11.15), is given by

$$\det(\bar{F}) = \det(c\_1 F\_1 + c\_2 F\_2) = \det(F\_2 - c\_1 \mathbf{a} \otimes \mathbf{N})$$

$$= \det(F\_2)(1 - c\_1 \mathbf{N} \cdot F\_2^{-1} \mathbf{a}).\tag{11.16}$$

Expanding and adding *c*<sup>1</sup> det(*F*2) − *c*<sup>1</sup> det(*F*2) = 0 to the right side results in

$$\det(\bar{F}) = \det(F\_2) - c\_1 \det(F\_2) \mathcal{N} \cdot F\_2^{-1} a + c\_1 \det(F\_2) - c\_1 \det(F\_2)$$

$$= c\_1 \det(F\_2) (1 - \mathcal{N} \cdot F\_2^{-1} a) + c\_2 \det(F\_2). \tag{11.17}$$

Using again identity (11.15) the volume change with respect to the laminate phase 1 yields

$$\det(F\_1) = \det(F\_2 - a \otimes \mathcal{N}) = \det(F\_2)(1 - \mathcal{N} \cdot F\_2^{-1}a) \tag{11.18}$$

and consequently, by plugging Eq. (11.18) into (11.17) one finds

$$\det(\bar{F}) = c\_1 \det(F\_1) + c\_2 \det(F\_2). \tag{11.19}$$

**Mapping of volume elements for a three phase rank 1 laminate** The change of volume with respect to the effective macroscopic deformation of the laminate, accounting for the geometrical compatibility condition and identity (11.15), is given by

$$\det(\bar{F}) = \det(c\_A \mathcal{F}\_A + c\_{M\_0} \mathcal{F}\_{M\_0} + c\_{M\_\triangle} \mathcal{F}\_{M\_\triangle})$$

$$= \det(\mathcal{F}\_A + c\_M(\mathcal{F}\_M - \mathcal{F}\_A)).\tag{11.20}$$

Accounting for the following identities

$$\begin{split} \mathbf{F}\_{M} &= \frac{1}{c\_{M}} (c\_{M\_{0}} \mathbf{F}\_{M} + c\_{M\_{\Delta}} \mathbf{F}\_{M\_{\Delta}}) \\ &= \bar{\mathbf{F}} + (1 - c\_{M}) \mathbf{a} \otimes \mathbf{N} - c\_{M\_{0}} (1 - \frac{1}{c\_{M}}) \mathbf{b} \otimes \mathbf{N}, \end{split} \tag{11.21}$$

$$F\_M - F\_A = (\mathbf{a} + \frac{c\_{M\_0}}{c\_M}\mathbf{b}) \otimes \mathbf{N} := \mathbf{d} \otimes \mathbf{N},\tag{11.22}$$

$$\det(F\_M) = \det(F\_M - F\_A + F\_A) = \det(F\_A)(1 + \mathcal{N} \cdot F\_A^{-1}d),\tag{11.23}$$

Eq. (11.20) yields

$$\det(\bar{F}) = \det(F\_A + c\_M d \otimes \mathcal{N}) = \det(F\_A)(1 + c\_M \mathcal{N} \cdot F\_A^{-1} d)$$

$$= (1 - c\_M) \det(F\_A) + c\_M (\det(F\_A) + \mathcal{N} \cdot F\_A^{-1} d)$$

$$= c\_A \det(F\_A) + c\_M \det(F\_M). \tag{11.24}$$

Applying a similar procedure accounting for

$$\det(F\_{M\_0}) = \det(F\_{M\_0} - F\_{M\_\triangle} + F\_{M\_\triangle})$$

$$= \det(\mathbf{b} \otimes \mathbf{N} + F\_{M\_\triangle})$$

$$= \det(F\_{M\_\triangle})(1 + \mathbf{N} \cdot F\_{M\_\triangle}^{-1}\mathbf{b})\tag{11.25}$$

*c<sup>M</sup>* det(*F <sup>M</sup>*) can be represented by

$$c\_M \det(\boldsymbol{F}\_M) = \det(\boldsymbol{c}\_M \boldsymbol{F}\_M) = \det(\boldsymbol{c}\_{M\_0} \boldsymbol{F}\_{M\_0} + (\boldsymbol{c}\_M - \boldsymbol{c}\_{M\_0}) \boldsymbol{F}\_{M\_\Delta})$$

$$= \det(\boldsymbol{c}\_M \boldsymbol{F}\_{M\_\Delta} + \boldsymbol{c}\_{M\_0} \boldsymbol{b} \otimes \boldsymbol{N}) = \boldsymbol{c}\_M \det(\boldsymbol{F}\_{M\_\Delta} + \frac{\boldsymbol{c}\_{M\_0}}{\boldsymbol{c}\_M} \boldsymbol{b} \otimes \boldsymbol{N})$$

$$= \boldsymbol{c}\_M \det(\boldsymbol{F}\_{M\_\Delta}) (1 + \frac{\boldsymbol{c}\_{M\_0}}{\boldsymbol{c}\_M} \boldsymbol{N} \boldsymbol{F}\_{M\_\Delta}^{-1} \boldsymbol{b}) + (\boldsymbol{c}\_M - \boldsymbol{c}\_{M\_0}) \det(\boldsymbol{F}\_{M\_\Delta})$$

$$= \boldsymbol{c}\_{M\_\Delta} \det(\boldsymbol{F}\_{M\_\Delta}) + \boldsymbol{c}\_{M\_0} \det(\boldsymbol{F}\_{M\_0}). \tag{11.26}$$

Thus, a correct mapping of the volume elements is obtained

$$\det(\bar{F}) = c\_A \det(F\_A) + c\_{M\_0} \det(F\_{M\_0}) + c\_{M\_\Delta} \det(F\_{M\_\Delta}).\tag{11.27}$$

**Mapping of area elements for a two phase rank 1 laminate** A similar calculation can be performed for the cofactor of *F*¯ , representing the change of an area element. The cofactor of a second order tensor *F* and a rank 1 connection *a* ⊗ *N* can be expressed, accounting for the definition of the cofactor cof(*A*) = det(*A*)*A* −T , identity Eq. (11.15) and the Sherman-Morrison-Woodbury formula<sup>1</sup>

$$\begin{split} \text{cof}(F + a \otimes N) &= \det(F + a \otimes N)(F + a \otimes N)^{\text{T}} \\ &= \det(F + a \otimes N)F^{\text{T}} \left( I - \frac{N \otimes F^{-1}a}{1 + N \cdot F^{-1}a} \right) \\ &= \det(F)F^{-T}(1 + N \cdot F^{-1}a) \left( I - \frac{N \otimes F^{-1}a}{1 + N \cdot F^{-1}a} \right) \\ &= \text{cof}(F) \left( (1 + N \cdot F^{-1}a)I - N \otimes F^{-1}a \right) . \end{split} \tag{11.29}$$

<sup>1</sup> Given an invertible second order tensor *A* and two rank-1 connected first order tensors *a* and *b* the inverse (*A* + *a* ⊗ *b*)−<sup>1</sup> exists if 1 + *b* · *A*−1*a* 6= 0 and the inverse is given by

$$(A + a \otimes b)^{-1} = A^{-1} \left( I - \frac{a \otimes A^{\cdot \mathbf{T}} b}{1 + b \cdot A^{-1} a} \right). \tag{11.28}$$

With Eq. (11.29) one finds that

$$\begin{split} \text{cof}(F\_1) &= \text{cof}(F\_2 - a \otimes N) \\ &= \text{cof}(F\_2) + \text{cof}(F\_2)(N \otimes F\_2^{-1}a - (N \cdot F\_2^{-1}a)I) . \end{split} \tag{11.30}$$

Considering now the macroscopic change of an area element cof(*F*¯ ), accounting for Eqs. (11.29) and (11.30) and in the same way as above adding *c*1cof(*F*2) − *c*1cof(*F*2) = **0**

$$\text{cof}(\bar{F}) = \text{cof}(c\_1 \mathbf{F}\_1 + c\_2 \mathbf{F}\_2) = \text{cof}(\mathbf{F}\_2 - c\_1 \mathbf{a} \otimes \mathbf{N})$$

$$= \text{cof}(\mathbf{F}\_2) \left( (1 - c\_1 \mathbf{N} \cdot \mathbf{F}\_2^{-1} \mathbf{a})I + c\_1 \mathbf{N} \otimes \mathbf{F}\_2^{-1} \mathbf{a} \right)$$

$$= (1 + c\_1 - c\_1) \text{cof}(\mathbf{F}\_2) + c\_1 \text{cof}(\mathbf{F}\_2) (\mathbf{N} \otimes \mathbf{F}\_2^{-1} \mathbf{a} - (\mathbf{N} \cdot \mathbf{F}\_2^{-1} \mathbf{a})I)$$

$$= c\_1 \text{cof}(\mathbf{F}\_1) + c\_2 \text{cof}(\mathbf{F}\_2). \tag{11.31}$$

This, as in the case of the volume elements, again gives a correct mapping of the area elements from the microscale to the macroscale.

**Mapping of area elements for a three phase rank 1 laminate** In a similar way as for the determinant the cofactor of the effective deformation yields

$$\text{cof}(\bar{F}) = \text{cof}(F\_A + c\_M(F\_M - F\_A)) = \text{cof}(F\_A + d \otimes N)$$

$$= \text{cof}(F\_A) \left( (1 + N F\_A^{-1} d)I - N \otimes F\_A^{-1} d \right). \tag{11.32}$$

Using the identity

$$\text{coef}(F\_M) = \text{coef}(F\_A) \left( (1 + \mathcal{N} F\_A^{-1} d)I - \mathcal{N} \otimes F\_A^{-1} d \right) \tag{11.33}$$

Eq. (11.32) results in

$$\text{cof}(\bar{F}) = (1 - c\_M)\text{cof}(F\_A) + c\_M \text{cof}(F\_M). \tag{11.34}$$

Furthermore,

$$\begin{split} c\_{M}\text{cof}(\mathbf{F}\_{M}) &= c\_{M}\text{cof}(\mathbf{F}\_{M\_{\triangle}} + \frac{c\_{M\_{0}}}{c\_{M}}(\mathbf{F}\_{M\_{0}} - \mathbf{F}\_{M\_{\triangle}})) \\ &= c\_{M}\text{cof}(\mathbf{F}\_{M\_{\triangle}} + \frac{c\_{M\_{0}}}{c\_{M}}(\mathbf{b} \otimes \mathbf{N})) \\ &= c\_{M}\text{cof}(\mathbf{F}\_{M\_{\triangle}}) \left( (1 + \frac{c\_{M\_{0}}}{c\_{M}}\mathbf{N} \cdot \mathbf{F}\_{M\_{\triangle}}^{-1}\mathbf{b})\mathbf{I} - \frac{c\_{M\_{0}}}{c\_{M}}\mathbf{N} \otimes \mathbf{F}\_{M\_{\triangle}}^{-1}\mathbf{b} \right) \\ &+ (c\_{M\_{0}} - c\_{M\_{0}})\text{cof}(\mathbf{F}\_{M\_{\triangle}}) \\ &= c\_{M\_{\triangle}}\text{cof}(\mathbf{F}\_{M\_{\triangle}}) + c\_{M\_{0}}\text{cof}(\mathbf{F}\_{M\_{0}}). \end{split} \tag{11.35}$$

This together with Eq. (11.34) shows the correct mapping of the area elements for the three phase rank 1 laminate

$$\text{cof}(\bar{F}) = (1 - c\_M)\text{cof}(F\_A) + c\_{M\_0}\text{cof}(F\_{M\_0}) + c\_{M\_\Delta}\text{cof}(F\_{M\_\Delta}).\tag{11.36}$$

## **Bibliography**


Krawietz, A., 1986. Materialtheorie. Springer-Verlag, Berlin Heidelberg.


## **List of figures**






## **List of tables**


#### **Schriftenreihe Kontinuumsmechanik im Maschinenbau Karlsruher Institut für Technologie (KIT) (ISSN 2192-693X)**


Die Bände sind unter www.ksp.kit.edu als PDF frei verfügbar oder als Druckausgabe bestellbar.


Die Bände sind unter www.ksp.kit.edu als PDF frei verfügbar oder als Druckausgabe bestellbar.


Dual phase steels are of significant importance for structural applications as they combine high strength and ductility. This is due to their two phase microstructure comprising ferrite and martensite. Martensite forms under rapid cooling of prior austenitic grains, accompanied by a change of the lattice symmetry. This is characterized by large directional deformations which induces dislocations in both, the austenitic and martensitic phase. On a continuum scale, the transformation of austenite into martensite can be described by the sharp interface theory which describes the motion of a singular surface separating two solids with different elastic properties and eigenstrains. In this work, a transformation model based on the sharp interface theory is developed. It is set in the finite strain context, and crystal plasticity effects of the participating phases as well as the dynamic of the singular surface are accounted for. Furthermore, a simple model to capture the inheritance of austenite dislocations by the forming martensite phase is included. Numerical applications of the model, using the Abaqus user subroutine UMAT, show its ability to qualitatively capture the experimentally observed thermal transformation behavior of low carbon austenite into martensite in dual phase steels.

J. Ruck

**Modeling martensite transformation based on a sharp interface**

ISSN 2192-693X ISBN 978-3-7315-1072-7

Gedruckt auf FSC-zertifiziertem Papier