A Numerical Model for the Prediction of Microstructure Distribution Across the Thickness of Quenched Steel Plates

In this work, we apply a numerical inverse analysis procedure based on the ICME framework to ensure a required microstructure and combination of properties in the quenched plate. The microstructure is decided first, and the cooling profile needed to obtain such microstructure is then calculated using the CALPHAD approach. Subsequently, an inverse analysis of the heat transfer problem provides the description of the convection mechanism in quenching that results in the demanded cooling profile. An additional constraint is set on the through-thickness thermal gradient to achieve a homogeneous microstructure. Finally, the resultant microstructure is computed. By means of the proposed model we are able to retrieve the necessary quenching process parameters and to quantify the influence of these parameters on the temperature and microstructure distribution within the plate after the quenching.


Introduction
Quenched high strength steel plates can offer an excellent combination of uniform deformation, yield strength, and toughness. It can be achieved through an appropriate and homogeneous microstructure that results from a precise control over the quenching process conditions. Accelerated cooling via water quenching is one of the most commonly used techniques. However, the governing mechanisms of the heat transfer in the process are difficult to quantify experimentally, especially for thicker plates. This becomes particularly complex when the experiment is conducted in industrial settings, i.e. on an actual quenching line. Furthermore, the impact of the quenching process settings on the microstructure is often difficult to correlate as other components like composition and geometry also influence the microstructural constituents of the steel. The microstructure formed during the quenching process determines the final properties of the quenched steel. In the industrial practice, the target range of the mechanical properties (hardness and yield strength in particular) in quenching are optimized for toughness with the tempering process. These mechanical properties would vary across the thickness of an as-quenched steel plate due to a gradient in the microstructure (i.e. different final distributions of each phase across the thickness). Given the uniform chemical composition prior to rolling and heating, this microstructure gradient is a product of a through-thickness thermal gradient during the quenching process.
In this paper, we explore the concept of integrated property engineering and process engineering. The approach is to set the target in terms of the desired combination of properties in the final product, and then to design a thermal cycle that results in these properties. From a practical point of view, a sequence of important questions needs to be answered to achieve that goal. First, is it thermodynamically possible to achieve a desired combination of properties? In the context of the plate quenching process, the question is if there is a cooling profile that for a given annealed microstructure of a certain chemical composition can result in the required properties. Secondly, we need to identify if the properties are strongly sensitive to the cooling profile, and what variation to the cooling profile still provides an acceptable range of properties. Next, it needs to be verified if the required cooling profile can be obtained for a specific plate thickness and composition, and also if there are no unacceptable deviations from the desired cooling rates everywhere in the material. Lastly, it is to be determined whether the cooling profile is within the operational range of the quenching line. At this point, we leave out other considerations, such as possible productivity issues and financial costs, and limit ourselves to the material and processing aspects.
In the recent years, Finite Element Method (FEM) models are used to assess the properties that are not always experimentally measurable. Rajeev [1] and Jin [2] developed FE-based models to predict the microstructural distortion and hardness distribution in an as quenched sample. Xiaohu [3] constructed a quantitative model for simulating microstructural evolution during the quenching and tempering process of steels. In the present paper, we try to combine both the microstructural evolution and the consequent mechanical properties related to the microstructure across the thickness of construction standard EN10025-6 [4] steel plates. JMatPro software [5] is used, which follows CALPHAD method to predict the major phase transformations, including the formation of ferrite, bainite and martensite, and to calculate the properties, namely strength and hardness, for an arbitrary cooling profile.
In order to get reliable results in the numerical analysis of the physical model an accurate description of the material properties and heat exchanged between the solid body and the fluid, defined through the heat transfer coefficients (HTC), should be provided. If on one hand, the thermalphysical properties can be evaluated quite easily by using the CCT diagrams and cooling curves, the same cannot be stated concerning the definition of the HTC values. For a given thermal treatment it is difficult to have a 'a priori' estimation of these coefficients since they are strictly dependent on the cooling capacity of the quenching machine, the surface conditions and the material of the steel component to be quenched. The solution strategies in literature include approaches based on trialand-errors and industrial experience, and inverse heat transport problems [7]. Due to the limitations of the aforementioned first class of methods, the solution of inverse and optimization problems became during the last decades a successful engineering tool for thermal processes design. References [7] and [8] provide a comprehensive overview of the solution techniques for inverse heat transfer problems.

Method
Metallurgical problem. The plate quenching process involves different interlinked processes, namely heat transfer and phase transformations. The resultant phase fractions defining the microstructure vary across the plate thickness due to different cooling profiles at each point. In the context of answering whether a cooling profile can result in required combination of properties, the important factors that influence the final properties after quenching are the cooling profile and the prior austenite grain size. Thermal problem. Let us consider the body ℬ that occupies a region of the three-dimensional Euclidean space ℇ with a regular boundary in its reference configuration. A variation of temperature T(x,t) within the body ℬ is governed by the heat equation with initial and boundary conditions: where , , are density, specific heat and heat conductivity, respectively, ̇ is the volumetric heat source and ( , ) is the heat transfer coefficient. Following the standard FEM approach [9], the weak form of Eq. 1 is obtained by employing the Galerkin method and is written in spatial configuration, adopting the Lagrangian framework. In the

Materials Science Forum Vol. 949
analyses of the transient simulations, performed by using the commercial software Ansys Mechanical [10], non-linear thermal properties ( , , ) and moving boundary conditions are considered. With specific regard to the coefficient ( , ), in this study the heat transfer coefficient is defined as a dot product of space -time functions where N is the vector size, ( , ) the vector of functions and the vector of parameters used for the parametrization of the heat transfer coefficient. Eqs 1 and 2 define the direct problem. The purpose of the inverse analysis presented in the later part of the paper is to solve the boundary inverse problem, which in our case is to identify the vector of parameters . Optimization problem. Let ( ) be the reference temperature at a material point, given as a function of time, i.e. ( ) = ( ) ( ). The thermal model is parametrized with a set of parameters p, and calculates the temperature as a function of time. The relative differences between the reference temperature ( ) and the temperature calculated at the corresponding spatial point in the model ( ) ( ) is given by where both the reference and calculated temperatures (•) = (•) ( ) are evaluated at fixed time points . The error function per material point is then defined as where is the number of time points. If the reference temperature is known at M material points, the error functions can be conveniently grouped into vector The objective of optimization is to find an ideal set of parameters * for which all components of ( ) are minimal. This essentially defines a multi-objective minimization problem. Since the ideal solution may not exist, it is sufficient to find a set of Pareto-optimal solutions, and select the most convenient one using an additional performance indicator, for instance the weighted average where the relative importance of the objectives is given by vector .
Since the evaluation of the error function in Eq. 4 involves numerical solution of a non-linear transient thermal problem (Eqs. 1-2), it is computationally expensive. Therefore, a direct iterative search for the parameters * is time consuming and may require a good initial guess within the basin of attraction of the global minimum. To circumvent this problem, the minimization of Eq. 5 is solved in three steps. In the first step, a surrogate model of the thermal problem is constructed using the Genetic Agglomeration response surface [11]. To this end, the direct thermal model needs to be evaluated on a number of points in the parameter space to obtain the calibration points. The calibration points are chosen according to the Optimal Space-Filling design of experiment, which is a variant of Latin Hypercube sampling. Once the calibration data are acquired, the Genetic Agglomeration constructs an ensemble of various metamodels (incl. 2 nd order polynomials, non-parametric

34
Simulation-Based Technology Development for Material Forming regression, Kriging, and moving least squares), and then it minimizes the root mean square error between the ensemble output and the calibration points. The details of the implementation can be found elsewhere [12]. The second step is to exploit the surrogate model to conduct the minimization of the approximated vector given by Eq. 5. The multi-objective genetic algorithm based on NSGA-II [13] and implemented in Ansys DesignXplorer [12] is used. The algorithm provides Pareto frontier as its output, from which the solution with minimal performance indicator (Eq. 6) is chosen. That solution is then used as the starting point in the third step, where a gradient-based minimization of Eq. 6 is done. The minimization employs the original direct model and not the surrogate model.

Results
Industrial experiment. The method described above is used in the analysis of plate quenching operation conducted on an industrial line at NLMK Clabecq. In that quenching process, the annealed plate passes five major zones of different cooling efficiency as illustrated in Figure 1. The most efficient water cooling can be achieved in the zones where high-pressure water is applied on the plate. In general, convection is the dominant heat transfer mechanism in this process. Due to accumulation of water on the plate top surface, the convective heat exchange at the top of the plate is systematically higher than at the plate bottom. The resulting asymmetry in cooling rate may be detrimental to the homogeneity of properties, and in general should be avoided. An industrial scale experiment was conducted on a 60 mm thick, 2.5 m wide and 12 m long plate of EN10025-6 conformant steel [14]. The plate was instrumented with a temperature measurement system. The system was equipped with three thermocouples, positioned at ¼, ½ and ¾ of the plate thickness, respectively, and located at the center of the RD-TD plane of the plate. The plate was then placed on the manufacturing line and processed according to the normal operational conditions of the line. The plate was annealed in the pre-heating furnace to the austenitization temperature (921 ºC) and directly afterwards quenched in the quenching unit (see Figure 1). The temperature profiles recorded by the instrumentation were time-correlated with the events of entering and exiting the quenching device. In order to identify the target range of hardness values, an analysis is performed which attempts to study the dependence of this property on the grain size and the cooling rate. A reference cooling profile at a point close to the surface (¼ th of the thickness) is taken to assess the influence of grain size on the hardness value. Figure 2(a) presents the relationship between the grain size and the calculated hardness. The hardness shows a clear minimum around grain size of 20 µm. It can be steeply increased by reducing the grain size, but also with the increase in grain size, albeit not that rapidly. The influential factors in achieving higher mechanical properties are grain size and martensite forming capability. For the smaller grain size though the martensite forming capability is hindered, it is compensated by finer grain size in obtaining higher hardness, while in case of large grain size the martensite forming capability is high enough to obtain the high hardness independent of the final grain size. The industrially quenched steels are further tempered in order to relieve the stresses developed due to fast cooling and increase the toughness at the expense of yield strength and hardness. The steels of higher grain sizes would lose the hardness significantly upon tempering. Moreover, large grains are generally detrimental to toughness. Therefore, smaller grain sizes are preferred for the quenching process of steel plate.
An average cooling rate is calculated for the reference cooling profile taken between the temperature range 750-300 °C where the bainitic and martensitic transition takes place. The hardness values are calculated for different cooling profiles by gradually increasing the average cooling rates. This is repeated for different grain sizes. As shown in Figure 2(b), the hardness gradually increases with increasing the average cooling rate. From the process engineering perspective, if the resultant mechanical properties do not change considerably with variation in the process parameters (cooling rate and the grain size here) then the process is stable and controllable. The main intent of optimizing the cooling profiles is to maximize the properties, while they are not excessively dispersed along the thickness. From an industrial point of view, it is also important to decrease the disparity in the properties from one plate to another plate of same composition. For this, we select the range of hardness values which are almost the same within a range of cooling rates and grain sizes.  Figure 1) with parametrized THC (Eq. 2) has been set up. The parametrization allows one to set the HTC independently for the plate top and bottom at each zone of the plate. In total, vector in Eq. 2 includes eight independent parameters. Since the temperature measurement is available at three points in the quenched plate, it is possible to construct three objective functions (Eq. 3) and formulate the multi-objective problem using Eqs. 4-5. The three-step approach described in section Optimization problem was followed. The FE thermal model was employed to compute time-dependent field of temperature. The temperature profile at spatial locations corresponding to the position of thermocouples in the experiment were then retrieved from the FE solution. Figure 3 compares the temperature evolution at the thermocouples and the temperature calculated by the model at the corresponding location by using the HTC obtained as the result of the optimization procedure. As can be seen, the cooling at the ¼ and ¾ thickness are very accurately reproduced, while at ½ thickness the cooling rate from the model is overestimated. However, the experimental temperature curves exhibit a spurious discontinuity at approx. 400 s, which is most likely an artifact of the measurement system. The largest temperature differences between the experiment and the model are restricted to that short period, as more clearly visible in Figure 3b.  Figure 6a. Subsequently, the cooling rates required to reach the hardness were calculated. This results in three target cooling curves, calculated at ¼, ½ and ¾ of the plate thickness, respectively. As the next step, inverse thermal problem was solved to determine the HTC needed in the quenching unit to approach the target cooling profile. The comparison between the target temperature and the thermal model prediction for the optimized HTC is shown in Figure 4. It is seen that at the outer control points (¼ and ¾ thickness) the desired target cooling profile is reached, while at the plate mid-thickness the material is overcooled in the lower temperature range. It is not merely a coincidental local minimum of the inverse thermal problem, though. Since the optimization procedure followed in solving the inverse problem includes sampling of the parameter space, it is possible to assess the relationship between the objective functions in Eq. 5. As seen in Figure 5, the objective function at ½ thickness is clearly in conflict with the other two criteria. From the wide Pareto fronts seen in the figure, one can conclude that only compromise solutions are possible, since there is no unique solution that minimizes all terms in Eq. 5.
The phase fractions and the essential mechanical properties (hardness, yield strength and tensile strength) as calculated for the cooling curves are presented in Figure 6. As expected, the maximum values of mechanical properties and volume fraction of martensite are obtained in the points exposed to the highest cooling rates. It is also seen that the differences between the target properties and the properties obtained in the optimization process described above are much limited. Interestingly, the deviations of properties at ½ of the thickness, where the discrepancy in terms of the temperature (cf Fig. 4) was the largest, are almost negligible. It can be explained by the fact that the target temperature in the phase transformation regime is very well matched by the optimized thermal model. This might also be related to the averaging of the cooling rates applied in designing the target cooling curves.

Conclusions
This paper discussed an integrated computational method for property-driven optimization of thermal treatment in thick plates. In this method, the target properties and cooling profiles required to reach that target are defined first. Subsequently, the inverse analysis is carried out to determine the HTC coefficients needed in the quenching operation to follow the cooling profile. The method was successfully tested on an industrial scale problem.
• The inverse analysis of thermal problem allows identification of the heat transfer coefficients in the manufacturing operation of the considered industrial quenching unit. • The selection of the cooling profile towards required target properties provided achievable goal for the inverse analysis. The discrepancies between the target properties and the properties obtained in optimization were minor. In the future work, an experimental validation is needed to verify the calculated properties in plates quenched according to the obtained cooling profiles. Accurate link between the heat transfer coefficients in the zones of the quenching unit and the parameters that control the actual quenching process needs to be established as well.