3D phase-field simulations to machine-learn 3D information from 2D micrographs

A novel approach is developed to support retrieval of 3D information from 2D experimental micrographs. The approach utilizes 3D phase-field simulations to train an artificial intelligence machine. In a first step, the phase-field simulations have to be validated to reproduce microstructural features which characterize elementary processes which govern processing and high temperature service exposure. The qualified 3D simulation setup is then applied to produce a high number of 2D simulated micrographs by automated sectioning. These simulated micrographs are then used to train a gradient boosting regression model together with the 3D information from the simulations. In the final step, the model is applied to 2D experimental micrographs to retrieve the hidden 3D features. The approach is generally applicable to all kinds of metallic materials, minerals or ceramics which can be treated quantitatively by phase-field simulations. In this paper we concentrate on the process of directional coarsening, referred to as ‘rafting’, in the field of creep of single crystal Ni-base superalloys. The experimental and modeling aspects of the evolution of the volume fraction of the γ′ phase during long term creep are discussed.

A novel approach is developed to support retrieval of 3D information from 2D experimental micrographs. The approach utilizes 3D phase-field simulations to train an artificial intelligence machine. In a first step, the phase-field simulations have to be validated to reproduce microstructural features which characterize elementary processes which govern processing and high temperature service exposure. The qualified 3D simulation setup is then applied to produce a high number of 2D simulated micrographs by automated sectioning. These simulated micrographs are then used to train a gradient boosting regression model together with the 3D information from the simulations. In the final step, the model is applied to 2D experimental micrographs to retrieve the hidden 3D features. The approach is generally applicable to all kinds of metallic materials, minerals or ceramics which can be treated quantitatively by phase-field simulations. In this paper we concentrate on the process of directional coarsening, referred to as 'rafting', in the field of creep of single crystal Ni-base superalloys. The experimental and modeling aspects of the evolution of the volume fraction of the γ ′ phase during long term creep are discussed. * Author to whom any correspondence should be addressed.
Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
Ni-base single crystal superalloys (SXs) are widely utilized in aircraft industries and power plants due to their outstanding creep resistance under harsh working environments [1,2]. The typical microstructure of this type of alloy is characterize by cuboidal L1 2 ordered γ ′ precipitates embedded in a fcc γ matrix, which will subject to changes and get topologically inverted during creep [3]. Aside from the base Ni-Al alloys, significant amount of alloying elements (such as Co, Cr, Mo, Re, Ru, W, Nb, Ta, Ti, etc) are added to improve their performance in many aspects, which also leverages the complexity of the alloy systems and thereby increases the difficulty to design new superalloys. Give that it is impossible to determine the optimal set of alloys through actual experimental study, creep models were introduced by integrating empirical equations or physical laws to predict the creep properties, among with those models that considered microstructural features such as γ ′ precipitate size, particle spacing and γ ′ volume fraction were proved to be a success [4][5][6][7][8][9][10][11][12].
While the distribution and evolution of γ ′ precipitate size and particle spacing had been broadly coupled into those creep models, the γ ′ volume fraction was assumed to be constant. However, an experimental investigation by Serin et al [13] revealed that during the creep deformation, volume fraction of γ ′ -phase would decrease with a variance that could reach up to ∼ 15%. A similar trend had also been reported by Gaubert et al [14] through phase-filed simulation of an AM1 type Ni-base superalloy, where 10% to 20% decrease of γ ′ was observed. Considering the significant influence of varying γ ′ volume fraction on the creep behavior confirmed by both experimental [15] and theoretical [11] studies, it is of great importance to investigate the time dependence of γ ′ -phase volume fraction.
Although there exists potential improvement of aforementioned models' performance by taking the evolution of f v γ ′ into consideration, explicit and accurate measurements remain scarce. The conventional ways of evaluating three-dimensional (3D) features are based on statistical analysis of two-dimensional (2D) micrographs, out of which the area fraction is always considered as an approximation of volume fraction. However, such methods require unbiased sampling of the specimen, which means either the microstructure is homogenous or enough samples are provided [16]. Apparently, the periodical microstructure possesses by Ni-base SXs is not able to fulfill the former criterion, especially when the precipitates coarsen heterogeneously during creep; the latter criterion can be reached by serial sectioning [17], nevertheless it is extremely expensive and labor-intensive. Modification such as considering the geometry relationship had also been made [18] but only applies to initial states of creep because a cuboidal shape of precipitates is required. In general, estimating 3D stereological information from 2D micrographs has always been a challenge.
Recently, applications of machine learning methods have changed the paradigm of materials researches [19]. However, lacking of sufficient amount of high quality data impedes the deployment of machine learning algorithms. To eliminate the limitation of data, Nwachukwu et al [20] investigated the feasibility to predict the strain level of micrographs obtained from creep tests of Ni-base superalloys using deep learning, where the model was trained using digital images generated by phase-field simulations. A similar strategy is also employed in this study, where thousands of 2D sections taken from 3D phase-field simulations were evaluated and served as training sets to retrieve the 3D information. In principle, the 2D sections contain full information of characteristic microstructural features as well as information about their evolution in dependence of processing and testing conditions. Consequently, the γ ′ volume fraction evolution during creep is investigated through 3D phase field simulation in the present work. The simulation and experimental results are compared through statistical analysis on quantitative microstructural descriptors extracted from 2D sections. A machine learning model is established to learn the underlying relationship between 2D and 3D information, by which the volume fraction can be predicted directly from experimentally measured 2D microstructures.

Theoretical and modeling approach
2.1. Simulation approach 2.1.1. Phase-field functional. In this study, a multi phase-field method [21][22][23] is employed to investigate the relationship between 2D microstructures and 3D volume fraction of a Ni-base SX during high temperature-low stress creep. The minimization of total energy of the whole system under certain conditions is formulated to be the driving force of phase transformation process [21], of which comprises contributions from interfacial, chemical and elastic free energies in the present work. Moreover, a crystal plasticity model is implemented as a modification of elastic free energy density to account for the non-homogeneous plastic deformation. As a consequence, the free energy functional at simulation domain Ω is described as: As for multi-component system, each phase (or domain of interest such as grains and precipitate) is assigned with a phase-field parameter ϕ α ∈ [0, 1] and diffuse interfaces are assumed to separate different phases (ϕ α and ϕ β ). The interfacial energies σ αβ are reformulated in the phases-fields and corresponding gradients: where η denotes the diffuse interface width. The chemical free energy is weighted linearly with the phase fractions: where µ is the chemical potential, c and c α are the chemical concentration of the system and phase α, f α (c α ) is the bulk free energy of an arbitrary phase α with respect to its concentration. The elastic free energy density is modeled as a function of elastic strain E α , and stiffness tensor C α : In such a manner, the microstructure evolution of the system could be simulated by tracking the interface with respect to space x and time t: in which M αβ represents the mobility of the interface between two phases. Considering the feasibility of simulations, a binary Ni-Al system is selected as target system together with adjusted parameters, which has been deliberated in details in the literatures [3,24,25]. Forty-eight 3D microstructures covering the strain level range of 0 to 1.2 % were constructed in the framework of Open Phase [26] software. Each of the cubic boxes comprises 256 × 256 × 256 grid points and represents a system size of (3.81 µm) 3 . The simulated temperature and tensile stress were adjusted to 950 • C and 350 MPa respectively, in order to provide an adequate comparison with the interrupted creep test results from Ruttert et al [27] and Horst et al [28].

2.1.2.
Dislocation density based crystal plasticity model. For the system's elastic and plastic response, deformation gradient description is employed.
where F is total deformation gradient and F e , F tr , F pl are elastic, transformation and plastic deformation gradient. The transformation deformation gradient F tr considers the lattice misfit between gamma and gamma prime phase in the system as following: and α γ and α ′ γ are the lattice parameters of the gamma and gamma prime. Plastic deformation gradient F pl captures the plastic behavior of the system, calculating by the physically based non-local crystal plasticity model [29,30]. The evolution of plastic deformation gradient is: where L pl consider the shear rates on all the slip systems: whereγ s is the shear rate, m s is slip direction and n s is normal to the slip plane of the slip system s. The shear rate is calculated by the Orowan law [31], which considers the dislocation density ρ s and their velocity v s on slip system s: The dislocation velocity takes into account local resolve shear stress τ s and dislocation density dependent critical resolved shear stress τ s where α is Taylor constant, G and b are shear modulus and burger vector and v 0 is the initial velocity and τ 0 is the initial critical resolve shear stress. ρ s is dislocation density i.e. sum of statistically stored dislocation density ρ s SSD and geometrically necessary dislocation density ρ s GND on slip systems. The evolution of ρ s SSD is based upon Kocks-Mecking law [32] and ρ s GND is evolved by taking the gradient of the rate of plastic strain rate. The detailed description is given in the [29,30].

Quantitative microstructure investigation
ParaView [33] is applied to perform the automated slicing on the 3D cubes. For each cube (or strain level), 254 slices that parallel to (1 0 0) plane in gray scale are conducted to mimic experimental behavior, which is deliberated in appendix A, figure A1. Subsequently, the slices are binarized to obtain a more convenient processing for the subsequent analysis. Four features that are mostly used in describing morphology of Ni-base superalloys, i.e. area fraction of γ ′ , width and height of γ ′ , horizontal γ channel width, are measured automatically out of these binary images. Illustration of measuring these features are summarized in appendix B.
To quantify the microstructure evolution, an absolute moment invariant is also adopted as shape factor in the present study [34]: where A is the area of individual γ ′ precipitate,μs are the second order moment invariants.

XGBoost algorithm
EXtreme Gradient Boosting (XGBoost ) is one of the ensemble algorithms that build up the model additively [35]. For a given sample i from the data set D, the output of the model is a collective construction of multiple base learner functions: where K is the number of base learner functions, and f k (x) is the function for the kth learner. The last base learner, f K (x i ), is created based on the evaluation of the output from the previous learnersŷ . Therefore, the loss function L for the data set D can be defined as: where n is the total number of samples. In this work, the classification and regression tree is used as the base learner, and a regularization term Ω(f) is introduced to evaluate the complexity of each tree. The objective function to create a tree f K can be written as: where g i , h i and ω(f K ) are defined as: , where γ and λ are the hyper-parameters, t is the number of leaves of the tree f K , ω j is the output from jth leave which is assigned to data points by I j = {i|q(x i ) = j}. Equation (16) can thus be rewritten as: where G j = ∑ i∈Ij g i and H j = ∑ i∈Ij h i . Therefore, the optimized leave score ω * j for a fixed tree structure can be determined as: with corresponding optimal structure core: The overall workflow of the proposed method is illustrated in figure 1.

Correlation between simulation results and experiments
It is always necessary to evaluate the deviation between the simulated results against the experiments, of which is even more critical in the present work, since both the training set and the test set are derived from phase-field simulations. Figure 2 compares the simulated and experimental creep curves. As can be seen from the comparison, the simulations are in good agreement with the experimental results, which indicates the established phase-field model captures the microstructure evolution of the target system pretty well.
However, since the current study aims to investigate the geometric relationship between 2D sections and 3D volume, it is more persuasive to compare the consistency between simulated and experimental microstructures directly. The 2D sections cut from simulated 3D microstructures and experimental results from literatures [27,28] at corresponding strain levels are displayed in figure 3. Experimental micrographs are cropped to the same size as the 2D cuts from simulations to eliminate the influence of observation window size. The initial microstructure in figure 3 is derived from simulating the precipitation hardening heat treatment process of a Ni-base SXs [36], and exhibits a typical morphology of Ni-base SXs where the cuboidal γ ′ precipitates with similar size and shape are distributed homogeneously in the γ matrix. Comparing to the experiment, the simulated γ ′ precipitates hold rounder corners. The size and shape of γ ′ obtained from simulations remain nearly unchanged at the onset of temperature and external stress, with a slight increment of the horizontal γ channel width. Subsequently, the microstructure undergoes directional coarsening and forms rafts (ε = 0.6%). At 1.0 % strain, the microstructure is considered to be topologically inverted, in which the γ phase is segmented by γ ′ precipitates and therefore the latter becomes the matrix. The morphology for the experiment is less regular, with noticeable disparity of γ ′ precipitate size. Furthermore, the samples experienced an integrated heat treatment through HIP (Hot Isostatic Pressing) [27], which can be the reason for elongated γ ′ precipitates at initial state, and presumably P-type rafting observed at 0.4% strain. Notwithstanding these flaws, the morphology of the simulation and experiment are quite comparable. Figure 4 presents the quantitative descriptors in this study. For each measured feature, the solid line represents the average values of individual 2D cuts at relative strain levels, the dashed lines enclose the intervals between 1st and 3rd quantile, the dotted-dash lines indicate the max and min values, and the points refer to results extracted from cropped experimental micrographs. As seen in figure 4(a), the γ ′ area fraction exhibits large variance for both simulated and experimental results, which suggests a potential relationship between f γ ′ a and the cutting positions. The area fractions of individual cuts at initial state is thus plotted in figure 5 as a function of cutting positions. It should be noted that the cutting position is actually referring to the indices of grid points over the whole box, not the physical scale. Obviously, f γ ′ a varies according to different cutting positions: a cut may largely intersect with precipitates, of mainly lie in within a channel. This effect will always appear in such kind of ordered precipitate arrangements, in experiment or simulation. It is probably quite strong in the simulations due to the relatively small simulation box. See section 4 for more details. The standard deviation σ a of γ ′ area fraction are plotted in figure 6 as a function of strain. The value of σ a is as high as 7% at the early stages, and deceasing to a relative small value till 1% where topological inversion takes place. Such a phenomenon is actually related with the stereographical problem introduced in section 1. The divergences of width and height of γ ′ precipitates are attributed to the size of observation window. For the shape factor s in figure 4(e), the value of it is close to 1 at early stages, where most of the precipitates keep cuboidal, and close to smaller values exhibits a mixture of P-type and N-type rafting, while only N-type rafting is considered in the present work; (2) the artifacts introduced  during preprocess and segmentation of micrographs. However, it could be concluded that the extracted features are in remarkable agreement between experimental and simulated microstructures. Figure 7(a) presents the volume fraction of γ ′ calculated from three-dimensional phase-field simulations. The starting value is set as 70% to keep in line with the experimental result as well as thermodynamic calculations [37]. The calculated result reveals that the total volume fraction of γ ′ , f γ ′ v , decreases fast at the early stage of creep by 4.5%; at ∼0.55% strain it increases rapidly by nearly 2%. f γ ′ v decreases again till ∼0.80% strain, where the amount of γ-phase increases slightly. Finally, after that f γ ′ v decreases monotonically till 1.2% strain. This behavior can be explained by the competition between chemical, elastic and interface energy contributions and the sharp increase of f γ ′ v at ∼0.55% strain is directly related to the coherency loss between precipitates and matrix, see section 4 for more details. The corresponding microstructures around non-monotonic change points, which are marked as hollow squares in figure 7(a), are presented in figure 7(b). At ε = 0.51%, most of the γ ′ precipitates are coarsened perpendicular to the loading stress, but still keep cuboidal shape; when ε = 0.53%, a small group of the coarsened γ ′ phase starts to merge with each other, where f γ ′ bursts. At strain level between 0.81% and 0.89%, the microstructure does not change much but follows an analogous correlation: the remaining isolated γ ′ precipitates start to merge, result in topologically inverted morphology.

Regression model trained by phase-field simulations
To eliminate the influence by outliers, percentile-based transformation [38] is utilized during preprocess of the data. The extracted feature set is then divided into two parts, where 80% of the data is utilized as train set and the left as test set. The hyperparameters tuning in the present work is executed via scikit-learn (version 1.1.2) [39], and organized as the following: at the first step, grid search cross validation is employed to roughly estimate the possible hyperparameter space for three critical hyperparameters: n_estimators, max_depth and learning_rate. Figure 8 presents the effects of different hyperparameters on the model performance, which is evaluated by mean squared error (MSE, the lower, the better) and coefficient of determination (R 2 , the higher, the better). It can be easily concluded from the plots that the optimal parameter intervals are 600 ± 50, 9 ± 1, 0.15 ± 0.05 for n_estimators, max_depth and learning_rate, respectively. The model is then re-trained and fine-tuned based on the estimated parameter intervals using randomized search cross validation. The derived hyperparameters are listed in table 1. The absolute error, which is based on volume fraction counts as percentage, on the test set is plotted in figure 9 where each scatter point represents an individual sample. The circumstances that show large discrepancies mainly focus at ε < 0.4%, which might be attributed to the large σ a as seen in figure 6. Overall, 96% of the test samples can be predicted successfully using the trained model with an absolute deviation less than 0.5%.

Discussion
According to the present phase field simulations, the volume fraction of γ ′ varies nonmonotonically during long-term exposure under high temperature and tensile stress, with a significant reduction of amount of γ ′ .
According to Goerler et al [40], the elastic energy stored in the channel maintains an equilibrium between the chemical driving force and the interfacial energy, resulting in the distinctive microstructure of Ni-base superalloys with uniform precipitate size and channel width. As a result of lattice mismatch, it is well known that the matrix endures compressive stress while the γ ′ precipices bear tensile stress [41][42][43], producing a higher volume fraction of γ ′ than its value under chemical equilibrium [44]. However, the introducing of external loading establishes heterogeneous stress fields and promotes the directional coarsening of γ ′ particles normal to the applied stress, during which the onset of creep strain in the channel relaxes local misfit stress [24], thereby promoting the growth of the phase and reducing f γ ′ v . As shown in figure 7(b), the rebound of f γ ′ v after passing the creep minimum(∼0.53%), is closely related to the merging of precipitates. At this stage, the loss of coherence permits the coalescence of the γ ′ phase and disrupts the equilibrium by decreasing the interfacial energy, which raises the volume fraction of the γ ′ phase. Taking the area fraction as an approximate estimate of the volume fraction, a similar pattern can be observed in experimental data (hollow circles in figure 4(a)), proving the validity of the current conclusion. In the final stage, where the topological inversion occurs under tension, the volume fraction of γ ′ again decreases to an optimal level under chemical equilibrium. We admit that there are artifacts stemming from the model's merging criteria, and a future quantitative examination of merging effects should be done, which is outside the scope of this article. In any scenario, treating the volume fraction of γ ′ phase as constant during creep would be inadequate.
As mentioned in the introduction, the conventional way to get the volume fraction is taken the area fraction as an approximation. Unfortunately, the prerequisite of large amounts of   On the one hand, the area fraction exhibits large deviation at early stages of creep, and makes it unreliable to approximate the value of volume fraction; on the other hand, as the microstructure gets rafted or topologically inverted, the deviation from area fraction and volume fraction becomes acceptable. While area fraction falls to be the representative value for volume fraction during the whole process, it is still considered as an important indicator as it is the most feasible measurement out of 2D images. It is evident that the probabilities of getting different types of cutting planes depend on the precipitates size, γ channel width and morphology of the microstructure. Such that the machine learning could be utilized to build up the implicit correlation among features and target value. As shown in figure 9, the trained machine gained excellent accuracy on the test data. The trained model is then applied to experimental data for validation. It should  be mentioned that there are no direct 3D experimental results so far and the prediction are compared with values from empirical method [37] as well as measured area fractions. The error bars are estimated from the discrepancies among different measuring point. The results are shown in figure 11, from which the predicted values are in good agreement with empirical methods.

Conclusion
3D phase field simulations show that the volume fraction of the γ ′ phase in single crystal Nibase superalloys evolves during creep. It is affected by an interaction of atomistic (interfacial energy), thermodynamic (evolving solubility with temperature) and micromechanic (elastic strain energy) aspects and should not be considered as being constant.
The present study explores the feasibility of predicting the volume fraction from 2D experimental images by training a machine that gains knowledge from physic-based simulations. The results are in good agreement with empirical observations regarding the evolution of microstructural features during high temperature creep. Further work is required to interpret aspects associated with the interpretation of experimental micrographs and the processing of simulated data.

Data availability statement
The dataset used to support the findings of this study is available from the (http://dx.doi.org/ 10.17632/2tn2yxj4dm.1) Figure A1. Schematic diagram of slicing process.

Appendix B. Automated feature extraction
γ ′ area fraction f γ ′ a : where N γ ′ and N total denote the number of pixels that represent γ ′ phase and the total image, respectively.

Width and height of precipitates:
OpenCV [45] is employed to detect the contours of γ ′ precipitates. The width and height are the bounding rectangle width and length, respectively. γ channel width: Ahmed et al [46]. proposed an automated quantification method to measure the horizontal γ channel width by means of counting the pixel lengths of vertical line profiles.