- Email: [email protected]
Heat conduction across an array of cylinders
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Contents lists available at ScienceDirect
Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound
Heat conduction across an array of cylinders C. Pozrikidis Department of Chemical Engineering, University of Massachusetts, Amherst, MA 01003, USA
a r t i c l e
i n f o
Keywords: Heat conduction Effective conductivity Integral equations Boundary-element method
a b s t r a c t Heat transport through an infinite medium across a one-dimensional array of cylinders consisting of another conductive material is considered with the objective of assessing the effective thermal conductivity of the cylinder band. The formulation accounts for the thermal contact resistance of the interface between two materials in imperfect contact. The effective conductivity is expressed in terms of the displacement of the temperature profile far from the cylinder band with respect to that prevailing in a homogeneous medium. Integral representations involving single-layer potentials are developed to describe the temperature field in the presence of the cylinders, pertinent integral equations are derived, and numerical results based on a boundary-element method are presented to illustrate the effect of the cylinder size and conductivity ratios on the effective conductivity of the array. © 2017 Elsevier Ltd. All rights reserved.
1. Introduction The physical and mechanical properties of particulate, composite, and multi-phase materials continues to be the subject of active research (e.g., [1]). A central goal of theoretical and laboratory investigations has been to predict the effective thermal, electrical, and other physical or mechanical properties of these heterogeneous media in terms of corresponding properties of the constituent materials. For example, asymptotic analyses, numerical investigations, and experimental studies have provided us with rigorous expressions and phenomenological correlations for the effective conductivity of dilute or concentrated, ordered or random dispersions of particles with different shapes. In some applications, the dispersed phase has the form of undesirable voids or inclusions embedded in manufacturing. In other applications, particles with high conductivity or permeability are intentionally dispersed in a matrix to enhance the transport properties of the dispersion. For example, in mixed-matrix membranes, small particles of a microporous material, identified as a filler, are dispersed in a dense nonporous polymer material, identified as a matrix, and then processed into a thin composite layer, identified as a membrane. The hope is that the filler, chosen for its high adsorption affinity or transport rate for a molecular species of interest, will improve the efficacy of the matrix in membranemediated separations (e.g., [2,3]). The effective thermal conductivity of a medium containing ordered arrays or random collections of spherical, spheroidal, or cylindrical particles has received considerable attention following Maxwell’s seminal analysis [4]. Effective medium theories have provided us with appropriate frameworks and suggested expressions for the effective diffusivity in
terms of fundamental solutions (e.g., [5]). Numerical solutions based on finite-element methods have been presented to extend and confirm the theoretical predictions. In this paper, we study the effective thermal conductivity of an infinite homogeneous medium hosting a one-dimensional array of parallel cylinders consisting of another conductive material. In engineering applications involving reinforced or engineered materials, the cylinders can be regarded as elongated fibers. The infinite medium, termed the matrix, and arrays of cylinders, identified as fillers, have different thermal conductivities. The interfaces between the two materials exhibit contact resistance due to imperfect mechanical contact. In applications, the actual contact area can be as low as one percent of the nominal contact area. The effect of dispersed particle contact resistance was considered by previous authors in the limit of infinite dilution and for other geometrical arrangements, such as regular dispersions of spheres (e.g., [6–9]). The configuration and method of analysis presented in this paper have not been studied or discussed by previous authors. There are similarities between the problem considered in this work and that of elastic wave propagation through a medium with an array of cracks (e.g., [10,11]). However, the two problems differ in the governing equations, for example, Helmholtz v. Laplace, and more important, in the boundary conditions imposed at the crack, void, or conductive cylinder surface. Differences are also noted with regard to the physical objectives. In the case of cracks, typically resembling flat plates, an interior solution is not required. Although formulations in terms of integral equations are possible in both cases, the derivation and final form of these equations are different. Fourier expansions for periodic configuration are typically employed in wave propagation problems (e.g., [12]).
E-mail address: [email protected] https://doi.org/10.1016/j.enganabound.2017.09.007 Received 16 January 2017; Received in revised form 18 August 2017; Accepted 16 September 2017 0955-7997/© 2017 Elsevier Ltd. All rights reserved.
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
contact resistance, and (b) continuity of flux, q, given by Fourier’s law of conduction, 𝑘c (𝑇layer − 𝑇matrix ) = ± 𝑞,
𝑞 ≡ −𝑘matrix
d𝑇layer d𝑇matrix = −𝑘layer d𝑧 d𝑧 (3)
at 𝑧 = ±𝛿, where kmatrix is the thermal conductivity of the matrix, klayer is the thermal conductivity of the layer, and kc is the thermal contact conductivity ranging from zero to infinity. As kc tends to infinity, the interfacial temperature drop tends to disappear at the interface. For future reference, we introduce the ratio of the layer-to-matrix thermal conductivity, a dimensionless thermal contact conductivity, and a dimensionless thermal contact resistance. 𝛼=
𝑘layer 𝑘matrix
,
𝛼c =
𝛿𝑘c , 𝑘matrix
𝑟c =
𝑘layer 𝛼 = . 𝛼c 𝛿𝑘c
(4)
For the configuration shown in Fig. 1, 𝛼 > 1. The solution of the boundary-value problem can be found by elementary methods and is given by ( ( ) ) 1−𝛼 1 𝛿 𝑇matrix (𝑧) = 𝛾𝑧 1 + (5) + 𝛼 𝛼c |𝑧| above or below the layer, |z| > 𝛿, and 𝑇layer (𝑧) =
Fig. 1. Illustration of one-dimensional diffusive transport through a matrix across a layer of another material with thickness 2𝛿 confined between two planes located at 𝑧 = ±𝛿.
1 − 𝛼 + 𝑟c 𝑑 1−𝛼 1 = = + . 𝛿 𝛼 𝛼c 𝛼
(7)
We observe that d becomes negative when 𝛼 > 1 + 𝑟c or, equivalently, 𝛼 > 𝛼c ∕(𝛼c − 1). We may select an arbitrary elevation, h ≥ 𝛿, and define an effective conductivity based on the equation 𝑘ef f (ℎ) ℎ ℎ ≡𝛾 = . 𝑘matrix 𝑇matrix (ℎ) ℎ+𝑑
(8)
This would be the thermal conductivity recorded when the heat transport is driven by a temperature difference between two isothermal plates located at 𝑦 = ±ℎ. Substituting the expression for d given in (7), we obtain
2. Transport across a conductive layer
𝑘ef f (ℎ) 1 = ( 𝑘matrix 1 + 1−𝛼 𝛼 +
To establish a point a reference, we consider steady, one-dimensional heat transport through an infinite homogeneous conductive medium, termed the matrix, across a layer of another conductive medium with thickness 2𝛿, as shown in Fig. 1. The temperature, T(z), satisfies Laplace’s equation in the matrix and in the layer,
1 𝛼c
)
𝛿 ℎ
=
𝛼 𝛼 + (1 − 𝛼 + 𝑟c )
𝛿 ℎ
.
(9)
As h/𝛿 tends to unity, the right-hand side tends to 𝛼∕(1 + 𝑟c ); correspondingly, the effective diffusivity tends to 𝑘layer ∕(1 + 𝑟c ). As the ratio h/𝛿 tends to infinity, the fraction on the right-hand side tends to unity and the effective conductivity reduces to the matrix conductivity.
(1)
3. Periodic arrangement of cylinders
Far above and below the layer, as z → ± ∞, the temperature assumes the asymptotic form 𝑇 ≃ 𝛾 (𝑧 ± 𝑑),
(6)
inside the layer, |z| < 𝛿. An inconsequential constant can be added to the right-hand sides of (5) and (6). Comparing (5) with (2), we conclude that the profile displacement is given by
In Section 2, a framework for evaluating the effective conductivity is established by discussing an elementary one-dimensional model involving heat transport across a compact layer. The effective diffusivity is quantified in terms of the displacement of the temperature profile with respect to that prevailing in the absence of the layer. Having established this framework, in Section 3 we consider the temperature distribution and profile displacement due to an array of cylinders. The Laplace equation governing the temperature distribution will be solved in terms of single-layer representations involving the periodic Green’s function. The formulation results in a coupled system of two integrals equations of the second kind over the cylinder contours. Numerical results by standard boundary-element methods will be presented and discussed with reference to the compact layer configuration and with respect to analytical expressions for a dilute dispersion.
d2 𝑇 = 0. d𝑧 2
1 𝛾𝑧 𝛼
In the main part of this work, we consider heat transport in an infinite medium across a periodic array of cylinders with arbitrary crosssectional shapes separated by distance L along the x axis, as shown in Fig. 2. The temperature satisfies Laplace’s equation in the xy plane,
(2)
where 𝛾 is a specified gradient and d is a profile displacement to be determined as part of the solution. An inconsequential constant can be added to the right-hand side of (2) to render the temperature positive within a domain of interest, including the matrix and the layer. At the interfaces between the layer and the matrix, we require two conditions: (a) a jump in the temperature determined by the thermal
∇2 𝑇 = 0,
(10)
subject to a periodicity condition in the direction of the x axis, 𝑇 (𝑥 + 𝐿, 𝑦) = 𝑇 (𝑥, 𝑦). Far above and below the cylinder array, as y → ± ∞, the diffusive field assumes the asymptotic form shown in (2), where 𝛾 is a specified gradient and the profile displacement, d, must be found as part of the 106
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 2. Illustration of heat transport in the xy plane across a periodic arrangement of cylinders with arbitrary cross-sectional shapes. A temperature profile is shown to illustrate the significance of the displacement length, d. The cylinder interface exhibits a thermal contact resistance due to imperfect contact.
solution. We recall that an inconsequential constant can be added to the right-hand side of (2). For future reference, we introduce the ratio of the cylinder-to-matrix thermal conductivity, a dimensionless thermal contact conductivity, and a dimensionless thermal contact resistance. 𝛼=
𝑘cylinder 𝑘matrix
,
𝛼c =
𝑎𝑘c 𝑘matrix
,
𝑘cylinder 𝛼 𝑟c = = , 𝛼c 𝑎𝑘c
must be satisfied. Simplifying and using (17), we obtain 𝑇matrix (𝐱) ≃ 𝛾 (𝑦 ± 𝑑), where 𝑑=
(11)
3.1. Integral representations
3.2. Integral equations
where 𝜙(𝓁) is a dimensionless strength density distribution, is the contour of any one cylinder, 𝓁 is the arc length around , (𝐱, 𝐱′ ) is the periodic Green’s function of Laplace’s equation in two dimensions representing the field due to a periodic collection of point sources,
Requiring continuity of flux at the particle surface, and employing the properties of the Green’s function, we obtain the equation 𝑛𝑦 (𝐱 ) −
(𝐱 , 𝐱 ′ ) ≃ −
|𝐱 − 1 ln 2𝜋 𝐿
.
1 |𝑦 − 𝑦′ | + ⋯ , 2𝐿
(14)
𝑛𝑦 (𝐱 ) = (15)
where
where the dots indicate exponentially decaying terms. Using the properties of the Green’s function, we find that, as y → ± ∞, the integral representation yields the far field ( ) 1 𝑇matrix (𝐱) ≃ 𝛾 𝑦 ± 𝜙(𝐱′ )(𝑦′ − 𝑦) d𝓁(𝐱′ ) , (16) 2𝐿 ∮
𝜒(𝐱) ≡
𝜙(𝐱) d𝓁(𝐱) = 0
1 2
PV
∮
𝜙(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ )
𝜓(𝐱) +
PV
∮
) 𝜓(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ) ,
(21)
1 2
1+𝛼 𝜒(𝐱) − ∮ 1−𝛼
𝜙(𝐱) + 𝛼𝜓(𝐱) , 1+𝛼
PV
𝜁 (𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ),
𝜁 (𝐱 ) ≡
𝜙(𝐱) − 𝛼𝜓(𝐱) 1−𝛼
(22)
(23)
are auxiliary functions. Conversely, 𝜙(𝐱) =
where the plus sign applies far above the array and the minus sign applies far below the array. In order for the integral on the right-hand side of (16) to tend to a constant far from the cylinder array, as y → ± ∞, the constraint ∮
𝜙(𝐱) +
where the point x lies at the cylinder contour, the gradient 𝛁 operates with respect to the evaluation point, x, and PV denotes the principalvalue integral. Rearranging, we obtain the integral equation
As |𝑦 − 𝑦′ | → ∞, we obtain a linear asymptotic behavior, (𝐱 , 𝐱 ′ ) ≃ −
1 2
( = 𝛼 𝑛𝑦 (𝐱 ) +
(13)
and 𝜚 = 2𝜋∕𝐿 is the wave number (e.g., [13]). As x → x′, we recover the free-space Green’s function given by 𝐱′ |
(19)
where 𝜓(𝓁) is a dimensionless strength density distribution. The two single-layer distributions, 𝜙 and 𝜓, must be adjusted to satisfy the cylinder-matrix interfacial conditions.
The temperature field in the matrix can be represented by a singlelayer harmonic potential, ( ) 𝑇matrix (𝐱) = 𝛾 𝑦 + 𝜙(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) , (12) ∮
( ) 1 (𝐱 , 𝐱 ) = − ln 2 {cosh[𝜚(𝑦 − 𝑦′ )] − cos[𝜚(𝑥 − 𝑥′ )]} , 4𝜋
1 𝜙(𝐱) 𝑦 d𝓁(𝐱) 2𝐿 ∮
is the profile displacement. Inside the cylinders the temperature field can be represented by another single-layer harmonic potential, ( ) 𝑇cylinder (𝐱) = 𝛾 𝑦 + 𝜓(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) , (20) ∮
where kc is the thermal contact conductivity and a is a designated cylinder radius.
′
(18)
(1 + 𝛼) 𝜒(𝐱) + (1 − 𝛼) 𝜁 (𝐱) , 2
𝜓(𝐱) =
(1 + 𝛼) 𝜒(𝐱) − (1 − 𝛼) 𝜁 (𝐱) , 2𝛼 (24)
and 𝜙(𝐱) − 𝜓(𝐱) =
) 𝛼2 − 1 ( 𝜒(𝐱) − 𝜁 (𝐱) . 2𝛼
(25)
If the pair of functions 𝜒 and 𝜁 is available, the corresponding pair 𝜙 and 𝜓 can be reconstructed, and vice versa.
(17) 107
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 3. Distribution of the single-layer density, 𝜙, for a periodic arrangement of circular cylinders with scaled radius (a) 𝑎∕𝐿 = 0.10 and (b) 0.49, for 𝛼c = 1.0 and 𝛼 = 0.05 (dots connected with dashed lines), 0.1, 0.2, 0.4, and 0.8. The dashed lines in (a) represent the analytical solution for isolated cylinders given in (31).
The interfacial condition for the jump in the temperature requires that
Half of the first fraction on the right-hand side is equal to the areal fraction of an interfacial band of thickness 2a occupied by the periodic array of cylinders, 𝑐𝑎 = 𝜋𝑎2 ∕(2𝑎𝐿). For an arbitrary collection of small cylinders inside an layer of thickness 2𝛿, we obtain
( ) 𝑘𝑐 𝜙(𝐱′ ) − 𝜓(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) 𝑘matrix ∮ = 𝑛𝑦 (𝐱 ) −
1 2
𝜙(𝐱) +
PV
∮
𝜙(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ),
1 − 𝛼 + 𝑟c 𝑑 𝑐 , ≃2 𝛿 1 + 𝛼 + 2𝑟c 𝑎
(26)
where the point x lies at the cylinder contour. In terms of the auxiliary functions 𝜒 and 𝜁 , we obtain the integral equation (1 + 𝛼) 𝜒(𝐱) + (1−𝛼) 𝜁 (𝐱) −
PV
−
(
∮
where ca is the cylinder areal fraction. The associated effective conductivity introduced in (9) is
( ) 1−𝛼 2 𝑘c 𝜒(𝐱′ )−𝜁 (𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) 𝛼 𝑘matrix ∮ (27)
𝐷ef f (ℎ) ℎ ℎ ≡𝛾 = ≃ 𝐷m 𝑇 (ℎ) ℎ+𝑑 1+2
∮
𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱) = − 1 , 2
𝐷ef f (𝛿) 1 − 𝛼 + 𝑟c 𝛿 = ≃1−2 𝑐 . 𝐷m 𝛿+𝑑 1 + 𝛼 + 2𝑟c 𝑎
∫
𝜒 d𝓁 = −(1 − 𝛼)
(35)
3.4. Numerical method
∫
𝜁 d𝓁,
The integral Eqs (22) and (28) were solved by a standard boundaryelement method involving particle contour discretization into straight segments and subtraction of the singularity of the single-layer potential in terms of the free-space Green’s function to ensure high accuracy (e.g., [13]). The double-layer potential is identically zero when evaluated at a point on a flat boundary element. The twelve-point Gauss–Legendre quadrature was employed to compute the non-singular and regularized element influence coefficients of the double-layer potential. Results of computations with 256 boundary elements generate results that are accurate to the third significant figure for nearly touching cylinders.
(30)
3.3. Small cylinders When the cylinder size is small compared to the separation, each cylinder can be studied as though it were in isolation. The solution of the integral equations for small circular cylinders of radius a in the limit as a/L → 0 is given by 𝛼c (1 − 𝛼) + 𝛼 1 − 𝛼 + 𝑟c sin 𝜃, sin 𝜃 = 2 𝛼c (1 + 𝛼) + 2𝛼 1 + 𝛼 + 2𝑟c
3.5. Results and discussion Numerical solutions of the integral equations for small circular cylinders with scaled radius 𝑎∕𝐿 = 0.10 are shown in Fig. 3(a) for several values of 𝛼 and dimensionless contact conductivity 𝛼c = 1.0. The numerical results are in excellent agreement with the analytical solution given in (31) for isolated cylinders, represented by the dashed lines. Corresponding solutions of the integral equations for nearly touching circular cylinders with scaled radius 𝑎∕𝐿 = 0.49 are shown in Fig. 3(b). These distributions differ markedly from the sinusoidal distributions for small
(31)
where 𝜃 is the polar angle measured around the cylinder center. Substituting this asymptotic solution into (19) and performing the integration, we obtain the scaled profile displacement 𝑑 𝜋𝑎 1 − 𝛼 + 𝑟c = . 𝑎 𝐿 1 + 𝛼 + 2𝑟c
(34)
(29)
which demonstrates that the integral constraint (17) is satisfied for any strength density distribution, 𝜙(𝓁), arising from the solution of the integral equation.
𝜙(𝜃) = 2
.
This formula generalizes results obtained by previous authors in the absence of contact resistance (e.g., [14]).
we find that (1 + 𝛼)
1 1−𝛼+𝑟c 𝛿 𝑐 1+𝛼+2𝑟c ℎ 𝑎
For ℎ = 𝛿, we obtain
) (1 + 𝛼) 𝜒(𝐱′ ) + (1−𝛼) 𝜁 (𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ) = 2 𝑛𝑦 (𝐱). (28)
Integrating both sides of the integral Eq. (22) with respect to arc length, 𝓁, around the cylinder contour, , and using the identity 𝑃𝑉
(33)
(32)
108
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 4. Scaled profile displacement, d/a, for a periodic array of circular cylinders with scaled radius 𝑎∕𝐿 = 0.01 ( × ), 0.10 (∗ ), 0.25 ( + ), and 0.49 (○), where L is the cylinder separation, and for (a) 𝛼c = 1 and (b) 10. The dotted lines represent the analytical solution for small cylinders given in (32). The solid line represents the predictions of the layer model discussed in Section 2 for layer half-thickness 𝛿 = 𝑎.
Fig. 5. Scaled profile displacement, d/a, and scaled effective conductivity, keff /kmatrix , for a periodic array of circular cylinders for scaled contact conductance (a, c) 𝛼c = 1 and (b, d) 10.
109
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
or well separated cylinders. It is worth noting that the distribution of 𝜙 over the main part of the cylinders is nearly flat. The scaled profile displacement, d/a, depends on three dimensionless parameters: the scaled cylinder radius, a/L, the cylinder to matrix conductivity ratio, 𝛼, and the contact conductivity ratio, 𝛼 c . Graphs of the scaled profile displacement against 𝛼 are shown in Fig. 4 for circular cylinders with scaled radius 𝑎∕𝐿 = 0.01, 0.10, 0.25, and 0.49, and contact conductance 𝛼c = 1 and 10. The numerical results for small cylinders are in excellent agreement with the analytical solution given in (32), represented by the dotted lines. In fact, the agreement is fair for cylinder radius as high as 𝑎∕𝐿 = 0.25. As the cylinders tend to touch, a/L → 0.5, the numerical results tend to those predicted by the layer model discussed in Section 2, represented by the solid line for layer half-thickness 𝛿 = 𝑎 . The layer model predicts that 𝛼 (1 − 𝛼) + 𝛼 𝑑 = c , 𝑎 𝛼𝛼c
conductive properties of the array. In practice, the contact conductance may be increased by interfacial pressure mediated by isotropic, uniaxial, or a more general type of imposed strain or stress. 4. Discussion We have presented a framework for computing the effective conductivity of a matrix hosting a layer of another material or a periodic array of cylinders in terms of the displacement of the temperature profile far from the layer or cylinders. The formulation accounts for interfacial contact resistance due to imperfect contact. The integral formulation and boundary-element methods for two-dimensional configurations discussed in this paper can be extended to describe periodic arrays of spheres in terms of the doubly periodic Green’s function of Laplace’s equation. A practical concern is that the computation of the Green’s function is expensive, even when expedited methods based in Ewald’s summation technique are employed. Typical formulations involving a periodic box where Dirichlet conditions are imposed at the top and bottom and periodicity conditions are imposed in two directions are a less elegant yet more tractable alternative.
(36)
independent of a/L. As 𝛼 tends to infinity, d/a tends to (1 − 𝛼c )∕𝛼c , which is zero when 𝛼c = 1 and negative when 𝛼c = 10, in agreement with the numerical results presented in Fig. 4. Expression (36) can be contrasted with that for small cylinders shown in (32). A linear combination of the layer and small-cylinder results can be used to describe the profile displacement in engineering applications where a high degree of accuracy is not required. A summary of numerical results is presented in Fig. 5(a) and (b) in the form of three-dimensional graphs of d/a against a/L and 𝛼 for 𝛼c = 1 and 10. Though the general behavior of the two top graphs presented in this figure is similar, significant quantitative differences are observed, as seen previously in Fig. 4. We note in particular, that the profile displacement is positive when 𝛼c = 1, indicating a reduction in the effective conductivity with respect to that of the matrix due to dominant effect of the contact resistance, even for large 𝛼. Physically, high-conductivity cylinders loosely embedded in the matrix act as near insulators. In contrast, the profile displacement may become negative when 𝛼c = 10, indicating an enhancement in the effective conductivity. This behavior is illustrated directly in corresponding graphs of the effective conductivity, computed as 𝑘ef f 𝑎 ≡ 𝑘matrix 𝑎+𝑑
References [1] Öchsner A, Murch GE. Heat transfer in multi-phase materials. New York: Springer; 2011. [2] Wang T-P, Kang D-Y. Predictions of effective diffusivity of mixed matrix membranes with tubular fillers. J Membr Sci 2015;485:123–31. [3] Pozrikidis C, Ford DM. Conductive transport through a mixed-matrix membrane. J Eng Math 2017. (In Press). [4] Maxwell JC. Electricity and magnetism. Oxford: Clarendon Press; 1873. [5] Jeffrey DJ. Conduction through a random suspension of spheres. Proc R Soc Lond A 1873;335:355–67. [6] Benveniste Y. Effective thermal conductivity of composites with a thermal contact resistance between the constituents: non-dilute case. J Appl Phys 1987;61:2840–3. [7] Cheng H, Torquato S. Effective conductivity of periodic arrays of spheres with interfacial resistance. Proc R Soc Lond A 1997;453:145–61. [8] Duschlbauer H, Pettermann E, Böhm HJ. Heat conduction of a spheroidal inhomogeneity with imperfectly bonded interface. J Appl Phys 2003;94. https://doi.org/10.1063/1.1587886. [9] Matt CF, Cruz ME. Heat conduction in two-phase composite materials with threedimensional microstructures and interfacial thermal resistance. In: Öchsent A, M̎urch GE, editors. Heat transfer in multi-phase materials, vol. 2; 2010. p. 63–97. [10] Achenbach JD, Kitahara M, Mikata Y, Sotiropoulos DA. Reflection and transmission of plane waves by a layer of compact inhomogeneities. Pure Appl Geophys 1988;128:101–18. [11] Achenbach JD. Use of elastodynamic reciprocity theorems for field calculations. In: S̎chiavone P, Constanda C, Mioduchowski A, editors. Integral methods in science and engineering. Birkhäuser; 2002. [12] Scarpetta E, Sumbatyan A. Wave propagation through a periodic array of inclined cracks. Eur J Mech A/Solids 2000;19:949–59. [13] Pozrikidis C. A practical guide to boundary-element methods with the software library BEMLIB. Boca Raton: Chapman & Hall/CRC; 2002. [14] Hashin Z, Shtrikman S. A variational approach to the theory of the effective magnetic permeability of multiphase materials. J Appl Phys A 1962;33:3125–31.
(37)
according to (8), as shown in Fig. 5(c) and (d). When 𝛼c = 1, the effective conductivity is less than unity due to the low contact conductance. When 𝛼c = 10, the effective conductivity may be higher than unity for sufficiently high 𝛼 due to the moderate contact resistance. These results emphasize the importance of interfacial asperity interlocking. Mechanical degradation may loosen the cylinders and significantly affect the
110
Contents lists available at ScienceDirect
Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound
Heat conduction across an array of cylinders C. Pozrikidis Department of Chemical Engineering, University of Massachusetts, Amherst, MA 01003, USA
a r t i c l e
i n f o
Keywords: Heat conduction Effective conductivity Integral equations Boundary-element method
a b s t r a c t Heat transport through an infinite medium across a one-dimensional array of cylinders consisting of another conductive material is considered with the objective of assessing the effective thermal conductivity of the cylinder band. The formulation accounts for the thermal contact resistance of the interface between two materials in imperfect contact. The effective conductivity is expressed in terms of the displacement of the temperature profile far from the cylinder band with respect to that prevailing in a homogeneous medium. Integral representations involving single-layer potentials are developed to describe the temperature field in the presence of the cylinders, pertinent integral equations are derived, and numerical results based on a boundary-element method are presented to illustrate the effect of the cylinder size and conductivity ratios on the effective conductivity of the array. © 2017 Elsevier Ltd. All rights reserved.
1. Introduction The physical and mechanical properties of particulate, composite, and multi-phase materials continues to be the subject of active research (e.g., [1]). A central goal of theoretical and laboratory investigations has been to predict the effective thermal, electrical, and other physical or mechanical properties of these heterogeneous media in terms of corresponding properties of the constituent materials. For example, asymptotic analyses, numerical investigations, and experimental studies have provided us with rigorous expressions and phenomenological correlations for the effective conductivity of dilute or concentrated, ordered or random dispersions of particles with different shapes. In some applications, the dispersed phase has the form of undesirable voids or inclusions embedded in manufacturing. In other applications, particles with high conductivity or permeability are intentionally dispersed in a matrix to enhance the transport properties of the dispersion. For example, in mixed-matrix membranes, small particles of a microporous material, identified as a filler, are dispersed in a dense nonporous polymer material, identified as a matrix, and then processed into a thin composite layer, identified as a membrane. The hope is that the filler, chosen for its high adsorption affinity or transport rate for a molecular species of interest, will improve the efficacy of the matrix in membranemediated separations (e.g., [2,3]). The effective thermal conductivity of a medium containing ordered arrays or random collections of spherical, spheroidal, or cylindrical particles has received considerable attention following Maxwell’s seminal analysis [4]. Effective medium theories have provided us with appropriate frameworks and suggested expressions for the effective diffusivity in
terms of fundamental solutions (e.g., [5]). Numerical solutions based on finite-element methods have been presented to extend and confirm the theoretical predictions. In this paper, we study the effective thermal conductivity of an infinite homogeneous medium hosting a one-dimensional array of parallel cylinders consisting of another conductive material. In engineering applications involving reinforced or engineered materials, the cylinders can be regarded as elongated fibers. The infinite medium, termed the matrix, and arrays of cylinders, identified as fillers, have different thermal conductivities. The interfaces between the two materials exhibit contact resistance due to imperfect mechanical contact. In applications, the actual contact area can be as low as one percent of the nominal contact area. The effect of dispersed particle contact resistance was considered by previous authors in the limit of infinite dilution and for other geometrical arrangements, such as regular dispersions of spheres (e.g., [6–9]). The configuration and method of analysis presented in this paper have not been studied or discussed by previous authors. There are similarities between the problem considered in this work and that of elastic wave propagation through a medium with an array of cracks (e.g., [10,11]). However, the two problems differ in the governing equations, for example, Helmholtz v. Laplace, and more important, in the boundary conditions imposed at the crack, void, or conductive cylinder surface. Differences are also noted with regard to the physical objectives. In the case of cracks, typically resembling flat plates, an interior solution is not required. Although formulations in terms of integral equations are possible in both cases, the derivation and final form of these equations are different. Fourier expansions for periodic configuration are typically employed in wave propagation problems (e.g., [12]).
E-mail address: [email protected] https://doi.org/10.1016/j.enganabound.2017.09.007 Received 16 January 2017; Received in revised form 18 August 2017; Accepted 16 September 2017 0955-7997/© 2017 Elsevier Ltd. All rights reserved.
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
contact resistance, and (b) continuity of flux, q, given by Fourier’s law of conduction, 𝑘c (𝑇layer − 𝑇matrix ) = ± 𝑞,
𝑞 ≡ −𝑘matrix
d𝑇layer d𝑇matrix = −𝑘layer d𝑧 d𝑧 (3)
at 𝑧 = ±𝛿, where kmatrix is the thermal conductivity of the matrix, klayer is the thermal conductivity of the layer, and kc is the thermal contact conductivity ranging from zero to infinity. As kc tends to infinity, the interfacial temperature drop tends to disappear at the interface. For future reference, we introduce the ratio of the layer-to-matrix thermal conductivity, a dimensionless thermal contact conductivity, and a dimensionless thermal contact resistance. 𝛼=
𝑘layer 𝑘matrix
,
𝛼c =
𝛿𝑘c , 𝑘matrix
𝑟c =
𝑘layer 𝛼 = . 𝛼c 𝛿𝑘c
(4)
For the configuration shown in Fig. 1, 𝛼 > 1. The solution of the boundary-value problem can be found by elementary methods and is given by ( ( ) ) 1−𝛼 1 𝛿 𝑇matrix (𝑧) = 𝛾𝑧 1 + (5) + 𝛼 𝛼c |𝑧| above or below the layer, |z| > 𝛿, and 𝑇layer (𝑧) =
Fig. 1. Illustration of one-dimensional diffusive transport through a matrix across a layer of another material with thickness 2𝛿 confined between two planes located at 𝑧 = ±𝛿.
1 − 𝛼 + 𝑟c 𝑑 1−𝛼 1 = = + . 𝛿 𝛼 𝛼c 𝛼
(7)
We observe that d becomes negative when 𝛼 > 1 + 𝑟c or, equivalently, 𝛼 > 𝛼c ∕(𝛼c − 1). We may select an arbitrary elevation, h ≥ 𝛿, and define an effective conductivity based on the equation 𝑘ef f (ℎ) ℎ ℎ ≡𝛾 = . 𝑘matrix 𝑇matrix (ℎ) ℎ+𝑑
(8)
This would be the thermal conductivity recorded when the heat transport is driven by a temperature difference between two isothermal plates located at 𝑦 = ±ℎ. Substituting the expression for d given in (7), we obtain
2. Transport across a conductive layer
𝑘ef f (ℎ) 1 = ( 𝑘matrix 1 + 1−𝛼 𝛼 +
To establish a point a reference, we consider steady, one-dimensional heat transport through an infinite homogeneous conductive medium, termed the matrix, across a layer of another conductive medium with thickness 2𝛿, as shown in Fig. 1. The temperature, T(z), satisfies Laplace’s equation in the matrix and in the layer,
1 𝛼c
)
𝛿 ℎ
=
𝛼 𝛼 + (1 − 𝛼 + 𝑟c )
𝛿 ℎ
.
(9)
As h/𝛿 tends to unity, the right-hand side tends to 𝛼∕(1 + 𝑟c ); correspondingly, the effective diffusivity tends to 𝑘layer ∕(1 + 𝑟c ). As the ratio h/𝛿 tends to infinity, the fraction on the right-hand side tends to unity and the effective conductivity reduces to the matrix conductivity.
(1)
3. Periodic arrangement of cylinders
Far above and below the layer, as z → ± ∞, the temperature assumes the asymptotic form 𝑇 ≃ 𝛾 (𝑧 ± 𝑑),
(6)
inside the layer, |z| < 𝛿. An inconsequential constant can be added to the right-hand sides of (5) and (6). Comparing (5) with (2), we conclude that the profile displacement is given by
In Section 2, a framework for evaluating the effective conductivity is established by discussing an elementary one-dimensional model involving heat transport across a compact layer. The effective diffusivity is quantified in terms of the displacement of the temperature profile with respect to that prevailing in the absence of the layer. Having established this framework, in Section 3 we consider the temperature distribution and profile displacement due to an array of cylinders. The Laplace equation governing the temperature distribution will be solved in terms of single-layer representations involving the periodic Green’s function. The formulation results in a coupled system of two integrals equations of the second kind over the cylinder contours. Numerical results by standard boundary-element methods will be presented and discussed with reference to the compact layer configuration and with respect to analytical expressions for a dilute dispersion.
d2 𝑇 = 0. d𝑧 2
1 𝛾𝑧 𝛼
In the main part of this work, we consider heat transport in an infinite medium across a periodic array of cylinders with arbitrary crosssectional shapes separated by distance L along the x axis, as shown in Fig. 2. The temperature satisfies Laplace’s equation in the xy plane,
(2)
where 𝛾 is a specified gradient and d is a profile displacement to be determined as part of the solution. An inconsequential constant can be added to the right-hand side of (2) to render the temperature positive within a domain of interest, including the matrix and the layer. At the interfaces between the layer and the matrix, we require two conditions: (a) a jump in the temperature determined by the thermal
∇2 𝑇 = 0,
(10)
subject to a periodicity condition in the direction of the x axis, 𝑇 (𝑥 + 𝐿, 𝑦) = 𝑇 (𝑥, 𝑦). Far above and below the cylinder array, as y → ± ∞, the diffusive field assumes the asymptotic form shown in (2), where 𝛾 is a specified gradient and the profile displacement, d, must be found as part of the 106
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 2. Illustration of heat transport in the xy plane across a periodic arrangement of cylinders with arbitrary cross-sectional shapes. A temperature profile is shown to illustrate the significance of the displacement length, d. The cylinder interface exhibits a thermal contact resistance due to imperfect contact.
solution. We recall that an inconsequential constant can be added to the right-hand side of (2). For future reference, we introduce the ratio of the cylinder-to-matrix thermal conductivity, a dimensionless thermal contact conductivity, and a dimensionless thermal contact resistance. 𝛼=
𝑘cylinder 𝑘matrix
,
𝛼c =
𝑎𝑘c 𝑘matrix
,
𝑘cylinder 𝛼 𝑟c = = , 𝛼c 𝑎𝑘c
must be satisfied. Simplifying and using (17), we obtain 𝑇matrix (𝐱) ≃ 𝛾 (𝑦 ± 𝑑), where 𝑑=
(11)
3.1. Integral representations
3.2. Integral equations
where 𝜙(𝓁) is a dimensionless strength density distribution, is the contour of any one cylinder, 𝓁 is the arc length around , (𝐱, 𝐱′ ) is the periodic Green’s function of Laplace’s equation in two dimensions representing the field due to a periodic collection of point sources,
Requiring continuity of flux at the particle surface, and employing the properties of the Green’s function, we obtain the equation 𝑛𝑦 (𝐱 ) −
(𝐱 , 𝐱 ′ ) ≃ −
|𝐱 − 1 ln 2𝜋 𝐿
.
1 |𝑦 − 𝑦′ | + ⋯ , 2𝐿
(14)
𝑛𝑦 (𝐱 ) = (15)
where
where the dots indicate exponentially decaying terms. Using the properties of the Green’s function, we find that, as y → ± ∞, the integral representation yields the far field ( ) 1 𝑇matrix (𝐱) ≃ 𝛾 𝑦 ± 𝜙(𝐱′ )(𝑦′ − 𝑦) d𝓁(𝐱′ ) , (16) 2𝐿 ∮
𝜒(𝐱) ≡
𝜙(𝐱) d𝓁(𝐱) = 0
1 2
PV
∮
𝜙(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ )
𝜓(𝐱) +
PV
∮
) 𝜓(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ) ,
(21)
1 2
1+𝛼 𝜒(𝐱) − ∮ 1−𝛼
𝜙(𝐱) + 𝛼𝜓(𝐱) , 1+𝛼
PV
𝜁 (𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ),
𝜁 (𝐱 ) ≡
𝜙(𝐱) − 𝛼𝜓(𝐱) 1−𝛼
(22)
(23)
are auxiliary functions. Conversely, 𝜙(𝐱) =
where the plus sign applies far above the array and the minus sign applies far below the array. In order for the integral on the right-hand side of (16) to tend to a constant far from the cylinder array, as y → ± ∞, the constraint ∮
𝜙(𝐱) +
where the point x lies at the cylinder contour, the gradient 𝛁 operates with respect to the evaluation point, x, and PV denotes the principalvalue integral. Rearranging, we obtain the integral equation
As |𝑦 − 𝑦′ | → ∞, we obtain a linear asymptotic behavior, (𝐱 , 𝐱 ′ ) ≃ −
1 2
( = 𝛼 𝑛𝑦 (𝐱 ) +
(13)
and 𝜚 = 2𝜋∕𝐿 is the wave number (e.g., [13]). As x → x′, we recover the free-space Green’s function given by 𝐱′ |
(19)
where 𝜓(𝓁) is a dimensionless strength density distribution. The two single-layer distributions, 𝜙 and 𝜓, must be adjusted to satisfy the cylinder-matrix interfacial conditions.
The temperature field in the matrix can be represented by a singlelayer harmonic potential, ( ) 𝑇matrix (𝐱) = 𝛾 𝑦 + 𝜙(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) , (12) ∮
( ) 1 (𝐱 , 𝐱 ) = − ln 2 {cosh[𝜚(𝑦 − 𝑦′ )] − cos[𝜚(𝑥 − 𝑥′ )]} , 4𝜋
1 𝜙(𝐱) 𝑦 d𝓁(𝐱) 2𝐿 ∮
is the profile displacement. Inside the cylinders the temperature field can be represented by another single-layer harmonic potential, ( ) 𝑇cylinder (𝐱) = 𝛾 𝑦 + 𝜓(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) , (20) ∮
where kc is the thermal contact conductivity and a is a designated cylinder radius.
′
(18)
(1 + 𝛼) 𝜒(𝐱) + (1 − 𝛼) 𝜁 (𝐱) , 2
𝜓(𝐱) =
(1 + 𝛼) 𝜒(𝐱) − (1 − 𝛼) 𝜁 (𝐱) , 2𝛼 (24)
and 𝜙(𝐱) − 𝜓(𝐱) =
) 𝛼2 − 1 ( 𝜒(𝐱) − 𝜁 (𝐱) . 2𝛼
(25)
If the pair of functions 𝜒 and 𝜁 is available, the corresponding pair 𝜙 and 𝜓 can be reconstructed, and vice versa.
(17) 107
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 3. Distribution of the single-layer density, 𝜙, for a periodic arrangement of circular cylinders with scaled radius (a) 𝑎∕𝐿 = 0.10 and (b) 0.49, for 𝛼c = 1.0 and 𝛼 = 0.05 (dots connected with dashed lines), 0.1, 0.2, 0.4, and 0.8. The dashed lines in (a) represent the analytical solution for isolated cylinders given in (31).
The interfacial condition for the jump in the temperature requires that
Half of the first fraction on the right-hand side is equal to the areal fraction of an interfacial band of thickness 2a occupied by the periodic array of cylinders, 𝑐𝑎 = 𝜋𝑎2 ∕(2𝑎𝐿). For an arbitrary collection of small cylinders inside an layer of thickness 2𝛿, we obtain
( ) 𝑘𝑐 𝜙(𝐱′ ) − 𝜓(𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) 𝑘matrix ∮ = 𝑛𝑦 (𝐱 ) −
1 2
𝜙(𝐱) +
PV
∮
𝜙(𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ),
1 − 𝛼 + 𝑟c 𝑑 𝑐 , ≃2 𝛿 1 + 𝛼 + 2𝑟c 𝑎
(26)
where the point x lies at the cylinder contour. In terms of the auxiliary functions 𝜒 and 𝜁 , we obtain the integral equation (1 + 𝛼) 𝜒(𝐱) + (1−𝛼) 𝜁 (𝐱) −
PV
−
(
∮
where ca is the cylinder areal fraction. The associated effective conductivity introduced in (9) is
( ) 1−𝛼 2 𝑘c 𝜒(𝐱′ )−𝜁 (𝐱′ ) (𝐱, 𝐱′ ) d𝓁(𝐱′ ) 𝛼 𝑘matrix ∮ (27)
𝐷ef f (ℎ) ℎ ℎ ≡𝛾 = ≃ 𝐷m 𝑇 (ℎ) ℎ+𝑑 1+2
∮
𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱) = − 1 , 2
𝐷ef f (𝛿) 1 − 𝛼 + 𝑟c 𝛿 = ≃1−2 𝑐 . 𝐷m 𝛿+𝑑 1 + 𝛼 + 2𝑟c 𝑎
∫
𝜒 d𝓁 = −(1 − 𝛼)
(35)
3.4. Numerical method
∫
𝜁 d𝓁,
The integral Eqs (22) and (28) were solved by a standard boundaryelement method involving particle contour discretization into straight segments and subtraction of the singularity of the single-layer potential in terms of the free-space Green’s function to ensure high accuracy (e.g., [13]). The double-layer potential is identically zero when evaluated at a point on a flat boundary element. The twelve-point Gauss–Legendre quadrature was employed to compute the non-singular and regularized element influence coefficients of the double-layer potential. Results of computations with 256 boundary elements generate results that are accurate to the third significant figure for nearly touching cylinders.
(30)
3.3. Small cylinders When the cylinder size is small compared to the separation, each cylinder can be studied as though it were in isolation. The solution of the integral equations for small circular cylinders of radius a in the limit as a/L → 0 is given by 𝛼c (1 − 𝛼) + 𝛼 1 − 𝛼 + 𝑟c sin 𝜃, sin 𝜃 = 2 𝛼c (1 + 𝛼) + 2𝛼 1 + 𝛼 + 2𝑟c
3.5. Results and discussion Numerical solutions of the integral equations for small circular cylinders with scaled radius 𝑎∕𝐿 = 0.10 are shown in Fig. 3(a) for several values of 𝛼 and dimensionless contact conductivity 𝛼c = 1.0. The numerical results are in excellent agreement with the analytical solution given in (31) for isolated cylinders, represented by the dashed lines. Corresponding solutions of the integral equations for nearly touching circular cylinders with scaled radius 𝑎∕𝐿 = 0.49 are shown in Fig. 3(b). These distributions differ markedly from the sinusoidal distributions for small
(31)
where 𝜃 is the polar angle measured around the cylinder center. Substituting this asymptotic solution into (19) and performing the integration, we obtain the scaled profile displacement 𝑑 𝜋𝑎 1 − 𝛼 + 𝑟c = . 𝑎 𝐿 1 + 𝛼 + 2𝑟c
(34)
(29)
which demonstrates that the integral constraint (17) is satisfied for any strength density distribution, 𝜙(𝓁), arising from the solution of the integral equation.
𝜙(𝜃) = 2
.
This formula generalizes results obtained by previous authors in the absence of contact resistance (e.g., [14]).
we find that (1 + 𝛼)
1 1−𝛼+𝑟c 𝛿 𝑐 1+𝛼+2𝑟c ℎ 𝑎
For ℎ = 𝛿, we obtain
) (1 + 𝛼) 𝜒(𝐱′ ) + (1−𝛼) 𝜁 (𝐱′ ) 𝐧(𝐱) ⋅ 𝛁(𝐱, 𝐱′ ) d𝓁(𝐱′ ) = 2 𝑛𝑦 (𝐱). (28)
Integrating both sides of the integral Eq. (22) with respect to arc length, 𝓁, around the cylinder contour, , and using the identity 𝑃𝑉
(33)
(32)
108
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
Fig. 4. Scaled profile displacement, d/a, for a periodic array of circular cylinders with scaled radius 𝑎∕𝐿 = 0.01 ( × ), 0.10 (∗ ), 0.25 ( + ), and 0.49 (○), where L is the cylinder separation, and for (a) 𝛼c = 1 and (b) 10. The dotted lines represent the analytical solution for small cylinders given in (32). The solid line represents the predictions of the layer model discussed in Section 2 for layer half-thickness 𝛿 = 𝑎.
Fig. 5. Scaled profile displacement, d/a, and scaled effective conductivity, keff /kmatrix , for a periodic array of circular cylinders for scaled contact conductance (a, c) 𝛼c = 1 and (b, d) 10.
109
C. Pozrikidis
Engineering Analysis with Boundary Elements 85 (2017) 105–110
or well separated cylinders. It is worth noting that the distribution of 𝜙 over the main part of the cylinders is nearly flat. The scaled profile displacement, d/a, depends on three dimensionless parameters: the scaled cylinder radius, a/L, the cylinder to matrix conductivity ratio, 𝛼, and the contact conductivity ratio, 𝛼 c . Graphs of the scaled profile displacement against 𝛼 are shown in Fig. 4 for circular cylinders with scaled radius 𝑎∕𝐿 = 0.01, 0.10, 0.25, and 0.49, and contact conductance 𝛼c = 1 and 10. The numerical results for small cylinders are in excellent agreement with the analytical solution given in (32), represented by the dotted lines. In fact, the agreement is fair for cylinder radius as high as 𝑎∕𝐿 = 0.25. As the cylinders tend to touch, a/L → 0.5, the numerical results tend to those predicted by the layer model discussed in Section 2, represented by the solid line for layer half-thickness 𝛿 = 𝑎 . The layer model predicts that 𝛼 (1 − 𝛼) + 𝛼 𝑑 = c , 𝑎 𝛼𝛼c
conductive properties of the array. In practice, the contact conductance may be increased by interfacial pressure mediated by isotropic, uniaxial, or a more general type of imposed strain or stress. 4. Discussion We have presented a framework for computing the effective conductivity of a matrix hosting a layer of another material or a periodic array of cylinders in terms of the displacement of the temperature profile far from the layer or cylinders. The formulation accounts for interfacial contact resistance due to imperfect contact. The integral formulation and boundary-element methods for two-dimensional configurations discussed in this paper can be extended to describe periodic arrays of spheres in terms of the doubly periodic Green’s function of Laplace’s equation. A practical concern is that the computation of the Green’s function is expensive, even when expedited methods based in Ewald’s summation technique are employed. Typical formulations involving a periodic box where Dirichlet conditions are imposed at the top and bottom and periodicity conditions are imposed in two directions are a less elegant yet more tractable alternative.
(36)
independent of a/L. As 𝛼 tends to infinity, d/a tends to (1 − 𝛼c )∕𝛼c , which is zero when 𝛼c = 1 and negative when 𝛼c = 10, in agreement with the numerical results presented in Fig. 4. Expression (36) can be contrasted with that for small cylinders shown in (32). A linear combination of the layer and small-cylinder results can be used to describe the profile displacement in engineering applications where a high degree of accuracy is not required. A summary of numerical results is presented in Fig. 5(a) and (b) in the form of three-dimensional graphs of d/a against a/L and 𝛼 for 𝛼c = 1 and 10. Though the general behavior of the two top graphs presented in this figure is similar, significant quantitative differences are observed, as seen previously in Fig. 4. We note in particular, that the profile displacement is positive when 𝛼c = 1, indicating a reduction in the effective conductivity with respect to that of the matrix due to dominant effect of the contact resistance, even for large 𝛼. Physically, high-conductivity cylinders loosely embedded in the matrix act as near insulators. In contrast, the profile displacement may become negative when 𝛼c = 10, indicating an enhancement in the effective conductivity. This behavior is illustrated directly in corresponding graphs of the effective conductivity, computed as 𝑘ef f 𝑎 ≡ 𝑘matrix 𝑎+𝑑
References [1] Öchsner A, Murch GE. Heat transfer in multi-phase materials. New York: Springer; 2011. [2] Wang T-P, Kang D-Y. Predictions of effective diffusivity of mixed matrix membranes with tubular fillers. J Membr Sci 2015;485:123–31. [3] Pozrikidis C, Ford DM. Conductive transport through a mixed-matrix membrane. J Eng Math 2017. (In Press). [4] Maxwell JC. Electricity and magnetism. Oxford: Clarendon Press; 1873. [5] Jeffrey DJ. Conduction through a random suspension of spheres. Proc R Soc Lond A 1873;335:355–67. [6] Benveniste Y. Effective thermal conductivity of composites with a thermal contact resistance between the constituents: non-dilute case. J Appl Phys 1987;61:2840–3. [7] Cheng H, Torquato S. Effective conductivity of periodic arrays of spheres with interfacial resistance. Proc R Soc Lond A 1997;453:145–61. [8] Duschlbauer H, Pettermann E, Böhm HJ. Heat conduction of a spheroidal inhomogeneity with imperfectly bonded interface. J Appl Phys 2003;94. https://doi.org/10.1063/1.1587886. [9] Matt CF, Cruz ME. Heat conduction in two-phase composite materials with threedimensional microstructures and interfacial thermal resistance. In: Öchsent A, M̎urch GE, editors. Heat transfer in multi-phase materials, vol. 2; 2010. p. 63–97. [10] Achenbach JD, Kitahara M, Mikata Y, Sotiropoulos DA. Reflection and transmission of plane waves by a layer of compact inhomogeneities. Pure Appl Geophys 1988;128:101–18. [11] Achenbach JD. Use of elastodynamic reciprocity theorems for field calculations. In: S̎chiavone P, Constanda C, Mioduchowski A, editors. Integral methods in science and engineering. Birkhäuser; 2002. [12] Scarpetta E, Sumbatyan A. Wave propagation through a periodic array of inclined cracks. Eur J Mech A/Solids 2000;19:949–59. [13] Pozrikidis C. A practical guide to boundary-element methods with the software library BEMLIB. Boca Raton: Chapman & Hall/CRC; 2002. [14] Hashin Z, Shtrikman S. A variational approach to the theory of the effective magnetic permeability of multiphase materials. J Appl Phys A 1962;33:3125–31.
(37)
according to (8), as shown in Fig. 5(c) and (d). When 𝛼c = 1, the effective conductivity is less than unity due to the low contact conductance. When 𝛼c = 10, the effective conductivity may be higher than unity for sufficiently high 𝛼 due to the moderate contact resistance. These results emphasize the importance of interfacial asperity interlocking. Mechanical degradation may loosen the cylinders and significantly affect the
110