In this work, a phase-field model for recrystallization is developed based on the conservation laws. There has been no attempt to develop a phase-field model of recrystallization based on the conservation laws, even though various phase-field simulation models to reproduce the recrystallization phenomenon have been proposed. However, it is unclear what conservation laws are required for such a model. In the previous paper, toward solving this problem, we developed conservation laws of mass, momentum, angular momentum, and energy and a law of entropy at the lattice scale for the process of recrystallization. In this paper, first, two continuous variables, i.e., the order parameter and crystal orientation, are introduced into the balance equation of mass for a single phase and that of angular momentum for the lattice, respectively. Next, the fluxes of the order parameter and crystal orientation are derived from the law of entropy by the use of rational thermodynamics. Moreover, the diffusion coefficient and mass source are modeled to derive the evolution equations, i.e., phase-field equations of the order parameter and crystal orientation. Finally, for the phase-field equation of the crystal orientation, neglecting the conservative part and integrating the equation with respect to time under the first-order approximation, a phase-field model that is used for stable calculations is developed. This work aims to develop a phase-field theory on the basis of the change in crystal lattice during recrystallization. This paper gives a physical background to the methodological phase-field approach in the case of recrystallization.

## I. INTRODUCTION

The grain size and texture of a metal affect its mechanical and material properties and, thus, have attracted considerable attention. Recrystallization is a process by which texture is created, and many studies on recrystallization have been reported.^{1–9} In the recrystallization process, the recrystallization area and matrix area switch with each other (Fig. 1). At the same time, owing to the major advances in computational ability, recrystallization simulations^{10–13} have been carried out using the phase-field models that are continuum models.^{14–16}

There are many well-known phase-field models that enable us to analyze material dynamics.^{17} The complicated morphology of dendrites can be reproduced by a phase-field model with an anisotropic diffusion term.^{14,18–23} The martensitic transformation can also be expressed by a phase-field model that considers elastic energy based on micromechanics.^{24–27} There are phase-field models which reproduce the phase transition in a binary alloy,^{28,29} as well as eutectic and peritectic growth and a magnetic phase.^{30–34} The phase-field crystal model is useful for expressing the periodic structure of atoms in terms of phase-field parameters.^{35,36} In addition, phenomena during solidification and grain growth and around dislocations have been studied using phase-field models^{37–42}

These phase-field models are usually derived through the energy principle^{43} rather than being based on conservation laws in which the physical quantities are defined clearly. While energy functionals usually have equivalent balance equations, the balance laws of the phase-field models have scarcely been discussed. There are some important formulations of phase-field models based on mixture theory,^{44} and models of the phase transformation phenomenon are modeled with an interface within the thermodynamical balance laws.^{45–48} Referring to those literature studies, the physical aspect of the order parameter of the phase-field model can be discussed based on mixture theory.^{49,50} Specifically, it is still unclear what conservation law is equivalent to the governing equation of the crystal orientation of the Kobayashi–Warren–Carter (KWC)-type phase-field model.^{51–54} Here, the KWC-type phase-field model has two order parameters. One is the order parameter and the other is the crystal orientation. However, it is quite important to understand what conservation law the phase-field model is based on.

In the phase-field crystal model,^{35,36} there is no parameter related to crystal orientation. It only has the order parameter indicating the probability of the existence of atoms, and its energy functional includes higher-order gradients in order to achieve the profile of the order parameter.

If models are based on the energy principle, we cannot introduce material properties through fluxes. Although there are some reports in which materials with phase transformation are thermodynamically discussed,^{55} they do not treat the continuum evolution of crystal orientation, which is observed in static recrystallization. Considering the above, formulating conservation laws for recrystallization and redeveloping the phase-field model based on conservation laws are meaningful.

Generally, the basic conservation laws of the continuum are not introduced assuming that the material is homogeneous; thus, nm-scale-phenomena, such as atoms or molecules, or mm-scale-phenomena, such as crystal grain behavior, are not considered. Modeling of recrystallization in which microscale crystal lattice scale behavior affects macroscopic phenomena requires making the behaviors coarse in the *μ*m or nm scale. The framework for such a continuum has been developed as the generalized continuum mechanics was developed.^{56} Well-known material models include strain gradient theory,^{57} n-order theory,^{58} micromorphic theory,^{59–61} micropolar theory,^{62–64} director theory,^{64} and multipolar theory.^{65,66} Specifically, the micropolar theory has been used for various materials so far.^{67–72} However, these models have a problem that the symmetry of the Cauchy stress disappears and asymmetric stress appears, which has never been observed in experiments. On the basis of these theories, a lot of studies attempt to bridge the discrete quantities to the continuum ones^{73,74} or to bridge microscale quantities with microstructures to the macroscale ones^{75,76} for various phenomena. As same as the micropolar theory, mixture theory has also been discussed for various phenomena.^{77–81} Additionally, the aforementioned theories have not been extended to mixture forms.

On the other hand, the balance equations for a continuum are derived assuming that a material including an internal structure can be seen as a continuum; thus, the microscopic physical quantities need not be considered. Microscopic quantities should be considered in the case that the behavior of atoms affects macroscopic phenomena such as recrystallization. Therefore, generalized continuum mechanics considering the microscopic structure in continuum theory has been formulated by Cosserat and others.

The authors previously modeled the behavior of materials during recrystallization by considering generalized continuum mechanics.^{82} On the basis of this model, summing and averaging the discrete conservation laws all over the Representative Volume Elements (RVEs) making up a mesoscale domain, the RVE-averaged conservation laws of mass, momentum, angular momentum, and energy were formulated in which, using the characteristic quantities at the bulk and lattice scales, the conservation law of angular momentum was divided into bulk and lattice parts. By converging the RVEs to a material point of the continuum model for the averaged conservation laws for phase *p*, the continuum conservation laws of mass, momentum, angular momentum, and energy for phase *p* were formulated. The obtained global forms of conservation laws were localized in the usual manner. Furthermore, by summing the localized conservation laws all over the phases and averaging them, the conservation laws of a mixture were derived. In addition, based on the conservation laws, the second law of thermodynamics was formulated. However, an order parameter that expresses the phases was not introduced and the constitutive equations of fluxes were not derived. The governing equations that numerically reproduce recrystallization must be developed from the basic equation.

To the above problems, in this study, continuous phase-field functions, i.e., the order parameter and crystal orientation, are introduced. The constitutive equations of fluxes are derived from the second law of thermodynamics developed in the previous paper through a thermodynamic discussion. Furthermore, by modeling the basic equations, phase-field equations are developed, with which stable calculations are possible. In addition, the result of simple simulations using the present model is shown. Additionally, the body force and radiant heat are not considered as conditions, unlike in the previous work.

Originally, the KWC-type model is a methodological model for practical applications such as recrystallization. Traditionally, such a methodological model is developed by the energy principle. Hence, the meaning of one of the order parameters of the KWC-type phase-field model, i.e., crystal orientation, has not been discussed. In order to reconsider what the order parameter is, it is important to develop the model from the viewpoint of a balance law and thermodynamic discussion. If we reformulate the phase-field model based on the balance law and thermodynamic discussion, these parameters should be driven from the four conservation laws and we can clarify the meaning of the parameters. Actually, one of the order parameters of the KWC-type phase-field model, i.e., crystal orientation, can be explained as the change in crystal lattice during recrystallization. The emphasis of this paper is on giving a physical background to the methodological phase-field approach in the case of recrystallization.

## II. BALANCE EQUATIONS

During the static recrystallization process, the atomic force is in unstable equilibrium owing to the thermal activation of atoms and stored dislocations on the boundary between the recrystallized phase and the matrix, as previously shown. Hence, crystal lattices rotate and undergo a transition to another stable state, as shown in Fig. 2. In this study, a crystal lattice is modeled as a rigid bar, i.e., a lattice element (LE), with the equivalent atomic mass shown in Fig. 3. The balance equations are as follows:

Equations (16) and (23) are the local forms of the variables, and there is no expression for the global form in this manuscript because it has already been formulated in a previous report.^{82} Configuration is assumed to be the current one.

Equations (1) and (2) were derived from the conservation laws of mass for a single phase and angular momentum for the lattice, respectively, in the previous work. Here, *ρ*_{p}, *v*_{p}, and $b\u2323p$ are the mass density, velocity, and mass source of phase *p*, and *ρ*, **s**, and **M** are the mass density, angular momentum per unit mass, and couple stress of the mixture, respectively. On the other hand, the energy equation can be written as

where *ε*, ** q**,

**, and**

*T***are the energy per unit mass, heat flux, Cauchy stress, and deformation rate of the mixture, respectively, and**

*D***≡ grad**

*Λ***. The quantity of phase**

*θ**p*is expressed as

with its averaged quantity *A* and relative quantity $A\u2322p$. In addition, $\xc0p$ is the material time derivative, which is defined as

Using Eqs. (4) and (5), Eq. (1) can be rewritten in the form of the usual material time derivative as

Furthermore, the usual continuous equation is derived as follows from the mass conservation for a mixture:

### A. Balance equation of the order parameter

A microscopic volume element *dv* can be expressed as the sum of the microscopic volume elements of the recrystallized phase and matrix as

Setting the volume fraction of the recrystallized phase as *ϕ*, the microvolume elements of the recrystallized phase and matrix are written as *dv*_{r} = *ϕdv* and *dv*_{m} = (1 − *ϕ*)*dv*, respectively. Moreover, considering the constant mass density all over the phases, a microscopic mass element can be calculated as

From the above, the following relations are obtained:

Defining the relative mass flux for phase *p* as

and changing the index to *r*, the following equation is obtained:

Using Eq. (10) so as to introduce the order parameter into Eq. (12), the next equation is derived as follows:

is obtained. Writing *ρ* = *const*. and defining the order flux $j\u2322r$ and mass source $b\u2323r$ as

the balance equation of the order parameter is obtained as

### B. Balance equation of crystal orientation

In the previous paper, the relations between the angular momentum per unit mass *s*_{p}^{(k)}, the crystal orientation *θ*_{p}^{(k)}, and the angular velocity *ω*_{p}^{(k)} were formulated as

where *J*^{(k)} is the inertia moment tensor of a lattice element. In the same manner as in the previous work, by averaging Eq. (17), the corresponding relations for a mixture are obtained as

where ** J** is the average inertia moment tensor, which was expressed as follows in the previous paper:

where $I\u2261\xi (k)\u2297\xi (k)m$ and $\Omega $ denotes the averaged spin angular velocity tensor of the lattice. Substituting Eq. (18)_{1} into Eq. (2) and dividing both sides by *ρ*, the following equation is obtained:

where the orientation flux is defined as

Substituting Eq. (19) into Eq. (20) and applying the previously obtained relations $J\u0307\omega =I\Omega \omega \u2212\Omega I\omega =\u2212\Omega I\omega $, we obtain

Equation (23) is the balance equation of the crystal orientation.

## III. CONSTITUTIVE EQUATIONS

### A. Clausius–Duhem inequality

The derivation procedure of the KWC-type phase-field model through a thermodynamic discussion is the same as that of the classical Clausius–Duhem relation, but the consideration of the rotation of the crystal lattice is different.

The second law of thermodynamics was formulated as follows in the previous work:

where *η* is the entropy density and ** s** is the local entropy flux considering the entropy flux due to the phase transition, which is expressed as

Writing *ρ* = *const*. and ** v** = −

*v*_{m}in Eq. (11), we obtain

where *μ* is the chemical potential, which is defined as

Substituting Eqs. (25) and (28) into inequality (24) and performing the divergence operation, the following inequality is derived:

where ** P** =

*ρ*

**≡**

*p**ρ*grad

*μ*. Substituting Eqs. (3) and (16) into inequality (30) and multiplying

*θ*

_{K}/

*ρ*to both sides, the Clausius–Duhem inequality is derived as

where *φ* is the dissipation fraction per unit mass. In addition, Eq. (21) is incorporated into inequality (31). On the other hand, Helmholtz’s free energy per unit mass *ψ* is written as

Taking the material time derivative of both sides of Eq. (32), we obtain

where $\theta \u0307K=0$ is considered owing to the assumption of an isothermal field. Therefore, by applying the Legendre transformation with Helmholtz’s free energy to inequality (31), the following inequality is derived:

### B. Conservative part of the constitutive equation

Dividing the orientation flux into conservative and dissipative parts using the relation $I\theta =I\theta (c)+I\theta (d)$, inequality (34) can be rewritten as

From inequality (35), the arguments of *ψ* can be chosen as

In addition, on the basis of the set of axioms, the arguments of the constitutive quantities can be chosen as

From Eq. (36), expanding $\psi \u0307$ using the chain rule, we obtain

The sixth, seventh, and eighth terms of the middle part of inequality (39) are linear for each rate. Considering the possible arguments of *ψ*, the conditions in the case that inequality (39) is true are given as

Therefore, the arguments of *ψ* can be reduced to

Note that the case of strong nonequilibrium is neglected and that local balance is assumed in which the state quantity at the material point is dependent on the thermodynamic parameters. The first and second terms of the middle part in inequality (39) are also linear in each velocity; thus, the following equations are derived:

From Eqs. (41) and (42), the arguments of the conservative part of the orientation flux and chemical potential are reduced to

where it is assumed that the chemical potential is a state quantity. In this study, *ψ* is to be expressed as a polynomial of the scalar invariants that are constituted by its arguments. Here, the invariants of *ϕ* and $\Lambda $ are given as

where the scalar invariant of ∇*ϕ* is included so as to consider the effect of the gradient of the order parameter. On the other hand, the vector potential of the angular velocity does not exist. Therefore, $div\u2009\theta \u0307\u22600$ can be assumed, whereas $div\u2009\omega =div\u2009\theta \u0307=0$ in a usual solenoidal field. Because $div\u2009\theta \u0307\u22600$, div ** θ** ≠ 0 can be assumed; thus,

*I*

_{θ1}= tr

**= tr(grad**

*Λ***) = div**

*θ***is included in the arguments of**

*θ**ψ*in this theory. On the other hand, the invariants whose order is higher than three are not considered so as not to generate high-order nonlinear terms. From the above, the fraction of free energy is given as

Calculating Eq. (46) assuming the independence of the arguments, the explicit constitutive equation of $I\theta (c)$ can be written as

where the coefficients are defined as

Note that *a*_{θc0}, *a*_{θc}, and $a\theta c\u2032$ in Eq. (48) are the functions of *ϕ* and ** Λ**.

### C. Dissipative part of the constitutive equation

Considering the quasi-compound process, *φ* can be divided into pure dissipation functions relating to fluxes *φ*_{J} and the mass source *φ*_{b} as

Furthermore, from Eq. (49), we obtain

Considering Eq. (43) and expanding ** p** = grad

*μ*using the chain rule, inequality (51) can be expressed as

where the region in which ** Λ** is fixed is nearly equal to the region where

**=**

*Λ**O*. The region where grad

**=**

*θ**O*is the region inside phases, i.e.,

*ϕ*= 0 or

*ϕ*= 1. Therefore,

*ψ*is constant and the third term on the right-hand side of Eq. (52) takes a value of zero. From the above, inequality (51) is reduced to

where *c*_{π} = *∂μ*/*∂ϕ*. Moreover, applying the normality rule to the curved surface of the dissipation function, the constitutive equations of the fluxes are obtained as

where *v* is constant. From inequality (51), the arguments of *φ*_{J} are reduced to

where ** p** is expanded with the arguments in Eq. (43)

_{2}. Applying Eq. (56) to Eq. (55), the arguments of the constitutive quantities are reduced to the same as those of

*φ*

_{J},

Moreover, representing *φ*_{J} as a polynomial of its scalar invariant, the invariants $\Lambda \u0307$ and ∇** θ** are given as

In the same way as for the conservative part of the dissipation function, neglecting the invariants whose order is higher than three, the dissipation function takes the following form:

Simplifying Eq. (61) further, the following equation is obtained:

where *a*_{ϕ} is defined as

where the constants in Eq. (65) are defined as

Note that *a*_{θd0}, *a*_{θd}, and $a\theta d\u2032$ in Eq. (66) are the functions of *ϕ*, ** Λ**, and $\Lambda \u0307$, respectively.

## IV. PHASE-FIELD EQUATIONS

### A. Phase-field equation of the order parameter

Assuming that the diffusion coefficient *a*_{ϕ} (m^{2}/s) is constant and *b*_{ϕ} is a function of *ϕ*, i.e., *b*_{ϕ} = *b*_{ϕ}(*ϕ*), Eq. (67) is rewritten as

where the term in the parentheses of the second term on the right-hand side of Eq. (68) denotes the argument. Next, *b*_{ϕ}(*ϕ*) is formulated as the mass source term due to the phase transition. First, *b*_{ϕ}(*ϕ*) is represented as the sum of the mass sources due to the autocatalytic reaction *b*_{ca} and grain boundary migration *b*_{gb},

The equation of the autocatalytic reaction can be written as

where *R* and *M* denote the recrystallized phase and matrix, respectively. From Eq. (70), *b*_{gb} can be represented as [see Figs. 4(a) and 4(b)]

where *k*_{a} is a coefficient and *τ* (s) is a representative interval of the autocatalytic reaction. Moreover, the autocatalytic reaction during recrystallization is due to the energy barrier and the stored dislocation energy. Therefore, *k*_{a} can be represented as the sum of the effects of both energies as

where *k*_{ag} is a coefficient that represents the effect of the stored dislocation energy *E*_{s} and *k*_{ab} is a coefficient that represents the effect of the energy barrier *W*_{b}. Furthermore, *E*_{s} is the driving force of nucleus growth, i.e., the matrix region undergoes a transition to the recrystallized phase owing to this driving force. According to the above, *k*_{ag} can be modeled as

where *a*_{ag} (m^{3}/J) is a constant. In addition, assuming that the energy barrier exists in the region between the matrix and recrystallized phase, i.e., *ϕ* = 1/2 and *ϕ* is to transit to the recrystallized phase when *ϕ* ≥ 1/2 and the matrix when *ϕ* ≥ 1/2, *k*_{ag} is modeled as

where *a*_{ab} (m^{3}/J) is a constant. Substituting Eqs. (72)–(74) into Eq. (71), *b*_{ca} is expressed as

Here, Fig. 5 shows the illustration of the effect of *E*_{s} and *W*_{b}. On the other hand, grain boundary migration without crystal lattice rotation is supposed as the nucleus shrinks to the curvature center of the nucleus [see Fig. 4(c)]. Therefore, the reaction equation is given as

From Eq. (76), the velocity equation can be expressed as

where *k*_{gb} is a coefficient. Moreover, the grain boundary energy is introduced for the modeling of *k*_{gb}. On the basis of the Read–Shockley concept,^{83} *E*_{gb} is dependent on the difference in the crystal orientation of grains. In a three-dimensional field, *E*_{gb} is dependent on the invariants that are conjugate with the crystal orientation gradient grad ** θ**. The invariants whose order is lower than two are written as

where ** Λ** = grad

**. In the invariants in Eq. (78), div**

*θ***does not form energy by itself; thus, the second-order quantities such as**

*θ*are adopted, where, considering that *E*_{gb} = 0 when *ϕ* = 0, i.e., a matrix, *E*_{gb} can be written as a polynomial of the invariants of Eq. (79) as

Here, *E*_{r} (J/m^{3}) is the reference grain boundary energy and $Egb*$ is a nondimensional quantity expressed as

where *ε*_{1} (m^{2}), *s*_{2} (m), *ε*_{2} (m^{2}), *s*_{3} (m), and *ε*_{3} (m^{2}) are constants. In Eq. (81), in accordance with the general energy form, *ϕ*^{2}/2 is employed. The constant *k*_{gb} in Eq. (77) must be synchronized with the temporal change in grain boundary migration when grad ** θ** is fixed, i.e., only due to the phase change. From Eq. (80), (

*∂E*

_{gb}/

*∂t*)

_{$\Lambda $}is obtained as

where, from Eq. (82), *k*_{gb} can be expressed as

Substituting Eqs. (75) and (84) into Eq. (69), the reaction term due to the mass source is written as

Substituting Eq. (85) into Eq. (68) and defining the temporal diffusion *α* (m^{2}) and mobility *M*_{ϕ} (1/s) as *α* = *a*_{ϕ}*τ* and *M*_{ϕ} ≡ 1/*τ*, respectively, the equation of the order parameter is obtained as

Here, we can see the same form as the phase-field equation in Eq. (86).

### B. Phase-field equation of crystal orientation

The orientation flux is expressed as the sum of the conservative and dissipative parts. Therefore, using Eqs. (47) and (65), the orientation flux **I**_{θ} is written as

Considering $\Lambda $ = grad ** θ**, Eq. (88) is also expressed as

From Eq. (89), it should be noted that the governing equation of the crystal orientation based on the conservation laws is distinctly different from that of the conventional KWC-type phase-field model based on the energy principle. The former takes the form of a sum of parabolic and hyperbolic functions, whereas the latter is expressed as a parabolic equation. Therefore, the proposed equation has the advantage that it can avoid the physical contradictions of a parabolic partial differential equation such as an infinite propagation speed of a material. However, Eq. (89) includes a second-order temporal partial differential term and the equation is given as a sum of parabolic and hyperbolic functions, making it difficult to derive a stable calculation. Therefore, we next use Eq. (89) to derive a model with which a relatively stable simulation can be carried out. Neglecting the conservative part of Eq. (89), only the dissipative part remains,

Next, both sides of Eq. (90) are integrated with respect to time. First, we consider the following equation:

Integrating both sides of Eq. (91) with respect to time, we obtain

where *c*(**x**) = 0 is assumed as the initial condition. Consider Eq. (92). Approximating the average velocity as $\u0227\theta d\u2248\delta a\theta d/\delta t$, the second term on the right-hand side of Eq. (92) is calculated as

where *δt* is an infinitesimal temporal interval. Moreover, according to Fig. 6, Eq. (93) can be expressed as the area of a trapezoid as

In addition, because grad *θ*_{t} = grad ** θ**−

*δ*(grad

**), the following equation is obtained:**

*θ*where small quantities whose order is higher than two are neglected. In addition, $\u0101\theta d=\u0101\theta d(\varphi ,\u2009\Lambda ,\u2207\varphi ,\u2009\Lambda \u0307)$ is assumed because $a\theta d=a\theta d(\varphi ,\Lambda ,\u2207\varphi ,\Lambda \u0307)$. Selecting the arguments *ϕ* and ||**Λ**||, $\u0101\theta d=\u0101\theta d(\varphi ,\u2009\Lambda )$ is obtained. In addition, $\u0101\theta d=a\theta d\u2212\delta a\theta d=a\theta dt+\delta a\theta d\u2212\delta a\theta d=a\theta dt$ is equal to *a*_{θd} at a known time. From the above, integrating both sides of Eq. (90) with respect to time, we obtain

Next, the diffusion coefficient of the orientation flux is modeled. Because Eq. (99) governs the change in orientation when grains impinge on each other, *a*_{θd0}, *a*_{θd}, and $a\theta d\u2032$ must be synchronized with the temporal change in *E*_{gb} when *ϕ* is fixed. Therefore, (*∂E*_{gb}/*∂t*)_{ϕ} can be expressed as

where $(grad\u2009\theta )\u2248grad\u2009\theta \u0307$ and $(div\u2009\theta )\u2248div\u2009\theta \u0307$ are assumed because of the extremely small time interval. In addition, ${grad\u2009\theta \u2009\u22c5\u2009(grad\u2009\theta )T}=2grad\u2009\theta \u2009\u22c5\u2009(grad\u2009\theta \u0307)T$ is assumed. Comparing Eq. (99) with Eq. (100), *a*_{θd0}, *a*_{θd}, and $a\theta d\u2032$ can be modeled as

where *a*_{d0} (m^{2}/s) is a constant. Substituting Eq. (101) into Eq. (99), the three-dimensional phase-field equation of the crystal orientation is obtained as

where *M*_{θ} ≡ *J*^{−1}*a*_{d0} (1/s) is the mobility tensor expressing the development of the crystal orientation. Moreover, considering that crystal lattice rotation occurs at the boundary region, the inertia moment per unit mass ** J** can be expressed as a function of the order parameter

**=**

*J***(**

*J**ϕ*) (m

^{2}). Therefore,

**can also be written as**

*M***(**

*M**ϕ*). In addition, Eq. (102) can be regarded as the KWC-type-phase-field equation of the crystal orientation extended to three dimensions.

### C. Two-dimensional phase-field equations

The three-dimensional quantities of the crystal orientation ** θ** and inertia moment tensor

**are given as**

*J*Now, reducing the problem to the two-dimensional problem in the *xy* plane, Eq. (103) is represented as

where ** θ** and

*J*^{−1}each have one component allowing them to be reduced to be written as the quantities

*θ*

_{z}≡

*θ*and

*J*

_{zi}≡

*J*

_{θi}, respectively. Considering that

*θ*

_{i}=

*θ*

_{3}≡

*θ*is the only component in the two-dimensional problem, the following calculation becomes possible:

Equation (107) is the conventional KWC-type phase-field equation of the order parameter. Equation (107) is represented with an index as

Since there is no rotation outside of the plane, the second and third terms on the right-hand side of Eq. (109) can be deleted. As a result, Eq. (109) is reduced to

where *M*_{θ} ≡ *M*_{θ33}, *s* ≡ *s*_{2}, and *ε* ≡ *ε*_{2}. Equation (110) coincides with the governing equation of the crystal orientation of the KWC-type phase-field model.

We would like to propose the use of the partial difference method for the final equations (107) and (110). However, in the case that Eq. (90) is solved instead of Eq. (110), the constrained interpolation profile scheme, which is one of the numerical schemes used to solve the hyperbolic partial differential equations, can be employed. The finite element method is generally used in the case that the deformation field is considered via the motion equation.

## V. CONCLUSIONS

We presented the mass conservation law of single phase and the conservation law of spin angular momentum, which have been developed in the previous work, by the use of the order parameter and crystal orientation. The constitutive equations of fluxes of the order parameter and crystal orientation were thermodynamically derived. Additionally, we modeled the reaction term of the basic equation of the order parameter and diffusion coefficient of the basic equation of the crystal orientation by synchronizing with the temporal change in the grain boundary energy. The conclusions we obtained from this study are summarized.

The balance equation of the order parameter is the same as the mass conservation law in which the volume fraction of the recrystallized phase is expressed by the order parameter. The mass source of the governing equation of the order parameter can be modeled by the reaction rate equation of grain growth and grain boundary migration.

The balance law of the crystal orientation is derived by the conservation law of the spin angular momentum for the mixture of the lattice. The governing equation in this study was derived in the form of the sum of the hyperbolic equation and parabolic equation. If we assume that the recrystallization is limited to the dissipative process and integrating the above equation with respect to time, the three-dimensional parabolic-type equation of the crystal orientation is obtained.

When the phase-field model for recrystallization is derived based on the conservation laws, the conservation laws should be developed for the generalized continuum in which the lattice is modeled as a kind of directa.

Currently, we do not consider the elastic energy in the model. However, we can take the elastic energy into account. In that case, the elastic energy can be considered in the source term based on reaction velocity theory. If we develop the energy of a grain boundary not only with the gradient of the crystal orientation but also with the elastic energy, this can be achieved. This model is expected to be interesting since it is reasonable to consider the elastic energy around dislocations when simulating recrystallization phenomena. The development of this model is a future work.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

### APPENDIX: THE TENSOR OF THE MOMENT OF INERTIA

We derive the material derivative of the moment of inertia averaged in the RVE $J\u0307$. First, we define *J*^{(k)}. When a lattice element is modeled as a bar (Fig. 7) and its center of gravity is fixed, the spin angular momentum of the lattice element is expressed as follows:

where *ξ*^{′(k)} is the position vector to the arbitrary point from the center of the gravity of the lattice element, *ρ*′ is the mass density of the lattice element, *Ω*^{(k)} is the spin angular velocity, defining *Ω*^{(k)} ≡ dual **ω**^{(k)} [ ], and *v*_{l} is the volume of one lattice. Defining *I*′^{(k)} ≡ *ξ*′^{(k)} ⊗ *ξ*′^{(k)}, the integrated function of Eq. (A1) can be expanded as

Here, *J*′^{(k)} is given as

Substituting Eq. (A2) into Eq. (A1) and considering *ω*^{(k)} is constant in the region of *v*_{l}, the left-hand side of Eq. (A1) is reduced to

Here, *J*^{(k)} is expressed as

Considering Eq. (A5) and writing *J*^{(k)} using *I*^{(k)},

where $ml$ means the mass averaging in the lattice element. Hence, taking a material derivation of both sides of Eq. (A6), we have

where the material derivative of *I*^{(k)}, $\u0130(k)$, can be written as

Using orthogonal tensor *R*^{(k)} and initial position vector $\xi \u20320(k)$, *ξ*^{′(k)} can be expressed as the follows:

Taking the material derivative of Eq. (A9), we obtain

where we define $\Omega (k)\u2261R\u0307(k)R(k)T$ [ ]. Substituting Eq. (A10) into Eq. (A8), the following equation is derived:

Here, ** Ω** and $\Omega \u0303(k)$ are the mass average and fluctuation of

*Ω*^{(k)}, respectively, and $m$ is the mass average in the RVE. Next, dividing

*I*^{(k)}into the average and fluctuation in the RVE as $I(k)=I(k)m+\u0128(k)\u2261I+\u0128(k)$, Eq. (A11) is reduced as follows:

Generally, the correlation of the fluctuation is not always infinitesimal. In this case, **I** has the order of square of the length of the lattice element, whereas $\u0128(k)$ has the order of the square of the amplitude of atom vibration around equilibrium position, which meaning that ⟨ ⟩ in right-hand side is negligible since sufficient small order to ( ) in right-hand side. Therefore, Eq. (A12) can be rewritten as follows:

Moreover, multiplying **ω** to Eq. (A13), and using the relationships (** Ωω**)

_{i}=

*Ω*

_{ij}

*ω*

_{j}=

*e*

_{ikj}

*ω*

_{k}

*ω*

_{j}= 0, we obtain the following equation: