Abstract

Current additive manufacturing (AM) technologies are typically limited by the minimum feature sizes of the parts they can produce. This issue is addressed by the microscale selective laser sintering system (μ-SLS), which is capable of building parts with single micrometer resolutions. Despite the resolution of the system, the minimum feature sizes producible using the μ-SLS tool are limited by unwanted heat dissipation through the particle bed during the sintering process. To address this unwanted heat flow, a particle scale thermal model is needed to characterize the thermal conductivity of the nanoparticle bed during sintering and facilitate the prediction of heat affected zones. This would allow for the optimization of process parameters and a reduction in error for the final part. This paper presents a method for the determination of the effective thermal conductivity of copper nanoparticle beds in a μ-SLS system using finite element simulations performed in ansys. A phase field model (PFM) is used to track the geometric evolution of the particle groups within the particle bed during sintering. Computer aided design (CAD) models are extracted from the PFM output data at various time-steps, and steady-state thermal simulations are performed on each particle group. The full simulation developed in this work is scalable to particle groups with variable sizes and geometric arrangements. The particle thermal model results from this work are used to calculate the thermal conductivity of the copper nanoparticles as a function of the density of the particle group.

Introduction

Commercially available additive manufacturing (am) systems are typically capable of producing parts with feature sizes of the order of hundreds of micrometers [1,2]. Although these manufacturing processes are versatile, this size limitation prevents their effective use in the micro-electronics industry, where parts require manufacturing resolutions of the order of a few micrometers. The leading process used to manufacture two and a half dimensional micro-electronics structures consists of a mixture of lithography, etching, and chemical deposition steps, which often require intricate setups and suffer from size limitations. The development of a microscale selective laser sintering (μ-SLS) system addresses the need for a more robust, flexible, and easy to use process while targeting a minimum feature size of the order of a single micrometer [36].

Selective laser sintering (SLS) is an additive manufacturing (AM) process that directly applies laser energy to a deposited powder bed. The laser energy increases the temperature within the bed in the areas of laser application which facilitates the fusion of the powder particles into a solid part. This sintering process is repeated for each subsequently deposited layer until a full part is formed [7]. One of the major advantages of SLS over other commercial AM systems is its ability to build true 3D parts, such as overhanging structures, due to support provided by the unsintered material adjacent to the sintered part. The μ-SLS system functions similarly to traditional SLS processes in that laser irradiation is applied to a particle bed and the process is repeated layer-by-layer. The primary difference between the two processes is the use of nanoparticles as the powdered material in the μ-SLS system, which facilitates the production of parts with single micrometer features. Given the scale of these particles and the required feature sizes, accurate characterization of the nanoparticle beds is needed to optimize the μ-SLS process.

Existing simulation work for SLS systems provides insight into the general process of heat transfer in powdered systems [810], although they focus on larger particles and typically model a melting and solidification process, as opposed to the solid-state sintering process experienced by the nanoparticles used in the μ-SLS system. These particles are of the order of 10–100 nanometers in diameter. The nanoparticles remain solid as they sinter, with grain boundary and surface diffusion dominating the underlying solid-state sintering mechanisms [11]. The effective thermal conductivity of metal nanoparticles within a powder bed system is found to be dominated by the evolving geometries of the nanoparticles as they diffuse into one another, and the grain boundary thermal contact resistance at the interfaces between adjacent nanoparticles [12]. The work presented in this paper focuses on these two dominating factors, while also considering the conductivities of the individual copper nanoparticles and the surrounding air medium.

In the μ-SLS process, the morphology of the nanoparticle powder bed and subsequently the effective thermal conductivity within the bed evolves as sintering takes place. As densification and thermal conductivity increase during sintering, the ability for heat to spread through the system is improved. A functional relationship between density and effective thermal conductivity is therefore needed to quantify heat spread and predict part shapes for the μ-SLS process.

Background

The modeling and simulation of the heat transfer characteristics and material properties of particle systems have been actively researched for years. Two-particle models have been developed to predict heat transfer between adjacent particles in thermal contact using analytical and computational methods [1316]. These models make key assumptions regarding particle–particle contact (Hertzian contact theory with semi-infinite solid approximation [14], and thermal pipe contact approximation [16]) that allow for numerical estimations of heat transfer in the two-particle system. This two-particle conduction can be expanded to particle systems with more than two particles given the assumed independence of particle–particle contacts. However, this typically requires extensive computational resources when particle counts become large. These models also fail to account for the time-varying, nonspherical, nanoscale geometries present during the sintering of metal nanoparticles. These additional considerations must be accurately modeled to effectively predict heat transfer and resulting temperature distributions during the μ-SLS process.

Models capable of estimating the effective material properties and heat transfer characteristics of aerogels consisting of packed nanoparticles have recently been developed [17,18]. These models analyze sets of particles and consider the relative importance of particle–particle conduction and the heat flows between particles and the surrounding gas medium. Although these conductivity models provide insight into relevant heat transfer parameters in nanoparticle systems, they analyze insulative materials and packing structures unlike those found in μ-SLS systems. They also do not consider the geometric evolution of the nanoparticles that occurs during a temperature driven solid state diffusion process.

Work has also been developed that considers the varying sizes and packing structures of nanoparticles of various materials [1921]. These research efforts determine the relationship between the porosity of the nanoparticle structure and its effective thermal conductivity using both experimental and analytical techniques. Although these research efforts developed relationships between particle packing porosity and effective thermal conductivity, they focus on porosity values that are much smaller than those experienced by the nanoparticles in the μ-SLS process. These porosity-thermal conductivity relationships also lack a full representation of the geometric evolution of particle geometries during the early sintering and neck formation stage, where porosities are high and effective thermal conductivities continue to be dominated by particle–particle contact resistances.

Expanding upon more general conduction models, particle heat transfer simulations for SLS systems have been developed for error minimization and part prediction [2224]. These models analyze sets of discrete particles that represent sections of larger particle beds. Although these SLS models effectively simulating the heat transfer process within powder beads for the initial spherical arrangement of particles, they model a particle melting and solidification process. In the μ-SLS system, sintering is a solid-state diffusion process, where the particles do not undergo phase transformations. Therefore, as the particles diffuse into each other, the resulting geometries of the particle groups change significantly despite the lack of a phase change. This work builds upon previous research efforts [12,25] and expands their scope to incorporate the complex geometric changes experienced by the nanoparticle powders undergoing laser sintering and estimates the effective thermal conductivity of these particle packings as a function of sintering duration.

Modeling Methods

A finite element model (FEM) was developed in this work to determine the effective thermal conductivity of copper nanoparticles as a function of densification in a μ-SLS system. A steady-state heat transfer simulation was performed in ansys mechanical to obtain temperature information within a heated nanoparticle packing. The finite element method was selected due to its proven accuracy, setup simplicity, and its widespread use in AM applications [26,27]. Nanoparticle packings, composed of approximately 20–30 nanoparticles (each of the order of 10–100 nanometers in diameter), are analyzed as they transition from discrete spherical particles to a coalesced, sintered unit [24]. Computer aided design (CAD) models are used to represent the packed nanoparticle geometries during the finite element analysis. Surface and grain boundary diffusion is simulated using a phase field model (PFM) to track the morphological evolution of the particles and obtain accurate solid models of the nanoparticles throughout the sintering process [28]. The PFM tracks the evolution of a fixed number of equally spaced points within a defined bounding box. These points can either represent solid particles or the surrounding air. The initial number of particles in the packing, the radii of each particle, and the dimensions of the bounding box are variable inputs to the PFM. For a given particle set, the PFM simulates the entire sintering process and records snapshots of the intermediate geometries that occur during sintering. The captured snapshots provide data files containing a set of phase field variables for each point within the model's bounding box. For each geometric snapshot, the phase field variables determine which points correspond to solid particles, and which points correspond to the surrounding air. This allows for visualizations of the geometry of the nanoparticle packings at different intermediate points during sintering. More specifically, the data files provide information regarding the particle density, ρ, of the region corresponding to the point in space, as well as a number of η variables equal to the total number of initial particles in the system. Both variables range from 0 to 1 and have cutoffs that determine the geometric arrangement of the nanoparticles. The ρ cutoff determines the threshold between the solid and vapor phases for each point in space (0 representing vapor, and 1 representing solid), and the η variable cutoff determines which initial particle each point in space belongs to. This allows the user to track the evolution of grain boundaries between particles during sintering. This dataset is the basis for the subsequent conversion and FEM processes discussed in this paper.

Conversion From Phase Field Model Data to Computer Aided Design Models.

Although the PFM data files provide a comprehensive description of the particle geometries throughout the sintering process, the raw data is not a usable geometric input for a finite element solver. FEMs accept solid models and surface bodies in several forms. To use the PFM data in a Finite Element simulation, a conversion process must be established to transform the PFM output data files into usable solid models for FEM input. To accomplish this, a three step conversion process is developed. To maintain the accuracy of the PFM output data, this conversion process must provide input CAD files that match the original PFM particle detail. The conversion process must also work for all possible initial particle arrangements, variable particle numbers, and the different geometries that occur as particles diffuse into one another. The steps of the conversion process are as follows.

The first step of the conversion process converts raw PFM data into Cartesian point clouds. The PFM output data is filtered using the ρ cutoff to determine the solid phase, and the η variable cutoffs are used to divide the domain into individual particles. The resulting collection of filtered points is then run through a surface extraction algorithm to isolate the exterior points and convert the data into a point cloud format. The algorithm evaluates each solid point to see if any of the surrounding points are in the vapor phase and writes a new data file with the filtered surface set. An example of a resulting point cloud is shown in Fig. 1(a).

Fig. 1
The full conversion process: (a) nanoparticles represented by point clouds, (b) nanoparticles represented by a triangular surface mesh derived from point clouds, (c) nanoparticles represented by solid models derived from triangular surface meshes, (d) nanoparticle solid models resting on a glass substrate, (e) hollowed out air block, and (f) full assembly of particles, air, and glass
Fig. 1
The full conversion process: (a) nanoparticles represented by point clouds, (b) nanoparticles represented by a triangular surface mesh derived from point clouds, (c) nanoparticles represented by solid models derived from triangular surface meshes, (d) nanoparticle solid models resting on a glass substrate, (e) hollowed out air block, and (f) full assembly of particles, air, and glass
Close modal

The second conversion step takes the point cloud coordinate files as input and outputs reconstructed surface meshes for each nanoparticle. These surface meshes are hollow shells composed of small surface triangles and are constructed based on the arrangement of points in the input point cloud. This step in the conversion process is performed in meshlab, an open-source 3D meshing, and reconstruction software, where parameters are adjusted to optimize a tradeoff between surface mesh accuracy and the number and quality of the resulting triangle surface meshes. The major substeps performed in meshlab include computation of surface normals, surface reconstruction, surface decimation (a reduction in the number of node points that make up the triangular surface mesh), mesh repair, and mesh smoothing. This surface meshing step is often the most critical, as it determines the usability, complexity, and surface quality of the downstream solid parts and subsequent finite element meshes. An example of the resulting surface mesh is shown in Fig. 1(b).

The third step of the conversion process takes triangular surface meshes as inputs and outputs solid CAD models with defined grain boundaries between adjacent particles. Particle surface meshes are first filled out to form solid CAD models in freecad, an open-source cad software. The particle models are then checked for holes that result from initial flaws in the surface mesh reconstruction process. Once the particle models are free of defects, a Python script is executed to detect intersections between adjacent nanoparticles in the particle set. This script checks to see if any two particles have overlapping points in their original Cartesian point cloud. The information for each pair of overlapping particles is added to an intersection matrix for further processing. The intersection matrix controls a separate Python script with freecad integration that performs particle–particle Boolean cuts. If a small and large particle overlap, the overlap region is removed from the smaller particle to simulate the process of smaller particles diffusing into larger particles. This occurs for each set of intersecting particles recorded in the intersection matrix. These Boolean cuts ensure that the grain boundaries between adjacent nanoparticles in the particle set are correctly modeled and that the interface between any set of two particles in contact is a perfect fit. This step is critical, as the thermal contact resistance between adjacent nanoparticles in a particle set is one of the most important input parameters in subsequent finite element simulations. This step concludes with two additional operations using freecad Python commands. First, a rectangular block of air with dimensions equal to the PFM bounding box is inserted into the modeling space. This block represents the air surrounding the nanoparticles and is included to allow the FEM to correctly model particle–air–particle conduction. This model neglects natural convection effects due to the extremely small Raleigh numbers, and correspondingly small heat transfer coefficients present in this nanoscale system. In the absence of external forced convection, the surrounding air can be modeled as a solid, thermally conductive medium that remains stationary for the duration of the simulation. To ensure perfect contact between the solid air block and the nanoparticle CAD models, another Python script is used to subtract each particle in the set from the surrounding air block. This forms a hollowed-out block that perfectly envelops the nanoparticle packing. Finally, a second block is inserted beneath the nanoparticles to model the glass substrate. Examples of the completed CAD models are shown in Fig. 1.

Conversion From Phase Field Model Data to Computer Aided Design Models.

The purpose of converting the PFM output data files to solid CAD models is to prepare the nanoparticle models for FEM. With proper setup, the FEM will solve for a 3D temperature distribution within the nanoparticle models and the surrounding air block, which will be subsequently used to extract the effective thermal conductivity of the particle set. The FEM solver used in this work is configured to solve the heat conduction equation
ρcpTt=·(kT)+Q
(1)

where k, ρ, cp, and Q represent the thermal conductivity, the density, the specific heat capacity at constant pressure, and the internal volumetric heat generation of a given solid object. The terms on the right-hand side of Eq. (1) describe the second-order spatial temperature gradient in the object summed with the heat generation occurring within that object, Q. The term on the left-hand side of Eq. (1) represents the transient thermal response of the object given its material properties. To solve this second-order partial differential equation for the spatial temperature distribution within an object, the nanoparticle CAD models must first be transferred to the finite element workspace, where the finite element meshing algorithm will divide each of the bodies into discrete elements and nodes on which an approximate solution can be solved. Then, heat transfer boundary conditions must be specified. Finally, values for the thermal contact conductance between adjacent nanoparticles and at particle–air interfaces are determined and set. Mesh discretization, boundary condition application, and thermal contact conductance determination were the key focuses of the finite element setup.

Mesh Generation.

Once the nanoparticle models, the enveloping air block, and the glass substrate are brought into the finite element workspace, a multistep meshing procedure is implemented. Due to the nature of the conversion process, the resulting particle models have surface profiles that match their corresponding triangular surface meshes. These profiles are composed of many small triangular faces of varying sizes, adding complexity to the geometry of each nanoparticle. The surrounding air block envelops these particles with matching faces at the contact points and thus has a similar level of geometric complexity. To account for this, heavily refined meshes are generated both within solid volumes and on planar surfaces. The surface and body meshes are generated separately with unique values to ensure an optimal mesh. Within the particle models and the surrounding air block, a fine tetrahedral mesh is applied. This mesh is restricted to a region of influence near the contact surfaces between the particle models and the solid air block. Hollow spheres of influence are generated and imported into the model to accomplish this volumetric meshing procedure. Each of the hollow spheres of influence encompasses the surface of one of the nanoparticles in the particle set. Both the air and particle models are meshed with this method, allowing for a fine mesh near contact regions, and a coarser mesh near the particle centers. This reduces mesh complexity and reduces overall computational time while emphasizing the regions with fine geometric detail. Additionally, a fine surface mesh is applied to all particle- air contact surfaces. This allows the meshing algorithm to handle a large number of small triangular faces present in the model. Finally, the global mesh is instructed to automatically refine the mesh further in tight spaces to improve the success rate of the meshing algorithm and facilitate subsequent Finite Element modeling.

Boundary Condition Application.

After the mesh settings are finalized for all surfaces and solid bodies, boundary conditions representing heat transfer processes are applied to the particle model. Before importing the solid models into ansys, the edges of the nanoparticle group are removed by a vertical cutting plane to allow for the consistent application of boundary conditions to smooth, vertical faces irrespective of initial conditions or time-steps. The locations of the vertical cuts are chosen such that the fraction of the resulting flat surface area that contains particle faces approximates the volume fraction of the nanoparticles with respect to the bounding box. This boundary modification helps maintain a uniform heat flow across different particle geometries. A coordinate direction is chosen, and a constant temperature of 22 °C is applied to the flat particle faces corresponding to the minimum value of this coordinate direction. A constant heat flux of 104MWm2 is then applied to the flat particle faces corresponding to the maximum value of the chosen coordinate direction. All remaining boundary faces in the model are insulated to simulate one-dimensional heat flow through the nanoparticle set. Radiation was neglected in the model, as it has been found to be negligible compared to particle–particle conduction and particle–air–particle conduction. This is largely due to the relatively low average sintering temperature of roughly 700 K [29]. Additionally, due to small Rayleigh numbers at the single micrometer scale, natural convection is neglected in this model as well. Applied boundary conditions are depicted in Figs. 2(c) and 2(d).

Fig. 2
Three-dimensional finalized solid models (note the 180-deg rotation about the z-axis between images c and d): (a) model with enveloping air block, (b) model without air block, (c) constant surface heat flux boundary condition, and (d) constant surface temperature boundary condition
Fig. 2
Three-dimensional finalized solid models (note the 180-deg rotation about the z-axis between images c and d): (a) model with enveloping air block, (b) model without air block, (c) constant surface heat flux boundary condition, and (d) constant surface temperature boundary condition
Close modal

Thermal Contact Conductance.

The thermal contact conductance values used in the simulation determine the temperature jumps at the interfaces between adjacent bodies in particle–particle contact and in particle- air contact. The thermal contact conductance for particle–air contact is set to 105MWm2K to accurately model perfect wetting of the particle surfaces by the surrounding air. The resulting temperature jump at these interfaces should be close to zero to approximate ideal contact between an object and a surrounding fluid. The thermal contact conductance for particle–particle contact is more difficult to estimate. Multiple thermal contacts occur from one end of the nanoparticle packing to the other along the specified coordinate axis. This is especially true for earlier time-steps, where more particle–particle contacts exist along the primary axis of heat flow. A value of 20 MWm2K was chosen for the particle–particle thermal contact conductance based on estimations from Refs. [29] and [30]. The former used the discrete element method (DEM) to compare relevant heat flows between spherical metal nanoparticles and selected the value of thermal contact conductance that corresponded to the transition point at which the temperature gradients within individual particles became negligible compared to the temperature drop across particle–particle contacts. The latter measured the thermal contact conductance between deposited metal films (Al, Ag, Sn, Zn, and In) and an underlying copper film using transient thermoreflectance methods. The values determined from Ref. [30] ranged from 8 to 95 MWm2K for the set of tested metals. The value obtained from Ref. [29], 20 MWm2K was determined for the same nanoparticle sintering system as covered in this work. Therefore, the value from Ref. [29] was used in this paper for consistency. Given this chosen value for thermal contact conductance between particles, particle–particle contact conductance, and particle–particle contact area are the most important heat transfer parameters in the system. Although the bulk thermal conductivity of copper changes as a function of the temperature and the size of the individual nanoparticles, these changes are negligible when compared to changes in the contact regions between adjacent nanoparticles. Therefore, the bulk thermal conductivity of copper (400 Wm*K) is used in this work.

A summary of the complete finite element analysis (FEA) process is shown in Fig. 3. Figure 3 shows the full finite element setup for three unique time-steps during a thermal simulation for a given particle bed. For each time-step, the heat flux and constant temperature boundary conditions are shown, along with the full 3D temperature distribution that results from the steady-state simulation. The first image in Fig. 3 shows the simulation just after the onset of sintering (and subsequent neck formation). The second image shows the particle group at a more advanced sintering stage. Finally, the third image shows the particle group at a late sintering stage, where simulated sintering has progressed well beyond the degree of sintering experienced by the particle beds in the real system. The full simulation is run the same way for all depicted time-steps, and the resulting temperature distributions are used to calculate the effective thermal conductivity at that time-step.

Fig. 3
Summary of finite element setup and sample results. Temperature scales are in units of Degrees Celsius.
Fig. 3
Summary of finite element setup and sample results. Temperature scales are in units of Degrees Celsius.
Close modal

Results and Discussion

The densification of ten one-by-one micrometer beds consisting of anywhere between 21 and 43 particles was simulated using the PFM setup described earlier. The simulation boxes containing the particle sets were 104 by 104 pixels in the xy plane and ranged from 56 to 86 pixels along the z-axis. PFM output data files representing nanoparticle geometries for each of the ten sintered beds were saved at 35 time-steps along the sintering timeline. Each output file was converted into a set of particle CAD models with well-defined contact regions. Each of these particle sets was then imported into an ansys mechanical FEM workspace where mesh generation and boundary condition application were performed before ultimately solving each model for the resulting steady-state temperature distribution within the nanoparticle system.

Upon completion of the FEM simulation for each time-step, the temperature distributions were analyzed to calculate the effective conductivity of the particle set. To calculate the effective thermal conductivity of a particle set, the steady-state, one-dimensional approximation to Fourier's law of heat conduction was rearranged and solved for the thermal conductivity parameter. This rearranged thermal conduction equation is shown in the following equation:
keff=LT2T1·q
(2)

where keff is the effective thermal conductivity, T2 and T1 are two independent temperature values located at different locations along the chosen axis, L is the distance between T2 and T1 along the same axis, and q is the heat flux magnitude at the boundary. The temperature distributions produced by the FEM solver provide temperature data corresponding to every node on every element within the finite element model. To apply Eq. (2) across the full particle set and obtain an effective thermal conductivity value for the given set, the output temperature information must be reduced. The temperature distributions were first filtered to only include the external faces corresponding to the minimum and maximum values of the x coordinates. The temperature datasets were again filtered to only include temperature data from particle nodes. This filtration only included particle faces corresponding to noninsulative boundary conditions and allowed for a more accurate temperature difference calculation across the set of particles. Finally, the filtered temperature data were averaged at both the maximum and minimum x-coordinate locations to obtain a single representative temperature value for each location. These average temperatures correspond to T2 and T1 from Eq. (2).

Once the output data was simplified to averaged temperature values, the effective thermal conductivity of the particle group could be determined. Equation (2) provides a first-pass estimate for this value, but it is unable to account for nonuniform surface areas in the primary direction of heat flow, and thus requires modification. The PFM data input is inherently stochastic, meaning the nanoparticle geometries present in each particle set are different. Therefore, the surface areas that are exposed to the constant surface temperature and the constant heat flux boundary conditions can vary considerably between beds and between time-steps for a given bed. Additionally, the cross-sectional area perpendicular to the axis of heat flow changes with location in a unique way for each bed. A correction factor is added to Eq. (2) to account for these sources of variation and provide more comparable values for the effective thermal conductivity of the particle set. The modified formulation is shown in the following equation:
keff=(AsVf)·LT2T1·q
(3)

Vf represents the volume fraction of the PFM bounding box taken up by the nanoparticle models. As represents the ratio of nanoparticle surface area at the heat flux boundary to overall surface area on the same boundary face. T2 and T1Correspond to the average temperatures on the filtered external particle faces, q represents the applied heat flux of 104MWm2K and L represents the distance between T2 and T1. The fraction shown in Eq. (3) corrects for differences between the surface area where the heat flux and temperature boundary conditions were applied and the average surface area along the axis of heat flow through the particle set. The full calculation was performed on ten beds for 35 sampled time-steps each. In addition to the thermal conductivity calculation, the density of each nanoparticle packing was calculated using a Python script. This script samples each bed and calculates the ratio of particle volume to the volume of a bounding box. Given the small number of nanoparticles, the bounding box used to determine the total volume was chosen to be a subset of the full PFM bounding box. This reduces the impact of boundary effects on the densification results by focusing on the interior particles that better represent the majority of particles in a much larger bed. The results of the simulations for the ten randomly generated particle beds were aggregated into a single dataset, shown in Fig. 4.

Fig. 4
Thermal conductivity versus density with linear fit line and 95% prediction bounds
Fig. 4
Thermal conductivity versus density with linear fit line and 95% prediction bounds
Close modal

Figure 4 includes the output data from the thermal simulations and density calculations along with a linear best fit line and prediction bounds. All simulations were run for a fixed number of time-steps, although each bed reached a different density value by the end of the PFM simulation. Many of these density values were substantially larger than the density values achievable during laser sintering (this does not include the annealing process that occurs postsintering). Therefore, thermal conductivity data points were excluded if their density values exceeded 5000 kgm3, as this is a reasonable upper bound for achievable densification using the μ-SLS system without postprocess annealing. The equation of the linear fit line is also shown in Fig. 4. Within this equation, keff represents the effective thermal conductivity of the particle packing and ρ represents the density of the particle packing. As inferred from the fit line, a 1000 kgm3 increase in the packing density increases the effective thermal conductivity of the particle packing by 7.1 Wm*K. The 95% prediction interval provides a confidence band for future predictions using the linear fit shown in Fig. 4. The large range for prediction results from the considerable uncertainty present for single-micron particle beds. In smaller particle systems, boundary effects have considerably more impact on the simulation results. To reduce this uncertainty, a section of a larger particle bed can be used at the cost of a significant increase in computational expense. The prediction interval in Fig. 4 will be used to quantify the uncertainty in downstream sintering models that incorporate the results from this work.

Due to the random nature of the particle bed generation process, beds have different initial densities. To account for this, the results shown in Fig. 4 are plotted as a function of densification percentage in Fig. 5. Densification percentage refers to the percent increase in density from the density value calculated before the onset of sintering. The equation of the linear fit line is also shown in Fig. 5. Within this equation, keff represents the effective thermal conductivity of the particle packing and ρ% represents the densification percentage experienced by the particle packing. For every 10% increase in the densification percentage, the effective thermal conductivity is increased by 2.6 Wm*K. Additionally, the fit equation predicts a value of 10.74 Wm*K for the effective thermal conductivity when densification percentage is near 0 (just after the onset of sintering.

Fig. 5
Thermal conductivity versus densification percentage with linear fit line and 95% prediction bounds
Fig. 5
Thermal conductivity versus densification percentage with linear fit line and 95% prediction bounds
Close modal

It is noted that these results describe early necking and initial sintering, but do not include the initial thermal conductivity of the nanoparticle packing before sintering occurs. The values for the initial packing condition were recorded separately. This allows for future simulations to have access to a thermal conductivity value for unsintered copper nanoparticles, in addition to a range of thermal conductivity values for copper nanoparticles as they sinter. The average effective thermal conductivities of the particle group before and after the onset of sintering are 1.05±0.26Wm*K and 10.45±1.37Wm*K, respectively. The associated densities before and after sintering are 3646.87±138.52kgm3 and 3648.15±137.73 kgm3, respectively. The increases in average effective conductivities and densities that result from the onset of sintering are 9.42 Wm*K (798.35% change from initial thermal conductivity) and 1.28 kgm3 (.035% change from initial density), respectively. Despite a small average increase in the density during initial sintering, the thermal conductivity increases nearly by a factor of 10. This is the result of the early neck formation that occurs when adjacent nanoparticles surpass the sintering threshold. During this initial necking phase, the contact regions between nanoparticles grow considerably, opening larger paths for heat to conduct between the particles. At the onset of this process, the formation of these heat pathways results in a very small change in the overall geometry, and therefore the density, of the packing. The laser sintering that occurs during the μ-SLS process is heavily impacted by the initial necking phase, so it is critical that the rapid increase in thermal conductivity during this phase is accurately captured.

The results developed in this work describe the relationship between copper nanoparticle densification and the corresponding increase in the effective thermal conductivity of a nanoparticle group. This relationship has been integrated into a preliminary part-scale, transient finite element model capable of predicting the thermal evolution of a powder bed during sintering. Each copper element of the part-scale model has its thermal conductivity updated at a series of time-steps based on the temperature history of the element.

Moving forward, this part-scale model will be validated against experimental results using a combination of IR camera images, profilometer profiles, and experimentally measured electrical properties. This part scale model validation will also allow for validation of the model input parameters, most notably the thermal conductivity curves that result from the nanoparticle model discussed in this work. While it is not feasible to directly measure the effective thermal conductivity of a 1 μm by 1 μm nanoparticle group, the heat spread in the part-scale model is primarily determined by the evolving thermal conductivity in the particle bed. The degree of alignment between the part-scale model results and experimental studies will determine the validity of the conductivity-densification relationship as a predictor of the true thermal conductivity evolution during laser sintering.

Additionally, this multiscale modeling effort will be used to input laser parameters and predict the resulting sintered part shapes. This predictive capability will allow for model-based control of laser parameters to achieve optimized part shapes with minimized heat affected zones. Additionally, larger, 2 μm by 2 μm beds can be used to reduce the uncertainty of the resulting linear fit line by reducing the impact of boundary effects on the results. By expanding the size of the bed, the thermal model can take a section from the larger bed and run the heat transfer analysis on an interior section, avoiding particles along the edges.

Conclusion

In this paper, a FEM and a fully scalable data conversion framework are developed to estimate the effective thermal conductivity of smalls sets of copper nanoparticles. This will facilitate the prediction of heat affected zones during μ-SLS. PFM data was converted into solid particle CAD models, and an automated steady-state heat transfer simulation was used to calculate temperature profiles and effective thermal conductivities of nanoparticle packings. A steady-state, 1D thermal conduction approximation was used alongside an area/volume fraction correction factor extracting thermal conductivity information out of resulting temperature data. The results show a rapid increase in the thermal conductivity during early particle necking and a reduced rate of change for later time-steps. The results also show that there is significant uncertainty in the thermal conductivity of the beds due to the variances in nanoparticle configurations throughout the bed but that there is a clear trend that as be particle bed starts to sinter together and densify, the thermal conductivity of the bed increases. Overall, the thermal conductivity of the bed almost doubles from initial sintering to the maximum densification that is typically achieved in the μ-SLS process.

Funding Data

  • NXP Semiconductors.

  • National Science Foundation (Grant No. 1728313; Funder ID: 10.13039/100000001).

Nomenclature

As =

surface area corresponding to the constant heat flux boundary condition

cp =

nanoparticle specific heat capacity

k =

nanoparticle thermal conductivity

keff =

effective thermal conductivity of the nanoparticle packing

L =

distance between locations corresponding to T2 and T1

q =

heat flux magnitude at the boundary

Q =

volumetric heat generation

T1 =

average temperatures of the external particle faces aligned with the positive y-axis

T2 =

average temperature of the external particle faces aligned with the negative y-axis

Vf =

volume fraction of the PFM bounding box taken up by the nanoparticle models

η =

variable detailing particle contribution of each point in space (PFM)

ρ =

density corresponding to each point in space (PFM)

References

1.
Frazier
,
W. E.
,
2014
, “
Metal Additive Manufacturing: A Review
,”
J. Mater. Eng. Perform.
,
23
, pp.
1917
1928
.10.1007/s11665-014-0958-z
2.
Sager
,
B.
, and
Rosen
,
D.
,
2002
,
Stereolithography Process Resolution
,
Georgia Institute of Technology
,
Atlanta, GA
.
3.
Roy
,
N. K.
,
Behera
,
D.
,
Dibua
,
O. G.
,
Foong
,
C. S.
, and
Cullinan
,
M.
,
2019
, “
A Novel Microscale Selective Laser Sintering (μ-SLS) Process for the Fabrication of Microelectronic Parts
,”
Microsyst. Nanoeng. 5
,
64
.
4.
Roy
,
N. K.
,
Foong
,
C. S.
, and
Cullinan
,
M. A.
,
2016
, “
Design of a Micro-Scale Selective Laser Sintering System
,”
Annual International Solid Freeform Fabrication Symposium
, Austin, TX, Aug. 9, pp.
1495
1508
.https://www.researchgate.net/publication/318528295_Design_of_a_Microscale_Selective_Laser_Sintering_System
5.
Roy
,
N.
,
Yuksel
,
A.
, and
Cullinan
,
M.
,
2016
, “
Design and Modeling of a Microscale Selective Laser Sintering System
,”
ASME
Paper No. MSEC2016-8569.10.1115/MSEC2016-8569
6.
Roy
,
N.
,
Dibua
,
O.
,
Foong
,
C. S.
, and
Cullinan
,
M.
,
2017
, “
Preliminary Results on the Fabrication of Interconnect Structures Using Microscale Selective Laser Sintering
,”
ASME
Paper No. IPACK2017-74173.10.1115/IPACK2017-74173
7.
Nelson
,
C.
,
McAlea
,
K.
, and
Gray
,
D.
,
1995
, “
Improvements in SLS Part Accuracy
,”
Annual International Solid Freeform Fabrication Symposium
, Austin, TX, 1995, pp.
159
169
.
8.
Dong
,
L.
,
Makradi
,
A.
,
Ahzi
,
S.
, and
Remond
,
Y.
,
2009
, “
Three Dimensional Transient Finite Element Analysis of the Selective Laser Sintering Process
,”
J. Mater. Process. Technol.
,
209
(
2
), pp.
700
706
.10.1016/j.jmatprotec.2008.02.040
9.
Moser
,
D.
,
Cullinan
,
M.
, and
Murthy
,
J.
,
2016
, “
Particle-Scale Melt Modeling of the Selective Laser Melting Process
,”
Annual International Solid Freeform Fabrication Symposium
, Austin, TX, 2016, pp.
247
256
.
10.
Moser
,
D.
,
Fish
,
S.
,
Beaman
,
J.
, and
Murthy
,
J.
,
2014
, “
Multi-Layer Computational Modeling of Selective Laser Sintering Processes
,”
ASME
Paper No. IMECE2014-37535.10.1115/IMECE2014-37535
11.
Dibua
,
O. G.
,
Yuksela
,
A.
,
Roya
,
N. K.
,
Foong
,
C. S.
, and
Cullinan
,
M.
,
2018
, “
Nanoparticle Sintering Model, Simulation and Calibration Against Experimental Data
,”
ASME
Paper No. MSEC2018-6383.10.1115/MSEC2018-6383
12.
Yuksel
,
A.
,
Yu
,
E. T.
,
Cullinan
,
M.
, and
Murthy
,
J.
,
2020
, “
Investigation of Heat Transfer Modes in Plasmonic Nanoparticles
,”
Int. J. Heat Mass Transfer
,
156
, p.
119869
.10.1016/j.ijheatmasstransfer.2020.119869
13.
Batchelor
,
G. K.
, and
O'Brien
,
R. W.
,
1977
, “
Thermal or Electrical Conduction Through a Granular Material
,”
Proc. R. Soc. London Ser. A Math. Phys. Sci.
,
355
, pp.
313
333
.https://doi.org/10.1098/rspa.1977.0100
14.
Sun
,
J.
, and
Chen
,
M. M.
,
1988
, “
A Theoretical Analysis of Heat Transfer Due to Particle Impact
,”
Int. J. Heat Mass Transfer
,
31
(
5
), pp.
969
975
.10.1016/0017-9310(88)90085-3
15.
Zhou
,
J. H.
,
Yu
,
A. B.
, and
Horio
,
M.
,
2008
, “
Finite Element Modeling of the Transient Heat 150 Conduction Between Colliding Particles
,”
Chem. Eng. J.
,
139
(
3
), pp.
510
516
.10.1016/j.cej.2007.08.024
16.
Shimizu
,
Y.
,
2006
, “
Three-Dimensional Simulation Using Fixed Coarse-Grid Thermal-Fluid Scheme and Conduction Heat Transfer Scheme in Distinct Element Method
,”
Powder Technol.
,
165
(
3
), pp.
140
152
.10.1016/j.powtec.2006.04.003
17.
Zhao
,
J. J.
,
Duan
,
Y. Y.
,
Wang
,
X. D.
, and
Wang
,
B. X.
,
2012
, “
A 3-D Numerical Heat Transfer Model for Silica Aerogels Based on the Porous Secondary Nanoparticle Aggregate Structure
,”
J. Non-Cryst. Solids
,
358
(
10
), pp.
1287
1297
.10.1016/j.jnoncrysol.2012.02.035
18.
Guo
,
J. F.
, and
Tang
,
G. H.
,
2019
, “
A Theoretical Model for Gas- Contributed Thermal Conductivity in Nanoporous Aerogels
,”
Int. J. Heat Mass Transfer
,
137
, pp.
64
73
.10.1016/j.ijheatmasstransfer.2019.03.106
19.
Wu
,
D.
, and
Huang
,
C.
,
2020
, “
Thermal Conductivity Study of SiC Nanoparticle Beds for Thermal Insulation Applications
,”
Phys. E
,
118
, p.
113970
.10.1016/j.physe.2020.113970
20.
Lin
,
Z.
,
Huang
,
C. L.
,
Zhen
,
W. K.
, and
Huang
,
Z.
,
2017
, “
Enhanced Thermal Conductivity of Metallic Nanoparticle Packed Bed by Sintering Treatment
,”
Appl. Therm. Eng.
,
119
, pp.
425
429
.10.1016/j.applthermaleng.2017.03.087
21.
Qin
,
F.
,
Hu
,
Y.
,
Dai
,
Y.
,
An
,
T.
, and
Chen
,
P.
,
2020
, “
Evaluation of Thermal Conductivity for Sintered Silver Considering Aging Effect With Microstructure Based Model
,”
Microelectron. Reliab.
,
108
, p.
113633
.10.1016/j.microrel.2020.113633
22.
Ganeriwala
,
R.
, and
Zohdi
,
T. I.
,
2016
, “
A Coupled Discrete Element- Finite Difference Model of Selective Laser Sintering
,”
Gran. Matter
,
18
(
2
), p.
21
.10.1007/s10035-016-0626-0
23.
Moser
,
D.
,
Pannala
,
S.
, and
Murthy
,
J.
,
2016
, “
Computation of Effective Thermal Conductivity of Powders for Selective Laser Sintering Simulation
,”
ASME J. Heat Transfer-Trans. ASME
,
138
(
8
), p.
082002
.10.1115/1.4033351
24.
Moser
,
D.
,
Yuksel
,
A.
,
Cullinan
,
M.
, and
Murthy
,
J.
,
2018
, “
Use of Detailed Particle Melt Modeling to Calculate Effective Melt Properties for Powders
,”
ASME J. Heat Transfer-Trans. ASME
,
140
(
5
), p.
052301
.10.1115/1.4038423
25.
Yuksel
,
A.
, and
Cullinan
,
M.
,
2016
, “
Modeling of Nanoparticle Agglomeration and Powder Bed Formation in Microscale Selective Laser Sintering Systems
,”
J. Addit. Manuf.
,
12
, pp.
204
215
.10.1016/j.addma.2016.07.002
26.
Li
,
Y.
, and
Gu
,
D.
,
2014
, “
Parametric Analysis of Thermal Behavior During Selective Laser Melting Additive Manufacturing of Aluminum Powder
,”
Mater. Des.
,
63
, pp.
856
867
.10.1016/j.matdes.2014.07.006
27.
Ancellotti
,
S.
,
Fontanari
,
V.
,
Molinari
,
A.
,
Iacob
,
E.
,
Bellutti
,
P.
,
Luchin
,
V.
,
Zappini
,
G.
, and
Benedetti
,
M.
,
2019
, “
Numerical/Experimental Strategies to Infer Enhanced Liquid Thermal Conductivity and Roughness in Laser Powder-Bed Fusion Processes
,”
J. Addit. Manuf.
,
27
, pp.
552
564
.10.1016/j.addma.2019.04.007
28.
Dibua
,
O.
,
Yuksel
,
A.
,
Roy
,
N.
,
Foong
,
C. S.
, and
Cullinan
,
M.
,
2018
, “
Nanoparticle Sintering Model, Simulation and Calibration Against Experimental Data
,”
ASME J. Micro Nanomanuf.
,
6
, p.
041004
.
29.
Yuksel
,
A.
,
Yu
,
E. T.
,
Cullinan
,
M.
, and
Murthy
,
J.
,
2020
, “
Thermal Transport in Nanoparticle Packings Under Laser Irradiation
,”
ASME J. Heat Transfer-Trans. ASME
,
142
, p.
032501
.10.1115/1.4045731
30.
Zheng
,
H.
, and
Jagannadham
,
K.
,
2014
, “
Interface Thermal Conductance Between Metal Films and Copper
,”
Miner., Met. Mater. Soc. ASM Int.
,
58
(
6
), pp.
67
74
.