- Email: [email protected]
Fluid flow characteristics of vascularized channel networks
Chemical Engineering Science 65 (2010) 6270–6281
Contents lists available at ScienceDirect
Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces
Fluid flow characteristics of vascularized channel networks Kee-Hyeon Cho a,b, Moo-Hwan Kim a,n a b
Department of Mechanical Engineering, Pohang University of Science and Technology, San 31, Hyoja-dong, Namgu, Pohang, Kyungbuk 790-784, Republic of Korea RIST (Research Institute of Industrial Science and Technology), Energy and Resources Research Department, Namgu, Pohang, Kyungbuk 790-600, Republic of Korea
a r t i c l e in f o
a b s t r a c t
Article history: Received 31 May 2010 Received in revised form 9 September 2010 Accepted 10 September 2010 Available online 17 September 2010
This paper reports the fluid flow characteristics of vascularized channel networks. To validate our vascular designs by an analytical approach, three-dimensional numerical works were performed. The numerical work covered the Reynolds number range of 2–1000, cooling channels volume fraction of 0.02, pressure drop range of 10–10 000 Pa, and six flow configurations: first, second, and third constructal structures with optimized hydraulic diameters (D1 and D2) and non-optimized hydraulic diameter (D) for each system size 10 10, 20 20 and 50 50, respectively. In these cases, the objective was to compare global flow resistance and mass flow rate distribution of the analytical solutions with those of numerical solutions subject to a fixed volume and a fixed pressure drop. This paper shows that the fluid flow performance of the second constructs is superior to that of the first and third constructs when the system size exceeds 20 20. The difference in flow resistance performance between the optimized and non-optimized structures was found to increase and manifests itself clearly as the system size increases. Results also reveal that flow uniformity become desirable with increasing the system size, and that the third construct configurations have better flow uniformity than the other architectures among the optimized and non-optimized channel configurations. The analytical results are also compared with numerical data and good agreement between the numerical and analytical results is found. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Mass transfer Optimization Laminar flow Hydrodynamics Constructal Self-healing
1. Introduction In a polymer electrolyte membrane fuel cell (PEMFC) system, during a high power density operation, up to 50% of the energy produced is heat. Without proper thermal management, fuel cell performance may not be maintained. A cell temperature that is too high may lead to membrane dehydration, and a cell temperature that is too low may result in water condensation or flooding (Larminie and Dicks, 2000). Thus, it is important to minimize the pressure drop in the flow of coolant through the cooling plates, as the pumping power reduces the overall efficiency of the system and results in a more non-uniform temperature, both of which are crucial to PEMFC performance (Lasbet et al., 2007). Numerous investigators have studied thermal and flow characteristics in micro-channels and mini-channels. Gosselin and Bejan (2005) proposed the optimization of fluid networks based on the minimization of the pumping power requirement. This work used the concept of constructal theory to develop a design tool for minimizing the pumping power required to join arbitrary points. Pence (2002) developed the fractal-like channel network to yield a 60% lower pressure drop for the same total flow
n
Corresponding author. Tel.: + 82 54 279 2165; fax: + 82 54 279 3199. E-mail address: [email protected] (M.-H. Kim).
0009-2509/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2010.09.020
rate and a 30 1C lower wall temperature under identical pumping power conditions with identical total convective surface areas. Ryu et al. (2003) found that channel width and depth were more crucial than others in relation to heat-sink performance among various design variables. They also discovered that optimal dimensions and corresponding thermal resistance have a powerlaw dependence on the pumping power. In general, the use of micro-channels may improve thermal performance due to its compactness, whereas, it may also lead to disadvantages of increased pumping power and uneven temperature distribution. Therefore, it is desirable to increase the temperature uniformity of the wall while decreasing the pressure drop to overcome the problems mentioned above. Note, however, that achieving better temperature uniformity and lower pumping power across the compact devices, such as the cooling plates of PEMFC or micro-electromechanical systems (MEMS), simultaneously is a difficult challenge. In recent years, the optimization technology in multi-scale heat sinks or cooling plates has attracted increasing interest. Several types of effective flow architectures are being proposed, because the design of flow architecture has a significant effect on the operating performance of engineering equipment. An important segment of this new literature is based on a general principle: the constructal law (Bejan, 2000; Bejan and Lorente, 2004, 2008), which was stated by Bejan in 1996 as follows: ’’For a finite-size system to persist in time (to live), it must evolve in such a
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
way that it provides easier access to the imposed currents that flow through it.’’ Constructal law focuses attention on the relationship between the architecture of the flow system and its global performance. The constructal law and the use of tree architecture were first proposed as a problem of pure heat conduction, and later were extended to structures for convective fins, fluid flow, and heat transfer. The constructal theory is a mental mindset wherein the generation of flow structures that occur everywhere in nature (e.g. river basins, lungs, atmospheric circulation, and vascular tissues) can be reasoned based on the evolutionary principle of the increase in flow access over time. The constructal law is the time direction of the ‘‘movie’’ of successive configurations (Bejan, 1997, 2000; Bejan and Lorente, 2004, 2008). When the designs evolve toward configurations that flow more easily, every design (an image) is replaced by an image that flows more easily. The succession of such images can be seen as the frames of a movie. The movie tape runs toward images that flow more easily. There is a growing volume of research on the development of flow architectures using the constructal law (e.g., Brod, 2003; Hernandez et al., 2003; Kraus, 2003; Bejan et al., 2004; Jones and Ghassemi, 2004; Lundell et al., 2004; Senn and Poulikakos, 2004a, 2004b; Tondeur and Luo, 2004; Muzychka, 2005; Pramanick and Das, 2005). One active direction in engineering is the vascularization of smart materials, so that they may offer new or improved volumetric functions: self-cooling, self-healing, and variable transport properties (Wang et al., 2009). A body with embedded tree-shaped flows that bathe the entire volume is a vascularized body. Recently, Kim et al. (2006), Lee et al. (2008), and Wang et al. (2006) illustrated how the method can be used to develop vascular structures for self-healing materials, which are widely encountered in natural systems, and can be inspired from these systems for design. Biomimetic design approaches are also being used in the development of self-healing systems for polymer composites (White et al., 2001; Rudraiah and Ng, 2004; Toohey et al., 2007; Trask et al., 2007, 2008). An important new direction is the development of smart materials with self-healing functionality. A continuous supply of healing agents throughout the volume of the self-healing material calls for the vascularization of the entire volume, so that the entire volume of the structural
6271
composite is protected against volumetric cracking. Otherwise, several cracks may form randomly, and at different sites, simultaneously and repeatedly (Kim et al., 2006; Lee et al., 2008). The recent progress on self-healing techniques points toward composite structures with embedded vasculatures of channels filled with healing fluids. When cracks occur in the solid structure, the channels are ruptured and the vasculature delivers the healing fluid to the crack sites. The vasculature is more effective when it delivers the fluid faster and more uniformly regardless of the (random) position of the crack site. The present paper is about the flow characteristics of vasculatures of this kind, with minimal global flow resistance and minimal non-uniformity in how the fluid invades the volume. Tree-shaped flow configurations offer maximum access between one point (e.g. inlet or outlet) and an infinite number of points (e.g. area or volume). Many superimposed tree flows are accommodated by a grid of channels resembling the grid of city traffic. For volumes that must be bathed uniformly by a single stream flowing in and out, the recommended dendritic flow architecture consists of two trees matched canopy-to-canopy (Bejan, 2000) (Fig. 1(a)). Saber et al. (2009, 2010) discussed the geometrical design of various channel networks based on a multiscale approach and proposed rapid design of channel multi-scale networks with minimum flow maldistribution. New classes of multi-scale trees matched canopy to canopy are illustrated in Fig. 1(b) and (c). The second construct in Fig. 1(b) has blocks of parallel channels serving as canopies for two matched trees, and the third construct in Fig. 1(c) has four blocks of parallel channels. This higher order construct has greater complexity. In this paper, we consider three types of design on a square flat volume consisting of 10 10 volume elements: a first-level construct (Fig. 1(a)), a second-level construct (Fig. 1(b)), and a third-level construct (Fig. 1(c)). We analyze the hydrodynamic performance of the optimized constructal channel configurations, evaluate their performance as vascular network systems for selfhealing or self-cooling, and compare them with non-optimized channel configurations. We also validate their performance by comparisons between the analytical and numerical results. To achieve this, in addition to the structure with 10 10 elements, we compare the structures with 20 20 and 50 50 volume elements.
Fig. 1. Schematic diagram of the three-dimensional geometry of cooling plates with tree-shaped cooling channel with 10 10 elements: (a) first construct; (b) second construct; and (c) third construct.
6272
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
2. Geometry Consider the vascularized self-healing body consisting of a square slab measuring X Y and having the thickness W, where W is the dimension of the solid body in the direction perpendicular to the plane X Y, as illustrated in Fig. 1. The size of the square domain is measured in terms of N N, where N is the number of small square elements counted along one side (Cho et al., in press). When different flow structures are compared, X and Y are fixed, but W (i.e., elemental length d) varies by increasing the system size, N2 ¼ 10 10, 20 20, or 50 50. For example, Fig. 1 represents a square slab with the system size, N2 ¼10 10. The path to higher performance involves increasing the freedom to morph the flow architecture (Bejan, 2000). This entails adding more length scales that can be varied; in this case, two diameter sizes instead of one. The ratio of diameters (D1/D2) is the additional degree of freedom. Optimized multiple scales D1 and D2 will be distributed non-uniformly through the available flow volume, as illustrated in Fig. 1. We built the first, second and third constructs with two diameter sizes on a square domain with 10 10, 20 20, and 50 50 elements from our earlier work (Cho et al., in press). According to constructal theory (Bejan, 2000; Bejan and Lorente, 2004, 2008), the thinner channels (D1) should be placed in the canopy and the thicker (D2) in the stem and main branches, as illustrated in Fig. 1. The shape of the channel cross-section is assumed to be fixed from one channel to the next. Note that the effect of local pressure losses can be further reduced by smoothing the transitions between the subsequent channels of different sizes. However, they are perpendicular at a channel intersection, while maintaining the original cross-section of subsequent channels for simplicity (Fig. 1). We also define the corresponding non-optimized configurations, where channels have only one size, on a square domain composed of 10 10, 20 20, and 50 50 elements. There is one channel size D, with channel volume Vc, fixed on the same basis as the results of the optimized configurations. The detailed geometrical dimensions for each configuration are summarized in Table 1. The volume fraction occupied by all channels is held constant:
f¼
total channel volume Vc ¼ total volume V
3.1. Analytical model We previously explored new vascular designs for the volumetric bathing of smart structures with volumetric functionalities (self-healing, cooling). One stream bathes the volume, while the vasculature is configured as two trees matched canopy to canopy. Several architecture types are optimized: one channel size vs. two channel sizes, increasing complexity (first, second, and third constructs) and increasing size (up to 50 50 elemental volumes). A square domain composed of 10 10 elements for an analytical approach is illustrated in Cho et al. (in press). The centers of the elements are indicated with black circles. The only degree of freedom in the present study is the channel diameter ratio D1/D2. The assumptions on which the analytical analysis is based are fully developed laminar flow, negligible pressure losses at junctions or bends and fluid with constant properties. The flow resistance for the Poiseuille flow through a straight duct with a polygonal cross-section can be written as DP nL p2 Po ¼ ð3Þ 2 _ m A 8V where n is the kinematic viscosity, p is the perimeter of the crosssection, Po is the Poiseuille constant (for example, Po ¼16 for round cross-sections), and A is the area of the cross-section that is fixed, because the total duct volume V and the duct length L are fixed, namely A¼V/L (Bejan, 2000; Bejan and Lorente, 2008). The group p2Po/A depends only on n, which is the number of sides of the polygonal cross-section, and accounts for how this last geometric degree of freedom influences global performance (Bejan and Lorente, 2008). In summary, the pressure drop along a channel with length Li, _ i is estimated as diameter Di, and mass flow rate m
DPi ¼ C
external length scale ðXYÞ1=2 ¼ internal length scale Vc 1=3
ð2Þ
Table 1 Geometric dimensions for the constructal configurations (X ¼Y ¼100 mm). System size
Complexity Vc (mm3)
d (mm)
D (D1 ¼ D2)
D2 (mm)
D1 (mm)
10 10
1st 2nd 3rd
2000 2000 2000
10 10 10
1.529 1.476 1.543
2.338 1.979 1.843
1.295 1.045 1.065
7.9 7.9 7.9
20 20
1st 2nd 3rd
1000 1000 1000
5 5 5
0.780 0.764 0.782
1.500 1.282 1.160
0.663 0.537 0.531
10.0 10.0 10.0
1st 2nd 3rd
400 400 400
2 2 2
0.316 0.313 0.316
0.842 0.737 0.656
0.274 0.228 0.222
13.6 13.6 13.6
50 50
3. Analytical model and method
ð1Þ
where V is the total volume and Vc is the total channel flow volume. The void fraction f may be viewed as the average porosity of the entire structure, with the observation that, in these designs, the channel volume is not distributed uniformly. Another geometric property of the flow structure is svelteness Sv, the ratio of the external length scale ((XY)1/2) divided by the internal length scale: Sv ¼
Sv is a global property of the flow architecture, playing an important role in the evolution of the flow architecture towards the best, or near-best, architecture in a fixed space (near the ‘‘equilibrium flow configuration’’ performance level (Bejan, 2000)).
Sv
_ i Li m D4i
ð4Þ
where C is a constant factor, for example, C ¼128n/p, if the duct cross-section is round and i indicates the channel rank, where DPi denotes the corresponding pressure drop through a duct that can either be a channel, a distributor or a collector portion (Bejan, 2000; Kim et al., 2006; Bejan and Lorente, 2008; Lee et al., 2008). 3.2. Analytical method _ consists of invoking Eq. (4) for all The calculation of DP=m channels, accounting for the mass continuity at every junction, and adding all the DPi’s along one flow path, from the inlet to the outlet of the flow structure. The details of deriving the global optimized channel configurations are illustrated in an earlier work (Cho et al., in press) and not repeated here. The dimensionless global flow resistance is defined as
c¼
DP 2 3 f d _ Cm
ð5Þ
where d is the elemental length and f is the porosity of a slab (Kim et al., 2006; Bejan and Lorente, 2008; Lee et al., 2008). The dimensionless flow resistance parameter c, which can be related to the well-known Kleiber’s rule (Singer et al., 2004) because the specific (mass-related) pressure drop across the optimized
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
6273
Table 2 Comparison of dimensionless pressure drop and dimensionless mass flow rate between analytical and numerical solutions for the first, second and third construct, with two diameters. System size
Complexity
WP (Pa)
Be
c¼
DP _ Cm
f2 d 3
Analytical (a) 10 10
1st 2nd 3rd
20 20
1st 2nd 3rd
50 50
1st 2nd 3rd
8
min
~ M
(D1 a D2) Numerical (b)
ab (%) a
Analytical (c)
Numerical (d)
5
5
cd (%) c
1.4 10 1.3 1010 1.4 108 1.3 1010 1.4 108 1.3 1010
10 1,000 10 1,000 10 1,000
2.170 2.170 2.438 2.438 2.729 2.729
2.248 6.340 2.705 7.034 3.025 8.242
3.6 192.2 11.0 188.5 10.8 202.0
1.546 10 1.437 107 1.372 105 1.231 107 1.224 105 1.103 107
1.493 10 4.920 106 1.236 105 4.266 106 1.105 105 3.652 106
3.4 65.8 9.9 65.3 9.7 66.9
4.1 108 3.9 1010 4.1 109 3.9 1011 4.1 109 3.9 1011
30 3,000 30 3,000 30 3,000
1.884 1.884 1.848 1.848 2.184 2.184
2.088 4.248 2.071 4.371 2.298 5.099
10.8 125.5 12.1 136.5 5.2 133.5
2.685 105 2.556 107 2.731 105 2.518 107 2.309 105 2.115 107
2.423 105 1.134 107 2.437 105 1.064 107 2.194 105 9.060 106
9.8 55.6 10.8 57.7 5.0 57.2
1.4 109 1.3 1011 1.4 109 1.3 1011 1.4 109 1.3 1011
100 10,000 100 10,000 100 10,000
1.527 1.527 1.230 1.230 1.477 1.477
1.555 2.365 1.361 1.909 1.630 2.348
1.8 54.9 10.7 55.2 10.4 59.0
4.426 105 4.325 107 5.487 105 5.184 107 4.569 105 4.304 107
4.346 105 2.792 107 4.958 105 3.339 107 4.140 105 2.707 107
1.8 35.4 9.6 35.6 9.4 37.1
channel networks increases with decreasing the system size (N2), as illustrated in Table 2, represents the ratio of the overall pressure drop to the mass flow rate across the slab in the case of the same size and porosity. The c value calculated for each assumed configuration has a single value. For example, the c value for Fig. 1(a) reaches the minimum value 2.17, see Table 2. Here the corresponding mass flow rate for each configuration can be determined from Eq. (5), see Table 2. The pressure loss DP is calculated by summing up all flow imperfections (Bejan and Lorente, 2008) X X 4L 1 1 DP ¼ f rU 2 þ K rU 2 ð6Þ Dh 2 2 d l d l where Dh and r are the hydraulic diameter and fluid density, respectively. The first summation refers to friction losses caused by wall friction due to the viscosity (both molecular and turbulent) of water in motion, and results from a momentum transfer between the molecules of adjacent fluid layers moving at different velocities. The geometrical factor f only depends on the shape of the channel. In this work, the roughness of the channel walls is ignored. The second summation in Eq. (6) refers to local losses, such as pressure drops caused by junctions, fittings, valves, inlets, outlets, enlargements, and contractions. Each local loss contributes to the total loss in proportion to KrU2/2, where K is the respective local loss coefficient. In general, local losses can be non-negligible in the functioning of complex tree-shaped networks. In other words, the friction losses gain in relative importance to the local (junctions) losses, as the svelteness of the flow architecture increases.
4. Numerical model and method 4.1. Numerical model The fluid flow performance of the constructal channel architectures are numerically simulated using a model for 3-dimensional fluid flow for each configuration. To simplify the numerical simulation, only the channel part of the vascular network body is included in the computational domain.
Self-healing agent or coolant is provided by an embedded three-dimensional channel network. It is pumped into the selfhealing network with a specified pressure at the inlet, as illustrated in Figs. 1 and 2. This pressure must be high enough to overcome the pressure losses when the self-healing network system is under operation. The porosity f is fixed at 0.02 for all configurations. The numerical work covers the overall pressure drop range DP¼10–10 000 Pa, which corresponds to the mass flow rate range 4.38 10 7–1.63 10 3 kg/s. The material properties of water used in this study are: density r ¼998.2 kg/m3 and dynamic viscosity m ¼8.51 10 5 kg/m s. To focus on the effect of the optimized and non-optimized channel configurations on the vascular network performance, the following assumptions are made:
The fluid flow is steady-state, isothermal, and 3-dimensional. The fluid is Newtonian and single phase, while the flow is laminar. All the fluid properties are constant. Based on the above assumptions, the governing equations for mass and momentum were solved numerically in the fluid domains of the cooling channels, as follows:
Mass conservation: @u @v @w þ þ ¼0 @x @y @z
ð7Þ
Momentum equation: ! @u @u @u @p @2 u @2 u @2 u ¼ þm r u þv þw þ þ 2 @x @y @z @x @x2 @y2 @z
r u
! @v @v @v @p @2 v @2 v @2 v þv þw ¼ þm þ þ @x @y @z @y @x2 @y2 @z2
! @w @w @w @p @2 w @2 w @2 w þv þw ¼þm r u þ 2 þ 2 @x @y @z @z @x2 @y @z
ð8Þ
ð9Þ
ð10Þ
Velocity components in the x, y, and z directions in the coordinate system were designated by u, v, and w, respectively.
6274
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 2. Schematic diagram showing a flow path for pressure distribution (10 10 elements): (a) first construct; (b) second construct; and (c) third construct.
The variables are defined in the Nomenclature. The work of numerically solving Eqs. (7)–(10) was based on a dimensionless formulation using the variables: x~ ¼
x , L
u~ ¼
u m , L DP
y , L
z~ ¼
v~ ¼
v m , L DP
y~ ¼
z , L
n~ ¼ ~ ¼ w
n L
w m L DP
ðPPout Þ P~ ¼ DP
ð11Þ
ð12Þ
ð13Þ
~ v~ and w ~ are where L is a reference length (L¼X¼Y) and u, ~ y, ~ z~ directions, while DP¼ the velocity components in the x, Pin Pout. The resulting dimensionless equations include:
Mass conservation: ~ @u~ @v~ @w þ þ ¼0 @x~ @y~ @z~
ð14Þ
Momentum equation: ! @u~ @u~ @u~ 1 @p~ 1 @2 u~ @2 u~ @2 u~ ~ þ v~ þw ¼ þ u~ þ þ @x~ @y~ @z~ Be @x~ Be @x~ 2 @y~ 2 @z~ 2
ð15Þ
! @v~ @v~ @v~ 1 @p~ 1 @2 v~ @2 v~ @2 v~ ~ þ v~ þw ¼ þ u~ þ þ 2 @x~ @y~ @z~ Be @y~ Be @x~ 2 @y~ 2 @z~
ð16Þ
! ~ ~ ~ ~ ~ ~ @w @w @w 1 @p~ 1 @2 w @2 w @2 w ~ þ v~ þw ¼ þ þ þ u~ 2 @x~ @y~ @z~ Be @z~ Be @x~ 2 @y~ 2 @z~ ð17Þ
natural convection. The Be value is fixed, because DP is fixed in this study. No-slip boundary conditions are applied at the channel surfaces. The channel surfaces are impermeable with no slip. The overall pressure drop across the self-healing networks is defined as DP¼Pin Pout, where Pin and Pout are the inlet and outlet pressures, respectively. The lowest pressure (Pout ¼0) is maintained as constant and uniform over the outlet planes of all channels. The resulting pressure drops across the cooling channels were such that the Reynolds number based on the inlet channel diameter was in the 2–1000 range of Reynolds numbers.
4.2. Numerical method and grid independence Computations were performed using a finite-volume package (Fluent Inc., 2006) with a pressure-based solver, node-based gradient evaluation, Semi-Implicit Method for Pressure-Linked Equations (SIMPLE) algorithm for pressure–velocity coupling and second order upwind scheme for momentum equations. In the grid generation, tetra grids were adopted. The independence of the solution, with respect to the grid size, was checked by examining the values of the mass flow rate and pressure drop between the inlet and outlet, and the skin friction coefficients for each geometrical configuration. Grid independence tests were carried out for the six configurations presented in the preceding sections. The number of cells varied from case to case; for example, the smallest number was 5 106 for the optimized first construct with 10 10 elements. Convergence is achieved when the residuals for the mass and momentum equation were smaller than 10 6. In addition, we closely examined the difference between the incoming and outgoing mass flow. For all of the cases simulated in this study, the relative difference was kept within 0.01%.
The dimensionless pressure drop number (Be) is defined by Be ¼
DPL2
mn
ð18Þ
where m and n are the fluid dynamic viscosity and fluid kinematic viscosity of the fluid, water, respectively. The dimensionless pressure drop number (Be) that Bhattacharjee and Grosshandler (1988) and Petrescu (1994) termed the Bejan number is clearly a measure of the relative magnitude of the heat transfer and fluid friction irreversibility. The pressure drop number plays the same role in forced convection that the Rayleigh number plays in
5. Results and discussion 5.1. Hydrodynamic characteristics Three-dimensional computational works were carried out to predict the fluid flow characteristics of the new vascular designs for the volumetric bathing of the smart structures. To simplify the problems, only the fluid flow performance of vascular channels, rather than the self-healing system or body, is discussed.
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
The analytical work illustrated that the pressure drop through the optimized vascular self-healing system was relatively lower than that through the non-optimized channels, given the identical channel volume and mass flow inlet. To verify the validity of our analytical methodology, the analytical results of the optimized vascular designs were compared with the numerical results. The total pressure drop versus the mass flow rate was found to be one of the most important details in reference to the design of the compact devices or micro-electromechanical systems (MEMS). The pressure drop was defined as the difference of the area-weighted average static pressure of the inlet and the outlets as follows: RR RR P dA P dA RRoutlet ð19Þ DP ¼ RRinlet inlet dA outlet dA The local mean velocity at the inlet (U) was the mass-weighted average velocity and defined as RR rvin dA R ð20Þ U ¼ Rinlet inlet r dA where vn is the velocity vector through surface A at the inlet. The Reynolds number based on the inlet hydraulic diameter (Din) and the local mean velocity at the inlet (U) is defined as
rUDin m
Re ¼
ð21Þ
The mass flow rate was also non-dimensionalized by writing in the order of the magnitude sense that for Poiseuille flow in a tube with length Y and hydraulic diameter D the pressure drop scales as
DP
_ nY m D4
ð22Þ
The channel diameter scale is related to the porosity of the structure:
f
D2 Y XYW
ð23Þ
where X Y. From Eqs. (22) and (23) and the Be definition Eq. (18), the proper dimensionless group for the mass flow rate emerges as ~ ¼ M
_ mY
nrf2 W 2
ð24Þ
6275
where n is the kinematic viscosity, r is the density, f is the porosity of cooling channel in a slab, W is the thickness of vascularized body and Y is the length of the vascularized body. Fig. 3 demonstrates the results of the comparison of the dimensionless global flow resistance c between the analytical and numerical results. For the optimized designs, as illustrated in Fig. 3(a), the global flow resistance c decreases slowly as the system size (N2) increases. Note, however, that the global flow resistance c increases steeply as the system size (N2) increases for the non-optimized designs, as illustrated in Fig. 3(b). Larger architectures flowed more easily if they were configured according to the developed method. The dimensionless global flow resistance c for the optimized vascular designs also increased as the overall pressure drop DP increased when the local losses were considered in the numerical analysis. In particular, they increased sharply for the system size, 10 10, because of the small Sv (7.6) and resulting vortices. The analytical and numerical results were very close to each other when the pressure drop number (Be) was relatively small. For example, the minimal relative difference between the analytical and numerical results was approximately 3.6% in the vicinity of Be ffi 1.4 108 for the first construct with 10 10 elements, as illustrated in Table 2. However, the third construct with 10 10 elements led to a large discrepancy (202.0%) in the vicinity of Be ffi 1.3 1010. This phenomenon was also due to the high local pressure losses, when the pressure drop number (Be) increased. These results are summarized in Table 2, which details the performance of the optimized architecture of the present study on the square domains N N that increased in size all the way to 50 50. Table 3 illustrates the relative percentage difference between the dimensionless mass flow rates of the optimized and nonoptimized structures. In comparing the dimensionless mass flow rate of optimized and non-optimized constructs, we see that the design of the optimized constructs is superior to that of the non-optimized constructs. For example, the dimensionless mass flow rate of the optimized design is larger by 137.6% based on the non-optimized construct in the vicinity of Be ffi 1.4 108. The maximum relative percentage difference between the dimensionless mass flow rates of the optimized and non-optimized structures is also about 370.4% based on the non-optimized construct for 20 20 elements. This trend is manifested more clearly when the system size increases: the optimized second construct designs for 50 50 elements provide greater flow
Fig. 3. Comparison of the global flow resistances between the analytical and numerical results: (a) D1 a D2 and (b) D1 ¼ D2.
6276
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Table 3 Comparison of dimensionless mass flow rate between the optimized and non-optimized constructs from numerical solutions for the first, second and third construct, with two diameters and one diameter. System size
Complexity
ba b (%)
~ M
Be
Optimized (D1 aD2) (a) 10 10
1st 2nd 3rd
20 20
1st 2nd 3rd
50 50
1st 2nd 3rd
8
Non-optimized (D1 ¼D2) (b)
1.4 10 1.3 1010 1.4 108 1.3 1010 1.4 108 1.3 1010
5
1.493 10 4.920 106 1.236 105 4.266 106 1.105 105 3.652 106
6.284 104 2.200 106 5.812 104 2.195 106 6.740 104 2.416 106
137.6 123.6 112.7 94.4 63.9 51.2
4.1 108 3.9 1010 4.1 108 3.9 1010 4.1 108 3.9 1010
2.423 105 1.134 107 2.437 105 1.064 107 2.194 105 9.060 106
5.137 104 3.035 106 5.181 104 3.260 106 6.122 104 3.565 106
371.7 273.6 370.4 226.4 258.4 154.1
1.4 109 1.3 1011 1.4 109 1.3 1011 1.4 109 1.3 1011
4.346 105 2.792 107 4.958 105 3.339 107 4.140 105 2.707 107
3.223 104 2.685 106 3.346 104 2.919 106 3.876 104 3.260 106
1,248.4 939.9 1,381.8 1,043.9 968.1 730.4
access, i.e., larger by a factor of about 13 based on the nonoptimized construct in the vicinity of Be ffi1.4 109. These results are summarized in Table 3. Fig. 4 illustrates how the local static pressure along the centerline of the system (z¼ 0) varied under the identical flow rate conditions for the three system sizes (10 10, 20 20 and _ of 10 10 elements was 50 50), where the mass flow rate m fixed at 5 10 4 kg/s, 20 20 elements at 3 10-4 kg/s and 50 50 elements at 1 10-4 kg/s, respectively. For comparison, the flow path was defined in Fig. 2 by presenting a number along one path from the lower-left to the upper-right corner of the diagram in a clockwise direction. As illustrated in Fig. 5, ‘‘0’’ is the start point at the inlet, while the final point is ‘‘4’’ at the outlet for the first and second constructs and ‘‘6’’ for the third constructs. The non-dimensionalized length w~ ðw~ ¼ w=Lw Þwas also introduced to plot the curve. For example, the dimensional length Lw is 0.19 m for the non-optimized first construct with 10 10 elements. The constant pressure gradients observed in the referenced channel (i.e. range ‘‘1–2’’ for the first, second, and third constructs, range ‘‘2–3’’ for the second and third construct) exhibit convincing evidence of the fully developed flow, considering that the left part of the curve corresponds to the inlet position. The non-linear pressure distribution in the inlet part (position ‘‘1’’) was caused by the entrance effect and significantly contributes to the total pressure drop. In comparing the pressure distribution of the optimized and non-optimized constructs, we found that the design of the optimized constructs was superior to that of the non-optimized constructs. For example, the pressure drop through the optimized first constructs with 10 10 elements was approximately 261.6 Pa lower than that through the non-optimized first constructs, whose value was around 131.4 Pa. Considering these results, a lower pressure drop across each channel for the optimized constructs was anticipated to result in a better flow uniformity, while a higher pressure drop across each channel for the non-optimized constructs was anticipated to yield a worse flow uniformity. As observed in Fig. 4, there were clear pressure recoveries at the first bifurcation position (1) for the non-optimized first construct with 10 10 elements, as expected; 24.6 Pa after decreasing approximately 35.2 Pa from the inlet. The pressure
recovery at position (1) was caused by the sharp change and recovery of the flow direction. A similar phenomenon occurred at the outlet, due to the same reason, and the pressure gradient at the outlet was the highest along the flow path, which resulted in a relatively large pressure drop. We also found that the difference of pressure drops between the optimized and non-optimized vascular designs increased significantly as the system size (N2) increased. Later, a similar phenomenon in hydrodynamic characteristics was also found. The skin-friction coefficient Cf, based on the area-averaged wall shear stress, a non-dimensional parameter defined as the ratio of the wall shear stress and the reference dynamic pressure is defined as
tw 2 r 2 vb
Cf ¼ 1
ð25Þ
where tw is the area-averaged wall shear stress, and r and vb are the reference density and velocity, respectively. Notice that Eq. (25) is in terms of bulk conditions. These conditions imply performing a volume-weighted average, e.g. R Li ui dAi vb ¼ A ð26Þ Vc where i indicates the channel rank and ui denotes the corresponding fluid velocity through surface Ai and length Li of straight channels. Fig. 5 illustrates the skin friction coefficient versus the Reynolds number (Re), based on the inlet hydraulic diameter for the two system sizes (10 10 and 50 50). From Fig. 5, it can be observed that the skin friction coefficient Cf decreased with increasing Re number increases. In addition, the drag on the wall of the optimized constructs was much lower than that of the nonoptimized constructs, except for the vascular designs with 10 10 elements (Fig. 5). Furthermore, the skin friction coefficient Cf increased gradually at the high Reynolds number, due to the much higher mass flow rate, as compared to those of the nonoptimized designs. Again, with increasing flow rates, the pressure drop caused by the bend, or junction, becomes more relevant. Figs. 6 and 7 illustrate how the non-dimensional mass flow ~ varied with the imposed pressure difference (WP) for the rate (M) two system sizes (10 10 and 50 50). As seen in Fig. 6, it
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
6277
Fig. 4. Pressure distributions along one path (‘‘0–1–2–3–4–5’’ for the first and second constructs, and ‘‘0–1–2–3–4–5–6’’ for the third constructs) subject to a fixed mass _ ¼ 5 104 kg/s; (b) 20 20 elements at m _ ¼ 3 104 kg/s; and (c) 50 50 elements at m _ ¼ 1 104 kg/s. flow rate: (a) 10 10 elements at m
~ increases with the becomes evident that the global mass flow rate M pressure drop number (Be) and that the mass flow rate variation is non-linear with its gradient increasing with the pressure drop. The non-linear relationship of the pressure drop and mass flow rate can be attributed to the junctions and the developing of the fluid flow at each junction or bend. At junctions or bends, the flow is disturbed and secondary flow motions are initiated, which then decays as it passes along the straight channel, such that the flow tends to develop again before reaching the next junction. In other words, with increasing mass flow rates, the pressure drop caused by the bend or junction becomes more relevant, due to the higher mass flow rates for the configurations with 10 10 elements. In fact, when Sv was smaller than 10; the effect of the local pressure losses (e.g. junctions, entrances) was non-negligible (Lorente and Bejan, 2005). These Poiseuille-type losses result in the discrepancy between the analytical and numerical results. On the contrary, the pressure drop caused by the bend, or junction becomes negligible for system size with 50 50 elements. Consequently, the results of analytical and numerical approaches are in relatively good agreement (see Fig. 7).
Note that the third constructal designs for the non-optimized channel configurations provide smaller global flow resistance among all system sizes (Figs. 6(b) and 7(b)). This means that the design with more branching levels is a good design. Figs. 6(b) and 7(b) show that although with the additional pressure drop at bifurcations the overall pressure drop of the third constructs is smaller than that of the first and second constructal structures, the third constructs distribute the flow better than the first and second constructs. This is especially true when the svelteness number is greater than 10 in the cases with 20 20 and 50 50 elements. 5.2. Flow non-uniformity Flow maldistribution was evaluated by using the flow nonuniformity ratio in the mass flow rate distribution in the multiscale channels. The non-uniformity ratio r was defined as r¼
_i m _t m
ð27Þ
6278
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 5. Variation of skin friction coefficients with Reynolds number: (a) 10 10 elements and (b) 50 50 elements.
Fig. 6. Comparison of the dimensionless mass flow rate versus the pressure drop number (Be) between the analytical and numerical results for 10 10 elements: (a) D1 aD2 and (b) D1 ¼ D2.
Fig. 7. Comparison of the dimensionless mass flow rate versus the pressure drop number (Be) between the analytical and numerical results for 50 50 elements: (a) D1 aD2 and (b) D1 ¼ D2.
_ i is the mass flow rate, only along the thin channels (D1), where m _ t was the total value of mass flow rates along the thin and m channels (D1), with a given number of channels (No) on the square domains of Fig. 2. The variable i indicates the channel, where channel number (No) is 10, 9, and 4 for the first, second, and third _i constructs with 10 10 elements, respectively. The flow rates m belong to the channels that cross the transverse ascending from the left to the right on the square domains of each construct, as
illustrated in Fig. 2. The values of the mass flow rate in each channel are calculated from the integration of the flow velocity inside each tube as follows: Z mi ¼ rui dAi ð28Þ A
where ui is the velocity vector through surface Ai of thin channels _ i is the mass flow rate along the thin channels (D1). The (D1) and m
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
higher non-uniformity ratio r indicates a poor distribution of fluid flow along the thin channels (D1). The results of the flow maldistribution in the thin channels with increasing Be are illustrated in Figs. 8 and 9 for the two system sizes (10 10 and 50 50). In comparing the maldistribution of optimized and non-optimized constructs, we see that the design of the optimized constructs is also superior to that of the non-optimized constructs, as we discussed in the previous hydrodynamic characteristics. As the pressure drop number (Be) increases, for both the optimized and non-optimized configurations, the mass flow rate distribution of the first, second and third constructs becomes generally undesirable, due to the minor losses when Be ( or DP) increases (see Figs. 8 and 9); thus, resulting in a more nonuniform flow in the thin channels. On the other hand, the third construct configurations have a better flow uniformity than the other architectures among the optimized and non-optimized channel configurations, whereas the first construct architecture is not desirable, because the uniformity effect decreases rapidly as Be increases and the system size N2 increases, as illustrated in Figs. 8 and 9. Similarly, the comparison between the analytical and numerical results indicates that there is a positive agreement for both the optimized and non-optimized configurations when the
6279
pressure drop number (Be) was relatively small. On the contrary, the difference between the analytical and numerical results started to increase when increasing the overall pressure drop (DP), because of the increase in local losses.
6. Conclusions In the present study, we investigated the fluid flow performance of constructal architecture based on the mass flow rates _ for the given pressure drop (DP). Six flow configurations were (m) described in this paper, and include the first, second and third constructal structures with optimized hydraulic diameters (D1 and D2) and non-optimized hydraulic diameters (D) for each system size, 10 10, 20 20, and 50 50. The validations were conducted by comparing the analytical results with the numerical results using a three-dimensional CFD approach. Numerical works confirmed the validity of the optimized designs by the analytical approach from our earlier work. A good agreement between the analytical and numerical results was demonstrated; more specifically, the relative devia~ were approximately tions of the dimensionless mass flow rate (M) 1.8% for the first construct with 50 50 elements.
Fig. 8. Comparison of the flow uniformity between the analytical and numerical results for 10 10 elements: (a) first construct; (b) second construct; and (c) third construct.
6280
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 9. Comparison of the flow uniformity between the analytical and numerical results for 50 50 elements: (a) first construct; (b) second construct; and (c) third construct.
The results imply that the fluid flow performance of the first constructal structure was superior to that of the second and third constructal structures when N ¼10. The best structure is the second constructal structure, when N ¼20 and 50 in the optimized constructal configurations. We also found that the difference in flow resistance performance between the optimized and nonoptimized structures increases and manifests itself clearly as N increases. On the other hand, the best architecture in nonoptimized configurations is the third construct, across all working conditions. These results are attributed to the hydrodynamic advantage documented in Figs. 4, 6 and 7. In general, the flow uniformity becomes desirable as N increases. The third construct configurations have better flow uniformity than the other architectures among the optimized and non-optimized channel configurations, whereas the first construct architecture is not desirable, because the uniformity effect decreases rapidly as Be increases. The optimized constructal self-healing channels offer a much better flow uniformity and a lower pressure drop than the nonoptimized self-healing channels (conventional types with one diameter). Therefore, the new vascular designs can be used to design the vascular channel networks of self-healing or selfcooling, microvascular lab-on-a-chip systems or microelectronic cooling applications resulting in uniform velocity distributions
between the channels. Further work is also required to validate the flow characteristics of new vascular designs and to quantify the effectiveness of using them on the self-healing system. An experimental study is in progress to confirm the thermal performances of new vascular designs.
Nomenclature A Be d Cf Dh Di D1 D2 f K Li Lt _ m _i m
area, m2 pressure drop number, Eq. (18) elemental length scale, m skin friction coefficient, Eq. (25) hydraulic diameter, m, Eq. (6) channel diameter, m hydraulic diameter of thin channels, m hydraulic diameter of thick channels, m friction factor, Eq. (6) local loss coefficient, Eq. (6) length of ith channel, m total channel length, m, Table 1 mass flow rate, kg/s mass flow rate of the ith channel only along the thin channels (D1), kg/s
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
_t m ~ M N No P Pin Pout Po Re Sv u, v, w U V vb Vc x, y, z X, Y, W
total mass flow rate only along the thin channels (D1), kg/s dimensionless mass flow rate, Eq. (24) number of elemental sequences d d in one direction number of thin channels, Fig. 2 local pressure, Pa inlet pressure, Pa outlet pressure, Pa Poiseuille constant, Eq. (3) Reynolds number, Eq. (21) svelteness number, Eq. (2) velocity components, m s 1, Eqs. (7)–(10) local mean velocity at inlet, m s 1, Eq. (21) total volume, m3 volume-averaged bulk velocity, m s 1, Eq. (26) total channel flow volume, m3 Cartesian coordinates dimensions of vascularized unit, m, Fig. 1
Greek symbols
DP
r f
m n c
tw
pressure difference, Pa fluid density, kg m 3 porosity, equation (1) fluid dynamic viscosity, kg/s m fluid kinematic viscosity, m2 s 1 non-dimensional global flow resistance, Eq. (5) area-averaged wall shear stress, N m 2, Eq. (25)
Subscripts b i in min out
bulk channel rank inlet minimum outlet
Superscript
dimensionless
Acknowledgment The authors would like to gratefully thank Dr. A. Bejan for his fruitful discussions and helpful suggestions. References Bejan, A., 1997. Advanced Engineering Thermodynamics 2nd ed. Wiley, New York. Bejan, A., 2000. Shape and Structure: from Engineering to Nature. Cambridge University Press, UK. Bejan, A., Dincer, I., Lorente, S., Miguel, A.F., Reis, A.H., 2004. Porous and Complex Flow Structures in Modern Technologies. Springer, New York. Bejan, A., Lorente, S., 2004. The constructal law and the thermodynamics of flow systems with configuration. International Journal of Heat and Mass Transfer 47, 3203–3214. Bejan, A., Lorente, S., 2008. Design with Constructal Theory. Hoboken, NJ. Bhattacharjee, S., Grosshandler, W.L., 1988. The formation of a wall jet near a high temperature wall under microgravity environment. ASME Heat Transfer Division 96, 711–716. Brod, H., 2003. Residence time optimized choice of tube diameters and slit heights in distribution systems for non-Newtonian liquids. Journal of Non-Newtonian Fluid Mechanics 111, 107–125.
6281
Cho, K.H., Lee, J., Kim, M.H., Bejan, A.. Vascular design of constructal structures with low flow resistance and nonuniformity. International Journal of Thermal Sciences, in press, doi:10.1016/j.ijthermalsci.2010.07.009. Fluent Inc., 2006. FLUENT Version 6.3 User’s Guide. Fluent Inc., Lebanon, NH. Gosselin, L., Bejan, A., 2005. Tree networks for minimal pumping power. International Journal of Thermal Sciences 44, 53–63. Hernandez, G., Allen, J.K., Mistree, F., 2003. Platform design for customizable products as a problem of access in a geometric space. Engineering Optimization 35, 229–254. Jones, G.F., Ghassemi, S., 2004. Thermal optimization of a composite heat spreader. In: Proceedings of the ASME Heat Transfer/Fluids Engineering Summer Conference, Charlotte, NC. Kraus, A.D., 2003. Constructal theory and the optimization of fin arrays. In: Proceedings of the ASME International Mechanical Engineering Congress and Exposition, Washington, DC. Kim, S.W., Lorente, S., Bejan, A., 2006. Vascularized materials: tree-shaped flow architectures matched canopy to canopy. Journal of Applied Physics 100, 063525. Larminie, J., Dicks, A., 2000. Fuel Cell Systems Explained. John Wiley & Sons, New York. Lundell, F., Thonon, B., Gruss, J.A., 2004. Constructal networks for efficient cooling/ heating. In: Proceedings of the Second Conference on Microchannels and Minichannels, Rochester, NY. Lorente, S., Bejan, A., 2005. Svelteness, freedom to morph, and constructal multi-scale flow structures. International Journal of Thermal Sciences 44, 1123–1130. Lasbet, Y., Auvity, B., Castelain, C., Peerhossaini, H., 2007. Thermal and hydrodynamic performances of chaotic mini-channel: application to the fuel cell cooling. Heat Transfer Engineering 28 (8–9), 795–803. Lee, J., Kim, S., Lorente, S., Bejan, A., 2008. Vascularization with trees matched canopy to canopy: diagonal channels with multiple sizes. International Journal of Heat and Mass Transfer 51, 2029–2040. Muzychka, Y.S., 2005. Constructal design of forced convection cooled microchannel cooling plates and heat exchangers. International Journal of Heat and Mass Transfer 48, 3119–3127. Petrescu, S., 1994. Comments on the optimal spacing of parallel plates cooled by forced convection. International Journal of Heat and Mass Transfer 37, 1283. Pence, D.V., 2002. Reduced pumping power and wall temperature in microchannel heat sinks with fractal-like branching channel networks. Microscale Thermophysical Engineering 6, 319–330. Pramanick, A.K., Das, P.K., 2005. Heuristics as an alternative to variational calculus for the optimization of a class of thermal insulation systems. International Journal of Heat and Mass Transfer 48, 1851–1857. Ryu, J.H., Choi, D.H., Kim, S.J., 2003. Three-dimensional numerical optimization of a manifold microchannel heat sink. International Journal of Heat and Mass Transfer 46, 1553–1562. Rudraiah, N., Ng, C., 2004. A model for manufacture of nano-sized smart materials free from impurities. Current Science 86, 1076–1091. Senn, S.M., Poulikakos, D., 2004a. Laminar mixing, heat transfer and pressure drop in tree-like microchannel nets and their application for thermal management in polymer electrolyte fuel cells. Journal of Power Sources 130, 178–191. Senn, S.M., Poulikakos, D., 2004b. Tree network channels as fluid distributors constructing double-staircase polymer electrolyte fuel cells. Journal of Applied Physics 96, 842–852. Singer, D., Engelhardt, N., Ince, A., Hallmann, B., 2004. Body size effects on tissue metabolic rate and ischemia tolerance in neonatal rat and mouse hearts. Thermochimica Acta 422, 19–23. Saber, M., Commenge, J.M., Falk, L., 2009. Rapid design of channel multi-scale networks with minimum flow maldistribution. Chemical Engineering Science 48, 723–733. Saber, M., Commenge, J.M., Falk, L., 2010. Microreactor numbering-up in multiscale networks for industrial scale applications: impact of flow maldistribution on the reactor performances. Chemical Engineering Science 65, 372–379. Tondeur, D., Luo, L., 2004. Design and scaling laws of ramified fluid distributors by the constructal approach. Chemical Engineering Science 59, 1799–1813. Toohey, K.S., Sottos, N.R., Lewis, J.A., Moore, J.S., White, S.R., 2007. Self-healing materials with microvascular networks. Nature Materials 6, 581–585. Trask, R.S., Williams, G.J., Bond, I.P., 2007. Bioinspired self-healing of advanced composite structures using hollow glass fibres. Journal of Royal Society Interface 13, 363–371. Trask, R.S., Williams, H.R., Bond, I.P., 2008. Self-healing polymer composites: mimicking nature to enhance performance. Express Polymer Letters 2, 238–250. White, S.R., Sottos, N.R., Geubelle, P.H., Moore, J.S., Kessler, M.R., Sriram, S.R., Brown, E.N., Viswanathan, S., 2001. Autonomic healing of polymer composites. Nature 409, 794–797. Wang, K.M., Lorente, S., Bejan, A., 2006. Vascularized networks with two optimized channel sizes. Journal of Phyics D: Applied Physics 39, 3086–3096. Wang, K.M., Lorente, S., Bejan, A., 2009. Vascular materials cooled with grids and radial channels. International Journal of Heat and Mass Transfer 52, 1230–1239.
Contents lists available at ScienceDirect
Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces
Fluid flow characteristics of vascularized channel networks Kee-Hyeon Cho a,b, Moo-Hwan Kim a,n a b
Department of Mechanical Engineering, Pohang University of Science and Technology, San 31, Hyoja-dong, Namgu, Pohang, Kyungbuk 790-784, Republic of Korea RIST (Research Institute of Industrial Science and Technology), Energy and Resources Research Department, Namgu, Pohang, Kyungbuk 790-600, Republic of Korea
a r t i c l e in f o
a b s t r a c t
Article history: Received 31 May 2010 Received in revised form 9 September 2010 Accepted 10 September 2010 Available online 17 September 2010
This paper reports the fluid flow characteristics of vascularized channel networks. To validate our vascular designs by an analytical approach, three-dimensional numerical works were performed. The numerical work covered the Reynolds number range of 2–1000, cooling channels volume fraction of 0.02, pressure drop range of 10–10 000 Pa, and six flow configurations: first, second, and third constructal structures with optimized hydraulic diameters (D1 and D2) and non-optimized hydraulic diameter (D) for each system size 10 10, 20 20 and 50 50, respectively. In these cases, the objective was to compare global flow resistance and mass flow rate distribution of the analytical solutions with those of numerical solutions subject to a fixed volume and a fixed pressure drop. This paper shows that the fluid flow performance of the second constructs is superior to that of the first and third constructs when the system size exceeds 20 20. The difference in flow resistance performance between the optimized and non-optimized structures was found to increase and manifests itself clearly as the system size increases. Results also reveal that flow uniformity become desirable with increasing the system size, and that the third construct configurations have better flow uniformity than the other architectures among the optimized and non-optimized channel configurations. The analytical results are also compared with numerical data and good agreement between the numerical and analytical results is found. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Mass transfer Optimization Laminar flow Hydrodynamics Constructal Self-healing
1. Introduction In a polymer electrolyte membrane fuel cell (PEMFC) system, during a high power density operation, up to 50% of the energy produced is heat. Without proper thermal management, fuel cell performance may not be maintained. A cell temperature that is too high may lead to membrane dehydration, and a cell temperature that is too low may result in water condensation or flooding (Larminie and Dicks, 2000). Thus, it is important to minimize the pressure drop in the flow of coolant through the cooling plates, as the pumping power reduces the overall efficiency of the system and results in a more non-uniform temperature, both of which are crucial to PEMFC performance (Lasbet et al., 2007). Numerous investigators have studied thermal and flow characteristics in micro-channels and mini-channels. Gosselin and Bejan (2005) proposed the optimization of fluid networks based on the minimization of the pumping power requirement. This work used the concept of constructal theory to develop a design tool for minimizing the pumping power required to join arbitrary points. Pence (2002) developed the fractal-like channel network to yield a 60% lower pressure drop for the same total flow
n
Corresponding author. Tel.: + 82 54 279 2165; fax: + 82 54 279 3199. E-mail address: [email protected] (M.-H. Kim).
0009-2509/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2010.09.020
rate and a 30 1C lower wall temperature under identical pumping power conditions with identical total convective surface areas. Ryu et al. (2003) found that channel width and depth were more crucial than others in relation to heat-sink performance among various design variables. They also discovered that optimal dimensions and corresponding thermal resistance have a powerlaw dependence on the pumping power. In general, the use of micro-channels may improve thermal performance due to its compactness, whereas, it may also lead to disadvantages of increased pumping power and uneven temperature distribution. Therefore, it is desirable to increase the temperature uniformity of the wall while decreasing the pressure drop to overcome the problems mentioned above. Note, however, that achieving better temperature uniformity and lower pumping power across the compact devices, such as the cooling plates of PEMFC or micro-electromechanical systems (MEMS), simultaneously is a difficult challenge. In recent years, the optimization technology in multi-scale heat sinks or cooling plates has attracted increasing interest. Several types of effective flow architectures are being proposed, because the design of flow architecture has a significant effect on the operating performance of engineering equipment. An important segment of this new literature is based on a general principle: the constructal law (Bejan, 2000; Bejan and Lorente, 2004, 2008), which was stated by Bejan in 1996 as follows: ’’For a finite-size system to persist in time (to live), it must evolve in such a
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
way that it provides easier access to the imposed currents that flow through it.’’ Constructal law focuses attention on the relationship between the architecture of the flow system and its global performance. The constructal law and the use of tree architecture were first proposed as a problem of pure heat conduction, and later were extended to structures for convective fins, fluid flow, and heat transfer. The constructal theory is a mental mindset wherein the generation of flow structures that occur everywhere in nature (e.g. river basins, lungs, atmospheric circulation, and vascular tissues) can be reasoned based on the evolutionary principle of the increase in flow access over time. The constructal law is the time direction of the ‘‘movie’’ of successive configurations (Bejan, 1997, 2000; Bejan and Lorente, 2004, 2008). When the designs evolve toward configurations that flow more easily, every design (an image) is replaced by an image that flows more easily. The succession of such images can be seen as the frames of a movie. The movie tape runs toward images that flow more easily. There is a growing volume of research on the development of flow architectures using the constructal law (e.g., Brod, 2003; Hernandez et al., 2003; Kraus, 2003; Bejan et al., 2004; Jones and Ghassemi, 2004; Lundell et al., 2004; Senn and Poulikakos, 2004a, 2004b; Tondeur and Luo, 2004; Muzychka, 2005; Pramanick and Das, 2005). One active direction in engineering is the vascularization of smart materials, so that they may offer new or improved volumetric functions: self-cooling, self-healing, and variable transport properties (Wang et al., 2009). A body with embedded tree-shaped flows that bathe the entire volume is a vascularized body. Recently, Kim et al. (2006), Lee et al. (2008), and Wang et al. (2006) illustrated how the method can be used to develop vascular structures for self-healing materials, which are widely encountered in natural systems, and can be inspired from these systems for design. Biomimetic design approaches are also being used in the development of self-healing systems for polymer composites (White et al., 2001; Rudraiah and Ng, 2004; Toohey et al., 2007; Trask et al., 2007, 2008). An important new direction is the development of smart materials with self-healing functionality. A continuous supply of healing agents throughout the volume of the self-healing material calls for the vascularization of the entire volume, so that the entire volume of the structural
6271
composite is protected against volumetric cracking. Otherwise, several cracks may form randomly, and at different sites, simultaneously and repeatedly (Kim et al., 2006; Lee et al., 2008). The recent progress on self-healing techniques points toward composite structures with embedded vasculatures of channels filled with healing fluids. When cracks occur in the solid structure, the channels are ruptured and the vasculature delivers the healing fluid to the crack sites. The vasculature is more effective when it delivers the fluid faster and more uniformly regardless of the (random) position of the crack site. The present paper is about the flow characteristics of vasculatures of this kind, with minimal global flow resistance and minimal non-uniformity in how the fluid invades the volume. Tree-shaped flow configurations offer maximum access between one point (e.g. inlet or outlet) and an infinite number of points (e.g. area or volume). Many superimposed tree flows are accommodated by a grid of channels resembling the grid of city traffic. For volumes that must be bathed uniformly by a single stream flowing in and out, the recommended dendritic flow architecture consists of two trees matched canopy-to-canopy (Bejan, 2000) (Fig. 1(a)). Saber et al. (2009, 2010) discussed the geometrical design of various channel networks based on a multiscale approach and proposed rapid design of channel multi-scale networks with minimum flow maldistribution. New classes of multi-scale trees matched canopy to canopy are illustrated in Fig. 1(b) and (c). The second construct in Fig. 1(b) has blocks of parallel channels serving as canopies for two matched trees, and the third construct in Fig. 1(c) has four blocks of parallel channels. This higher order construct has greater complexity. In this paper, we consider three types of design on a square flat volume consisting of 10 10 volume elements: a first-level construct (Fig. 1(a)), a second-level construct (Fig. 1(b)), and a third-level construct (Fig. 1(c)). We analyze the hydrodynamic performance of the optimized constructal channel configurations, evaluate their performance as vascular network systems for selfhealing or self-cooling, and compare them with non-optimized channel configurations. We also validate their performance by comparisons between the analytical and numerical results. To achieve this, in addition to the structure with 10 10 elements, we compare the structures with 20 20 and 50 50 volume elements.
Fig. 1. Schematic diagram of the three-dimensional geometry of cooling plates with tree-shaped cooling channel with 10 10 elements: (a) first construct; (b) second construct; and (c) third construct.
6272
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
2. Geometry Consider the vascularized self-healing body consisting of a square slab measuring X Y and having the thickness W, where W is the dimension of the solid body in the direction perpendicular to the plane X Y, as illustrated in Fig. 1. The size of the square domain is measured in terms of N N, where N is the number of small square elements counted along one side (Cho et al., in press). When different flow structures are compared, X and Y are fixed, but W (i.e., elemental length d) varies by increasing the system size, N2 ¼ 10 10, 20 20, or 50 50. For example, Fig. 1 represents a square slab with the system size, N2 ¼10 10. The path to higher performance involves increasing the freedom to morph the flow architecture (Bejan, 2000). This entails adding more length scales that can be varied; in this case, two diameter sizes instead of one. The ratio of diameters (D1/D2) is the additional degree of freedom. Optimized multiple scales D1 and D2 will be distributed non-uniformly through the available flow volume, as illustrated in Fig. 1. We built the first, second and third constructs with two diameter sizes on a square domain with 10 10, 20 20, and 50 50 elements from our earlier work (Cho et al., in press). According to constructal theory (Bejan, 2000; Bejan and Lorente, 2004, 2008), the thinner channels (D1) should be placed in the canopy and the thicker (D2) in the stem and main branches, as illustrated in Fig. 1. The shape of the channel cross-section is assumed to be fixed from one channel to the next. Note that the effect of local pressure losses can be further reduced by smoothing the transitions between the subsequent channels of different sizes. However, they are perpendicular at a channel intersection, while maintaining the original cross-section of subsequent channels for simplicity (Fig. 1). We also define the corresponding non-optimized configurations, where channels have only one size, on a square domain composed of 10 10, 20 20, and 50 50 elements. There is one channel size D, with channel volume Vc, fixed on the same basis as the results of the optimized configurations. The detailed geometrical dimensions for each configuration are summarized in Table 1. The volume fraction occupied by all channels is held constant:
f¼
total channel volume Vc ¼ total volume V
3.1. Analytical model We previously explored new vascular designs for the volumetric bathing of smart structures with volumetric functionalities (self-healing, cooling). One stream bathes the volume, while the vasculature is configured as two trees matched canopy to canopy. Several architecture types are optimized: one channel size vs. two channel sizes, increasing complexity (first, second, and third constructs) and increasing size (up to 50 50 elemental volumes). A square domain composed of 10 10 elements for an analytical approach is illustrated in Cho et al. (in press). The centers of the elements are indicated with black circles. The only degree of freedom in the present study is the channel diameter ratio D1/D2. The assumptions on which the analytical analysis is based are fully developed laminar flow, negligible pressure losses at junctions or bends and fluid with constant properties. The flow resistance for the Poiseuille flow through a straight duct with a polygonal cross-section can be written as DP nL p2 Po ¼ ð3Þ 2 _ m A 8V where n is the kinematic viscosity, p is the perimeter of the crosssection, Po is the Poiseuille constant (for example, Po ¼16 for round cross-sections), and A is the area of the cross-section that is fixed, because the total duct volume V and the duct length L are fixed, namely A¼V/L (Bejan, 2000; Bejan and Lorente, 2008). The group p2Po/A depends only on n, which is the number of sides of the polygonal cross-section, and accounts for how this last geometric degree of freedom influences global performance (Bejan and Lorente, 2008). In summary, the pressure drop along a channel with length Li, _ i is estimated as diameter Di, and mass flow rate m
DPi ¼ C
external length scale ðXYÞ1=2 ¼ internal length scale Vc 1=3
ð2Þ
Table 1 Geometric dimensions for the constructal configurations (X ¼Y ¼100 mm). System size
Complexity Vc (mm3)
d (mm)
D (D1 ¼ D2)
D2 (mm)
D1 (mm)
10 10
1st 2nd 3rd
2000 2000 2000
10 10 10
1.529 1.476 1.543
2.338 1.979 1.843
1.295 1.045 1.065
7.9 7.9 7.9
20 20
1st 2nd 3rd
1000 1000 1000
5 5 5
0.780 0.764 0.782
1.500 1.282 1.160
0.663 0.537 0.531
10.0 10.0 10.0
1st 2nd 3rd
400 400 400
2 2 2
0.316 0.313 0.316
0.842 0.737 0.656
0.274 0.228 0.222
13.6 13.6 13.6
50 50
3. Analytical model and method
ð1Þ
where V is the total volume and Vc is the total channel flow volume. The void fraction f may be viewed as the average porosity of the entire structure, with the observation that, in these designs, the channel volume is not distributed uniformly. Another geometric property of the flow structure is svelteness Sv, the ratio of the external length scale ((XY)1/2) divided by the internal length scale: Sv ¼
Sv is a global property of the flow architecture, playing an important role in the evolution of the flow architecture towards the best, or near-best, architecture in a fixed space (near the ‘‘equilibrium flow configuration’’ performance level (Bejan, 2000)).
Sv
_ i Li m D4i
ð4Þ
where C is a constant factor, for example, C ¼128n/p, if the duct cross-section is round and i indicates the channel rank, where DPi denotes the corresponding pressure drop through a duct that can either be a channel, a distributor or a collector portion (Bejan, 2000; Kim et al., 2006; Bejan and Lorente, 2008; Lee et al., 2008). 3.2. Analytical method _ consists of invoking Eq. (4) for all The calculation of DP=m channels, accounting for the mass continuity at every junction, and adding all the DPi’s along one flow path, from the inlet to the outlet of the flow structure. The details of deriving the global optimized channel configurations are illustrated in an earlier work (Cho et al., in press) and not repeated here. The dimensionless global flow resistance is defined as
c¼
DP 2 3 f d _ Cm
ð5Þ
where d is the elemental length and f is the porosity of a slab (Kim et al., 2006; Bejan and Lorente, 2008; Lee et al., 2008). The dimensionless flow resistance parameter c, which can be related to the well-known Kleiber’s rule (Singer et al., 2004) because the specific (mass-related) pressure drop across the optimized
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
6273
Table 2 Comparison of dimensionless pressure drop and dimensionless mass flow rate between analytical and numerical solutions for the first, second and third construct, with two diameters. System size
Complexity
WP (Pa)
Be
c¼
DP _ Cm
f2 d 3
Analytical (a) 10 10
1st 2nd 3rd
20 20
1st 2nd 3rd
50 50
1st 2nd 3rd
8
min
~ M
(D1 a D2) Numerical (b)
ab (%) a
Analytical (c)
Numerical (d)
5
5
cd (%) c
1.4 10 1.3 1010 1.4 108 1.3 1010 1.4 108 1.3 1010
10 1,000 10 1,000 10 1,000
2.170 2.170 2.438 2.438 2.729 2.729
2.248 6.340 2.705 7.034 3.025 8.242
3.6 192.2 11.0 188.5 10.8 202.0
1.546 10 1.437 107 1.372 105 1.231 107 1.224 105 1.103 107
1.493 10 4.920 106 1.236 105 4.266 106 1.105 105 3.652 106
3.4 65.8 9.9 65.3 9.7 66.9
4.1 108 3.9 1010 4.1 109 3.9 1011 4.1 109 3.9 1011
30 3,000 30 3,000 30 3,000
1.884 1.884 1.848 1.848 2.184 2.184
2.088 4.248 2.071 4.371 2.298 5.099
10.8 125.5 12.1 136.5 5.2 133.5
2.685 105 2.556 107 2.731 105 2.518 107 2.309 105 2.115 107
2.423 105 1.134 107 2.437 105 1.064 107 2.194 105 9.060 106
9.8 55.6 10.8 57.7 5.0 57.2
1.4 109 1.3 1011 1.4 109 1.3 1011 1.4 109 1.3 1011
100 10,000 100 10,000 100 10,000
1.527 1.527 1.230 1.230 1.477 1.477
1.555 2.365 1.361 1.909 1.630 2.348
1.8 54.9 10.7 55.2 10.4 59.0
4.426 105 4.325 107 5.487 105 5.184 107 4.569 105 4.304 107
4.346 105 2.792 107 4.958 105 3.339 107 4.140 105 2.707 107
1.8 35.4 9.6 35.6 9.4 37.1
channel networks increases with decreasing the system size (N2), as illustrated in Table 2, represents the ratio of the overall pressure drop to the mass flow rate across the slab in the case of the same size and porosity. The c value calculated for each assumed configuration has a single value. For example, the c value for Fig. 1(a) reaches the minimum value 2.17, see Table 2. Here the corresponding mass flow rate for each configuration can be determined from Eq. (5), see Table 2. The pressure loss DP is calculated by summing up all flow imperfections (Bejan and Lorente, 2008) X X 4L 1 1 DP ¼ f rU 2 þ K rU 2 ð6Þ Dh 2 2 d l d l where Dh and r are the hydraulic diameter and fluid density, respectively. The first summation refers to friction losses caused by wall friction due to the viscosity (both molecular and turbulent) of water in motion, and results from a momentum transfer between the molecules of adjacent fluid layers moving at different velocities. The geometrical factor f only depends on the shape of the channel. In this work, the roughness of the channel walls is ignored. The second summation in Eq. (6) refers to local losses, such as pressure drops caused by junctions, fittings, valves, inlets, outlets, enlargements, and contractions. Each local loss contributes to the total loss in proportion to KrU2/2, where K is the respective local loss coefficient. In general, local losses can be non-negligible in the functioning of complex tree-shaped networks. In other words, the friction losses gain in relative importance to the local (junctions) losses, as the svelteness of the flow architecture increases.
4. Numerical model and method 4.1. Numerical model The fluid flow performance of the constructal channel architectures are numerically simulated using a model for 3-dimensional fluid flow for each configuration. To simplify the numerical simulation, only the channel part of the vascular network body is included in the computational domain.
Self-healing agent or coolant is provided by an embedded three-dimensional channel network. It is pumped into the selfhealing network with a specified pressure at the inlet, as illustrated in Figs. 1 and 2. This pressure must be high enough to overcome the pressure losses when the self-healing network system is under operation. The porosity f is fixed at 0.02 for all configurations. The numerical work covers the overall pressure drop range DP¼10–10 000 Pa, which corresponds to the mass flow rate range 4.38 10 7–1.63 10 3 kg/s. The material properties of water used in this study are: density r ¼998.2 kg/m3 and dynamic viscosity m ¼8.51 10 5 kg/m s. To focus on the effect of the optimized and non-optimized channel configurations on the vascular network performance, the following assumptions are made:
The fluid flow is steady-state, isothermal, and 3-dimensional. The fluid is Newtonian and single phase, while the flow is laminar. All the fluid properties are constant. Based on the above assumptions, the governing equations for mass and momentum were solved numerically in the fluid domains of the cooling channels, as follows:
Mass conservation: @u @v @w þ þ ¼0 @x @y @z
ð7Þ
Momentum equation: ! @u @u @u @p @2 u @2 u @2 u ¼ þm r u þv þw þ þ 2 @x @y @z @x @x2 @y2 @z
r u
! @v @v @v @p @2 v @2 v @2 v þv þw ¼ þm þ þ @x @y @z @y @x2 @y2 @z2
! @w @w @w @p @2 w @2 w @2 w þv þw ¼þm r u þ 2 þ 2 @x @y @z @z @x2 @y @z
ð8Þ
ð9Þ
ð10Þ
Velocity components in the x, y, and z directions in the coordinate system were designated by u, v, and w, respectively.
6274
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 2. Schematic diagram showing a flow path for pressure distribution (10 10 elements): (a) first construct; (b) second construct; and (c) third construct.
The variables are defined in the Nomenclature. The work of numerically solving Eqs. (7)–(10) was based on a dimensionless formulation using the variables: x~ ¼
x , L
u~ ¼
u m , L DP
y , L
z~ ¼
v~ ¼
v m , L DP
y~ ¼
z , L
n~ ¼ ~ ¼ w
n L
w m L DP
ðPPout Þ P~ ¼ DP
ð11Þ
ð12Þ
ð13Þ
~ v~ and w ~ are where L is a reference length (L¼X¼Y) and u, ~ y, ~ z~ directions, while DP¼ the velocity components in the x, Pin Pout. The resulting dimensionless equations include:
Mass conservation: ~ @u~ @v~ @w þ þ ¼0 @x~ @y~ @z~
ð14Þ
Momentum equation: ! @u~ @u~ @u~ 1 @p~ 1 @2 u~ @2 u~ @2 u~ ~ þ v~ þw ¼ þ u~ þ þ @x~ @y~ @z~ Be @x~ Be @x~ 2 @y~ 2 @z~ 2
ð15Þ
! @v~ @v~ @v~ 1 @p~ 1 @2 v~ @2 v~ @2 v~ ~ þ v~ þw ¼ þ u~ þ þ 2 @x~ @y~ @z~ Be @y~ Be @x~ 2 @y~ 2 @z~
ð16Þ
! ~ ~ ~ ~ ~ ~ @w @w @w 1 @p~ 1 @2 w @2 w @2 w ~ þ v~ þw ¼ þ þ þ u~ 2 @x~ @y~ @z~ Be @z~ Be @x~ 2 @y~ 2 @z~ ð17Þ
natural convection. The Be value is fixed, because DP is fixed in this study. No-slip boundary conditions are applied at the channel surfaces. The channel surfaces are impermeable with no slip. The overall pressure drop across the self-healing networks is defined as DP¼Pin Pout, where Pin and Pout are the inlet and outlet pressures, respectively. The lowest pressure (Pout ¼0) is maintained as constant and uniform over the outlet planes of all channels. The resulting pressure drops across the cooling channels were such that the Reynolds number based on the inlet channel diameter was in the 2–1000 range of Reynolds numbers.
4.2. Numerical method and grid independence Computations were performed using a finite-volume package (Fluent Inc., 2006) with a pressure-based solver, node-based gradient evaluation, Semi-Implicit Method for Pressure-Linked Equations (SIMPLE) algorithm for pressure–velocity coupling and second order upwind scheme for momentum equations. In the grid generation, tetra grids were adopted. The independence of the solution, with respect to the grid size, was checked by examining the values of the mass flow rate and pressure drop between the inlet and outlet, and the skin friction coefficients for each geometrical configuration. Grid independence tests were carried out for the six configurations presented in the preceding sections. The number of cells varied from case to case; for example, the smallest number was 5 106 for the optimized first construct with 10 10 elements. Convergence is achieved when the residuals for the mass and momentum equation were smaller than 10 6. In addition, we closely examined the difference between the incoming and outgoing mass flow. For all of the cases simulated in this study, the relative difference was kept within 0.01%.
The dimensionless pressure drop number (Be) is defined by Be ¼
DPL2
mn
ð18Þ
where m and n are the fluid dynamic viscosity and fluid kinematic viscosity of the fluid, water, respectively. The dimensionless pressure drop number (Be) that Bhattacharjee and Grosshandler (1988) and Petrescu (1994) termed the Bejan number is clearly a measure of the relative magnitude of the heat transfer and fluid friction irreversibility. The pressure drop number plays the same role in forced convection that the Rayleigh number plays in
5. Results and discussion 5.1. Hydrodynamic characteristics Three-dimensional computational works were carried out to predict the fluid flow characteristics of the new vascular designs for the volumetric bathing of the smart structures. To simplify the problems, only the fluid flow performance of vascular channels, rather than the self-healing system or body, is discussed.
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
The analytical work illustrated that the pressure drop through the optimized vascular self-healing system was relatively lower than that through the non-optimized channels, given the identical channel volume and mass flow inlet. To verify the validity of our analytical methodology, the analytical results of the optimized vascular designs were compared with the numerical results. The total pressure drop versus the mass flow rate was found to be one of the most important details in reference to the design of the compact devices or micro-electromechanical systems (MEMS). The pressure drop was defined as the difference of the area-weighted average static pressure of the inlet and the outlets as follows: RR RR P dA P dA RRoutlet ð19Þ DP ¼ RRinlet inlet dA outlet dA The local mean velocity at the inlet (U) was the mass-weighted average velocity and defined as RR rvin dA R ð20Þ U ¼ Rinlet inlet r dA where vn is the velocity vector through surface A at the inlet. The Reynolds number based on the inlet hydraulic diameter (Din) and the local mean velocity at the inlet (U) is defined as
rUDin m
Re ¼
ð21Þ
The mass flow rate was also non-dimensionalized by writing in the order of the magnitude sense that for Poiseuille flow in a tube with length Y and hydraulic diameter D the pressure drop scales as
DP
_ nY m D4
ð22Þ
The channel diameter scale is related to the porosity of the structure:
f
D2 Y XYW
ð23Þ
where X Y. From Eqs. (22) and (23) and the Be definition Eq. (18), the proper dimensionless group for the mass flow rate emerges as ~ ¼ M
_ mY
nrf2 W 2
ð24Þ
6275
where n is the kinematic viscosity, r is the density, f is the porosity of cooling channel in a slab, W is the thickness of vascularized body and Y is the length of the vascularized body. Fig. 3 demonstrates the results of the comparison of the dimensionless global flow resistance c between the analytical and numerical results. For the optimized designs, as illustrated in Fig. 3(a), the global flow resistance c decreases slowly as the system size (N2) increases. Note, however, that the global flow resistance c increases steeply as the system size (N2) increases for the non-optimized designs, as illustrated in Fig. 3(b). Larger architectures flowed more easily if they were configured according to the developed method. The dimensionless global flow resistance c for the optimized vascular designs also increased as the overall pressure drop DP increased when the local losses were considered in the numerical analysis. In particular, they increased sharply for the system size, 10 10, because of the small Sv (7.6) and resulting vortices. The analytical and numerical results were very close to each other when the pressure drop number (Be) was relatively small. For example, the minimal relative difference between the analytical and numerical results was approximately 3.6% in the vicinity of Be ffi 1.4 108 for the first construct with 10 10 elements, as illustrated in Table 2. However, the third construct with 10 10 elements led to a large discrepancy (202.0%) in the vicinity of Be ffi 1.3 1010. This phenomenon was also due to the high local pressure losses, when the pressure drop number (Be) increased. These results are summarized in Table 2, which details the performance of the optimized architecture of the present study on the square domains N N that increased in size all the way to 50 50. Table 3 illustrates the relative percentage difference between the dimensionless mass flow rates of the optimized and nonoptimized structures. In comparing the dimensionless mass flow rate of optimized and non-optimized constructs, we see that the design of the optimized constructs is superior to that of the non-optimized constructs. For example, the dimensionless mass flow rate of the optimized design is larger by 137.6% based on the non-optimized construct in the vicinity of Be ffi 1.4 108. The maximum relative percentage difference between the dimensionless mass flow rates of the optimized and non-optimized structures is also about 370.4% based on the non-optimized construct for 20 20 elements. This trend is manifested more clearly when the system size increases: the optimized second construct designs for 50 50 elements provide greater flow
Fig. 3. Comparison of the global flow resistances between the analytical and numerical results: (a) D1 a D2 and (b) D1 ¼ D2.
6276
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Table 3 Comparison of dimensionless mass flow rate between the optimized and non-optimized constructs from numerical solutions for the first, second and third construct, with two diameters and one diameter. System size
Complexity
ba b (%)
~ M
Be
Optimized (D1 aD2) (a) 10 10
1st 2nd 3rd
20 20
1st 2nd 3rd
50 50
1st 2nd 3rd
8
Non-optimized (D1 ¼D2) (b)
1.4 10 1.3 1010 1.4 108 1.3 1010 1.4 108 1.3 1010
5
1.493 10 4.920 106 1.236 105 4.266 106 1.105 105 3.652 106
6.284 104 2.200 106 5.812 104 2.195 106 6.740 104 2.416 106
137.6 123.6 112.7 94.4 63.9 51.2
4.1 108 3.9 1010 4.1 108 3.9 1010 4.1 108 3.9 1010
2.423 105 1.134 107 2.437 105 1.064 107 2.194 105 9.060 106
5.137 104 3.035 106 5.181 104 3.260 106 6.122 104 3.565 106
371.7 273.6 370.4 226.4 258.4 154.1
1.4 109 1.3 1011 1.4 109 1.3 1011 1.4 109 1.3 1011
4.346 105 2.792 107 4.958 105 3.339 107 4.140 105 2.707 107
3.223 104 2.685 106 3.346 104 2.919 106 3.876 104 3.260 106
1,248.4 939.9 1,381.8 1,043.9 968.1 730.4
access, i.e., larger by a factor of about 13 based on the nonoptimized construct in the vicinity of Be ffi1.4 109. These results are summarized in Table 3. Fig. 4 illustrates how the local static pressure along the centerline of the system (z¼ 0) varied under the identical flow rate conditions for the three system sizes (10 10, 20 20 and _ of 10 10 elements was 50 50), where the mass flow rate m fixed at 5 10 4 kg/s, 20 20 elements at 3 10-4 kg/s and 50 50 elements at 1 10-4 kg/s, respectively. For comparison, the flow path was defined in Fig. 2 by presenting a number along one path from the lower-left to the upper-right corner of the diagram in a clockwise direction. As illustrated in Fig. 5, ‘‘0’’ is the start point at the inlet, while the final point is ‘‘4’’ at the outlet for the first and second constructs and ‘‘6’’ for the third constructs. The non-dimensionalized length w~ ðw~ ¼ w=Lw Þwas also introduced to plot the curve. For example, the dimensional length Lw is 0.19 m for the non-optimized first construct with 10 10 elements. The constant pressure gradients observed in the referenced channel (i.e. range ‘‘1–2’’ for the first, second, and third constructs, range ‘‘2–3’’ for the second and third construct) exhibit convincing evidence of the fully developed flow, considering that the left part of the curve corresponds to the inlet position. The non-linear pressure distribution in the inlet part (position ‘‘1’’) was caused by the entrance effect and significantly contributes to the total pressure drop. In comparing the pressure distribution of the optimized and non-optimized constructs, we found that the design of the optimized constructs was superior to that of the non-optimized constructs. For example, the pressure drop through the optimized first constructs with 10 10 elements was approximately 261.6 Pa lower than that through the non-optimized first constructs, whose value was around 131.4 Pa. Considering these results, a lower pressure drop across each channel for the optimized constructs was anticipated to result in a better flow uniformity, while a higher pressure drop across each channel for the non-optimized constructs was anticipated to yield a worse flow uniformity. As observed in Fig. 4, there were clear pressure recoveries at the first bifurcation position (1) for the non-optimized first construct with 10 10 elements, as expected; 24.6 Pa after decreasing approximately 35.2 Pa from the inlet. The pressure
recovery at position (1) was caused by the sharp change and recovery of the flow direction. A similar phenomenon occurred at the outlet, due to the same reason, and the pressure gradient at the outlet was the highest along the flow path, which resulted in a relatively large pressure drop. We also found that the difference of pressure drops between the optimized and non-optimized vascular designs increased significantly as the system size (N2) increased. Later, a similar phenomenon in hydrodynamic characteristics was also found. The skin-friction coefficient Cf, based on the area-averaged wall shear stress, a non-dimensional parameter defined as the ratio of the wall shear stress and the reference dynamic pressure is defined as
tw 2 r 2 vb
Cf ¼ 1
ð25Þ
where tw is the area-averaged wall shear stress, and r and vb are the reference density and velocity, respectively. Notice that Eq. (25) is in terms of bulk conditions. These conditions imply performing a volume-weighted average, e.g. R Li ui dAi vb ¼ A ð26Þ Vc where i indicates the channel rank and ui denotes the corresponding fluid velocity through surface Ai and length Li of straight channels. Fig. 5 illustrates the skin friction coefficient versus the Reynolds number (Re), based on the inlet hydraulic diameter for the two system sizes (10 10 and 50 50). From Fig. 5, it can be observed that the skin friction coefficient Cf decreased with increasing Re number increases. In addition, the drag on the wall of the optimized constructs was much lower than that of the nonoptimized constructs, except for the vascular designs with 10 10 elements (Fig. 5). Furthermore, the skin friction coefficient Cf increased gradually at the high Reynolds number, due to the much higher mass flow rate, as compared to those of the nonoptimized designs. Again, with increasing flow rates, the pressure drop caused by the bend, or junction, becomes more relevant. Figs. 6 and 7 illustrate how the non-dimensional mass flow ~ varied with the imposed pressure difference (WP) for the rate (M) two system sizes (10 10 and 50 50). As seen in Fig. 6, it
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
6277
Fig. 4. Pressure distributions along one path (‘‘0–1–2–3–4–5’’ for the first and second constructs, and ‘‘0–1–2–3–4–5–6’’ for the third constructs) subject to a fixed mass _ ¼ 5 104 kg/s; (b) 20 20 elements at m _ ¼ 3 104 kg/s; and (c) 50 50 elements at m _ ¼ 1 104 kg/s. flow rate: (a) 10 10 elements at m
~ increases with the becomes evident that the global mass flow rate M pressure drop number (Be) and that the mass flow rate variation is non-linear with its gradient increasing with the pressure drop. The non-linear relationship of the pressure drop and mass flow rate can be attributed to the junctions and the developing of the fluid flow at each junction or bend. At junctions or bends, the flow is disturbed and secondary flow motions are initiated, which then decays as it passes along the straight channel, such that the flow tends to develop again before reaching the next junction. In other words, with increasing mass flow rates, the pressure drop caused by the bend or junction becomes more relevant, due to the higher mass flow rates for the configurations with 10 10 elements. In fact, when Sv was smaller than 10; the effect of the local pressure losses (e.g. junctions, entrances) was non-negligible (Lorente and Bejan, 2005). These Poiseuille-type losses result in the discrepancy between the analytical and numerical results. On the contrary, the pressure drop caused by the bend, or junction becomes negligible for system size with 50 50 elements. Consequently, the results of analytical and numerical approaches are in relatively good agreement (see Fig. 7).
Note that the third constructal designs for the non-optimized channel configurations provide smaller global flow resistance among all system sizes (Figs. 6(b) and 7(b)). This means that the design with more branching levels is a good design. Figs. 6(b) and 7(b) show that although with the additional pressure drop at bifurcations the overall pressure drop of the third constructs is smaller than that of the first and second constructal structures, the third constructs distribute the flow better than the first and second constructs. This is especially true when the svelteness number is greater than 10 in the cases with 20 20 and 50 50 elements. 5.2. Flow non-uniformity Flow maldistribution was evaluated by using the flow nonuniformity ratio in the mass flow rate distribution in the multiscale channels. The non-uniformity ratio r was defined as r¼
_i m _t m
ð27Þ
6278
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 5. Variation of skin friction coefficients with Reynolds number: (a) 10 10 elements and (b) 50 50 elements.
Fig. 6. Comparison of the dimensionless mass flow rate versus the pressure drop number (Be) between the analytical and numerical results for 10 10 elements: (a) D1 aD2 and (b) D1 ¼ D2.
Fig. 7. Comparison of the dimensionless mass flow rate versus the pressure drop number (Be) between the analytical and numerical results for 50 50 elements: (a) D1 aD2 and (b) D1 ¼ D2.
_ i is the mass flow rate, only along the thin channels (D1), where m _ t was the total value of mass flow rates along the thin and m channels (D1), with a given number of channels (No) on the square domains of Fig. 2. The variable i indicates the channel, where channel number (No) is 10, 9, and 4 for the first, second, and third _i constructs with 10 10 elements, respectively. The flow rates m belong to the channels that cross the transverse ascending from the left to the right on the square domains of each construct, as
illustrated in Fig. 2. The values of the mass flow rate in each channel are calculated from the integration of the flow velocity inside each tube as follows: Z mi ¼ rui dAi ð28Þ A
where ui is the velocity vector through surface Ai of thin channels _ i is the mass flow rate along the thin channels (D1). The (D1) and m
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
higher non-uniformity ratio r indicates a poor distribution of fluid flow along the thin channels (D1). The results of the flow maldistribution in the thin channels with increasing Be are illustrated in Figs. 8 and 9 for the two system sizes (10 10 and 50 50). In comparing the maldistribution of optimized and non-optimized constructs, we see that the design of the optimized constructs is also superior to that of the non-optimized constructs, as we discussed in the previous hydrodynamic characteristics. As the pressure drop number (Be) increases, for both the optimized and non-optimized configurations, the mass flow rate distribution of the first, second and third constructs becomes generally undesirable, due to the minor losses when Be ( or DP) increases (see Figs. 8 and 9); thus, resulting in a more nonuniform flow in the thin channels. On the other hand, the third construct configurations have a better flow uniformity than the other architectures among the optimized and non-optimized channel configurations, whereas the first construct architecture is not desirable, because the uniformity effect decreases rapidly as Be increases and the system size N2 increases, as illustrated in Figs. 8 and 9. Similarly, the comparison between the analytical and numerical results indicates that there is a positive agreement for both the optimized and non-optimized configurations when the
6279
pressure drop number (Be) was relatively small. On the contrary, the difference between the analytical and numerical results started to increase when increasing the overall pressure drop (DP), because of the increase in local losses.
6. Conclusions In the present study, we investigated the fluid flow performance of constructal architecture based on the mass flow rates _ for the given pressure drop (DP). Six flow configurations were (m) described in this paper, and include the first, second and third constructal structures with optimized hydraulic diameters (D1 and D2) and non-optimized hydraulic diameters (D) for each system size, 10 10, 20 20, and 50 50. The validations were conducted by comparing the analytical results with the numerical results using a three-dimensional CFD approach. Numerical works confirmed the validity of the optimized designs by the analytical approach from our earlier work. A good agreement between the analytical and numerical results was demonstrated; more specifically, the relative devia~ were approximately tions of the dimensionless mass flow rate (M) 1.8% for the first construct with 50 50 elements.
Fig. 8. Comparison of the flow uniformity between the analytical and numerical results for 10 10 elements: (a) first construct; (b) second construct; and (c) third construct.
6280
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
Fig. 9. Comparison of the flow uniformity between the analytical and numerical results for 50 50 elements: (a) first construct; (b) second construct; and (c) third construct.
The results imply that the fluid flow performance of the first constructal structure was superior to that of the second and third constructal structures when N ¼10. The best structure is the second constructal structure, when N ¼20 and 50 in the optimized constructal configurations. We also found that the difference in flow resistance performance between the optimized and nonoptimized structures increases and manifests itself clearly as N increases. On the other hand, the best architecture in nonoptimized configurations is the third construct, across all working conditions. These results are attributed to the hydrodynamic advantage documented in Figs. 4, 6 and 7. In general, the flow uniformity becomes desirable as N increases. The third construct configurations have better flow uniformity than the other architectures among the optimized and non-optimized channel configurations, whereas the first construct architecture is not desirable, because the uniformity effect decreases rapidly as Be increases. The optimized constructal self-healing channels offer a much better flow uniformity and a lower pressure drop than the nonoptimized self-healing channels (conventional types with one diameter). Therefore, the new vascular designs can be used to design the vascular channel networks of self-healing or selfcooling, microvascular lab-on-a-chip systems or microelectronic cooling applications resulting in uniform velocity distributions
between the channels. Further work is also required to validate the flow characteristics of new vascular designs and to quantify the effectiveness of using them on the self-healing system. An experimental study is in progress to confirm the thermal performances of new vascular designs.
Nomenclature A Be d Cf Dh Di D1 D2 f K Li Lt _ m _i m
area, m2 pressure drop number, Eq. (18) elemental length scale, m skin friction coefficient, Eq. (25) hydraulic diameter, m, Eq. (6) channel diameter, m hydraulic diameter of thin channels, m hydraulic diameter of thick channels, m friction factor, Eq. (6) local loss coefficient, Eq. (6) length of ith channel, m total channel length, m, Table 1 mass flow rate, kg/s mass flow rate of the ith channel only along the thin channels (D1), kg/s
K.-H. Cho, M.-H. Kim / Chemical Engineering Science 65 (2010) 6270–6281
_t m ~ M N No P Pin Pout Po Re Sv u, v, w U V vb Vc x, y, z X, Y, W
total mass flow rate only along the thin channels (D1), kg/s dimensionless mass flow rate, Eq. (24) number of elemental sequences d d in one direction number of thin channels, Fig. 2 local pressure, Pa inlet pressure, Pa outlet pressure, Pa Poiseuille constant, Eq. (3) Reynolds number, Eq. (21) svelteness number, Eq. (2) velocity components, m s 1, Eqs. (7)–(10) local mean velocity at inlet, m s 1, Eq. (21) total volume, m3 volume-averaged bulk velocity, m s 1, Eq. (26) total channel flow volume, m3 Cartesian coordinates dimensions of vascularized unit, m, Fig. 1
Greek symbols
DP
r f
m n c
tw
pressure difference, Pa fluid density, kg m 3 porosity, equation (1) fluid dynamic viscosity, kg/s m fluid kinematic viscosity, m2 s 1 non-dimensional global flow resistance, Eq. (5) area-averaged wall shear stress, N m 2, Eq. (25)
Subscripts b i in min out
bulk channel rank inlet minimum outlet
Superscript
dimensionless
Acknowledgment The authors would like to gratefully thank Dr. A. Bejan for his fruitful discussions and helpful suggestions. References Bejan, A., 1997. Advanced Engineering Thermodynamics 2nd ed. Wiley, New York. Bejan, A., 2000. Shape and Structure: from Engineering to Nature. Cambridge University Press, UK. Bejan, A., Dincer, I., Lorente, S., Miguel, A.F., Reis, A.H., 2004. Porous and Complex Flow Structures in Modern Technologies. Springer, New York. Bejan, A., Lorente, S., 2004. The constructal law and the thermodynamics of flow systems with configuration. International Journal of Heat and Mass Transfer 47, 3203–3214. Bejan, A., Lorente, S., 2008. Design with Constructal Theory. Hoboken, NJ. Bhattacharjee, S., Grosshandler, W.L., 1988. The formation of a wall jet near a high temperature wall under microgravity environment. ASME Heat Transfer Division 96, 711–716. Brod, H., 2003. Residence time optimized choice of tube diameters and slit heights in distribution systems for non-Newtonian liquids. Journal of Non-Newtonian Fluid Mechanics 111, 107–125.
6281
Cho, K.H., Lee, J., Kim, M.H., Bejan, A.. Vascular design of constructal structures with low flow resistance and nonuniformity. International Journal of Thermal Sciences, in press, doi:10.1016/j.ijthermalsci.2010.07.009. Fluent Inc., 2006. FLUENT Version 6.3 User’s Guide. Fluent Inc., Lebanon, NH. Gosselin, L., Bejan, A., 2005. Tree networks for minimal pumping power. International Journal of Thermal Sciences 44, 53–63. Hernandez, G., Allen, J.K., Mistree, F., 2003. Platform design for customizable products as a problem of access in a geometric space. Engineering Optimization 35, 229–254. Jones, G.F., Ghassemi, S., 2004. Thermal optimization of a composite heat spreader. In: Proceedings of the ASME Heat Transfer/Fluids Engineering Summer Conference, Charlotte, NC. Kraus, A.D., 2003. Constructal theory and the optimization of fin arrays. In: Proceedings of the ASME International Mechanical Engineering Congress and Exposition, Washington, DC. Kim, S.W., Lorente, S., Bejan, A., 2006. Vascularized materials: tree-shaped flow architectures matched canopy to canopy. Journal of Applied Physics 100, 063525. Larminie, J., Dicks, A., 2000. Fuel Cell Systems Explained. John Wiley & Sons, New York. Lundell, F., Thonon, B., Gruss, J.A., 2004. Constructal networks for efficient cooling/ heating. In: Proceedings of the Second Conference on Microchannels and Minichannels, Rochester, NY. Lorente, S., Bejan, A., 2005. Svelteness, freedom to morph, and constructal multi-scale flow structures. International Journal of Thermal Sciences 44, 1123–1130. Lasbet, Y., Auvity, B., Castelain, C., Peerhossaini, H., 2007. Thermal and hydrodynamic performances of chaotic mini-channel: application to the fuel cell cooling. Heat Transfer Engineering 28 (8–9), 795–803. Lee, J., Kim, S., Lorente, S., Bejan, A., 2008. Vascularization with trees matched canopy to canopy: diagonal channels with multiple sizes. International Journal of Heat and Mass Transfer 51, 2029–2040. Muzychka, Y.S., 2005. Constructal design of forced convection cooled microchannel cooling plates and heat exchangers. International Journal of Heat and Mass Transfer 48, 3119–3127. Petrescu, S., 1994. Comments on the optimal spacing of parallel plates cooled by forced convection. International Journal of Heat and Mass Transfer 37, 1283. Pence, D.V., 2002. Reduced pumping power and wall temperature in microchannel heat sinks with fractal-like branching channel networks. Microscale Thermophysical Engineering 6, 319–330. Pramanick, A.K., Das, P.K., 2005. Heuristics as an alternative to variational calculus for the optimization of a class of thermal insulation systems. International Journal of Heat and Mass Transfer 48, 1851–1857. Ryu, J.H., Choi, D.H., Kim, S.J., 2003. Three-dimensional numerical optimization of a manifold microchannel heat sink. International Journal of Heat and Mass Transfer 46, 1553–1562. Rudraiah, N., Ng, C., 2004. A model for manufacture of nano-sized smart materials free from impurities. Current Science 86, 1076–1091. Senn, S.M., Poulikakos, D., 2004a. Laminar mixing, heat transfer and pressure drop in tree-like microchannel nets and their application for thermal management in polymer electrolyte fuel cells. Journal of Power Sources 130, 178–191. Senn, S.M., Poulikakos, D., 2004b. Tree network channels as fluid distributors constructing double-staircase polymer electrolyte fuel cells. Journal of Applied Physics 96, 842–852. Singer, D., Engelhardt, N., Ince, A., Hallmann, B., 2004. Body size effects on tissue metabolic rate and ischemia tolerance in neonatal rat and mouse hearts. Thermochimica Acta 422, 19–23. Saber, M., Commenge, J.M., Falk, L., 2009. Rapid design of channel multi-scale networks with minimum flow maldistribution. Chemical Engineering Science 48, 723–733. Saber, M., Commenge, J.M., Falk, L., 2010. Microreactor numbering-up in multiscale networks for industrial scale applications: impact of flow maldistribution on the reactor performances. Chemical Engineering Science 65, 372–379. Tondeur, D., Luo, L., 2004. Design and scaling laws of ramified fluid distributors by the constructal approach. Chemical Engineering Science 59, 1799–1813. Toohey, K.S., Sottos, N.R., Lewis, J.A., Moore, J.S., White, S.R., 2007. Self-healing materials with microvascular networks. Nature Materials 6, 581–585. Trask, R.S., Williams, G.J., Bond, I.P., 2007. Bioinspired self-healing of advanced composite structures using hollow glass fibres. Journal of Royal Society Interface 13, 363–371. Trask, R.S., Williams, H.R., Bond, I.P., 2008. Self-healing polymer composites: mimicking nature to enhance performance. Express Polymer Letters 2, 238–250. White, S.R., Sottos, N.R., Geubelle, P.H., Moore, J.S., Kessler, M.R., Sriram, S.R., Brown, E.N., Viswanathan, S., 2001. Autonomic healing of polymer composites. Nature 409, 794–797. Wang, K.M., Lorente, S., Bejan, A., 2006. Vascularized networks with two optimized channel sizes. Journal of Phyics D: Applied Physics 39, 3086–3096. Wang, K.M., Lorente, S., Bejan, A., 2009. Vascular materials cooled with grids and radial channels. International Journal of Heat and Mass Transfer 52, 1230–1239.