Thermal Modeling in Metal Additive Manufacturing Using Graph Theory: Experimental Validation With Laser Powder Bed Fusion Using In Situ Infrared Thermography Data

Fast and accurate computational approaches to predict the temperature distribution (thermal history) in additive manufactured parts are valuable for the understanding and prevention of flaw formation, among other critical functions exemplified in Fig. 1 [1–4].

In a previous paper, we proposed a graph theory-based computational heat transfer approach for predicting the thermal history in additive manufactured parts in near real time [5]. In that paper, the graph theory-predicted temperature trends were verified with (a) exact analytical Green’s function-based solutions, (b) finite element (FE) implementation of Goldak’s double ellipsoid moving heat source model [6,7], and (c) a commercial software for thermal simulation in additive manufacturing (AM) (Autodesk Netfabb). Results from our prior work showed that the graph theory approach was about ten times faster than the benchmark Goldak’s model implemented in a commercial FE software (Abaqus) [5]. The mean absolute percentage error (MAPE) of the graph theory-derived predictions relative to FE analysis was less than 10%.

The objective of this paper is to validate the graph theory approach in the specific context of the laser powder bed fusion (LPBF) AM process using in situ infrared thermal measurements. In LPBF, metal powder is deposited on a bed (build plate) and selectively melted layer-upon-layer with a laser. The temperature gradients induced in the part during LPBF is one of the main causes for flaws such as cracking and distortion in shape [2,4].

To realize the foregoing objective, we frame the following four tasks:

1. Two test part geometries are made in two separate builds on a Renishaw AM250 production-grade LPBF system. The two test parts are described below.

   • A cylinder of diameter 8 mm and height 60 mm. The test cylinder is sintered alongside eight other identical cylindrical-shaped parts on the build plate. The total time for the build is 171 min (1200 layers).

     The build plan for the test cylinder has three phases. First, for the first 20 mm vertical height of the test cylinder, it is scanned simultaneously with the rest of the eight cylinders. In the second phase, the test cylinder is processed to a total height of 40 mm, while the scanning of the rest of the eight other cylinders is paused. The third and concluding phase is identical to the first phase—the test cylinder is processed to its final height of 60 mm along with the rest of the eight cylinders.

     The preceding build strategy, where some parts are intermittently scanned, leads to a variation in the time required by the laser to process a layer, because the laser takes longer to complete a layer in the first and third phases when more parts are scanned, compared with the second phase. Accordingly, the time elapsed between the processing of two successive layers—called the interlayer cooling time (ILCT)—varies across the build. The layer-to-layer variation in ILCT leads to microstructural heterogeneity in LPBF parts [8].

   • Inverted cone shape with a bottom diameter of 2 mm, top diameter of 20 mm, and vertical build height of 11 mm. The build time is nearly 51 min (220 layers). The gradual increase in the surface area of the cone as a function of its vertical build height causes a variation in the ILCT.

     Simple cylinder and cone-shaped test parts are used because the ILCT can be readily determined compared with complex-shaped parts; the ILCT is a critical input parameter for model validation.

2. Surface temperature measurements for the two test parts are acquired layer-by-layer using an in situ longwave infrared (LWIR) thermal camera. The surface temperature measurements is recorded over their entire build duration consisting of 1200 and 220 layers for the cylinder and inverted cone shapes, respectively. To the best of the authors’ knowledge, model validation efforts documented in the literature use in situ temperature measurements from at most 25 layers.

3. The steady-state surface temperature for the two test parts is predicted using the graph theory-based approach and subsequently validated layer-by-layer relative to the experimental measurements acquired using the thermal camera.

4. The layer-by-layer thermal history predictions from the graph theory predictions are compared with a previously published FE model [9]. The comparison of FE and graph theory approach includes predictions of the steady-state surface temperature, as well as the temperature at a point in the interior of each of the two test parts.

The rest of this paper is structured as follows. In Sec. 2, we review the prior literature and delineate the challenges involved in the acquisition of thermal history in LPBF. Section 3 describes the experimental methodology, and adaption of the graph theory and FE approaches to predict the surface temperature. Section 4 reports results concerning the validation of the graph theory with experimental measurements, and comparison with FE analysis. Lastly, conclusions and avenues for future work are summarized in Sec. 5.

Articles by Yan et al. [10] and Tapia and Elwany [11] review in situ thermal measurement approaches in AM. There are two approaches to obtain part-level in situ temperature in LPBF. The first approach is to embed thermocouple(s) inside the part or in the substrate, for measuring temperature trends at a few discrete locations. The second approach uses thermal imaging to measure the part surface temperature [12,13]. This work applies the second strategy. In Secs. 2.1 and 2.2, we highlight the key challenges in both the thermocouple and infrared thermography approaches, respectively, and justify the use of infrared thermography as a viable means to validate the graph theory-based thermal model.

The temperature profile at discrete points in LPBF parts is obtained by embedding thermocouple(s) (i): inside the substrate, (ii) incorporating thermocouples inside a prebuilt part(s) and then building the test part(s) over the prebuilt part(s). To the best of our knowledge, there are no examples in the literature that describe stopping the LPBF process to instrument thermocouples inside the test part.

Researchers have acquired the temperature trends at the underside of the part by brazing thermocouples on the surface of the build plate in a manner such that the head of the thermocouple is barely exposed [12,13]. However, it is observed that the temperature signals obtained by the thermocouple inside the substrate are considerably attenuated as the part grows in size.

For example, Dunbar et al. embedded both a thermocouple and strain gage array within the build plate of an EOS M270 machine to validate their predicted distortion trends [12,13]. In their setup, the sensor array is coupled to a battery-powered data acquisition system incorporated underneath the build plate. Thermocouples are drilled through the build plate, such that the tip of a thermocouple is exposed (≈0.25 mm) above the build plate. The build dimensions for the test coupons used by Dunbar et al. are 6.25 mm × 6.25 mm × 2.33 mm (vertical build height), and the experiment lasts 10 min, in which temperature data are acquired for a maximum of three layers.

Promoppatum et al. [14] used a setup similar to that of Dunbar to acquire temperature data for a large 165 mm × 60 mm × 70 mm (vertical build height) stainless steel part. Temperature trends at five discrete points on the underside of the part were tracked using thermocouples. The temperature readings recorded by the thermocouples at the underside of the part reduced to a steady-state temperature of 200 °C within 25 layers.

Similar attenuation of the temperature signature acquired by a thermocouple embedded in the substrate is also observed in a recent work of Wood et al. [15]. Experiments were conducted on a custom experimental open architecture LPBF setup at Edison Welding Institute. Prebuilt cuboid-shaped stainless steel coupons of ∼12 mm × 12 mm × 12 mm size were embedded with four thermocouples at different layer heights, and one thermocouple was embedded into the build plate. In one of their tests, Wood et al. [15] deposit a total of five layers (200 μm total thickness) on the prebuilt coupons. Thermocouple data were acquired for roughly 7 s per layer. The thermocouple embedded in the build plate did not register any temperature variation.

Researchers have obtained the temperature trends at discrete points inside the part and substrate in the directed energy deposition (DED) AM process. In DED, because the part is not surrounded by powder, a thermocouple can be readily embedded onto the substrate or spot welded on the surface of the part by stopping the process as demonstrated in the work of Heigel et al. [16]. In one of our prior works, we validated the graph theory approach for the DED process with temperature data obtained from thermocouples embedded inside the substrate for titanium alloy test parts [17].

The preceding literature review highlights the difficulty in obtaining the thermal history at discrete points via contact-based thermocouples embedded inside the LPBF substrate or part as the temperature signature attenuates after a few layers.

Given the challenge in measurement of temperature trends in the inside of the part with thermocouples in LPBF, researchers frequently use a thermal camera to obtain the relative temperature trends on the top surface of the part before the next layer of powder material is deposited [11,18–21]. The infrared camera is typically mounted either inside or outside the chamber at an angle to the powder bed—called staring configuration. However, the surface temperature data acquired by the thermal camera is a relative measurement and not the absolute surface temperature [18,19,22–24]. Hence, the temperature readings captured by the thermal camera must be calibrated a priori under practical LPBF conditions so that they can be scaled to absolute temperature measurements.

Researchers at the National Institutes of Standards and Technology report different approaches to calibrate thermographic measurements in LPBF [18,19,22]. One such approach uses the concept of a black body emitter to calibrate the thermal camera measurements [18]. This approach is exemplified by Rodriguez et al. [25]. The key idea is to embed a thermocouple inside a deep cavity drilled in an AM test part. When the part is heated, the cavity inside the object behaves as a black body emitter per Planck’s law. As the part is heated in a controlled manner to a steady-state temperature, the absolute temperature readings measured by the thermocouple are used to scale (calibrate) the surface temperature readings measured by a thermal camera. In this work, we use the black body emitter approach to calibrate the thermal camera. The calibration procedure is summarized in Sec. 3.2 and described in detail in Ref. [8].

The schematic of the experimental setup is shown in Fig. 2; further details are reported in Ref. [8]. A LWIR thermal camera (FLIR A35X) with a spectral range of 7.5 μm to 13 μm is incorporated within the build chamber of a Renishaw AM250 LPBF machine. The thermal camera is sealed inside a vacuum-tight box with a germanium window and focused onto the build plate inclined at an angle of 66 deg from the horizontal. The configuration of the infrared thermal camera allows the measurement of surface temperature of the entire top surface of the part.

Thermal images are captured at a resolution of 320 × 256 pixels, providing a pixel resolution of approximately 1 mm², and recorded at a rate of 60 frames per second. The response time for the sensor is approximately 12 μs. The calibration process used in our previous works is briefly described in Sec. 3.2 [8].

To calibrate the temperature trends captured by the thermal camera, a cylinder-shaped test artifact is made using LPBF. The calibration test artifact is identical in geometry, material (SAE 316L), and LPBF processing parameters used for the two experimental builds (Sec. 3.3). The calibration setup shown in Fig. 3 is adapted from Ref. [8].

The temperature of the calibration artifact is controlled using a 200 W cartridge heater embedded in a recess in the bottom. The calibration artifact is heated, and its resulting surface temperature is recorded using two thermocouples located in two respective recesses milled on its top surface. One of these thermocouples (TC1 in Fig. 3) is used as a feedback control for the cartridge heater, while the other (TC2) records the temperature trends used for calibration.

The thermal camera is calibrated in the range of 300 K–800 K because the maximum temperature to which the cartridge heater is operational is 800 K. A 9-pixel × 9-pixel (9 mm²) sample of thermal intensity values in the center of the top surface of the calibration artifact are extracted from the thermal camera data.

A calibration function (Fig. 3(c)) is obtained by fitting the average intensity over the 9-pixel × 9-pixel sample area recorded by the thermal to the mean temperature recorded by the thermocouple TC2. To ascertain the uncertainty in the thermal camera readings, the calibration procedure is repeated ten times (n = 10). The 95% confidence interval in temperature readings in the interval of 300 K to 800 K ranged from 0.1% to 1% of the mean temperature reading [8]. For temperature readings beyond 800 K, we expect the calibration function to remain valid, as it is derived from Planck’s law and the emissivity would not change significantly until melting occurs (viz., 1643 K).

The calibration procedure is repeated with a thin layer of unmelted powder deposited on top of the calibration artifact, and the test data are used to derive another calibration curve. Such a two-part calibration procedure, with a solid part, and with unmelted powder layered on top, ensures that the temperature readings account for the change in material emissivity in LPBF after a layer is fused (but before a new layer is recoated), and after a new layer is added (but before it is melted).

In this work, we make two LPBF test parts in two different builds that are designed to influence the surface temperature in the part through variation in the ILCT. The scan pattern, process parameters, and material properties for the two builds are reported in Table 1.

The test part is a cylinder of diameter 8 mm and height 60 mm in the center of the build plate. This cylindrical test part is built in three phases, as depicted in Fig. 4. The test part is built with a laser power of 200 W, while the rest of the eight other cylinders are built at 5 W; all parts are built without anchoring supports. As we will explain shortly, building the rest of the eight cylinders prevents their collapse during phase 3 of Build 1.

In phase 1, the test part is built along with eight other identical cylinders arranged in a grid pattern. The ILCT in phase 1 is roughly 10.5 s. After a build height of 20 mm is reached (400 layers, each layer is 50 μm), the processing of the rest of the cylinders is stopped, marking the end of phase 1, and start of phase 2. In phase 2, only the test sample, i.e., the center cylinder is processed until a total build height of 40 mm is reached (800 layers). Because only one cylinder is processed, the ILCT reduces to nearly 6.6 s from 10.5 s in phase 1. Lastly, in phase 3, all nine cylinders are again processed for a total build height of 60 mm (1200 layers). Accordingly, in Phase 3, the ILCT again increases from 6.6 s to approximately 10.5 s. The total build time is about 171 min.

In phase 3, because there is unmelted powder underneath the rest of the eight cylinders—there are no anchoring supports below the part, the parts will tend to move and cause a build failure. Hence, the laser power for melting of the eight cylinders around the periphery of the test sample is always set at a minimum of 5 W. In other words, the scanning of the rest of the eight cylinders at low power allows them to be built without supports.

The test part devised for this build is shown in Fig. 5; it is an inverted cone whose diameter gradually increases from 2 mm to a diameter 20 mm over a vertical build height of 11 mm (50 μm layer thickness, 220 layers). The build time is about 51 min.

In this test part, the ILCT increases almost linearly in proportion to the build height from 10 s at the start of the build to 16 s at the final 11 mm build height. Furthermore, the temperature of the top surface increases progressively with the deposition of new layers as the narrower cross section of the part in the preceding layers impedes the diffusion of heat.

We process the surface temperature data (T(t)) acquired by the thermal camera to obtain the steady-state surface temperature between two immediate layers. The steady-state temperature between layers k and k + 1 is represented as $T¯kk+1$⁠. The steady-state temperature $T¯kk+1$ is derived from the time-varying surface temperature T(t) acquired from the thermal camera using the following steps.

As explained in Sec. 3.2, the surface temperature T(t) is the temperature averaged over the 9 mm² area encompassing the center of the part, which corresponds to a 9-pixel × 9-pixel region of the thermal camera image (Fig. 3). Referring to Fig. 6, the surface temperature signal T(t) is distilled into three steps common to both Build 1 and Build 2. The y-axis of Fig. 6 is T(t) in Kelvin. The x-axis is time in seconds; each data point is processed from a frame of the thermal camera (frame rate 60 Hz).

• Step 1: Large upward spike denoting the beginning and end of melting

• In this stage of the process, the laser is active (ON) and is currently scanning the powder bed. A large upward peak is observed when the laser is directly sintering the 9 mm² area sampled from the thermal image. The large upward spike lasts less than 0.5 s (30 thermal image frames). The time from the end of the large upward spike to the start of the next upward spike is the ILCT.

• Precise quantification of the ILCT is critical for model validation purposes; the time t = ILCT in Eq. (1), Sec. 3.5. However, the ILCT is not constant, but can change between layers depending on the build plan and shape of the part. As we will show in Sec. 4.1, for Build 1 the ILCT for phase 1 and phase 3 is approximately 10.5 s, which is considerably longer than the ILCT for phase 2 which is 6.6 s. In Build 2, the ILCT increases continually over the build from 10 s to 16 s. The ILCT is tracked using a spike detection procedure using the findpeaks command in matlab command.

• Step 2: First downward spike due to the recoater blocking the field-of-view of the IR camera when it returns to the powder reservoir

• After the end of melting of a layer, the recoater returns to fetch fresh powder. During step 2, the bed is lowered so that the recoater can pass freely over the powder bed and avoid contact with the part. As the recoater returns to fetch fresh powder, the IR camera field-of-view is momentarily blocked leading to a large downward spike in temperature lasting less than a 1/50th of a second.

• Step 3: Second downward spike due to new powder being deposited on the powder bed when the recoater rakes a new layer of powder on the surface of build plate

• As the recoater makes another pass to deposit a fresh layer of powder, it again momentarily blocks the field-of-view of the IR camera causing a large downward spike in the signal. Because the recoater speed is considerably slower than in the previous step 2, hence the downward peak lasts for close to 1/5th of a second.

• When the powder is initially spread, it is at ambient temperature (300 K). Therefore, the powder will extract heat from the solidified part surface which is still at a higher temperature. The heat required to raise the powder temperature will cause the surface temperature of the part to decrease. This drop in the surface temperature of the part due to the deposition of fresh powder is accounted in the separate powder-related step in the calibration of the infrared camera as described in Sec. 3.2. The corrected temperature signature is overlaid on the surface temperature signal in Fig. 7(a). The insulating nature of the powder causes the corrected temperature to be lower than the bare part temperature. The steady surface temperature readings for layer k and k + 1, $T¯kk+1$⁠, is the temperature in the relatively flat portion of the curve in Fig. 7(b) before layer k + 1 is processed. Described another way $T¯kk+1$ is the minimum temperature recorded just before melting of the new layer.

The graph theory approach is illustrated schematically in Fig. 8 in the context of Build 1. These steps are discussed in detail in our previous works [5,26]. The graph theory approach, as explained in our previous work, converts the part geometry into a set number of discrete nodes [5].

To maintain consistency with the calibrated thermal data, the temperature distribution predictions for the graph theory model are validated against the same 9 pixel × 9 pixel sample region on the surface of the test part shown in Fig. 3. While validation of the graph theory predictions with temperature measurements nearer to the edges of the test part would be valuable, we are constrained by the limited 1 mm² resolution of the thermal camera used in this paper. The part-powder boundary involves complex phenomena encompassing convective and conductive heat transfer modes, compared with the dominance of conduction-based heat loss near to the center of the part, therefore, the measurement uncertainty at the edges of the part would become overwhelmingly large.

To mitigate the computation time, instead of simulating the deposition one individual layer at a time (layer height 50 μm), we adapt the graph theory approach to simulate the deposition of several layers at a time. Such a layer consolidated from two or more individual layers is called a super layer or meta-layer and is commonly used in coarse FE modeling of the AM process to reduce the computation time [9].

Using the super layer approach is particularly well suited to the graph theory method as the precision is independent of the simulated time-step. This is because the time t for which the heat is diffused in the part in Eq. (1) can be set to one large time-step without computing the temperature at intermediate discrete steps as in FE analysis. The time t is set to the ILCT accrued over super layers, and the super layer is varied from 3 mm (consisting of 60 individual layers of 50 μm each) to 0.3 mm (six individual layers).

The graph theory simulation studies require tuning of two types of factors.

1. Number of nodes (N)

   Selecting the total number of nodes (N) into which the part is discretized involves a tradeoff in computation time and accuracy [5]. In our previous work for a complex geometry part, selecting a higher number of nodes results in a smaller error in comparison with benchmark FE studies, while degrading the computational efficiency [5]. We evaluated the effect of varying the number of nodes from 1000 to 5000 in steps of 1000.

2. Model parameters related to heat diffusion

   In the graph theory approach, two model parameters related to the heat diffusion must be determined, namely, the gain factor (g) Eq. (1) and the neighborhood distance (ɛ) which governs the connectivity of the nodes [5]. There is an interaction between these two parameters. To mitigate this complexity, and need for extensive tuning, in this paper we have made one change to the graph theory model, instead of setting ε to an absolute distance in mm, we now connect the nearest 50 neighbors of a node with edges. The number 50 is selected based on extensive offline studies.

   We report the MAPE and root mean square error (RMSE, Kelvin) for each tested combination of super layer thickness (SLT) and number of nodes. To obtain the gain factor (g), we fix the total number of nodes at 1000 and conduct a grid search with respect to the steady-state infrared thermal measurements obtained for phase 1 of Build 1. To make the calibration more rigorous, the layer height for simulation is set at 50 μm, which is the same as the actual layer height of the build.

   The value of g is changed with the number of nodes (N) fixed at 1000 and layer thickness 50 μm. The graph theory approach is applied for the first 20 mm of the build height, i.e., the graph theoretic model is calibrated for the temperature readings from phase 1 of Build 1. The results from the model calibration procedure are shown in Fig. 9(a). The value of g that minimizes the MAPE and RMSE is selected. The results from the grid search are shown in Fig. 9(b). The value of g that minimizes MAPE and RMSE is 1.5 × 10⁴; this value is set constant for all subsequent simulation studies, including phase 2 and phase 3 of Build 1, and entirety of Build 2.

The rest of the material-related constants and simulation parameters are described in Table 2. The simulations are conducted in the matlab environment on a desktop personal computer with an Intel Core i7-6700 CPU, clocked at 3.40 GHz with 32 GB of onboard memory.

The FE approach used for predicting the thermal history in LPBF parts is detailed in our previous publication [9]. To maintain parity, the FE model uses the identical meta-layer or super layers implemented for the graph theory approach. In our prior work, we obtained both the temperature distribution and distortion in an LPBF part by simulating the deposition of super layers. The FE predicted thermal-induced distortions are within 10% of offline measurements [9].

To ensure equitable comparison of FE and graph theory approaches, the following steps are taken: (i) for both the FE and graph theory implementations, MAPE and RMSE are quantified in for the same number of nodes and resolution (super layer thickness); (ii) we compare the computation time required by the FE predictions to converge to approximately the same MAPE and RMSE of the graph theory predictions at an identical level of resolution. Lastly, we qualitatively compare the FE and graph theory predictions for a point in the interior of the two test parts; there is no additional computation cost associated with calculating the temperature distribution at an interior point.

The surface thermal signatures recorded for Build 1 are shown in Fig. 10. The y-axis in Fig. 10(a) is the surface temperature T(t). In Fig. 10(b), the steady-state surface temperature between two successive layers $(T¯kk+1)$ is tracked on the y-axis. Changes in the temperature trend across the three phases of Build 1 are more clearly evident in Fig. 10(b) on processing the temperature signatures T(t) from Fig. 10(a).

A gradual increase in the steady-state surface temperature is observed during phase 1, succeeded by a sharp increase observed at the start of phase 2, and finally followed by a drop at the start of phase 3. These changes in the temperature correspond to the change ILCT; the reason for the sharp increase in temperature in phase 2 is the decrease in ILCT to roughly 6.6 s (Fig. 10(d)), compared with 10.5 s in phase 1 and phase 3 (Fig. 10(c)).

To summarize these observations, in Build 1, the ILCT, and consequently, the surface temperature distribution of the cylindrical test part changes considerably from phase 1 through 3. These temperature trends from Build 1 have two practical implications, as shown in our previous work [8]. First, tasks that require stopping the build, e.g., replenishing the power, re-filling the chamber with inert gas, that entail a change in the ILCT are liable to cause microstructural heterogeneity. Second, it is not viable to optimize the process parameters for one type of geometry, and consider this knowledge as transferable to other situations—the process parameters must be demarcated through in silico thermal experiments for every build if there is any change in the part geometry, orientation, build layout, number of parts, and scanning strategy.

Fig. 11(a) maps the effect of changing the SLT on the steady-state surface temperature distribution predicted by the graph theory approach for number of nodes N = 3000. In Fig. 11(b), the converse case, i.e., the SLT is maintained constant (0.3 mm, six individual layers of 50 μm each) and the steady-state surface temperature distribution with varying N is predicted. A more detailed sensitivity analysis determining the effect of N and SLT on the MAPE, RMSE, and computation time is reported in-depth in Appendix  A.

In general, the prediction accuracy improves (MAPE and RMSE reduces) as the SLT is decreased, and N is increased. However, the relationship is not linear. An amicable balance in both accuracy and computation time is obtained by setting SLT = 0.3 mm (6 individual layers) and N = 3000. The error under these conditions (MAPE) is close to 13%, and the results are obtained in approximately 26 min (≈1/6th of the actual build time of 171 min). The MAPE reduces to ∼9% for N = 4000 and SLT = 0.3, however, the computational time increases to 65 min.

The predictions from the graph theory approach are compared with FE analysis in Fig. 12 and Table 3. As explained in Sec. 3.6, to ensure equitable comparison, we implemented the super layer approach in a FE in a commercially software (abaqus), the detailed implementation of the FE analysis is described in Ref. [9].

In Fig. 12(a), representative thermal trends for two super layer settings 0.3 mm (six individual layers) and 0.5 mm (eight individual layers), with N = 3000 for the graph theoretic approach are compared with FE analysis under identical conditions relative to the experimental temperature measurements. Likewise, Fig. 12(b) shows the analysis repeated for N = 5000. The comparison between FE and graph theory results is quantified in Table 3.

With N = 3000 and SLT 0.3 mm, the MAPE for the FE analysis is approximately 16%, and the results are obtained within 2048 s (34 min). Using the graph theory approach, the MAPE is 14%, and the trends are obtained in 1655 s (27 min) of computation.

Next, we fix the MAPE of ∼9%, and RMSE 16.5 ± 1 K and compare the computational time for graph theory and FE approach for an identical resolution (SLT = 0.3). For the graph theory approach, the MAPE and RMSE reduced to less than 9% and 16 K on increasing N = 4000 with corresponding computation time of 65 min. To achieve the same level of prediction error, it requires the FE approach 5000 nodes, and nearly 175 min. Effectively, the graph theory approach requires 40% of the computation time of FE to reach approximately the same level of MAPE and RMSE with an identical level of resolution (super layer thickness). The computational advantage of the graph theory approach is retained when the number of nodes N = 5000 for both FE and graph theory; the graph theory approach converges ∼30% faster than FE.

The graph theory approach is currently implemented in matlab (an interpreted computer language) which does not allow multicore processing, and the code is not optimized for parallelization. In contrast, the FE analysis is conducted in a commercial package (abaqus), which is threaded through multiple cores. Porting the graph theory approach to a compiled language, such as c++ with code optimization will further increase its computational efficiency.

The procedure described earlier in Sec. 3.4 is used to preprocess the thermal signatures obtained for Build 2. As shown in Fig. 13(a), in Build 2, a gradual increase in the steady-state surface temperature is observed. The ILCT increases with the build height (Fig. 13(b)), because, as the top surface area increases, the time required to scan successive layers also increases.

Consequently, while the temperature increases for layers near the top, there is also more time for the layer to cool because it takes longer time for the laser to scan a larger surface area. Hence, there is a gradual increasing trend in surface temperature as the accumulation of heat near the top surface occurs concurrently with an increase in ILCT. The steady-state surface temperature nearly reaches the liquidous temperature of 316L stainless steel (1600 K). The practical implication of Build 2 is that the process parameters must be adapted layer-by-layer as opposed to a fixed parameter set so that the part temperature remains consistent. This can be achieved by varying the scan pattern, laser power, and ILCT (by pausing the process between layers for a longer time).

The inverted cone shape is more complex than the cylinder as its cross-sectional area changes with the build height, and hence more number of nodes are required to capture the heat flux. As before, a smaller SLT and larger number of nodes both improve the accuracy of the solution.

Sensitivity analysis tracking the effect of number of nodes (N) and super layer thickness on MAPE, RMSE, and computation time is reported in Appendix  B. Representative results are shown in Fig. 14. For instance, in Fig. 14(a), when the number of nodes (N) is set at 4000, and the SLT is 0.2 mm the MAPE ∼6%, and computation time is close to 41 min. In Fig. 14(b), the SLT is held equal to 0.3 mm and N = 4000, which results in MAPE ∼7%, and computation time ∼33 min; the actual build time for Build 2 is close to 51 min.

The graph theory-derived predictions for Build 2 are compared with the FE analysis in Fig. 15 and Table 4. For equitable comparison, the FE analysis is set to a super layer thickness of 0.2 mm and 0.3 mm, and the number of nodes (N) is set at 4000. As apparent from Fig. 15(a), both the FE and graph theory approaches track the increasing surface temperature trends evident in the experimental data. Furthermore, for the results shown in Fig. 15(b), we increased number of nodes (N) for the FE analysis until it converged to a nearly identical accuracy level of accuracy in terms of MAPE and RMSE observed for graph theory.

As exemplified in Fig. 15, and quantitatively in Table 4, for a fixed resolution (SLT), and for an RMSE of 28.5 ± 1 K, and MAPE ∼6%, the FE analysis required N = 6800 and 96 min of computation time. By contrast, for the foregoing degree of prediction error, the graph theory approach required N = 4000, and converged in 41 min. In other words, the graph theory required 40% fewer nodes and converged within 40% of the time required by FE. These results affirm the computational advantages of the graph theory approach over FE.

We compare the graph theory and FE predictions for a point in the interior volume of the two builds. The results, shown in Fig. 16, substantiate that both the graph theory and FE approaches capture the cyclical heating and cooling characteristic of LPBF as the material is deposited and melted layer-upon-layer.

In the case of Build 1, Fig. 16(a), the measurement (sample) location for testing both FE and graph theory models is on the central axis cylinder, 4 mm from the bottom (viz., 1/5th of the height of Phase 1 build of 20 mm; the total build height over the three layers is 60 mm). In Build 2, we observe the temperature at a location on the central axis of the cone, and at a distance of 2 mm from the bottom (∼1/5th of the build height of 11 mm). The x-axis of Fig. 16 corresponds to the build height, in the context of time. Closer examination of the plots, particularly Fig. 16(b), shows a close alignment in the peaks and valleys of the temperature predictions from the graph theory approach and FE predictions.

Observation of Fig. 16(a), in the context of Build 1, reveals that there is a slight increase in the temperature corresponding to start of phase 2, and decrease at the start of phase 3. The change in temperature at the interior point shown in Fig. 16(a), corresponding to the phases of Build 1 is not as prominent as in the steady-state surface temperature plots in due to attenuation of the temperature signature as the part grows in size. The attenuation of the thermal signature harkens to the limitations in acquiring the temperature data by embedding a thermocouple inside of the substrate discussed in Sec. 2.1.

Concerning Build 2, in Fig. 16(b) the attenuation again affects the cyclical pattern as the build progresses. However, the temperature at the measurement point for Build 2 does not decrease as in Build 1, because, as explained in the context of Fig. 13, Sec. 4.2.1, the surface temperature for the conical test part increases continually throughout the process. Although a close agreement in the trends predicted by the FE and graph theory approaches is observed in Fig. 16, however, we note the difference in peak temperature predictions. In the absence of a temperature data in the interior of the part, this difference is attributed to the limited number of nodes used for the FE and graph theory simulations.

This work provides experimental evidence to substantiate the computational efficiency and accuracy of the graph theoretic approach proposed in our previous work [5]. We arrive at the following conclusions through two experimental builds conducted in the specific context of the LPBF AM process.

1. In Build 1, a cylindrical part 8 mm diameter, 60 mm vertical build height is processed with a phased build plan, such that other parts are intermittently added and removed during the build. Consequently, the interlayer cooling time varies over the 1200 layers of the build, which in turn influences the thermal history and microstructure.

   The graph theoretic approach predicts the resulting complex thermal history within 65 min with a MAPE less than 10% and 16 K root mean squared error (RMSE), which is substantially smaller than the actual build time of 171 min. For a comparable level of MAPE and RMSE, a coarse finite element approximation requires 177 min.

   From a practical perspective, Build 1 shows that the graph theory approach is capable of emulating a complex multi-part build plan with test parts being removed and added during the process.

2. In Build 2, a conical part of diameter 20 mm and vertical build height 11 mm is processed in a way such that the diameter of its circular end progressively increases with the build height (an inverted cone). The actual build time is 51 min. The steady-state surface temperature for this test part gradually increases with the addition of new layers, despite processing under constant LPBF parameters. The graph theory approach predicts this increasing steady-state temperature trends within 41 min the surface temperature distribution with MAPE less than 7% and RMSE of 32 K relative to the experimental observations. In contrast, for nearly identical level of MAPE and RMSE, the finite element approach requires 96 min.

The graph theory approach has the potential to facilitate physics-based optimization of process parameters (laser power, hatch pattern, etc.), build strategy, support placement, etc., to minimize warping and distortion. In our forthcoming works, we will validate the graph theory approach for complex geometry parts. However, the prediction of microstructural evolution with graph theory remains a formidable challenge. This is because the microstructure evolution is a function of both the part-level temperature and the meltpool-level thermal-fluid flow phenomena. The graph theory approach currently does not incorporate meltpool-level phenomena.

One of the authors (PKR) thanks the NSF for funding his research through the following grants CMMI-1719388, CMMI-1739696, CMMI-1752069, and OIA-1929172 at University of Nebraska-Lincoln. Specifically, the concept of using graph theory for modeling in metal additive manufacturing applications is funded through CMMI-1752069. The foregoing NSF grant also provided additional supplemental funding toward the lead author (Mr. Reza Yavari) to intern at Edison Welding Institute, at Columbus Ohio through the NSF INTERN program. The authors aprreciate the efforts of the journal associate editor Professor Steve Schmid in shepherding this manuscript through the review process. Sincere thanks are due to the three anonymous reviewers whose comments and suggestions have doubtlessly improved the rigor of this work.

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper. Data provided by a third party listed in Acknowledgment.

Table and accompanying Fig. 17 shows the effect of number of nodes and super layer thickness on the prediction accuracy and computation time for Build 1. Build 1—Cylinder built in multiple phases.

60 mm build height, 1200 layers, 50 μm layer thickness, 171 min (10,260 s) build time

Table and accompanying Fig. 18 shows the effect of number of nodes and super layer thickness on the prediction accuracy and computation time for Build 2. Build 2—Inverted cone.

11 mm build height, 220 layers, 50 μm layer thickness, 51 min (3060 sec) build time