Gaussian-process based modeling and optimal control of melt-pool geometry in laser powder bed fusion

Abstract

Studies have shown that melt-pool characteristics such as melt-pool size and shape are highly correlated with the formation of porosity and defects in parts built with the laser powder bed fusion (L-PBF) additive manufacturing (AM) processes. Hence, optimizing process parameters to maintain a constant melt-pool size during the build process could potentially improve the build quality of the final part. This paper considers the optimal control of laser power, while keeping other process parameters fixed, to achieve a constant melt-pool size during the laser scanning of a multi-track build under L-PBF. First, Gaussian process regression (GPR) is applied to model the dynamic evolution of the melt-pool size as a function of laser power and thermal history, which are defined as the input features of the GPR model. Then a constrained finite-horizon optimal control problem is formulated, with a quadratic cost function defined to minimize the difference between the controlled melt-pool size and its reference value. A projected gradient descent algorithm is applied to compute the optimal sequence of laser power in the proposed control problem. The GPR modeling is demonstrated using simulated data sets, a mix of simulated and experimental data sets, or pure experimental data sets. Numerical verification of the control design of laser power is performed on a commercial AM software, Autodesk’s Netfabb Simulation. Simulation results demonstrate the effectiveness of the proposed GPR modeling and model-based optimal control in regulating the melt-pool size during the scanning of multi-tracks using L-PBF.

Appendices

A. Computation of initial temperature

Consider the multi-track shown in Fig. 1. Define a fixed coordinate system with the origin at the starting point of the first track, x-axis pointing to the laser scanning direction of the first track, and y-axis pointing to the direction where subsequent tracks will be scanned. Without loss of generality, consider that the laser is melting the ith track and it is going to scan a point with coordinates $\mathbf {p} = [x; y]$ on this track. Noting that a multi-track single-layer is considered here, thus 2-dim coordinates are sufficient to describe any point in the build plane.

When the laser scans the first track (i = 1), the initial temperature $T_{init}$ equals to the ambient temperature $T_a$, i.e., $T_{init} = T_a$. When $i > 1$, to account for the temperature contribution of each past track j ($j = 1, \ldots , i-1$) to the point $\mathbf {p}$, a virtual heat source j is assigned to replace the original laser heat source as soon as it finishes scanning track j. The virtual heat source continues to travel along the same direction at the same scanning speed $v_j$ and with the same power value ($Q_j$) as the physical laser heat source at the end of track j. A pair of virtual heat source and heat sink could be defined to better characterize the temperature cooling after the physical laser heat source finishes scanning the track j, but the temperature difference between using a pair of virtual heat source/sink and using a single virtual heat source is negligible in this study, and hence the latter is used for simplicity.

Let $\mathbf {p}_j = [x_j; y_j]$ denote the coordinates of the virtual heat source j. Computation of such coordinates could easily account for the build plan and inter-hatch dwell time (including skywriting time and any additional dwell), during which the virtual heat source keeps traveling. Let $q_j$ denote the net power of the virtual heat source of track j, which equals to the product of laser power $Q_j$ and the laser absorption efficiency. Then the initial temperature at point $\mathbf {p}$ on the ith track is computed as the summation of the ambient temperature and temperature contributions from all virtual heat sources j, $j = 1, \ldots , i-1$, given as follows:

$$\begin{aligned} T_{init} (\mathbf {p}) = T_a+\sum _{j=1}^{i-1} \frac{q_j}{2\pi k \Vert \mathbf {p}-\mathbf {p}_j\Vert _2} e^{-\frac{v_j(w_j+\Vert \mathbf {p}-\mathbf {p}_j\Vert _2)}{2a}} \end{aligned}$$

(28)

with each term in the summation representing the Rosenthal’s solution (Rosenthal 1946) due to a virtual heat source j. In (28), $w_j = x - x_j$, $\Vert \cdot \Vert _2$ denotes the 2-norm. Constant values of the thermal conductivity k and thermal diffusivity a for Inconel 625 listed in Table 1 are used here, in contrast to temperature-dependent material properties used in FEA. For the ith track, to compute $T_{init}(s)$ using (28), the coordinates $\mathbf {p} = [x;y]$ can be easily computed as follows: $x = s$ for an odd-number track and $x = L-s$ for an even-number track; $y = (i-1)\cdot h$ with h denoting the hatch space. Our prior study showed that this analytical computation has a good agreement, to a certain extent, with FEA simulations (Li et al. 2017).

B. Computation of GPR estimated melt-pool volumes and partial derivatives

Consider that the design matrix $\mathbf {X}$ consists of n training data, i.e., $\mathbf {X} = [\mathbf {x}_1, \cdots , \mathbf {x}_j, \cdots , \mathbf {x}_n]$ where each training data $\mathbf {x}_j = [V_j, Q_j, T_{init,j}]^T$ consists of three elements $\mathbf {x}_j(1) = V_j$, $\mathbf {x}_j(2) = Q_j$, and $\mathbf {x}_j(3) = T_{init,j}$. By the squared exponential covariance function in (10), the matrix $K(\widehat{\varvec{\xi }}(s), \mathbf {X})$ of dimension $1 \times n$ is given in (29) below. For $\mathbf {w} = [w_1, \cdots , w_j, \cdots , w_n]^T$, the estimated states $\widehat{V}(s+1)$ in (78) - (9) and their partial derivatives are given in (30) - (32) below.

$$\begin{aligned}&K(\widehat{\varvec{\xi }}(s), \mathbf {X}) = \sigma ^2_f \cdot \left( \begin{array}{c} e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_1(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_1(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_1(3))^2/\sigma ^2_{l3}]} \\ \vdots \\ e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_j(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_j(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_j(3))^2/\sigma ^2_{l3}]} \\ \vdots \\ e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_n(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_n(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_n(3))^2/\sigma ^2_{l3}]} \end{array}\right) ^T \end{aligned}$$

(29)

$$\begin{aligned}&\widehat{V}(s+1) = \sum _{j=1}^n w_j \sigma ^2_f \cdot e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_j(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_j(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_j(3))^2/\sigma ^2_{l3}]} \end{aligned}$$

(30)

$$\begin{aligned}&\frac{\partial \widehat{V}(s+1)}{\partial \widehat{V}(s)} = - \sum _{j=1}^n w_j \sigma ^2_f \cdot \frac{\widehat{V}(s) - \mathbf {x}_j(1)}{\sigma ^2_{l1}} \cdot e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_j(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_j(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_j(3))^2/\sigma ^2_{l3}]} \end{aligned}$$

(31)

$$\begin{aligned}&\frac{\partial \widehat{V}(s+1)}{\partial Q(s)} = - \sum _{j=1}^n w_j \sigma ^2_f \cdot \frac{Q(s) - \mathbf {x}_j(2)}{\sigma ^2_{l2}} \cdot e^{-\frac{1}{2}[(\widehat{V}(s)-\mathbf {x}_j(1))^2/\sigma ^2_{l1} + (Q(s)-\mathbf {x}_j(2))^2/\sigma ^2_{l2} + (T_{init}(s)-\mathbf {x}_j(3))^2/\sigma ^2_{l3}]} \end{aligned}$$

(32)