Computational investigation of factors affecting thermal conductivity in a particulate filled composite using finite element method

International Journal of Engineering Science 56 (2012) 86–98

Contents lists available at SciVerse ScienceDirect

International Journal of Engineering Science journal homepage: www.elsevier.com/locate/ijengsci

Computational investigation of factors affecting thermal conductivity in a particulate filled composite using finite element method Muhammad Zain-ul-Abdein a,⇑, Sajjad Azeem b, Syed Mushtaq Shah c a

Faculty of Material Science Engineering, GIK Institute of Engineering Sciences & Technology, Topi, Pakistan Department of Chemical Engineering, UET Peshawar, Pakistan c Université de Lyon, CNRS, INSA-Lyon, LaMCoS UMR5259, F69621, France b

a r t i c l e

i n f o

Article history: Received 28 February 2012 Accepted 20 March 2012 Available online 20 April 2012 Keywords: Polymer matrix composite Bakelite–graphite powder Thermal conductivity Finite element method

a b s t r a c t Particulate filled polymeric composites with enhanced thermo-physical properties are highly demanded in electronic industry. This paper presents an experimental and computational investigation of the thermal conductivity enhancement in a bakelite–graphite composite material. The experimental work illustrates an effect of the graphite addition in different volume fractions upon the effective thermal conductivity of the composite. Computational investigation was performed in two parts. The first part explains a development of experimentally validated finite element models for the estimation of effective thermal conductivity, while the second part demonstrates a detailed analysis of the factors affecting thermal conductivity of the composite. The factors that were examined include particle size with individual constituent properties, and air gaps/voids and interface additions in terms of packing density. The findings showed that not only the finite element simulations may be exploited for the prediction of effective thermal conductivity in a composite material; they may also be helpful in suggesting the optimum particle size and packing density factors to suit the industrial design requirements. Ó 2012 Elsevier Ltd. All rights reserved.

1. Introduction Recent developments in composite materials have successfully overcome several limitations in the applications of conventional materials faced by manufacturing industries. The electronic industry, for instance, is in a permanent quest of materials with high corrosion resistance, high strength, less weight and high electrical and thermal conductivities. Polymer matrix composites (PMCs) with high electrical and thermal conductivities are a better alternate to metal matrix composites (MMCs), since they offer less weight and ease of fabrication in comparison to the latter. Significant research has lately been devoted, to the development of such PMCs (Azeem & Zain-ul-abdein, 2012; Boudenne, Ibos, Fois, Majeste, & Gehn, 2005; Kim, Choi, Lee, & Kang, 2008; Mamunya, Davydenko, Pissis, & Lebedev, 2002; Tekce, Kumlutas, & Tavman, 2007), where the particulate filled composites (PFCs) with improved thermal/electrical conductivity were proposed. An exhaustive study has also been dedicated over a period of time to the theoretical understanding of the variations in thermo-physical properties of the composite materials. Wiener (1912) and Hashin and Shtrikman (1962) bounds, for example, provide maxima and minima of these properties where little information is available for a two-phase medium. From the early works of Clausius–Mossotti (Clausius, 1879) and Maxwell (1884), a self-consistent scheme (SCS) evolved, which later on led to the development of an effective medium scheme (EMS) and a differential scheme (DS) of Bruggeman (1935). Further study yielded an effective field method (EFM) with two variants as Mori and Tanaka (1973) scheme and Kanaun and ⇑ Corresponding author. Tel.: +92 938 271 858; fax: +92 938 271 865. E-mail address: [email protected] (M. Zain-ul-Abdein). 0020-7225/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijengsci.2012.03.035

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

87

Levin (1994) scheme. Since these methods propose the non-empirical formulas, they are used extensively by researchers. However, a major limitation of these schemes is that they assume the presence of second phase particle as an isolated non-interacting inclusion. This implies that for a higher volume fraction of second phase particles improved numerical models are required. The reader is referred to the works of Markov (2000), Sevostianov and Kachanov (2007), for a comprehensive review and applications of these schemes. Various authors have used the computational methods for a direct / inverse identification of effective thermal conductivity (Ke) of the polymeric composites. In their work, Telejko and Malinowski (2004) proposed that the finite element method (FEM), as an inverse identification tool, can be successfully used for the estimation of thermal conductivity with only a few experimental results. Song and Youn (2006) evaluated Ke for carbon nanotube filled polymeric composite using the control volume FEM. They have concluded that homogenization technique employing FEM is a promising tool for the estimation of anisotropic material properties such as thermal conductivity, elastic modulus etc. Decarlis and Jaeger (2001) exploited a self-consistent FEM for the prediction of Ke of a heterogeneous two-phase material for spherical, ellipsoidal and cylindrical shaped inclusions. Yamada, Igawa, Taguchi, and Jitsukawa (2002) investigated numerically and experimentally the effect of SiC fiber upon thermal conductivity in SiC/SiC composite. They have found that the FEM calculations can reveal an improvement in Ke of a composite with different filler addition. Likewise, Bakker (1997) studied the influence of complex porosity and inclusion structures on the thermal and electrical conductivities using FEM and developed a simple relationship that could transfer the 2D conductivity into 3D conductivity. Dede (2010) has recently suggested an FE simulation and optimization technique for the heat flow within an anisotropic composite material. However, verification from experimental results is yet to be established. Yinping and Xingang (1995) used numerical analysis for the prediction of Ke of mixed solid materials. They suggested that the volume fraction, not the particle size / shape, of each phase influence the Ke. This is, however, in contradiction with the theoretical models proposed by Mori-Tanaka and Bruggeman (see Azeem & Zain-ul-abdein, 2012, for reference). Although much attention has been given to the investigation of effective thermal conductivity in the PFCs through FEM, nevertheless, how small the particle size of a PFC should be is rarely quantified in literature. As a general rule, it is assumed that the smaller the particle size, the higher the Ke. It is, however, shown in this work that Ke varies non-linearly with particle size. Moreover, in numerical estimation of Ke of a two-phase particulate medium, the significance of a third phase (generally air) is less highlighted. In other words, little thought has been given to estimating the effect of packing density upon the Ke value. This has also been addressed in this paper. From the design perspective, it is important to know beforehand that which factors influence strongly / weakly the material properties so that a better design can be reached with minimum effort involved. A computational investigation of these problems using properties of the individual constituents for a particulate filled bakelite–graphite composite has been the subject of present work. A simplified approach, where representative volume elements (RVEs) were used, is being proposed so as to examine the role of the ‘key’ factors affecting Ke. 2. Strategy Computational analysis of the factors that influence the effective thermal conductivity (Ke) of a particulate filled composite (PFC) is by no means an easy task. Difficulties arise when micro-sized powders are mixed with each other in various proportions. For instance, fineness of particle size, clustering of like particles, presence of voids / air gaps, interfacial contact etc. are the most important parameters that contribute toward the bulk Ke of a composite material. Owing to several complexities involved in the measurement of these factors, various simplifications were assumed during simulation, which are discussed in Section 4. The entire FE simulation scheme and its findings are being reported here for the first time. The strategy adopted, including experimental observations, is presented in the following steps: Step-I: Experimental measurement of Ke of the bakelite–graphite composite samples with varying volume percent of graphite (30%, 35%, 40%, 45%, 50% and 55%) using P.A Hilton Heat Conduction Unit H-940 and following ASTM standard E1225-99 ASTM (1999). Step-II: Development of an experimentally validated numerical model for a material of known thermal conductivity (brass in this case); a prototype simulation. Step-III: Inverse identification of Ke of bakelite–graphite composites through FE simulations and comparison with experimental results. Step-IV: Computation of the effect of increasing particle size or decreasing clusters with a fully random distribution of bakelite–graphite powder for 50 vol.% graphite sample. Step-V: Numerical calculation of the effect of packing density (range: 100–67%) and interface addition due to powdered particles with a random distribution of air gaps in a 50 vol.% graphite sample. The entire strategy is explained at length in the subsequent sections. 3. Experimental work A detailed description of the experimental work has been provided in Azeem and Zain-ul-abdein (2012). A brief synopsis, relevant to the computational investigation, is presented here. The bakelite–graphite powders were mixed in several volume

88

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

percent in dry state. The samples were molded in a cylindrical shape of diameter 25 mm through compression molding using PRONTOPRESS-2. A P.A. Hilton Heat Conduction Unit H-940 was used to provide the heat flow along the cylindrical axis of the composite specimens. This was achieved by placing the composite samples in between the two brass terminals termed as a ‘hot end’ and a ‘cold end’. The samples were wrapped in teflon tape to avoid heat losses to the surroundings. The unwrapped contact surfaces were coated with a conducting paste to reduce the thermal contact resistance. A schematic representation of the experimental setup, adopted from Azeem and Zain-ul-abdein (2012), is shown in Fig. 1. The ‘hot end’ refers to the terminal which was heated with a constant flux of 10 W, while the ‘cold end’ was equipped with a water-cooled brass bar having a fluid circulation rate of 400 ml/min. The heating and cooling setups were installed at the flat faces of the terminals only, while the remainder of the cylindrical bar was insulated so as to induce the directional heat flow. Temperatures were recorded at the core of the cylindrical bar with the help of several Pt-100 thermocouples (TCs) shown as TC1, TC2, . . ., TC9 in Fig. 1. The TC locations are given in Table 1. Note that TC4, TC5 and TC6 were used only when the thermal conductivity of brass, instead of composite sample, was measured. Also, Th and Tc refer to the temperature values extrapolated from TC1–TC3 and TC7–TC9, respectively. Since the slope of the temperature gradient is proportional to the thermal conductivity of the material, the Ke of the composite samples was calculated when a steady state (temperature difference < ±0.1 °C/min) had reached. The calculations were done in accordance with the standard ASTM: E1225-99 ASTM (1999) using the following equation:

Ke ¼

Q av g  DZ S  DT C

ð1Þ

where Ke = Effective thermal conductivity (W m1 K1) of composite sample DZ = Length (m) of composite sample (Table 2) S = 4.91  104 m2 = Cross-sectional area of the sample DTC = Th – Tc = Temperature (K) difference across the composite sample Qavg = Average heat input (W) as defined in Eq. (2)

Q av g ¼

qh þ qc 2

ð2Þ

where h i 3 = heat flow in hot zone (W) qh ¼ K brass S TZ11 T Z 3 qc ¼ K brass S

h

T 7 T 9 Z 7 Z 9

i

= heat flow in cold zone (W)

with Kbrass = 115 W m1 K1 (D.E. Tyler, Copper Alloy Products, Selection: Nonferrous Alloys, & vol. 2. ASM Hand book, 1998)

Brass terminal (hot end)

TC1

TC2

TC3 Th

Composite Tc TC7 TC8 TC9 sample

25mm

Heat flux

O

O’

ΔZ

30mm

Temperature (°C)

30mm

Distance (mm) Fig. 1. Schematics of experimental setup and temperature measurement.

Brass terminal (cold end)

Water cooled end Insulation

89

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98 Table 1 Thermocouple locations. Thermocouples

TC1

TC2

TC3

TC4

TC5

TC6

TC7

TC8

TC9

Distance (mm) from hot end

5

15

25

35

45

55

65

75

85

Table 2 Experimentally measured thermal conductivities for various samples. Sample type

Graphite (vol.%)

Length DZ (mm)

Effective Thermal conductivity Ke (W m1 K1)

Brass Composite-1 Composite-2 Composite-3 Composite-4 Composite-5 Composite-6

– 30 35 40 45 50 55

30 13 8 13 8 13 13

114.98 4.84 5.92 7.84 9.23 10.78 12.28

T1  T3 = temperature difference between TC1 and TC3 T7  T9 = temperature difference between TC7 and TC9 Z1  Z3 = distance between TC1 and TC3 Z7  Z9 = distance between TC7 and TC9 The measured thermal conductivities for various specimens are reported in Table 2. Although the conductivity of the brass sample was already known, yet it was measured to ensure the accuracy of the experimental procedure. 4. Computational investigation 4.1. Heat transfer analysis FE software AbaqusÒ was employed to identify the Ke value of composite samples. Several heat transfer analyses were performed for different cases. The heat equation which is solved to calculate the temperature field, T, in time (t) and 2D space (x, y) is written as:

    @ @T @ @T @T þ þ Q ¼ qC p Ke Ke @x @x @y @y @t

ð3Þ

where Ke is the effective thermal conductivity in W m1 K1, q is the density in kg m3, Cp is the specific heat in J kg1 K1 and Q is the internal heat generation rate in W m3. Since a steady state with a temperature difference of less than ±0.1 °C/ min must be reached, the specific heat term can be ignored. Owing to the axial symmetry of the specimen, 2D axisymmetric models were used throughout, where the element type was DCAX4. The model geometry is shown schematically in Fig. 2. The model is subdivided into three parts: the hot and cold ends with material as brass of known thermal conductivity (115 W m1 K1) and a central grey-colored part representing material whose conductivity is to be estimated. 4.2. Boundary conditions identification The success of a numerical model is critical to a precise identification of boundary conditions. There were several loading and boundary conditions used in the experimental work. These conditions are also shown schematically in Fig. 2. Although

Direction of heat flow Applied flux (10W)

Hot-End

12.5mm

q = 0 (Insulated surface)

Cold-End

Heat transfer due to thermal conductance

Axis of symmetry 30mm

ΔZ

30mm

q = 0 (Plane of symmetry) Fig. 2. Schematics of loading and boundary conditions and geometry of axisymmetric model.

90

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

the equipment was designed to avoid heat losses as much as possible; the effective flux, which heats up the hot-end, would always be less than the applied flux (10 W). Similarly, heat extraction from the cold-end by means of a water-cooled brass bar involves contact surfaces and hence gives rise to the thermal conductance, a surface roughness and pressure dependent interfacial property. The determination of effective heat flux and thermal conductance for different composite samples is not only difficult but it also complicates the numerical model unnecessarily, and hence increases the overall computation time. A simple solution to the problem is an identification of the boundary conditions as time dependent prescribed temperatures. This was achieved by an extrapolation of the temperatures at TC1–TC3 to point O and TC7–TC9 to point O’, as shown in Fig. 1. The thermal boundary conditions, so calculated, were then applied at all the respective nodes of the hot and cold end faces. Once the boundary conditions are defined, the temperature gradient within the numerical model will depend upon the thermal conductivities of the brass terminals and composite samples. Since the conductivity of brass is known, the only unknown left in the model is the Ke of composite, which is said to be identified when experimentally measured temperatures at TC1–TC3 and TC7–TC9 are reproduced through simulation.

Heat flow Axis of symmetry

Composite mesh

Particle size: 6.25 x 6.5 mm2

Particle size: 3.125 x 3.25 mm2

Particle size: 1.56 x 1.625 mm2

Particle size: 0.781 x 0.812 mm2

Particle size: 0.391 x 0.406 mm2

Particle size: 0.195 x 0.203 mm2

Bakelite

Graphite

Representative Volume Element (RVE)

Fig. 3. Composite mesh with material properties of bakelite (50 vol.%) and graphite (50 vol.%) as individual constituents for different particle sizes.

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

91

4.3. FE mesh and material properties A two-dimensional axisymmetric mesh was created for all the FE models. Four nodes quadrilateral elements with linear interpolation (type: DCAX4) were used where temperature was the only degree of freedom. Maximum number of nodes and elements in any model were 8845 and 8704, respectively, where the size of the smallest element, within the composite mesh, was 195 lm  203 lm. As far as thermal material properties of the composite samples are concerned, two different approaches were employed. In the first approach, the composite was assumed to have the properties of a single homogeneous bulk material. This is because the bakelite and graphite particles were mixed at micro-scale, while their response to thermal loading was determined at macro-scale and hence the term effective thermal conductivity (Ke) of the composite was used. Whether this assumption was safe will be depicted in Section 5.3. The second approach assumes thermal properties of the individual constituents viz. bakelite (KB = 1.4 W m1 K1 Perry, 1997), graphite (KG = 130 W m1 K1 Callister, 2000) and air (KA = 0.0257 W m1 K1 Engineering Toolbox, 2011). This approach is closer to the reality because, no matter how fine the particle size was, the behavior of bakelite and graphite powders is dependent upon their intrinsic properties. It should be noted that the thermal properties of air were also used. The reason is that the mixing and compression molding of powdered particles can never lead to 100% packing density; there is always some percentage of air gaps/voids involved due to mismatch in particles’ shape and interlocking phenomenon.

Fig. 4. Composite mesh with material properties of bakelite, graphite and air as individual constituents with equal vol.% of bakelite and graphite and different vol.% of air/voids.

92

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

Moreover, in dry state, the medium at the interface of graphite and bakelite powders is also air, which, in turn, gives rise to the thermal contact resistance. The inclusion of air gaps as a constituent of PFC can, therefore, not be ignored. Based upon different material properties and distribution of individual components within the mesh, a total of four types of FE models were developed, which are described as separate cases in the following: a. Case-1: Mesh generation with material properties as a single bulk material for the estimation of Ke of composite samples with varying volume percent of graphite viz. 30%, 35%, 40%, 50% and 55%. b. Case-2: Mesh generation (Fig. 3) with material properties of bakelite and graphite as individual constituents with 50 vol.% of graphite for the investigation of effects of heterogeneity in the composite sample upon temperature distribution and particle size / clustering upon Ke. c. Case-3: Mesh generation (Fig. 4) with material properties of bakelite, graphite and air as individual constituents with equal vol.% of bakelite and graphite and varying vol.% of air (i.e. 0%, 1%, 5%, 10%, 15%, 20%, 25%, 30% and 33%) for the investigation of effects of packing density and interface addition upon Ke. d. Case-4: Mesh generation with material properties as a single bulk material for the estimation of equivalent thermal conductivity of composite samples in Cases 2 and 3. One of the assumptions in generating meshes for cases 2 and 3 was that the particle size for both the bakelite and graphite powders was identical throughout and was square-shaped such that an RVE contains equal numbers of non-adjacent bakelite/graphite elements. Another simplification assumed random distribution of individual components (bakelite, graphite and air) within the composite mesh. This was done to maximize the interaction of unlike elements. An obvious outcome of these assumptions is a heterogeneous uniformly-sized randomly distributed composite material. 5. Results and discussion 5.1. Prototype simulation As the name suggests, a prototype simulation was meant to develop an FE model which could verify an appropriate integration of all the necessary inputs in the heat transfer analysis. Since a target simulation is meant to identify the Ke of composite, it is important to ensure that the identified value of Ke is reasonably accurate. This can be achieved by performing a simulation where nothing is unknown. Hence, an analysis was carried out where the entire FE model consisted of a single material brass with the thermal conductivity value of 115 W m1 K1. The thermal histories were first recorded at all the thermocouples (TC1–TC9) through experimentation as discussed in Section 3. Having applied the calculated thermal boundary conditions, the simulated time-temperature curves were compared with experimental curves (Fig. 5). It was found that a value of 115 W m1 K1 yields a close approximation to experimental results. Some discrepancies observed may be attributed to the imprecision in thermocouple locations. It may also observed from Fig. 5 that the temperature difference between TC3 and TC4 (DTTC3–TC4) is almost identical to the difference between TC2 and TC3 (DTTC2–TC3) for both the experimental and simulated results. Note that there is an interface between TC3 and TC4, while the material between TC2 and TC3 is continuous. It was mentioned in the experimental work that a conductive paste was used at the interfaces of the brass terminals and the brass/composite samples in order to reduce the thermal contact resistance. The effectiveness of the paste is obvious from the curves at TC3 and TC4 and also at TC6 and TC7 (involving another interface). Had the thermal contact resistance been significant, the term DTTC3–TC4 would have been much greater than DTTC2–TC3. A similar effect is noted for TC6, TC7 and TC8, where DTTC6–TC7 is approximately the same as DTTC7–TC8.

Temperature (oC)

60

TC1

55

TC2

50

TC3 TC4

45 TC5

40

TC6

35

TC7 TC8

30 TC9

25 20 0

2000

4000

6000

8000

10000

12000

14000

16000

Time (s) Fig. 5. Comparison between experimental (markers) and simulated (solid lines) thermal histories for thermal conductivity estimation of brass at TC1–TC9.

93

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

5.2. Inverse identification of Ke for composite samples Having developed a prototype numerical model, the FE simulations were performed for the inverse identification of Ke of composite samples with graphite volume percent of 30, 35, 40, 45, 50 and 55. Fig. 6 presents a comparison between experimental and simulated results for time-temperature curves at all the TCs except TC4, TC5 and TC6 for a composite sample with 50 vol.% graphite. This is because it was practically much difficult to install these TCs at the specified locations within core of the compaction molded composite samples. Nevertheless, the temperature histories at remaining TCs are sufficient to develop the comparison and estimate the Ke. For the inverse identification of Ke, two extreme values, conductivities of bakelite (KB = 1.4 W m1 K1) and graphite (KG = 130 W m1 K1), were initially assumed as lower and upper limits. It was found that the calculated temperature differences DTTC1–TC3 and DTTC7–TC9 became negligibly small when the minimum value of conductivity had been used. In addition, the simulated TC1–TC3 curves were lying close to the experimental TC1 curve; while the simulated TC7–TC9 curves were in close approximation to TC9 curve. This means that there was practically no or very little heat flow across the composite sample, or in other words, a correct thermal gradient with brass terminals could only be achieved by increasing the Ke value of composite. Likewise, the maximum value of conductivity had shown a relatively small thermal gradient than the one observed experimentally at the steady state. This, in turn, implies that the Ke value was too high to maintain an experimentally observed temperature gradient; rather it quickly allowed the heat flow through the entire model. A few more simulation runs by gradually limiting the Ke value finally resulted in a close estimation of the thermal histories recorded experimentally (Fig. 6). It must be mentioned here that an estimated Ke value is said to be identified if, at steady state, it yields a temperature difference of less than 1 °C at all the respective TCs i.e. DT = |TCiExp – TCiSim| < 1°C with i = 1, 2, 3, 7, 8 and 9. Note that Fig. 6 shows a total of 9 simulated curves (solid lines) and 6 experimental curves (markers). 70

Temperature (oC)

65

TC1 TC2 TC3

60

Brass Terminal (hot end)

TCS3

55 50

TCS7

45

Composite sample (length=13mm)

TCS9

40 TC7 TC8 TC9

35 30

Brass Terminal (cold end)

25 20 0

2000 4000 6000

8000 10000 12000 14000 16000 18000 20000

Time (s) Fig. 6. Comparison between experimental (markers) and simulated (solid lines) thermal histories for Ke estimation of 50 vol.% graphite composite at TC1– TC3 and TC7–TC9.

Table 3 Summary of the experimental and simulated results for all the composite samples. Composite (Graphite vol.%)

Temperatures at steady-state (°C) TC1

TC2

TC3

TC7

TC8

TC9

Effective thermal conductivity Ke (W m1 K1)

30

Experiment Simulation

86.59 87.14

84.78 85.42

83.04 83.71

29.15 29.18

28.41 27.48

26.72 25.77

4.84 4.84 (0.05)

35

Experiment Simulation

69.09 69.13

67.45 67.34

65.69 65.63

37.11 36.85

35.48 35.13

33.77 33.39

5.92 5.92 (+0.05)

40

Experiment Simulation

70.91 71.15

69.12 69.46

67.45 67.78

33.39 33.91

32.03 32.23

30.18 30.54

7.84 7.84 (0.05)

45

Experiment Simulation

62.93 63.14

61.28 61.41

59.49 59.69

40.41 40.79

38.71 39.07

37.09 37.35

9.23 9.23 (0.05)

50

Experiment Simulation

64.35 65.12

62.51 63.37

60.89 61.61

35.26 35.99

33.61 33.81

31.87 32.06

10.78 10.78 (0.05)

55

Experiment Simulation

67.63 68.12

65.93 66.35

64.11 64.59

41.29 41.27

39.58 39.51

37.84 37.75

12.28 12.28 (+0.05)

94

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

A comparison of Fig. 5 (brass sample) and Fig. 6 (50 vol.% graphite composite) shows a strong influence of identified Ke. Here, the former indicates a uniform temperature distribution throughout the numerical model due to a constant thermal conductivity value used; whereas the latter illustrates a clear difference in temperature distribution within the brass terminals and the composite sample due to the change in thermal conductivity along the cylindrical axis. In Fig. 6 TCS3, TCS7 and TCS9 refer to the location at which simulated curves were obtained within the composite sample. For example, TCS3 was a curve plotted within the composite mesh at a distance of 3 mm from the interface between hot-end brass terminal and composite sample. Likewise, TCS7 and TCS9 were drawn at 7 mm and 9 mm, respectively. The Ke values of other composite samples were identified in an identical manner and are reported in Table 3. Note that the simulated results also provide upper / lower bound for the estimated Ke value. In addition, Table 3 summarizes the measured and calculated temperatures at steady-state over all the TCs for all the composite samples used. 5.3. Effect of particle size variation The case-2 type simulations were performed using the composite mesh as shown in Fig. 3. The aim of these simulations was to examine the effect of particle size variation upon the temperature distribution and the effective thermal conductivity.

Fig. 7. Temperature contours with different particle sizes; a. Composite mesh, b. A = 6.25  6.5 mm2, c. A = 3.125  3.25 mm2, d. A = 1.56  1.625 mm2, e. A = 0.781  0.812 mm2, f. A = 0.391  0.406 mm2, g. A = 0.195  0.203 mm2, h. same as g with homogeneous bulk material properties.

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

95

1/Area (1/mm2)

Effective thermal conductivity (W/m-K)

0

5

10

15

20

25

30

80 Ke vs A-1

70 60 50 40 30 20

-0.14 y =K46.78x e vs A R² = 0.911

10 0 0

5

10

15

20

25

30

35

40

45

Particle area (mm2) Fig. 8. Effective thermal conductivity (Ke) as a function of particle size/area (A).

Such an influence could only be investigated if the properties of individual constituents (bakelite and graphite) were used. The volume percent of both the constituents were kept equal i.e. 50% bakelite + 50% graphite. Six different particle sizes were used in the order of their increasing fineness (Fig. 3), or conversely, decreasing clusters; since a large-sized particle may also be assumed as a cluster of fine particles. Fig. 7 demonstrates the contours of temperature distribution resulting from different particle sizes used. It may be noticed that with increasing particle fineness (or decreasing cross-sectional area, A), the temperature contours progressively flatten as in a homogeneous medium. While the heterogeneity of the medium is noted in case of large particle size where non-uniform contours are predominant. This implies that the more the particle size is fine, the lesser the heterogeneity of the medium is effective. It can, therefore, be established that the assumption of using bulk material properties as a homogeneous medium (the first approach mentioned in Section 4.3) is appropriate for the estimation of Ke when micro-sized particles are involved. As the temperature contours change with particle size variation, the heat conduction through the composite should also change. This further implies that the effective thermal conductivity, Ke, should depend upon the particle size/cross-sectional area (A). Since it was not possible to get the Ke value directly from the composite mesh using the individual constituent properties, an indirect way is to develop an equivalent model with composite mesh as a homogeneous medium and then compare the temperature histories at TC locations for both the heterogeneous and homogeneous composite meshes (Fig. 7. g. and h., for example). The effective thermal conductivities (Ke), so obtained from the equivalent models, are plotted against the cross-sectional area (A) of the particle in Fig. 8. It has been found that there is a logarithmic decay in Ke values with increasing particle size/area (A). The Ke values are significantly small for larger A; they, however, quickly rise to an almost constant value beyond a critical particle size, which in the present case is 1.56  1.625 mm2. It not only indicates that the Ke value of a PFC depends strongly upon the particle size, but also that the micro-sized particles are generally sufficient to yield maximum or at least close to maximum Ke and an increasingly smaller size beyond a critical limit would not affect the Ke value significantly. An approximate relation between Ke and A for bakelite–graphite PFC takes the form of power law i.e. K e ¼ 6:44ð1=A0:14 Þ . It is important to note at this stage that the experimentally measured Ke value of 10.78 W m1 K1 (Table 3) for the 50 vol.% graphite micro-sized composite, where an average particle size was 14.7 lm, is not even close to the minimum predicted value of 26.1 W m1 K1 for the macro-sized particles (6.5  6.25 = 40.625 mm2). Moreover, with decreasing particle size, the estimated Ke increases further. This shows that there must be some other strong factors which are responsible for the shift of Ke values toward the lower limits. Given that the composite samples were dry compacts, the packing density of the powdered particles and the thermal contact resistance should play an important role in suggesting the possible values of Ke. Effect of these factors is being discussed in the subsequent section. 5.4. Effect of change in packing density Air gaps are almost always present in the powder compacts due to the interlocking of particles with each other. Additionally, the surface roughness of particles at micro level gives rise to the thermal contact resistance, where the medium between the contact surfaces is, usually, air for dry powders. Collectively, the presence of these voids affects the ‘packing density’ of the PFCs. This means that the packing density (D), practically, may never reach a value of 100% with a simple compaction molding of micro-sized powders. There is always some inclusion of air, which may influence the properties of the composite considerably. From the perspective of numerical modeling, it may, therefore, be assumed that the air gaps and the medium at the contact surfaces jointly represent a third constituent, called air/voids, dispersed evenly in a matrix of bakelite and graphite.

96

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

Fig. 9. Temperature contours with different packing densities; a. Composite mesh, b. 0% voids (100% D), c. 1% voids (99% D), d. 5% voids (95% D), e. 10% voids (90% D), f. 15% voids (85% D), g. 20% voids (80% D), h. 33% voids (67% D).

Numerical simulations were performed for a composite sample of 50 vol.% graphite and 50 vol.% bakelite with the smallest element size of 0.195  0.203 mm2 and varying volume fraction of air i.e. 0%, 5%, 10%, 15%, 20%, 25%, 30% and 33%. The higher the volume fraction of air, the lesser the packing density. Excluding the elements representing air inclusions from the composite mesh the remaining elements were equally divided into bakelite and graphite (see Fig. 4) to maintain their equal proportions. Note that the composite mesh with 0% air is the same as in Fig. 3 with minimum particle size. Fig. 9 illustrates the temperature contours for various packing densities in the order of increasing void fraction. It may be noted that the lesser the packing density, the lesser the heat conduction through the composite sample. Fig. 10 presents the change in effective thermal conductivity (Ke) as a function of packing density (D) and percentage voids. It may be observed that Ke reaches its maximum value at 100% D, shows a sharp decrease at only 1% addition of voids (i.e. 99% D), while decreases linearly with a relatively shallow slope below 99% D. Extrapolation of the curve from 99–67% D values yields upper and lower limits of Ke. An empirical relationship may be written as:

K e ¼ 104:9D  62:08

ð4Þ

It is evident from the above equation that the Ke value drops to zero when D is about 60% (or voids  40%). From the practical point of view, such a high percentage of air suggests that the bakelite and graphite powders are loosely placed close to each

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

97

Fig. 10. Effective thermal conductivity as a function of % age packing density/air gaps.

other and even if there is some contact between the particles it is simply insufficient to allow the heat flow through the composite. Let, for the sake of argument, D be 100% in the above relation then the Ke value becomes 42.82 W m1 K1. However, this value, already known from Fig. 8, is 64.6 W m1 K1. This difference can be understood from the simplification assumed in the previous section where the bakelite and graphite particles had a perfect contact with each other i.e. no inclusion of air even at particles’ interfaces. This is, however, not possible in dry molded compacts where a well defined interface between particles exists almost inevitably. Ignoring interlocking and large air pockets, a very little addition of voids (say 1%) in numerical model should, therefore, be assumed equivalent to the effect of an interface medium between the bakelite and graphite powders and hence a strong influence of thermal contact resistance should arise. It may, thus, be suggested that an initial drop in Ke values is a consequence of interface addition within the composite sample with air as medium of interaction; while large air gaps are responsible for further linear decrease. Now using Eq. (4) for the 50 vol.% graphite composite having Ke value of 10.78 W m1 K1 (Table 3), the packing density, D, is estimated up to 70% approx. This implies that the composite contains almost 30 vol.% air as voids/air gaps or interface medium. This also suggests that the effective thermal conductivity of the composite may be increased merely by increasing the packing density through compaction, heat treatment etc. Moreover, it should be mentioned that Eq. (4) is valid only for equal volume fraction of bakelite and graphite. For other composite samples with unequal bakelite and graphite fractions, the same approach can be used to develop a relationship between Ke and D. 6. Conclusions From several experimental and computational findings presented above, following conclusions can be drawn: – Computational means may be employed successfully for the Ke estimation of PFCs. The Ke of PFC depends upon the volume fraction of individual constituents as in the present case of bakelite–graphite PFC; the higher the graphite content, the higher the Ke. – One of the most important factors that increases the Ke is particle size. Fine particle size yields better Ke. The computational approach, however, reveals that Ke does not depend linearly upon the particle size; rather there is a critical level of particle refinement beyond which the increase in Ke is trivial. This finding is of key importance from the design point-ofview. – Numerical models also reveal that a dry interface addition has a profound effect upon reducing the Ke value, while further decrease in Ke is dependent upon the presence of air pockets / voids. This means that a limited improvement in Ke may be expected with increased packing densities. However, a significant increase in Ke should result if the dry interface is eliminated/minimized, say, by heat treatment. This work may further be extended for the development of computational models with a mixture of different particles distributed within the composite. Yet another investigation may be carried out for spherical and irregular shapes of the powdered particles. References ASTM Standard E1225–99 (1999). Standard Test Method for thermal conductivity of solids by means of the guarded-comparative-longitudinal heat flow technique. Annual Book of ASTM Standards (Vol. 14.02). Azeem, S., & Zain-ul-abdein, M. (2012). Investigation of thermal conductivity enhancement in bakelite–graphite particulate filled polymeric composite. International Journal of Engineering Science, 52, 30–40.

98

M. Zain-ul-Abdein et al. / International Journal of Engineering Science 56 (2012) 86–98

Bakker, K. (1997). Using the finite element method to compute the influence of complex porosity and inclusion structures on the thermal and electrical conductivity. International Journal of Heat and Mass Transfer, 40(15), 3503–3511. Boudenne, A., Ibos, I., Fois, M., Majeste, J. C., & Gehn, E. (2005). Electrical and thermal behavior of polypropylene filled with copper particles. Composites Part A Applied Science and Manufacturing, 36, 1545–1554. Bruggeman, D. A. G. (1935). Calculation of various physical constants in heterogeneous substances. I. Dielectric constants and conductivity of composites from isotropic substances (German). Annalen Physics, 24(5), 636–679. Callister, W. D. Jr., (2000). Fundamentals of Materials Science and Engineering (fifth ed.). Wiley. ISBN:047139551X, p. 459. Clausius, R. (1879). Die mechanische Behandlung der Elektricitat. Braunshweig: Vieweg. Decarlis, A., & Jaeger, M. (2001). Effective thermal conductivity of heterogeneous two-phase material using the self-consistent finite element method. Scripta Materialia, 44, 1955–1958. Dede, E. M. (2010). Simulation and optimization of heat flow via anisotropic material thermal conductivity. Computational Materials Science, 50, 510–515. Engineering Toolbox, 2011, http://www.engineeringtoolbox.com/air-properties-d_156.html, Accessed on July 13, 2011. Hashin, Z., & Shtrikman, S. (1962). A variational approach to the theory of the effective magnetic permeability of multiphase materials. Journal of Applied Physics, 33, 3125–3131. Kanaun, S. K., & Levin, V. M. (1994). Effective field method in mechanics of matrix composite materials. In K. Z. Markov (Ed.), Recent advances in mathematical modelling of composite materials (pp. 1–58). World Sci.. Kim, S. W., Choi, B., Lee, S. H., & Kang, K. H. (2008). Measurement of thermophysical properties of particulate-filled polymer composites. High Temperatures – High Pressures, 37(1), 21–30. Mamunya, Y. P., Davydenko, V. V., Pissis, P., & Lebedev, E. V. (2002). Electrical and thermal conductivity of polymers filled with metal powders. European Polymer Journal, 38, 1887–1897. Markov, K. Z. (2000). Elementary micromechanics of heterogeneous media. In K. Z. Markov & L. Preziozi (Eds.), Heterogeneous media: Micromechanics modeling methods and simulations (pp. 1–162). Boston: Birkhauser. Maxwell, J. C. (1884). A treatise on electricity and magnetism (3rd ed.). USA: Dover Publications Inc., pp. 440–441. Mori, T., & Tanaka, K. (1973). Average stress in matrix and average elastic energy of materials with misfitting inclusions. Acta Materialia, 21, 571–574. Perry, R. H., & Green, D. W. (1997). Perry’s chemical engineers’ handbook (7th ed.). McGraw-Hill, p. 337. Sevostianov, I., & Kachanov, M. (2007). Effect of interphase layers on the overall elastic and conductive properties of matrix composites applications to nanosize inclusion. International Journal of Solids and Structures, 44, 1304–1315. Song, Y. S., & Youn, J. R. (2006). Evaluation of effective thermal conductivity for carbon nanotube/polymer composites using control volume finite element method. Carbon, 44, 710–717. Tekce, H. S., Kumlutas, D., & Tavman, I. H. (2007). Effect of particle shape on thermal conductivity of copper reinforced polymer composite. Journal of Reinforced Plastics and Composites, 26(1), 113–121. Telejko, T., & Malinowski, Z. (2004). Application of an inverse solution to the thermal conductivity identification using the finite element method. Journal of Materials Processing Technology, 146, 145–155. Tyler, D. E. (1998). Wrought copper and copper alloy products, properties and selection: nonferrous alloys and special purpose materials. ASM hand book (Vol. 2, pp. 241–264). ASM International. Wiener, O. (1912). Die Theorie des Mischkorpers fur das Feld des stationaaren Stromung, Erste Abhandlung die Mttelswertsatze fur Kraft. Polarisation und Energie Abh Math-Physichen Klasse Konigl Sacsh Gessel Wissen, 32(6), 509–604. Yamada, R., Igawa, N., Taguchi, T., & Jitsukawa, S. (2002). Highly thermal conductive, sintered SiC fiber-reinforced 3D-SiC/SiC composites: Experiments and finite-element analysis of the thermal diffusivity/conductivity. Journal of Nuclear Materials, 307–311, 1215–1220. Yinping, Z., & Xingang, L. (1995). Numerical analysis of effective thermal conductivity of mixed solid materials. Materials and Design, 16(2), 91–95.