International Journal of Thermal Sciences 148 (2020) 106135
Contents lists available at ScienceDirect
International Journal of Thermal Sciences journal homepage: http://www.elsevier.com/locate/ijts
A hybrid approach for predicting the effective thermal conductivity of sintered porous materials John Polansky a, *, Nick Jeffers b, Jeff Punch a a b
CONNECT Centre, Stokes Labs, Bernal Institute, University of Limerick, Limerick, Ireland Nokia Bell Labs, Blanchardstown, Dublin 15, Ireland
A R T I C L E I N F O
A B S T R A C T
Keywords: Effective thermal conductivity Porous media Sintered porous materials Granular assemblies
Estimating the effective conductivity of porous media has been widely based on homogenisation techniques with a limited number of discrete analyses. Though bulk estimates may be sufficient within prescribed bounds and applications, localised properties are not readily available. This work employs resistive network analytics in conjunction with 2D and 3D porous media modelling, to estimate the effective conductivity through local interrogation of porous sintered bronze disks. The resulting hybrid method uses rigid body simulations to generate the porous network and analytical mathematics to evaluate the effective conductivity. The analytics are herein verified and validated for both 2D and 3D cases. The 2D analysis used printed circuit boards to validate the analytics and produced an experimental mean and standard error of 308.7 � 12.6 μΩ for networks having resistors ranging from 10 to 1000 Ω. A 3D analysis of sintered porous bronze disks was performed for both electrical and thermal conductivities. A set of simulated samples consistent with the physical specimens were generated using an open source 3D creation software and tested in manners consistent with experiments. The electrical tests across the diameter of the bronze samples yielded equivalent experimental and theoretical re sistances of 1.93 � 0.42 mΩ and 2.03 � 0.17 mΩ respectively. A Xenon flash analysis across the thickness of the samples produced an experimental and theoretical Effective Thermal Conductivity (ETC) of 21.9 � 3.4 W/mK and 18.1 � 0.3 W/mK respectively. This method was able to successfully estimate the ETC within 17% of the measured mean, while traditional contact-based empirical correlations were found to underpredict the ETC by as much as 58%. The models align well with both experimental studies in both thickness and diametric analysis. Extension of the models to samples of differing particle sizes revealed the ETC to increase with particle size.
1. Introduction The conductive properties of porous media are important to a range of engineering applications such as thermal management [1], bio materials [2], fuel cells [3], composites [4], and pebble bed reactors [5] among others. How heat migrates through porous media is dictated by the complexity of the network and the interplay of heat transfer modes. This work is focused on the analysis of porous networks to estimate the Effective Thermal Conductivity (ETC) of sintered porous bronze, how it may be measured experimentally, and its relation to models and analytics. The conductive properties of porous media have been estimated using variety of techniques including but not limited to volume aver aging (based on the constituent thermal conductivities), network con nectivity and contact resistances [6]. Although versatile,
homogenisation techniques are not well suited to the development of engineered porous media. Rather, discrete approaches are sought to capture the local properties otherwise ignored in traditional ETC methods. Discrete approaches like FEM and FEA allow for the capture of local thermal properties for a given network architecture. Though robust, such methods can be computationally expensive due to meshing and iterative convergence, thereby making them potentially unsuitable for statistical or parametric analysis. Early porous media models characterised the upper and lower bounds of equivalent networks using resistance networks based on their respective series and parallel circuits [7]. These two cases clearly demonstrate the importance of the formed network as they have the same porosity through a 90� rotation, yet define the upper and lower ETC bounds. The parallel circuit model was enhanced by Liang et al. [8] and bridged the two limits using calibrated coefficients. Alternatively, a hybrid mapping approach was developed by Wang et al. [9] to analyse
* Corresponding author. E-mail address:
[email protected] (J. Polansky). https://doi.org/10.1016/j.ijthermalsci.2019.106135 Received 9 July 2019; Received in revised form 4 October 2019; Accepted 9 October 2019 Available online 11 November 2019 1290-0729/© 2019 Elsevier Masson SAS. All rights reserved.
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
contact resistances on each particle to make computations tractable. The simulation time was noted to take seconds while FEM took several hours. Of the simulated samples, the ETC was found to increase with compaction, particle size and solid volume fraction. These finding sug gest that the thermal performance of porous materials is, at least in part, driven by the formed network and may be a weak function of porosity if at all. Given the importance of the formed network, detailed studies of physical specimens are critical for validation. Grandjean et al. [20] analysed physical tin oxide porous structures using a scanning electron microscope and converted them into 2D mesh structures for use in finite element analysis. The simulation results were accurate for close packing arrangements but deviated for pore volume fractions greater than 40%. Similarly, Roussel et al. [3] used GeoDict software to discretize the sample slices into voxels and create a three-dimensional structure. The minimum voxel ratio per particle results in an increase in computations as a multiple of the particle count. This approach was noted to take several minutes of computational time depending on the voxel resolu tion yet proved inferior when compared to resistor-based approaches for small contact resistances [3]. Depending on how porous media is formed, it is liable to impact the overall thermal performance of the material. Though porous media may take many different physical forms, this work is focused on packed bed forms with the particles being spherical in nature. For such cases, simple porous media structures can be described as a variation of the circle packing problem. Fortunately, the circle packing problem has been widely studied for both mono- and unequally-sized circles. A brief introduction is provided here as it applies to sintered bronze; additional material can be found in Refs. [21,22]. In the case of two-dimensional packing problems, Huang et al. [23] developed deterministic algorithms to pack a prescribed rectangle with a predefined set of unequal circles. Alternatively, George et al. [24] used a heuristic approach to pack the space with the greatest number of cir cles, using a priori knowledge of the circle set to selectively place the circles for the greatest packing density. While this achieves maximal packing structures, the packs were noted to be geometrically stable but physically unstable [24]. The method used to mimic porous media manufacturing must be given some consideration as to its suitability to problem. A form more consistent with manufactured porous media was developed by Visscher and Bolsterli [25] who used Monte Carlo simu lations to fill a rectangular domain subject to gravity and domain shaking to promote settling. The packing density for 2D assemblies was limited by the lack of interstitial passages for smaller circles. The Visscher and Bolsterli [25] study included 3D simulations using mono-sized particles. In the more general case of irregularly shaped particles, Gan et al. [13] investigated the ETC of ellipsoid particles in a packed bed. The particles were noted to arrange themselves such that the major axis was normal to gravity as one might expect. The resulting packed beds produced a higher ETC than spherical packed beds, which was attributed to the change in coordination number and contact area [13]. Further investigation demonstrated that the interplay between the conduction network and the fluid network changed the ETC, while particle size had no significant effect for similar packing arrangements. With the formed samples, comes a need for local analytical interro gation of the samples. Of the methods utilised in prior works, point-topoint resistance methods show promise regarding accuracy and speed. Numerous point-to-point resistance calculation methods have been studied as summarised by Vos [26]. Many are limited to networks that are small, symmetric or infinite. While there have been advancements made for infinite lattices [27], they provided no tractable means of solving large finite networks like those of porous media. This gap was addressed by Wu [27], who developed a closed form analytical approach to calculating the two-point resistance for finite networks. The network is represented as a weighted, undirected graph with the nodes and edges all encapsulated in a weighted Laplace matrix. The derivation yielded a single expression, Eq. (1), which was verified for regular (structured)
Nomenclature Latin A d D h i k m n R V
Area [m2] Sample diameter [m] Particle diameter [m] Sample height [m] Summand index Thermal conductivity [W/mK] Sample mass [kg] Particle count Resistance [K/W or Ω] Volume [m3]
Subscripts eq Equivalent exp Exposed Greek
α β λ
π ρ ϕ
ψ
Reference node Reference node Eigenvalue Mathematical Const Material density [kg/m3] Porosity Eigenvector
random porous media samples as their equivalent structured networks. In contrast to the one-to-one mapping approach, bulk based approaches were developed to describe the ETC of randomly structured samples, namely Effective Medium Theory (EMT) and Maxwell-Eucken. These models consider stochastic networks of multi-component mixtures and account for percolation effects [7]. Such simplified network models are well suited to investigating materials such as carbon fibre reinforced polymers [4], patterned structures [10] and burn char residues [11]. The resistive network limits and interstitial models are graphically illustrated in Table 1, and more extensively described in Refs. [7,11]. The EMT and Maxwell-Eucken models can capture shifts in bulk conductivity due to changes in the network, however they are princi pally focused on the global properties rather than local. To address this, researchers looked to the Discrete Element Method (DEM) to replicate the discrete nature of the individual particles of porous media. The foundations of DEM are based on soil and rock mechanics which was extended to porous media analysis as they share packing structures, formed networks, interfacial junctions and surface contacts [9,12–17]. Such assemblies are common to a range of problems including bio ag gregates [2], soils [18], and pebble bed reactors [5] among others. Yun et al. [19] used DEM to simulate 3D particulate assemblies of Ottawa 20–30 sand with variable external forces. Changes in the applied forces were found to affect the network and contact areas through which heat was conducted [19]. Each sample simulated was found to have different conductance networks and corresponding thermal properties while maintaining the same porosity. Similarly, Zhang et al. [12] simulated granular assemblies using DEM with imposed isothermal Table 1 Range of ETC models for heat loads travelling from top to bottom. Parallel
MaxwellEucken
EMT
MaxwellEucken
Series
2
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
lattices with Neumann and Dirichlet boundary conditions in 1D [27], 2D €bius strip [27], cobweb [29], fan [30]) and (rectangular lattice [28], Mo 3D (cuboid [27], Klein bottle [27], embedded globe [31]). The weighting (thermal/electrical resistance) of the network was limited to a singular value in each of the principle cartesian directions. � �2 n � X ψ ia ψ iβ � 1 ¼ (1) Req λi i¼2
diameter of 22 � 0.05 mm and a peak height of 3 � 0.05 mm, as illus trated in Fig. 1 (A). The bottom and perimeter of the samples are consistent with the mould shape and dimensions; however, the top has an uneven structure as illustrated in Fig. 1 (B). The changes in surface height were measured using the Keyence 3D microscope and software using 4 sampling lines across the top of each of the 15 samples tested. The sample tops were found to have a height range with a mean and standard deviation of 2.51 � 0.23 mm respectively. The geometric properties of the samples limit experimental interrogation to perimeter contact measurements and optical for the top and bottom. As such, electrical probe-based measurements were limited to the perimeter and thermal flash method was performed across the top and bottom. An accurate analysis of the samples required knowledge of the ma terial alloy. For this, a Scanning Electron Microscope (SEM) equipped with an Energy Dispersive X-ray (EDX) was used to determine the composition of three specimens. Two sets of EDX measurements were performed for each of the three specimens, at the bulk and neck junc tions between particles. The material was determined to be phosphor bronze: UNS C52400 (88% Cu, 11% Sn) based on the EDX analysis, as exemplified in Fig. 2 (A). The phosphorous was noted to have a high concentration in the necks formed between the particles. Having iden tified the bronze composition, the electrical and thermal conductivity values taken to be 6.37 � 106 1/Ωm and 50 W/mK respectively [33]. In addition to the material characterisation, the SEM provided a means of estimating the size range of the particles and sintered junc tions, as illustrated in Fig. 2 (B–C). From the SEM analysis, the particles were estimated to have a mean diameter of 0.75 mm. The diameters were assumed to have a normal distribution with a standard deviation of 0.042 mm. Finally, the formed sinter necks between particles were estimated to have a mean diameter of 0.3D (0.225 mm) with a normal distribution having a standard deviation of 0.013 mm. With the properties of the porous sintered bronze samples charac terised, a set of simulations were created to estimate the conductive properties of the materials. The analyses begin with a 2D approach to verify and validate the analytics for non-orthogonal networks and demonstrate its suitability for 3D analysis.
The formulation of Eq. (1) has a singularity stemming from the Laplacian which was overcome by Izmailian et al. [29] through their use of the second minor of the networks Laplace matrix. This approach was noted to address the unique cases of the cobweb and fan networks, with a modified version capturing the globe network using a single summand. A later study of point-to-point resistances by Vos [26] was able to demonstrate that the resistance could be obtained using a pseudo-inverse Laplace. Through these derivations, the point-to-point resistance was found to be a function of the spanning trees through the two points and the total number of spanning trees [26]. Subsequent derivations were able to relate global changes in resistance to changes in resistance of a single edge. The causal link between the local and global resistance changes for non-symmetric Laplace matrices was studied by Cernanova et al. [32], in which they theorised its application to elec trical fault analysis of integrated circuits. From the outlined literature, it is apparent that the study of heter ogenous porous material properties is multifaceted problem spanning a variety of topics. Of the studies, some have demonstrated that tradi tional methods such as homogenisation are unable to capture the true nature of heterogeneous materials. Although discrete methods of anal ysis (FEA, FEM, DEM, etc.) have proven effective, they can be compu tationally intensive due to meshing requirements and further compounded by increasing particle counts. By comparison, resistive network methods are limited to sparse square matrices of the same size as the number of particles in the sample. Given the size and nature of the sintered bronze porous media (typically 3000 particles and 8500 sin tered junctions), an implementation of Wu’s [27] resistor analytics is applied to 2D and 3D porous media simulations. The following sections outline the physical sample properties, 2D electrical verification and validation, 3D electrical and thermal verification and validation with variance in sample composition.
3. Methodology 3.1. Verification and validation - 2D
2. Sample properties
To generate a suitable set of 2D porous samples, a constructive al gorithm was developed and implemented in MATLAB. The twodimensional particles were based on idealistic circles stacked within a prescribed rectangular set of bounds. The boundaries of the domain were specified as the bottom and side walls with the top being open. The domain reflects a vessel which is to be filled with particles in an upward filling manner consistent with the manufacturing method, with the noted exception regarding singular and sequential placement. The
The conductive properties of sintered materials are dictated by the material properties and the formed network in each sample. For the purposes of this work, sintered bronze disks were selected as a suitable test geometry given their simplicity for manufacture, consistent geo metric tolerancing and accessibility for testing. The sintered bronze disks as received from the supplier, had a
Fig. 1. Sintered bronze sample (A) physical sample, (B) 3D microscope scan of the sample. 3
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
Fig. 2. Sample composition analysis: (A) SEM-EDX analysis, (B) Particle diameter (C) Necks size.
bottom and side walls were taken to be rigid boundaries for which the particles can contact through a tangent condition. The top boundary condition was specified as a threshold condition where a particle can only be positioned such that it is wholly below the top or tangent to it. Filling of the domain was done sequentially using a pseudo-random particle, seeded from a pre-defined distribution. To promote diversity in the formed networks, a distribution [34] was used with the diameters ranging from with 10–100 μm. The bottom row was sequentially con structed from left to right until filled. The particle tops were then mapped to form a singular topology. All subsequent particles were placed sequentially with the topological map being continuously
updated. The particles were placed in a mechanically stable position (local minima), with the process continuing until the last attempted placement resulted in a placement exceeding the top of the domain. The resulting arrangement was analysed to identify the connectivity, finally giving the connectivity network. This process is summarised in Fig. 3. The resulting connectivity network is an unweighted and undirected graph. The particles (nodes) and sinter necks (edges) were used to create an adjacency matrix. The adjacency matrix was weighted based on the sample properties (i.e. particle pair radii, neck area, bulk conductivity). While the calculated weightings were representative of the physical specimen, the magnitude and range of the resistance values did not suit
Fig. 3. 2D sample creation: (A) Packing, (B) Network creation, (C) Resulting network (graph). 4
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
the manufacture of Printed Circuit Boards (PCBs) used for validation. To address this, E96 surface mount resistors were selected, each having a discrete value within the range of 10–1000 Ω inclusive, with a manu facturer specified tolerance of 0.1%. The weighting of the graph edges were assigned pseudo-randomly based on the standard E96 resistor values. Three random networks were selected for the validation study, with duplicates of each for repeatability. The PCBs were manufactured with the interconnecting traces having a manufacturer specified thickness of 35 μm, width of 500 μm and resistivity of 16.8 nΩ m. The PCBs were populated with the prescribed E96 resistors and soldered in a reflow oven. The PCBs were optically inspected for defects and then probed electrically using a Kiethley 2450 SMU in a 2-wire configuration. Losses attributed to the traces were calculated and applied to the weighted Laplace matrix to minimise error. Additionally, probe lead losses (251 � 0.17 mΩ) were applied as a bias to the measured resistance values. Probing of the PCBs was done using the first node as the refer ence for all measurements. Initial testing revealed several discrepancies between the experi mental and analytical values. The discrepancies presented an opportu nity to test the theory of Cernanova et al. [32] that the analytical solutions could be used to identify and locate faults in the circuit. The faults were sought out by measuring the point-to-point resistances and moving the reference node closer to the discrepancy until the resistances across each resistor were measured (in the network). The nodes with the largest discrepancies had the resistors removed and individually tested. This approach led to the identification of components that were out of specification or had faulty solder junctions. The faulty components were replaced with new components within specification. The PCBs were then tested with one node as a reference node with respect to all other nodes on the circuit. The tests required 434 individual measurements, yielding a mean and standard error of 308.7 � 12.6 μΩ. The study successfully validated the analytical method in addition to confirming the fault detection theory of Cernanova et al. [32].
mean and standard error of 54.5 � 7.0 μΩ respectively. The anisotropic case produced a mean and standard error of 67.8 � 4.6 μΩ respectively. The observed differences were attributed to rounding errors in the cal culations done in Simulink and MATLAB. Given the consistency of the 3D verification study, the analytical solution was deemed suitable for use in analysing more complex and irregular structures like those of porous media. 4. 3D simulated samples 4.1. Sample generation With the physical properties characterised, a set of 1000 simulated samples were created in an open source 3D creation software - Blender 2.79. The particles were simulated as ico-spheres with 320 triangular faces per particle. The particle diameters were assigned using Python’s numpy normally distributed random number generator within the range of 0.75 � 0.126 mm. Each simulated sample was initialised with a new pseudo-random distribution of particle sizes ensuring variance between each generated sample. The particles were arranged in an array above a rigid cylindrical container having a diameter of 22 mm and walls higher than the desired 3 mm sample height. The number of particles was calculated using equivalent volumes to ensure an excess of particles were used to fill the container. The particles were assigned a gravitation body force and rigid body collision properties, allowing the particles to interact with the container and each other as they fell in. The body collisions were assigned a damping factor of 0.1 to ensure the particles would reach a static condition. The simulations were run until the particles were stationary, upon which the positional data and particle size was saved to file and exported for subsequent analysis in MATLAB. In matching the physical sample geometry, particles found above the 3 mm sample height were removed using Boolean logic, yielding sam ples consistent with the bronze sinters, as exemplified in Fig. 5 (A). The resulting simulated samples were then analysed to identify how the particles were interconnected. All disconnected particle sets, if any, were removed from the analysis result, as is consistent with the manufacturing process. The result is a single connected network as exemplified in Fig. 5 (B). The physical attributes of the sintered particle network (diameter, neck diameter, material properties) were then applied to the Laplace matrix, yielding a weighted Laplace matrix. The resulting Laplace matrices form the basis for all subsequent electrical and thermal analysis.
3.2. Extension to 3D Having validated the 2D cases, it is expected that an explicitly defined 3D network should yield a similar result. This can be demon strated using a simple cube as illustrated in Fig. 4. The graphs illustrated in Fig. 4, are isomorphic and describe the same network in both two- and three-dimensions. Given the 2D cube to be valid, then the 3D isomorphic graph must also be valid. Therefore, the analytics are taken to be valid for both 2D and 3D networks. To support this, a pair of Simulink simulations were devised using 3D orthonormal networks with each principle direction having four re sistors. Two test cases (isotropic & anisotropic) were used to verify the analytics. The isotropic case used 1 Ω resistors and the anisotropic case used pseudo-random resistors ranging from 1 to 100 Ω inclusive. The resistance values were computed using Eq. (1) and compared to those produced by the Simulink simulation. The computed resistance values were differenced and normalised by the analytical value and then averaged over the 100 test points per case. The isotropic case produced a
4.2. Electrical probe analytics The perimeter of the disk was selected for probing as it offered a consistent surface with an increased length scale for a given set of dia metrically opposite points on the disk. Though probing the sample could be done through a point contact, a distributed line segment across the thickness of the disk was used to be consistent with the experiments. The probes were incorporated into the modelled network through the addition of two points to the network, each representing a single probe. The points on the sample contacted by the probes then add new con nections (edges) to the network as illustrated in Fig. 6. The probe points were diametrically opposite each other and were measured at 18 points around the perimeter with each having an angular separation of 10� . When the probing locations are moved around the sample, a new set of connections are formed resulting a modified form of the original porous network. This yields 18 different networks for each simulated sample, giving a total of 18 000 simulated electrical data points. Each probed network is therefore completely described using weighted Laplace matrices. The eigenvalues and eigenvectors of the Laplace matrices are then used in Eq. (1) to obtain the equivalent elec trical resistance between the probes.
Fig. 4. Isomorphic representations of a cube in 3D (left) and 2D (right). 5
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
Fig. 5. Simulated sample: (A) Particles (blue-interior, green-perimeter, yellow-probed), (B) Network. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
4.3. Thermal flash analytics
used in Eq. (1) to obtain the source-to-sensor thermal resistance.
As with the electrical analysis, the thermal flash model analytics were crafted to replicate the experimental conditions as they were performed across the top and bottom of the samples. The same 1000 simulated samples were re-used subject to a line of sight visibility analysis. This analysis discretised the individual particles into a square array of pixels. The individual pixel size was calculated for each sample based on the smallest particle diameter in the sample which would result in a minimum resolution of 100 pixels per particle diameter. The pro jection of each particle from top-to-bottom was then calculated giving the exposure fraction of each particle to the sensor. The same visibility analysis was done from bottom-to-top to obtain the exposure fractions for those particles able to see the light source. The visibility analysis was then used to create the network probe points and associated links be tween the individual particles and the source and sensor of the thermal flash. The extent of the visibility network was limited to the portion of the sample exposed to the source and sensor of the thermal flash as illustrated in Fig. 7. In each case, the simulated samples, like the phys ical ones, had a small fraction of area that was uninterrupted, thereby forming a small number of pathways from the source to the sensor. The uninterrupted pathways were excluded from the formed network as they would yield a trivial and non-physical result. Furthermore, the impin gent radiation is assumed to be completely absorbed with no losses and no radiative networks between the individual particles were considered. The resistances (R) associated with the weighting of the connections between the sample and the source/sensor were calculated based on the visibility fractions (Aexp ) and particles size (D) as given by Eq. (2). These values were then incorporated into the samples weighted Laplace matrix thereby giving a completely described resistance network. The eigen values and eigenvectors of the weighted Laplace matrices were then
R¼
Fig. 6. The porous network (Blue) with the added probes and connections (Red). (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
Fig. 7. Porous network (Blue) with probes and 5% of the formed connections visualised (Red). (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
D 2kAexp
(2)
5. Experimental methods 5.1. Electrical resistivity To perform the point-to-point electrical resistance measurements, an experimental test rig was constructed as illustrated in Fig. 8. The porous bronze sample (yellow) was fixed to a rotating pedestal (white) on an angular turntable (black) with a pair of copper probes (brown) placed at diametrically opposite sides of the disk in a clamshell style clamp (grey). The copper probes spanned the entire thickness of the bronze samples. When clamped on the disk, the copper probes contact multiple points along the thickness of the disk’s perimeter. Once clamped, the copper probes were connected to a Keithley 2450 SMU in a 4-wire configuration for measuring the electrical resistance. A total of 15 sintered bronze disk samples were tested at 18 points with 10� each, giving a total of 270 point-to-point resistance measurements. 5.2. Thermal flash The thermal flash analysis was performed on the sintered bronze samples using a Linseis Xenon Flash XFA 500 unit. This measurement approach is well suited to samples having varied topology like those of the sintered bronze samples. Though this measurement technique is typically applied to solid samples having a smooth and homogenous
6
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
volume (1.14 cm3). The sample volume can be obtained using either the maximal or mean dimensions of the irregular surfaces as noted previ ously. For this study, both were calculated, assuming maximal and mean heights of 3 mm and 2.45 mm respectively as is consistent with the 3D microscope analysis. The maximal and mean height cases yielded mean porosities of 0.47 and 0.35 respectively. The simulated samples fall within this range, though closer to the mean height case. Given that the particle diameters, junction sizes and porosities are all comparable, the simulation driven analytics were deemed comparable to the physical test specimens. 6.2. Electrical resistivity analysis Fig. 8. Experimental setup for measuring the electrical resistivity across the sintered bronze samples.
The electrical resistance measurements are illustrated in Fig. 9, with the radial direction indicating the resistance values and the theta di rection indicating the angular probe locations. The resistance values appear to be well distributed, having a mean and standard deviation of 1.93 � 0.42 mΩ respectively. The same probing method was applied to the 1000 simulated sam ples, with the same 18 probe points around the perimeter. The electrical resistance values calculated have a mean and standard deviation of 2.03 � 0.17 mΩ as illustrated in Fig. 10. Comparing the experimental and simulated results, the simulated cases are found to be consistent with the experimental results, with a narrower distribution. As such, the simulations and analytics can describe the conductive properties of porous sintered bronze disks.
surface, ASTM E1461-13 thermal flash standard [35], does provide exception to this and permits application to porous samples like the sintered bronze disks. The bronze samples were coated on the top and bottom with three coats of graphite to give the surfaces a homogeneous and low reflectivity surface finish. The analysers thermal sensor was cooled with liquid ni trogen for 2 h prior to the commencement of testing. The graphite coated bronze samples were loaded into a carousel and placed into the analyser. Each sample was subjected to a single flash analysis once every 5 min, with at least five measurements per sample. 6. Results & discussion
6.3. Thermal flash analysis
6.1. Sample properties
The sintered samples were visually inspected and observed to have a relatively low pass through area and deemed suitable for thermal flash analysis. The thermal flash experiments were performed for each of the 15 bronze samples, with a test result exemplified in Fig. 11. The thermal signal exhibits an early spike followed by a gradual increase toward steady state. The early thermal spike was attributed to the porosity of the samples which permits a small portion of the source radiation to impact the sensor directly. The 15 bronze samples were tested at least five times each, with the thermal rise profiles fitted using the Linseis Xeon Flash software. The software was given the following properties (Bulk density - 8.78 g/cm3, heat capacity - 0.380 J/gK and sample thickness - 0.3 cm) [33]. With the provided properties and measured thermal diffusivity the software calculated the ETC. Of the results, only those having a fit quality of 95% or greater were retained. The sample ETC values are illustrated in Fig. 12, with the samples having a calculated ETC with a mean and standard deviation of 21.9 � 3.4 W/mK. The simulated sam ples had an ETC with a mean and standard deviation of 18.1 � 0.30 W/mK. The ETC of the simulated samples is found to be consistent with the physical samples measured using the flash method.
The simulated samples appear to be consistent with the physical samples with some exceptions. The top surface of the physical samples had a greater variance in height as compared to the simulated samples. This variance is attributed to the differences between the manufacturing process and the particle removal by way of Boolean logic. For a given geometry, the porosity of the samples can be readily obtained using Eq. (3). How the volume of the void and total sample is computed can have a significant impact on the sample porosity. For the simulated samples, the top particles are relatively even while the physical samples have height variances up to a single particle. The simulated samples have a narrow porosity range with a mean and standard deviation of 0.38 � 0.002. φ¼
Vvoid ¼1 Vsample
4m
ρπd2 h
(3)
The porosity of the physical samples was calculated using the measured mass, bulk material density (8.89 g/cm3) and the sample
Fig. 9. Electrical resistance measurements of 15 sintered bronze samples over 18 probe locations.
Fig. 10. Simulated electrical resistance values and their mean values for 1000 simulated samples and 18 probe locations. 7
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
Drawing comparison with porosity based empirical correlations like those of Grootenhuis et al. [36] and Koh et al. [37], the ETC is under-estimated with values ranging from 9.3 W/mK and 11.5 W/mK respectively. Comparing the empirical correlation values and those of this work, the traditional contact-based empirics underpredict the ETC by as much as 58% and fail to achieve a value with the experimental error extremums. By implementing the radiative condition, this work better estimates the mean ETC with a 17% under-estimation, which is within the experimental error range. The empirical correlations of Grootenhuis et al. [36] are based on the use of a contact condition in the experimental setup. Using such a method is expected to give a different ETC as the boundary condition will yield a different network with fewer connections as exemplified in Fig. 13. The contact network is observed to be a subset of the radiation network connections. The linking re sistances are principally driven by the contact area, which is a function of the applied pressure in the case of contact experiments. The length scale used in determining the ETC is taken as the geometric height of the porous sample. This length scale favours the radiative ETC flash analysis due to the internally accessed paths, which have a shorter relative path length. Fig. 11. Thermal flash data and fitted curve giving the effective thermal conductivity.
6.4. Particle size analysis To compliment the 0.75 mm sintered sample analysis, a set of sup plementary simulations were generated to investigate the changes in ETC over a range of particle sizes. Four additional sample sets were created with mean diameters of 1, 0.875, 0.625 and 0.5 mm. Each nominal size had 1000 samples, each having the same Gaussian distri butions relative to the nominal diameter and sinter junctions (0.3D). The simulated samples were analysed using the radiative condition, with the ETC plotted in Fig. 14. The ETC is observed to increase in variance and magnitude with particle size. The increase in ETC is asymptotically bounded by the bulk conductivity. 7. Conclusion The effective conductivity of porous media is a complex and multi faceted problem and may be studied using a variety of methods. Ho mogenisation methods are generally limited to bulk estimations and unable to investigate local properties or changes in network architec ture. Though some finite methods such as FEA, FEM, and DEM can investigate local properties, they are computationally intensive due to meshing requirements. To overcome this, a hybrid method was devel oped using rigid body simulations to generate the porous network and
Fig. 12. Thermal conductivity value ranges for the 15 samples tested with fit qualities �95%.
Fig. 13. Formed heat flow networks for: (A) Contact and (B) Radiation.
Fig. 14. Effective thermal conductivity for radiation of modelled samples. 8
J. Polansky et al.
International Journal of Thermal Sciences 148 (2020) 106135
analytical mathematics to evaluate the effective conductivity. The formed network structures can be described by a weighted Laplacian matrix and solved using analytics in a fraction of the time. To verify and validate the analytics, a set of 2D and 3D studies were performed. The 2D analysis successfully verified and validated the analytics for regular and irregular weighted networks. Furthermore, the 2D study confirmed the diagnostic capability postulated by Cernanova et al. [32]. The extension to 3D was demonstrated to be consistent with the 2D approach in principle, as demonstrated by the iso-morphic equivalency. In further support, the analytics were verified using Simulink and found to be consistent and therefore suitable for estimating the ETC of all 3D networked structures. A hybrid approach to porous 3D networks was developed using open source 3D creation software and the extended 3D analytical methods. The electrical resistance of the bronze samples was measured experimentally and found to be consistent with the analytical estimation, having mean and standard deviation values of 1.93 � 0.42 mΩ and 2.03 � 0.17 mΩ respectively. The complementary thermal flash study yielded an ETC with a mean and standard deviation of 21.9 � 3.4 W/mK with the simulations having 18.1 � 0.3 W/mK. Again, the ETC simulations values fall within the experimental results. By contrast, traditional contact-based empirical correlations under predict the ETC by as much as 58% and fail to achieve a value with the experimental error extremums. A wider study of the particle size was performed through simulation and analytics. The ETC of the samples with differing particle size show an increase in ETC with diameter for the radiative case. In conclusion, resistive network analytics can be used in concert with physical models to accurately predict the ETC of porous materials quickly and efficiently.
[4] H. Yu, D. Heider, S. Advani, A 3D Microstructure Based Resistor Network Model for the Electrical Resistivity of Unidirectional Carbon Composites, 2015, pp. 740–749. Composite Structures. [5] W. van Antwerpen, C. du Toit, P. Rousseau, A Review of Correlations to Model the Packing Structure and Effective Thermal Conductivity in Packed Beds of MonoSized Spherical Particles, 2010, pp. 1803–1818. Nuclear Engineering and Design. [6] M. Kaviany, Principles of Heat Transfer in Porous Media, Springer, London, 1999. [7] J. Wang, J.K. Carson, M.F. North, D.J. Cleland, A new approach to modelling the effective thermal conductivity of heterogeneous materials, Int. J. Heat Mass Transf. (2006) 3075–3083. [8] Y. Liang, X. Li, A new model for heat transfer through the contact network of randomly packed granular material, Appl. Therm. Eng. (2014) 984–992. [9] F. Wang, X. Li, The stagnant thermal conductivity of porous media predicted by the random walk theory, Int. J. Heat Mass Transf. (2017) 520–533. [10] C. Hsu, P. Cheng, K. Wong, Modified Zehner-Schlunder models for stagnant thermal conductivity of porous media, Int. J. Heat Mass Transf. (1994) 2751–2759. [11] J. Staggs, Estimating the thermal conductivity of chars and porous residues using thermal resistor networks, Fire Saf. J. (2002) 107–119. [12] H. Zhang, Q. Zhou, H. Xing, H. Muhlhaus, A DEM study on the effective thermal conductivity of granular assemblies, Powder Technol. (2011) 172–183. [13] J. Gan, Z. Zhou, A. Yu, Effect of particle shape and size on effective thermal conductivity of packing beds, Powder Technol. (2017) 157–166. [14] C. Chueh, A. Berti, J. Pharoah, C. Nicolella, Effective conductivity in random porous media with convex and non-convex porosity, Int. J. Heat Mass Transf. (2014) 183–188. [15] C. Argento, D. Bouvard, Modeling the effective thermal conductivity of random packing of spheres through densification, Int. J. Heat Mass Transf. (1996) 1343–1350. [16] S. Masamune, M. Smith, Thermal Conductivity of Beds of Spherical Particles, 1963, pp. 136–143. I&EC Fundamentals. [17] W. Siu, S.-K. Lee, Transient temperature computation of spheres in threedimensional random packings, Int. J. Heat Mass Transf. (2004) 887–898. [18] C. Lee, L. Zhuang, D. Lee, S. Lee, I.-M. Lee, H. Choi, Evaluation of Effective Thermal Conductivity of Unsaturated Granular Materials Using Random Network Model, 2017, pp. 76–85. Geothermics. [19] T.S. Yun, T.M. Evans, Three-dimensional random network model for thermal conductivity in particulate materials, Comput. Geotech. (2010) 991–998. [20] S. Grandjean, J. Absi, D. Smith, Numerical calculations of the thermal conductivity of porous ceramics based on micrographs, J. Eur. Ceram. Soc. (2006) 2669–2676. [21] A. Lodi, S. Martello, M. Monaci, Two-dimensional packing problems: a survey, Eur. J. Oper. Res. 141 (2002) 241–252. [22] B. Korte, J. Vygen, Combinatorial Optimization: Theory and Algorithms, Springer, 2000. [23] W. Huang, Y. Li, H. Akeb, C. Li, Greedy algorithm for packing unequal circles into a rectangular container, J. Oper. Res. Soc. (2005) 539–548. [24] J.A. George, J.M. George, B.W. Lamar, Packing different-sized circles into a rectangular container, Eur. J. Oper. Res. (1995) 693–712. [25] W.M. Visscher, M. Bolsterli, Random packing of equal and unequal spheres in two and three dimensions, Nature 239 (1972) 504–507. [26] V.S.S. Vos, Methods for Determining the Effective Resistance, Universiteit Leiden, 2016. [27] F. Wu, Theory of resistor networks: the two-point resistance, J. Phys. A Math. Gen. (2004) 6653–6673. [28] J. Essam, F. Wu, The exact evaluation of the corner-to-corner resistance of and MxN network: asymptotic expansion, J. Phys. A Math. Theor. (2009), 025205, 110. [29] N.S. Izmailian, R. Kenna, F. Wu, The two-point resistance of a resistor network: a new formulation and application to the cobweb network, J. Phys. A Math. Theor. (2014), 035003, 1-10. [30] J. Essam, Z.-Z. Tan, F. Wu, Resistance between two nodes in general position on an mxn fan network, Phys. Review E (2014), 032130, 1-5. [31] Z.-Z. Tan, J. Essam, F. Wu, Two-point resistance of a resistor network embedded on a globe, Phys. Review E (2014), 012130, 1-7. [32] V. Cernanova, J. Brenkus, V. Stopjakova, Non-symmetric finite networks: the twopoint resistance, J. Electr. Eng. (2014) 283–288. [33] ASM Metals Handbook, in: Properties and Selection: Nonferrous Alloys and Special-Purpose Materials, vol. 2, ASM International, 1990. [34] H.A.M. GmbH, Heraeus Additive Manufacturing, June 2017 [Online]. Available: www.heraeus-additive-manufacturing.com. [35] ASTM E1461-13, Standard Test Method for Thermal Diffusivity by the Flash Method, ASTM International, 2013. [36] P. Grootenhuis, R. Powell, R. Tye, Thermal and electrical conductivity of porous metals made by powder metallurgy methods, Proc. Phys. Soc. 65 (1952) 502–511. [37] K J C Y, A. Fortini, Prediction of thermal conductivity and electrical resistivity of porous metallic materials, Int. J. Heat Mass Transf. 16 (1973) 2013–2021.
Declaration of competing interest In regards to our manuscript titled: A hybrid approach for predicting the effective thermal conductivity of sintered porous materials, we formally declare that we have no conflicts of interest associated with this manuscript or its findings. Acknowledgements This research is conducted with the financial support of Science Foundation Ireland (SFI) under Grant Number 13/RC/2077 and has been part funded by the European Regional Development Fund through the SFI Research Centres Programme. Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi. org/10.1016/j.ijthermalsci.2019.106135. References [1] W. Libeer, F. Ramos, C. Newton, M. Alipanahrostami, C. Depcik, X. Li, Two-phase heat and mass transfer of phase change materials in thermal management systems, Int. J. Heat Mass Transf. 100 (2016) 215–223. [2] F. Pennec, A. Alzina, N. Tessier-Doyen, B. Nait-ali, N. Mati-Baouche, H. De Baynast, D. Smith, A combined finite-discrete element method for calculating the effective thermal conductivity of bio-aggregates based materials, Int. J. Heat Mass Transf. (2013) 274–283. [3] D. Roussel, A. Lichtner, D. Jauffres, R.K. Bordia, C.L. Martin, Effective transport properties of 3D multi-component microstructures with interface resistance, Comput. Mater. Sci. (2015) 277–283.
9