A computationally efficient thermal model for selective laser melting

Selective laser melting (SLM) is a widely used additive manufacturing method for building metal parts in a layerby-layer manner thereby imposing almost no limitations on the geometrical layout of the part. The SLM process has a crucial impact on the microstructure, strength, surface quality and even the shape of the part, all of which depend on the thermal history of material points within the part. In this paper, we present a computationally tractable thermal model for the SLM process which accounts for individual laser scanning vectors. First, a closed form solution of a line heat source is calculated to represent the laser scanning vectors in a semi-infinite space. The thermal boundary conditions are accounted for by a complimentary correction field, which is computed numerically. The total temperature field is obtained by the superposition of the two. The proposed semi-analytical model can be used to simulate manufacturing geometrically complex parts and allows spatial discretisation to be much coarser than the characteristic length scale of the process: laser spot size, except in the vicinity of boundaries. The underlying assumption of linearity of the heat equation in the proposed model is justified by comparisons with a fully non-linear model and experiments. The accuracy of the proposed boundary correction scheme is demonstrated by a dedicated numerical example on a simple cubic part. The influence of the part design and scanning strategy on the temperature transients are subsequently analysed on a geometrically complex part. The results show that overhanging features of a part obstruct the heat flow towards the base-plate thereby creating local overheating which in turn decrease local cooling rate. Finally, a real SLM process for a part with an overhanging feature is modelled for validation of the proposed model. Reasonable agreement between the model predictions and the experimentally measured values can be observed.


Introduction
Selective laser melting (SLM) is an additive manufacturing (AM) method, where metal powder is consolidated into a solid part in a layerby-layer fashion. In SLM, first, a layer of powder is dispensed uniformly across the built area, enclosed within an atmosphere controlled build chamber. Subsequently, a high energy laser beam is scanned over the powder bed in order to melt and fuse the powder particles in a selective manner. The 2D cross-sectional layout of each layer is attained by slicing the CAD drawing of the 3D part. Each layer hence requires a laser scanning pattern in accordance with the corresponding 2D cross-sectional layout. Finally, the build platform is lowered by a distance equal to the thickness of the powder layer so that the build area is levelled with the powder dispensing system. These three steps are repeated until all the layers are deposited to fabricate a 3D part.
When the part has down-facing or overhanging surfaces, with respect to the building direction, sacrificial support structures are usually included. After the completion of the laser processing, a stress relieving heat treatment is often applied to reduce the residual stresses before the part is detached from the base-plate and support structures are removed. A final finishing process may be performed to reduce the surface roughness and/or induce compressive surface stresses [2] to improve the fatigue resistance.
The thickness of the powder layer in SLM ranges between 20-100 μm [3], whereas the laser spot radius is typically tens of micrometres. The implication of these two characteristic length scales are two-fold: (i) SLM is suitable for manufacturing a small batch of relatively small parts with dimensions ideally no larger than few centimetres since approximately 150 layers are needed for every cm of the part in the build direction, and (ii) geometrically complex, intricate parts can be built with a spatial resolution of one tenth of a mm. Moreover, because of the layer-by-layer realisation of the net shape, building time and manufacturing cost of a part are insensitive to the degree of geometrical complexity the part entails. In addition, almost no material is wasted besides building the supports. These are the distinguishing features of SLM when compared to conventional manufacturing technologies such as casting, forging, milling and CNC. Consequently, SLM is especially promising for manufacturing topology optimised metal parts that are typically geometrically complex. The realisation of topology optimized next generation engineering parts with the aforementioned conventional manufacturing techniques are either too costly or impossible. AM also becomes more attractive when the material of interest, e.g. Ti does not allow for high speed milling or machining.
For SLM to reach its full potential and further increase permitted design freedom, the realisation of defect free parts with high dimensional accuracy and microstructural control are required. Defects such as porosity, surface roughness and formation of dross are particularly detrimental for mechanical performance. Dimensional accuracy required for high precision parts is attained only when thermally induced distortions during the build are kept below tight tolerances. Moreover, residual stresses may also lead to part distortion when the part is detached from the base-plate and/or supports [4][5][6]. Finally, control over the microstructure, i.e. the nature and morphology of the constituting solid phases is also a key factor for determining mechanical properties. Microstructural development is dictated by the process induced (re) heating/cooling cycles and subsequent heat treatments experienced.
All these considerations described above depend on the temperature history of the part during the SLM process. In-situ temperature measurement that could act as a part of the feedback loop for temperature control is an active research topic but not yet available in the state-ofthe-art SLM machines. Visual methods based on pyrometry and infrared cameras [7,8] have shown good potential for in-situ inspection of the SLM process. However, so far these in-situ measurement methods only yield information from the surface, suitable for melt-pool monitoring. Real-time observation of the building process is hindered by limited fields of view, highly temporal load and large amounts of data to be handled [9]. Moreover, because of the harsh processing environment, the monitoring equipment cannot be conveniently integrated into the AM machines, and instead, externally mounted and less capable solutions are in place for temperature monitoring [9].
Thermal process modelling on the other hand, comprising computational prediction of the temperature transients of a part is of paramount importance. Using thermal process models one can quantify the effects of process parameters (laser speed, power and scanning strategy) in combination with the choice of material and the part topology, on the development of temperature. In light of these findings, optimal process parameters and design rules for additive manufacturing can be identified leading to high quality parts.
The physics of the SLM process manifests coupling of a variety of phenomena such as heat transfer, fluid dynamics and solid mechanics. Moreover, the transient nature of the heat transfer phenomena is further complicated by the mismatch of length and time scales that are associated with the spot radius of the scanning laser beam and the part dimensions. The latter is typically on the order of centimetres which implies a total manufacturing time of hours. The implication of the separation of scales is that a very fine spatial/temporal discretisation is required when a standard numerical method, such as finite elements (FE) is employed to solve for temperature transients. This in turn results in a high number of degrees of freedom and a very small step size for time integration which make the model computationally prohibitive.
In so-called continuum models, for the sake of simplicity, the powder is commonly assumed as a continuous medium with a constant laser absorption ratio. Furthermore, the effect of laser induced heating is simplified by a moving point, surface or volumetric heat source [10]. In this context, Solberg et al. [11] developed a general purpose implicit, nonlinear FE code which is able to utilise commodity parallel-processing platforms to model additive manufacturing. Hodge et al. [12] performed thermo-mechanical SLM process FE simulations for a 1 mm 3 cubic part with a layer thickness of 50 μm. In the latter since the layer thickness is representative of the actual SLM process, more accurate insight regarding the thermo-mechanical behaviour is obtained, compared to models where tens of real SLM layers are lumped into one simulation layer for the sake of computational efficiency [13]. However, the computational costs associated with high accuracy is unfortunately prohibitive.
In our previous paper [16], we developed an immensely computationally efficient thermal process model for SLM which describes the moving laser spot with a set of point heat sources introduced along the laser scanning vectors. The closed form analytical solution for a point heat source in a semi-infinite medium is utilised to capture the steep temperature gradients analytically in the vicinity of the laser spot in a transient manner. Correct boundary conditions of the problem are then enforced utilising the superposition principle, i.e. with the combination of image point heat sources situated outside the domain and a complementary smooth correction field. The former are also described analytically, while the latter is solved numerically. Since the steep temperature gradients in the vicinity of the laser spot are accounted for with an analytical description, a coarse spatial discretisation for the complementary numerical correction field becomes acceptable. Consequently, computational costs of the thermal process simulations were significantly reduced, enabling simulations of multiple layers for parts with dimensions in the cm length scale. However, image fields can be applied to bodies having convex surfaces only. Moreover, finding the necessary number of image sources and their locations becomes cumbersome even for relatively simple geometric layouts. Since SLM is intended for manufacturing parts having complex geometries, in this paper we will no longer use image fields to enforce the boundary conditions. In contrast, we investigate how to improve the accuracy of the numerical correction field by means of mesh refinement only in the close proximity of part surfaces. Moreover, in order to further improve both the computational efficiency and the accuracy, laser scanning vectors will be described by means of a line source description instead of an array of point sources. The line source solution will be derived analytically based on the point heat source equation that we used in [16] previously.
The outline of the paper is as follows. Section 2 details the proposed novel thermal model for SLM process. In Section 3, the assumed linearity of the heat equation and the accuracy of the proposed line source equation are first examined by comparisons with a fully non-linear model and experiments. Then two numerical examples are presented to investigate the required degree of mesh refinement near the boundaries and to demonstrate our ability to model the building process of relatively complex geometries. Finally, a real SLM process for a part with an overhanging feature is modelled to further validate the accuracy of the proposed model. The article concludes with a reiteration of the most salient points of the study.

Model description
Consider a three dimensional solid part V with an arbitrary shape that is perfectly bonded onto a base-plate B and being built in an ongoing SLM process as illustrated in Fig. 1. At time t = 0, the part is submerged into the powder bed and a thin layer of powder has been dispensed on its top surface ∂V top (see Fig. 1a). The lateral surface of the part ∂V lat , and the top surface of the base-plate which is not bonded to the part V, subsequently denoted with ∂B top , are also covered by metal powder. At time t = 0, the laser starts to scan over the uppermost powder layer with a laser scanning vector that makes an angle θ with the x 1 -axis. In response, the temperature of the body V and the baseplate B increase as governed by the heat equation where T is the temperature, Q v is the rate of volumetric heat generation, i.e. the heat source term, ρ, c p and k are the density, the constant- pressure specific heat and the thermal conductivity, respectively. It is assumed that the part V and the base-plate B are made up the same material. Since the mean conductivity of the powder (assumed to be a homogeneous continuum) covering ∂V lat ∪ ∂V top ∪ ∂B top is approximately hundred times smaller than that of the solid body [17], it is assumed that the powder has a negligible heat conductivity [16]. Moreover, the amount of heat lost due to radiation and convection from the top surface ∂V top of body V is also negligible in comparison to the amount of heat transferred by conduction within the solid [18,19]. This assumption is justified for low and moderate power densities that are of interest in this study but becomes questionable for instance in the keyhole regime. As a result, no heat-flux (insulating) boundary conditions are appropriate for the surfaces that are in contact with the powder. The base-plate is thick and the build area is much larger than the contact area between the part and the base-plate. Therefore, the temperature at the bottom and lateral surfaces of the base-plate, ∂B bot and ∂B lat respectively, are assumed to be maintained at a fixed value during the entire laser scanning process [17,20]. Consequently, constant temperature boundary conditions are prescribed for ∂B bot ∪ ∂B lat . As long as the conductivity of the powder and the heat losses due to convection and radiation are negligible [18,19], the energy present in the melt-pool is almost entirely transmitted to the body V through the surface ∂V top which is immediately beneath the uppermost layer of powder. Therefore, in our model the uppermost powder layer (hence its heat capacity) is discarded and instead heat sources representing laser energy are positioned directly on ∂V top .
Eq. (1) is nonlinear due to (i) the phase transitions between the solid, liquid and vapour phases, and (ii) the temperature dependence of thermal properties k and c p . The melting-solidification and evaporationcondensation phenomena occurring in body V are localised to the meltpool [21]. Moreover, a material point undergoing melting because of laser heating subsequently solidifies rapidly since the laser spot having a high velocity moves away from the point of interest quickly. Therefore, it can be argued that the latent heat initially absorbed during melting is subsequently released during solidification within a short duration. However, the same logic does not apply to the latent heat of vaporization, since the energy lost during evaporation is not gained back by condensation. Hodge et al. [12] employed a reduced effective laser power to implicitly account for the energy loss by evaporation, radiation, and mass ejection. Following Hodge et al. [12], an effective power of P e = 2P/3 (where P is the nominal power) is adopted in our study. Consequently, the nonlinearity of Eq. (1) is now solely associated with the temperature dependence of thermal parameters c p and k. However, Chidls et al. [22] reported that upon fixing the values of k and c p for an assumed value of temperature, a solution of Eq. (1) for the SLM process can still be well-approximated (see also [16]). Hence, by neglecting the temperature dependence of thermal properties, Eq. (1) can finally be linearised as where α = k/ρc p is the thermal diffusivity at a temperature remains to be determined.

Laser scanning vectors
It remains to account for the effect of laser scanning vectors, i.e. the Q v term in Eq. (2). A scanning vector such as the one depicted in Fig. 1 can be discretised into an array of instantaneous point heat sources, as shown in Fig. 1b. These point sources are activated sequentially, in accordance with the direction of the laser scanning vector and the scanning speed v. The analytical closed-form solution of Eq. (2) that exists for a point source in a semi-infinite space 1 is represented by T I ( ) , which reads where t is the time instance of interest and the index I denotes the Ith point heat source [23]. The parameter A in Eq. (3) is the fractional absorptivity of the powder for laser energy and is determined as described in [17]. Time t I the time shift for the point source to diffuse to a radius equal to the laser spot radius r l [16]. The distance between the material point of interest x i and the point source position x i I ( ) is represented by R. The energy Q v associated with an individual heat source can be determined as P e Δt where P e = 2P/3 [12] is the effective laser power and Δt is the time step equal to the duration between the activation of two consecutive sources, i.e.
with L being the length of the scanning vector, v the laser scanning speed and N the total number of point sources for the scanning vector of interest. Since Eq. (2) is linear, the temperature field due to a scanning vector in the semi-infinite space can be computed by the superposition of individual temperature fields due to point sources, i.e.
where I = 1, … M and M is the number of point sources that are active at a given time t and M ≤ N.
The accuracy of the temperature prediction via Eq. (5) undoubtedly depends on temporal resolution, i.e. the value of Δt. For instance, increasing the total number of sources N or, equivalently, decreasing the value of Δt for a fixed value of v increases the accuracy of Eq. (5). However, from a computational point of view, increasing N adds up to the memory requirements and increases the number of floating point operations. If the time step Δt → 0, the array of discretised point sources converges to a single continuous line source. Consequently, the summation in Eq. (5) can be converted into an integral as The line source starts at time instance t 0 and terminates at t f . If t ≤ t f , the upper bound of the integral t e = t, otherwise, t e = t f . Recall that the discrete line source is represented by a set of discrete point sources where individual activation instances of point sources are given by t I 0 ( ) .
In contrast, the variable t 0 is the local activation instance of the continuous line source given as , where ζ is the local coordinate of a line source as illustrated in Fig. 1a. It is important to note that Eq. (6) can also be cast into a closed form  The error function erfc(ϕ) is defined as For a total number of Θ laser scans, the corresponding T is expressed as where the subscript L denotes it is the expression for a continuous line source, and the superscript j represents the jth scan vector, where j = 1, 2, …, Θ.
The remarkable advantage of using Eq. (7a) instead of Eq. (5) is two-fold. Firstly, the solution of the linear heat equation is now temporally exact, since Eq. (7a) is valid in the limit of infinitesimal time step. Secondly, since the effect of a complete scanning vector is cast into a closed form as a function of t, the temperature due to a line source in a semi-infinite space can be obtained directly using Eq. (7). Using a continuous line source is also computationally more efficient compared to the discrete line heat source comprising an array of point sources that requires the summation of closed form solutions through Eq. (5).
In order to demonstrate the advantages of the line source solution, a 1 mm long single laser scanning vector with the scanning velocity of = v 0.2 m/s applied on a very large and thick base-plate as shown in Fig. 2 is considered. Since the dimensions of the base-plate are much larger than that of the laser scanning vector, it is safe to assume the base-plate to be a semi-infinite medium and the temperature evolution therein can be calculated as described above. Now consider the Point G, located at the end of the laser track (see Fig. 2). In order to compute temperature T at point G estimated with point sources, one can consider an array of point sources as in Eq. (5) each of which has a contribution to the resulting temperature field given through Eq. (3), and the total number of the point sources N depends on Δt c.f. Eq. (4). Alternatively, we can directly use the line source solution derived above, given in Eq. (7). We denote the results attained with the point source solution given in Eq. (5) by T P and the line source solution given in Eq. (7) by T L .
The ratio of T P /T L as a function of Δt is plotted in Fig. 3. Recall that T L is insensitive to the chosen value of Δt. The value of T P /T L thus indicates the amount of error introduced into the solution of T P due to the finite value of the time step Δt. When the value of T P /T L ≈ 1, the T P solution is well-approximated. For instance, for a Δt = 10 −4 s which   Additive Manufacturing 31 (2020) 100955 corresponds to N = 50, the temperature at point G estimated with point sources is approximately 11 times of that of the exact line source solution taken as the reference. This is clearly not a meaningful representation of the scanning vector. The total CPU time for computing T p and T L are represented by t P and t L respectively, 2 and the ratio t P /t L as a function of the step size Δt is depicted in Fig. 4. Recall that, as the value of Δt decreases or equivalently N increases, the accuracy of the solution T P is improved. Conversely, the CPU time required increases, manifesting a compromise between accuracy and computational efficiency. For instance, in order to improve the accuracy of T P from the T P /T L value from 1.744 (at Δt = 10 −5 s) to 1.001 (at Δt = 10 −8 s), the ratio t P /t L shown in Fig. 4 increases from 1 to 307.7, which means the computational time t P increases more than 300 fold. The computational time t L associated with the line source is not only independent from the temporal discretisation of the problem, but also computationally inexpensive with t L ≈0.0013 s in this example. For a regular SLM product, which is typically finished with hundreds of thousands of laser scanning tracks, the difference in computational tractability between line and point source solutions will be immense. Therefore, in the remainder of the paper, the laser scanning will be accounted for by the line source approach exclusively.

Boundary correction fields
In Section 2.1, the temperature prediction for a laser scanning vector on a semi-infinite domain mimicking a very large base-plate is performed using Eq. (2), i.e. the solution of the linear heat equation. However, in order to model the SLM process of a finite size part, a set of boundary conditions should be prescribed. For that purpose, we once again rely on the superposition principle and thus the linearity of the heat equation. The total temperature field T(x i , t) is decomposed into T , comprising the total contribution of the all line source solutions in a semi-infinite space and a complimentary field T to account for the actual boundary conditions as illustrated in Fig. 5. Therefore the total temperature field reads Since the source term Q v due to laser energy has been accounted for in T , the complimentary field T can be obtained by solving with the boundary/initial conditions detailed below. The insulating boundary conditions for the total temperature field T, prescribed on ∂B top ∪ ∂V lat are enforced as follows where m i and l i , illustrated in Fig. 1b Finally, the initial condition at t = t 0 is given as where the initial temperature T c is set to the room temperature.
Note that the temperature field T acts as a correction field when superimposed onto T and prescribe the desired initial/boundary conditions of the problem at hand. Provided T and T Grad(˜) are finite, Eq. (9) describes a smooth T field that can be solved by any standard numerical method such as FE analysis. However, when the laser scanning vectors approach to the lateral boundary ∂V lat , the steep gradients of the analytical field T in the vicinity of the line source causes T to be less smooth and hence a fine discretisation of T is necessary. In our previous work [16], in order to circumvent this problem, image sources too were introduced to account for the boundary correction in the vicinity of a boundary. The temperature field associated with the image sources is denoted by T , and the superposition scheme is then modified as However, the image sources can be used for convex boundaries only [16]. In order to extend the capabilities of our semi-analytical method, here we will adopt the decomposition of T as in Eq. (9) without introducing any image source. Instead a mesh refinement is employed in the vicinity of the boundaries so that the complimentary T field can be resolved numerically by conventional FE analysis. Meanwhile, a coarse spatial discretisation is still used for the rest of the part. As a result, the thermal modelling of building geometrically complex SLM products can be described while achieving considerable computational efficiency.

Numerical examples and validation of results
Four numerical examples are presented in this section. The aim of the first example is to justify the use of linear heat equation upon arriving to the closed form solution of Eq. (7a) for a line source on a semiinfinite space. For that purpose, once again a single laser scanning line on a very large base-plate is analysed, where the effect of boundary conditions can be neglected. Consequently, the temperature due to a single laser scan line can be well-approximated by the temperature field T alone. The melt-pool width estimated by Eq. (7a) is then compared with the predictions of a fully non-linear powder scale model reported in [24]. This model explicitly accounts for the temperature dependent properties, phase transitions and the 3D volumetric nature of a moving heat source. The second example is to demonstrate the accuracy of superposition scheme given in Eq. (9) in comparison with the more rigorous but impractical superposition scheme given in Eq. (11) including image sources. For that purpose a laser scanning vector applied near the edge of a cubic part will be investigated with both superposition schemes. Subsequently, thermal process modelling of a more intricate part having overhanging surfaces and a cylindrical hole and thereby concave surfaces are studied. Finally, a real SLM process of a part with overhanging features reported in [25] is simulated to evaluate the overall accuracy of the proposed model.
Two materials are considered in the subsequent numerical examples. Stainless steel (SS) 316L is selected for the first numerical example and discussed in Section 3.1. Ti-6Al-4V is used for the remaining three examples. The thermal properties of SS 316L and Ti-6Al-4V are given in Table 1. The values quoted therein are representative for a temperature close to the melting point of SS 316L and Ti-6Al-4V respectively. In our previous study [16] it has been shown that the linear heat equation Eq. (2) can accurately predict the melt-pool dimensions for Ti-6Al-4V when the thermal properties are chosen at a temperature approximately equal to the melting point [16].
For all the results presented except for Figs. 10 and 11, the temperature field T is calculated using our in-house MATLAB code and the temperature field T is solved using the commercial FE software ABAQUS. The user subroutines dflux and disp are used to apply the Neumann (Eqs. (10b) and (10c)) and Dirichlet (Eq. (10d)) boundary conditions, respectively. The user subroutine uexternaldb is used to input the laser scanning information comprising the start and the end locations and the scanning directions of every laser scanning vector. For ease of plotting, the contour plots in Figs. 10 and 11 make use of ABAQUS subroutine umatht where the T calculation is incorporated. However, the values of T with umatht can be calculated at the integration points of the elements only. Consequently, the nodal T values can only be extrapolated from the T values determined at the integration points. The in-house MATLAB code is used to directly calculate the temperature T at any point of interest without extrapolation, for all the results presented in this paper except for Figs. 10 and 11. 3

A line heat source in semi-infinite medium
A single laser scanning vector is applied on a very large base-plate, which can be considered as a semi-infinite medium and thus, in the absence of any boundary conditions, the temperature field is predicted by calculating T only. Melt-pool widths are predicted by monitoring the temperature field and are subsequently compared with the predictions of the non-linear powder scale model and experimental findings reported in [24]. The laser spot radius is taken 27 μm and an absorptivity of 0.35 is used in [24]. As stated in Section 2, an effective power of P e = 2P/3 is employed in the proposed model to compensate for the evaporation, convection and radiation associated heat losses. To make a fair comparison, an absorptivity of A = 0.525 (as 2A/3 = 0.35) is used in the proposed model. The simulation and the corresponding experimental results are tabulated in Table 2. Table 2 shows that the proposed model leads to slightly wider meltpools than both the powder scale non-linear model and the experiments [24]. We expect the discrepancy is due to the lack of 3D volumetric description of the laser source, which leads to higher predicted temperatures. Therefore, by utilising a set of thermal properties valid for a temperature around the melting point, melt-pool widths predicted by the linear heat equation agree reasonably well with those of the nonlinear model and the experiments. The biggest error between the proposed model and the experiments occurs at a power of 300 W with a laser speed of 1.8 m/s. For these findings, the error is 13% with respect to the experiments. For the other two sets of process parameters, the errors are both less than 3%.

Simple cubic part
It remains to investigate the accuracy of the superposition scheme of Eq. (9) and study what degree of mesh refinement is needed in the vicinity of ∂V in order to achieve an acceptable accuracy. For this purpose consider a simple cubic component with an edge length of 2 mm, as illustrated in Fig. 6. Now consider the scanning vector depicted in Fig. 6 that runs parallel to one of the edges of the lateral boundary. An effective laser power of P e = 82.5 W is used and the laser speed is 0.5 m/s while the laser spot radius is 20 μm. The absorptivity is chosen to be 0.77. The scanning vector is separated from the edge of the boundary with a distance of 100 μm only. Since the lateral surfaces of the cube are in contact with powder and hence assumed to be insulated, the total heat flux on ∂V lat is prescribed to be zero. The bottom surface which is bonded to the base-plate is fixed to a constant temperature of T c = 200°C for simplicity. Since the scanning vector is situated very close to the lateral boundary, the smoothness of the field T is questionable. In order to capture the anticipated non-smooth behaviour of T numerically, mesh refinement is performed near the boundary. A FE mesh comprising 8-node hexahedral linear elements is employed and the mesh is refined in the proximity of the boundary as shown in Fig. 7. Four different mesh densities with the minimum mesh size δ e ranging from 0.0125 mm to 0.1 mm are considered. We introduce the dimensionless parameters η = δ e /r l to quantify the degree of mesh refinement, where r l is the laser spot radius which dictates the wavelength of the temperature field T for heat sources close to the boundary. Since a single scanning vector is applied to the top surface of the cube, the element size is also gradually increased along the negative x 3 direction. The nodes that are located at the interface between the outer rim with refined mesh and the inner core with coarse mesh are constrained to have the same temperature by applying the tie constraints in ABAQUS.  [14,15] where the temperature field is decomposed into T , valid in semi-infinite space and known analytically and the complimentary correction field T , calculated numerically that imposes the correct boundary conditions of the problem at hand. The temperature field T of the cube is calculated using the two different superposition schemes Eqs. (9) and (11). The temperature history for points C1-C5, which are located on the boundary edge (see Fig. 6), are used to compare the prediction of the two superposition schemes. To have a fair comparison, T fields for both superposition schemes are calculated using the line source description derived in Section 2.1. In our previous paper [16] it was demonstrated that for the simple cubic geometry shown in Fig. 6, the image source method is convenient to apply and the superposition scheme given in Eq. (11) was shown to be extremely accurate for predicting the temperature even with coarse spatial discretisation. Hence, the temperature calculated by the superposition scheme of Eq. (11) is taken as the reference. A single line image source, 0.1 mm away from the lateral boundary as the mirror image of the scanning vector is present for this reference calculation. For both superposition schemes, the T field is solved using FE analysis. For the reference calculation a fixed mesh size of δ e = 0.5 mm, i.e. η = 25 is used. The laser scanning duration is 0.0036 s. A total time duration of 0.01 s is considered for calculating the temperature history.
The temperature for the points C1-C5 is plotted as a function of time in Fig. 8a with a time step of Δt = 10 −5 s. The black solid lines denote results calculated by the superposition scheme of Eq. (11) and constitute the reference, while the purple and yellow dashed lines present results calculated by Eq. (9) with the lowest (η = 5) and highest (η = 0.625) level of mesh refinement, respectively. The blue lines in Fig. 8a on the other hand corresponds to the temperature prediction of Eq. (9) with a uniform mesh size of η = 25 (see the inset in Fig. 8a).
The blue curves in Fig. 8a clearly demonstrate the need for mesh refinement to attain accurate results. For η = 5 (purple dashed lines), the degree of mesh refinement is insufficient to capture the boundary effect very accurately; the corresponding temperature histories of points C1-C5 still deviate from the reference black curves. For the finest mesh with η = 0.625, the corresponding results (yellow dashed lines) have a perfect agreement with the reference (black solid lines). The results obtained by η = 2.5 and η = 1.25 are not shown in Fig. 8a, since they are visually indistinguishable from the yellow dashed curves obtained for η = 0.625. These results thereby indicate that a mesh size of η = 2.5 is suitable for the desired level of accuracy.
Recall that the calculation of T is independent of the time step. However, the accuracy of the solution of the complementary field T is sensitive to the temporal discretisation. The time integration scheme employed makes use of the quantities T x t ( , ) i in order to calculate T(x i , t + Δt). Therefore finite size of the Δt still introduces an error to the prediction of T(x i , t + Δt). In order to asses the degree of this error due to temporal discretisation we evaluate the dimensionless parameter for points C1-C5 where t total = 0.01 s is the total duration for the scanning vector to be applied. In essence, ϕ indicates integral difference between the entire temperature histories, i.e. the difference between the area underneath the T versus t prediction of given mesh size and a time step with that of the reference calculation (illustrated in Fig. 8a by black curves). Recall that the reference calculation is performed with the superposition scheme Eq. (11) using a time step of Δt = 10 −5 s. When the time step is chosen as 10 −5 s, the value of error ϕ for points C1-C5 are plotted in Fig. 8b. The error ϕ for all the points C1-C5 is less than 5% when η = 2.5. This demonstrates η = 2.5 is sufficient for the desired level of accuracy. Then, for η = 2.5, the effect of the time step is depicted in Fig. 8c. When the time step is Δt = 10 −3 s, the error become as large as 30%, while when the time step is 10 −4 s, the error rapidly reduces to less than 5%.
The calculations of T are performed by the parallel computing capability of ABAQUS with 32 CPUs. The total CPU time and the wall clock time of the four mesh densities considered for the time step value of Δt = 10 −5 s are tabulated in Table 3. It can be seen that with the increasing mesh density, both the wall clock time and the total CPU time increases rapidly. According to Fig. 8, as considerable accuracy is achieved with the mesh density η = 2.5 and Δt = 10 −4 s. These two values are hence used for the remaining numerical examples.

Process modelling of a wedge-shaped part
We proceed with examining the SLM process of a wedge-shaped part with two overhanging lateral surfaces and a cylindrical hole as shown in Fig. 9a. The overhanging lateral surfaces both have an inclination of 45°with respect to the x 1 -axis. In this example we aim to demonstrate the process modelling capability of our semi-analytical model on a relatively complex geometry and investigate the effect of laser scanning strategy on the temperature transients during manufacturing. Note that the cylindrical hole introduces concave surfaces that rules out the possibility of using image sources for imposing boundary conditions. Therefore the superposition scheme given in Eq. (9) is used together with mesh refinement in the proximity of boundaries.
An effective laser power of P e = 40 W is applied and the laser speed is 0.5 m/s with the laser spot radius of 20 μm. The absorptivity is also chosen to be 0.77. The layer thickness is 50 μm. Recoating a new layer of powder usually takes around 10 s [27]. In this numerical example, the recoating time between each layer is assumed to be long enough so that the temperature of the part with high conductivity returns to the initial temperature before scanning each new layer for simplicity. Hence the thermal history for scanning each layer is independent. We note in passing that residual heat may build up in the part due to the processing of previous layers. This build up of heat can be easily accounted for by the proposed model, by choosing the initial T field to be equal to the total temperature field before scanning a new layer instead of using Eq. (10e).
The part illustrated in Fig. 9a is also divided into two regions as the outer rim and the inner core. The outer rim of the part is discretised with a fine FE mesh while the inner core and the base-plate are Table 2 The melt-pool width (μm) at different power P (W) and speed v (m/s).
The proposed model 112 106 106 Experiment [24] 109 104 94 Non-linear model [24] 89 ± 4 94 ± 12 96 ± 8 discretised with a much coarser mesh. The 8-node hexahedral elements with linear interpolation functions are used for the entire part and baseplate. The inner core of the part comprises elements with η ≈ 25, i.e. average element size is around 25 times of the laser spot radius r l . In Section 3.2 it was found that an element size of δ e = 0.05 mm (η = 2.5) was sufficient for the desired level of accuracy of the numerical solution of the temperature field T when time integration is performed with a time step of Δt = 10 −4 s. Numerical parameters identical to those in Section 3.2 are adopted here. The temperature values of nodes meeting along the interface between the outer rim and the inner core of the part are constrained to be equal using the tie constraint in ABAQUS. Consecutive scanning vectors having opposite directions (see Fig. 9b), which is a common practice in SLM [28], is utilised to build the part. The orientation of scanning vectors are kept constant for all layers but two different scanning orientations with θ = 0°and θ = 45°are investigated. The border offset in Fig. 9b is λ = 50 μm and the hatch spacing is h = 100 μm. The temperature distributions attained while building the 6th and 14th layers are shown in Figs. 10 and 11 at the instance when the depicted laser track is finished. In Fig. 10, θ = 0°while in Fig. 11, θ = 45°. Figs. 10a and 11 a give the snapshots of temperature distribution while scanning the 6th layer, whereas Figs. 10b and 11 b show the snapshots of temperature distribution while scanning the 14th layer. In Fig. 10c, T (the temperature field due to line sources in a semi-infinite domain only) is also illustrated while scanning the 14th layer with a scanning orientation of θ = 0°, in order to emphasise the contribution of boundary effects through T field to the resulting total temperature field. The contour levels above the melting point are designated with grey colour in order to illustrate the size of the melt-pool at the given instance.
As shown in Figs. 10a and 10 b, the surface of the cylindrical hole during the SLM process, initially constitutes an up-facing surface (see for example layer 6) which gradually becomes a down-facing (overhanging) surface as the build continues (see for example layer 14). It is well-known that [29] heat generated by the laser can be quickly transferred from up-facing surfaces towards the base-plate which acts as a heat sink, whereas the heat transfer from down-facing (overhanging) surfaces are obstructed by the poorly conducting powder hindering the path towards the base-plate. This effect is naturally captured in our semi-analytical model as observed in Figs. 10 and 11. When the laser is near the up-facing surface of the circular hole (Figs. 10a and 11 b) less heat accumulation is observed compared to that of when the laser is near the down-facing surface of the circular hole (Figs. 10b and 11 b). It is also observed that local overheating in the vicinity of overhangs is present irrespective of the scanning orientation θ. The effect of insulating boundaries on the rise of local temperature can also be seen by comparing the temperature distributions in Fig. 10b and c. In Fig. 10c since boundary effects are not accounted for, the size of the overheated zone is smaller compared to that of in Fig. 10b.
The influence of scanning orientation on the temperature distributions can also be observed from Figs. 10 and 11. The scanning orientation θ = 45°is able to heat up the part more effectively compared to θ = 0°. For instance, in Fig. 10a, it can be observed at the left side of the part away from the cylindrical hole, temperature can be as low as 300°C, while in Fig. 11a, higher temperatures are observed. This implies that scanning with θ = 45°would induce a smaller thermal gradient which in turn result in lower residual stresses.
The temperature histories of selected points are also investigated. Fig. 12 schematically shows the left half of the top-view of the uppermost layer whilst building the 6th and 14th layer. Fig. 12a and c represent the 6th layer with the scanning orientation of θ = 45°and θ = 0°, respectively. Fig. 12b and d represent the 14th layer with the scanning orientation of θ = 45°and θ = 0°, respectively. Point F1 is located at the boundary edge of the cylindrical hole, while point F2 is located at the lateral overhanging edges.
In Fig. 13 temperature is plotted against time for point F1 while building the 6th and 14th layers with the scanning orientation of θ = 0°. The blue dotted lines are the corresponding T , which assumes the laser scanning is performed in a semi-infinite medium. The solid lines are the temperature predictions, i.e. = + T T Tˆwhere boundary effects are accounted for through T .
The peak T value, for layer 6 is slightly lower than the temperature T, which is around 1300°C. This is due to the fact that boundary effects associated with the up-facing surface (see for instance Fig. 10a) whilst depositing the 6th layer is not a major impediment for heat flow towards the base-plate. In contrast, the boundary effects of an overhanging surface appearing in the case of the 14th layer significantly obstruct heat flow and, thus, in turn induce a higher temperature. Consequently, the peak value of the total temperature T at F1 for layer 14 is much higher than the associated T and also the peak temperature T whilst building the 6th layer.
It can also be observed from Fig. 13 that the average cooling rate for point F1 from its peak temperature is higher for the 6th layer in comparison to 14th layer which is a direct consequence of the presence of an overhanging feature in the part geometry which reduces the cooling rate locally. SLM imposed cooling rate is essential for the morphology and fraction of α and β phases of Ti-6Al-4V alloy [30]. Fig. 13 also indicates the boundary effects are essential to accurately calculate the cooling rate experienced by the part, especially for designs having overhanging features. It is thus demonstrated that the proposed model is able to pick up the variance of cooling rates with different overhanging angles. Fig. 14 shows the temperature history of point F1 upon depositing the 6th and 14th layer with the scanning orientation of θ = 45°and thereby demonstrating a strong dependence of temperature evolution to the scanning orientation c.f. Fig. 13. The peak temperature values depicted in Fig. 14 when θ = 45°are lower when compared to the peak values in Fig. 13. Moreover, the shape of the curves for layer 14 differ extensively between Figs. 13 and 14. The discrepancy between the peak temperatures can be explained by means of measuring the smallest distance between the nearest scanning vectors and the point F1. When the scanning orientation is chosen to be θ = 45°, the average distances between the point F1 and two nearest scanning vectors responsible for the peaks are shown to be larger in Fig. 12, both for layer 6 and 14. Moreover, for layer 14, when the scanning orientation is θ = 45°, the first one of the two nearest scanning vectors to the point F1 moves away from point F1 whereas the second one moves towards the point F1 and vice versa when θ = 0°(see Fig. 12). The former results with an increase in the temperature just after the global peak as can be seen in Fig. 14. It is therefore important to note that scanning orientation alone can modify the temperature transients which in turn have an impact on the resulting microstructure.
In Fig. 15, the temperature histories of point F2 while building the 14th layer with two different scanning orientations are plotted. Different temperature evolution for point F2 can be clearly observed. As  illustrated in Fig. 12b and d, before the peak temperature is obtained with θ = 45°, a total of four laser scanning vectors has been applied whereas with θ = 0°this number is six. Moreover, the individual laser scanning vectors are longer in the case of θ = 0°. Consequently, the part has been exposed to a higher energy level with θ = 0°resulting with a higher peak temperature compared to the scanning orientation of θ = 45°. Another factor that contributes to the difference between the peak temperatures is that point F2 is 0.0707 mm apart from the two nearest scanning vectors when θ = 0°while the distance between point F2 and the associated nearest scanning vectors becomes 0.0866 mm when θ = 45°. Finally one can also compare the peak temperature attained at point F1 (see Fig. 13) and F2 (see Fig. 15) at layer 14 with θ = 0°. The higher peak temperature attained at point F1 implies the overhanging feature due to the cylindrical hole is more problematic in terms of overheating compared to the overhang at the lateral surface. The complimentary numerical field T is calculated using the parallel computing capability of ABAQUS with 32 CPUs. The wall clock time for a single layer ranges from 7 to 15 minutes. The wall clock time increases with layer number as the number of DOFs also increases. Normally the computational time for a comprehensive fully non-linear model is very long and it may take hundreds of hours to compute a 3D model with several layers of real-time SLM process, see for example [10]. In comparison, the proposed model shows a great advantage in computational efficiency.

Validation of results with a real 3D part
In this section we compare our model predictions with the measurements performed during SLM of a real 3D part, as reported by Hopper [25]. The front view of the part is shown in Fig. 16a. The smallest angle of the overhang made with the x 2 -axis is 26.7°. The layer thickness is 30 μm and the laser speed is 1.25 m/s. The scanning angle is 16.6°with respect to the x 1 -axis and the hatch spacing is 65 μm. The laser powers for scanning the high power region and low power region are 170 W and 155 W, respectively (see Fig. 16b). The corresponding effective power P e = 2P/3 is considered in the simulation as suggested in [12]. The absorptivity is 0.77 and the initial temperature is 20°C. Tetrahedral elements are used to mesh the part and the average size of the elements are specified as 0.5 mm and 0.05 mm for coarse and fine mesh regions, respectively. The time step is 10 −4 s. In our simulation, the recoating time between each layer is again assumed to be long enough so that the temperature of the part cools down to the initial value before scanning each new layer. Since the thermal history for scanning each layer is independent, only scanning the last layer is simulated to make the comparison. Fig. 17 shows the temperature history of point G (see Fig. 16b). Time t′ is defined such that t′ = 0 when the laser is exposed at point G. Some oscillations in the experimental curve can be observed. The oscillations are mainly caused by the noise of the signals. As mentioned in Fig. 9. A wedge-shaped part with two overhanging surfaces with an overhang angle of 45°and containing a cylindrical hole, built on a base-plate. The built part is divided into the outer rim, which is discretised with a fine mesh, and the inner core, which is discretised with a coarse mesh. (b) An alternating scanning strategy with a scanning orientation of θ and a hatch spacing of h is applied for building each layer. The border offset is denoted as λ.
[25], the exposure time of the sensor cannot be too long to capture the rapid change of the temperature. On the other hand, shortening the sensor exposure time induces noise in the temperature measurement. Moreover, the melt pool plume and ejected particles during the laser scanning can also cause some errors in the temperature measurement. It can be observed that the oscillations mainly occur when the laser approaches point G, when t′ < 0. When t′ > 0 and without the interference of the laser scanning, the experimental curve becomes smoother and agrees well with the predictions of simulation.
There are several reasons for the discrepancy observed between the simulated and measured temperature. The proposed model uses the reduced effective power to implicitly account for any form of energy loss. However, neglecting the latent heat and the 3D volumetric nature of the laser source both lead to an overestimation of temperature. Besides, the simplification of neglecting the powder conductivity may also contribute to the overestimation of the temperature [31]. On the other hand, as mentioned in [25], saturation of the image sensor can cause underestimation of experimentally measured values. Moreover, since the measurement device can only work in a certain temperature range, the temperature below 1000°C cannot be measured. Consequently, the simulated curve agrees reasonably well with the experiment.

Conclusions
A semi-analytical thermal model for predicting the complete thermal history during SLM process for a part with arbitrary shape is presented. We have derived a closed form temperature solution for a continuous line heat source in a semi-infinite space which has been shown to be a computationally efficient and accurate representation of scanning vectors applied during the SLM process. The line source solution is easy to implement for arbitrary scanning orientations. In order to predict temperature evolution in finite bodies, the superposition principle is used. A complementary temperature field to impose the physically relevant boundary conditions is solved numerically. This complementary field is smooth away from the boundaries but requires a relatively fine discretisation in the vicinity of boundaries. Consequently, an mesh refinement near the surfaces becomes a pragmatic approach to address the multi-scale nature of the problem that limits the computational costs without compromising the accuracy. The accuracy of the line source equation and the superposition scheme is validated by a non-linear model and the experiments from the literature.
The proposed semi-analytical model has been demonstrated to successfully analyse complex geometries. Numerical examples also showed that it is essential to take into account the boundary effect for Fig. 10. Snapshots of the temperature distribution for (a) layer 6 and (b) layer 14 for a scanning orientation of θ = 0°. (c) Temperature field T due to a line source in semi-infinite space for the layer 14. Points labelled as F1 and F2 will be used to compare pointwise temperature histories.   predicting the temperature distribution of a SLM part. This is especially important to asses part designs containing overhanging features. With the predicted temperature, one may further develop a model which relates the temperature and lack of fusion between laser tracks, and thus the proposed model is potentially helpful to predict and avoid the lack of fusion behaviour by optimizing the laser power, speed and hatching distance. According to the simulation results, overhanging surfaces are prone to overheating and more so around a hole and thereby experience slower cooling rates. Moreover, different scanning orientations also affect the predicted temperature distribution and history henceforth both scanning strategy and part topology can be independently designed to control the microstructure.

Conflicts of interest
The authors declare no conflicts of interest. where L and U are the corresponding lower and upper bound of the integration, respectively. Since the laser scan is assumed to be terminated at t f , the upper bound U would be t r 1/2 when t ≤ t f while U = (t − t f + t r ) −1/2 if t > t f . The parameters of Y, C and L are defined in Eq. (7b) with