- Email: [email protected]
Theoretical and experimental studies of microflows in silicon microchannels
Available online at www.sciencedirect.com
Materials Science and Engineering C 28 (2008) 910 – 917 www.elsevier.com/locate/msec
Theoretical and experimental studies of microflows in silicon microchannels L. Renaud ⁎, C. Malhaire, P. Kleimann, D. Barbier, P. Morin Institut des Nanotechnologies de Lyon, INL, Université de Lyon and INSA de Lyon, Lyon, F-69003, France; Université Lyon 1, Lyon, F-69622, France Available online 18 October 2007
Abstract The determination of fluid flows in silicon microchannels is important for the design of microfluidic systems. In this paper, experimental investigations on the characteristics of low fluid flows (few μl h− 1) in silicon trapezoidal microchannels (21 μm in depth, length and width ranging from 200 to 440 mm and 58 to 267 μm, respectively) are presented. The test-devices have been fabricated using micromachining technologies. A double KOH etching process has been used to achieve microchannels in (100)-oriented silicon wafers as well as deep in-plane cavities used for capillary connections. Silicon has been finally anodically bonded on Pyrex substrates. The experimental set-up, based on the measurement of a differential pressure and a liquid–air interface displacement in a gauged tube, is fully detailed in terms of fluidic connections and measurement principle. The experimental results are in good agreement with the Navier–Stokes theory, solved by a simple iterative method. However, finite element modelling has been used to study complex 3D problems that were found in the devices and the experimental set-up. Finally, we propose abacuses for three different channel cross-sections that may be used to easily compute the flow in a microchannel. © 2007 Elsevier B.V. All rights reserved. Keywords: Microfluidic; Microflow; Navier–Stokes equations; Abacuses; Finite element modelling
1. Introduction The determination of flow velocity profile in microchannels is important for the design of microfluidic systems. In fact, micromachined channels are essential parts in many applications from heat pipes to micro-pumps, μ-TAS and Lab-On-AChip devices [1]. With the device size reduction, many investigations of the liquid flow in microchannels, with different experimental methods, have been presented in the literature to verify the validity of macrofluidic flow theory to the small scale [2]. In fact, the microchannels are characterised by very large surface to volume ratios and surface phenomena that are typically ignored in larger conduits might become important at the microscale. Some significant flow characteristic deviations from theoretical predictions have been experimentally observed but some contradictions among these results appear, according to the work of I.P.I. Papautsky et al. [3]. For example, J. Pfahler et al. [4] have found that experiences carried out with microchannels having hydraulic diameter (Dh) of 76.1 μm were in rough agreement with the prediction from the Navier–Stokes equations. Significant deviations from the theory concerning the ⁎ Corresponding author. Tel.: +33 4 72 17 73 99. E-mail address: [email protected] (L. Renaud). 0928-4931/$ - see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.msec.2007.10.070
friction coefficient (up to 3 times) were only found for smaller channels (Dh = 1.6 and 3.3 μm). B. Xu et al. [5] studied microchannels with Dh ranging from 30 to 344 μm and observed a maximum 20% overestimated relationship between friction factor and Reynolds number. Even for water flow in nanometer depth channels, J.T. Cheng and Giordano validated the Poiseuille flow theory with difference less than 10% between experimental and theoretical flows [6]. On the contrary, Q. Weilin et al. [7,8] observed higher pressure gradient and flow friction factor (up to 40%) in trapezoidal microchannels (hydraulic diameter ranging from 51 to 169 μm) than those given by the conventional laminar flow theory. They proposed a roughness-viscosity model to interpret the experimental results. Those reported deviations from theory could be explained by temperature effects, pressure dependant properties or wall roughness, as noticed by H. Herwig and Hausner [9] or by experimental artifacts: errors due to the fabrication of the samples, fluid leakage, bubbles trapping, and dusts. In this paper, we present experimental investigations on the characteristics of fluid flows in KOH-etched trapezoidal microchannels (Dh ranging from 45.7 to 74.3 μm) for very low Reynolds number value (Re b 1). The experimental set-up is fully detailed and we show some practical aspects in microfluidic experiments like an original capillary connection
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
911
flow is always laminar under low pressure drop (typically a few bar) [7], leading to a unidirectional flow and a uniform absolute pressure in the cross-section. For a permanent flow and a fixed pressure drop ΔP between the inlet and the outlet of the channel, Eq. (1) simplifies in the cross-section to: DP ∂2 w ∂2 w þ 2 þ 2 ¼ 0: gL ∂y ∂z Fig. 1. Geometrical definition of the Navier–Stokes numerical solving.
and microflow measurements. A simple and proven iterative method for solving the Navier–Stokes equations in the case of a 2D flow is presented for a non-specialist audience [10]. Nevertheless, finite element modelling has been used to study the complex three-dimensional problems that were found in the experimental set-up. Finally we give abacuses for three different channel cross-sections in the Appendix that may be used to quickly compute the flow in a microchannel and proposed to support a microfluidic design work. 2. Theory of pressure-driven flows The relationship between the fluid velocity and the absolute pressure for an incompressible viscous liquid is given by the classical fluid dynamics theory and the well-known Navier– Stokes equation: ∂Yv þ ∂t
Y Y
vj
Y
Y
v ¼ j
P Y þ mD v q
ð1Þ
Where Y v stands for the fluid velocity vector with components (u, v, w) each expressed for a set of Euler components (x, y, z, t). P and ρ are the absolute pressure and the voluminal mass, respectively (Fig. 1). The kinematic viscosity ν is the ratio between the dynamic viscosity η and the voluminal mass ρ. In the case of a microfluidic horizontal straight channel (x-direction), the
ð2Þ
The boundary condition at the wall surface is the no-slip condition, assumed for hydrophilic walls (w(y,z) = 0 at the wall). Eq. (2) is a classical elliptic equation that can be solved using a second-order equidistant central difference scheme and a Successive Over Relaxation (SOR) algorithm. The detailed general explanations are given in [11]. Iterations have been performed until convergence of the calculated mean flow. An example of velocity profile in a trapezoidal microchannel crosssection is given in Fig. 2. Since the flow is laminar, the flow rate Q (defined as the product of the mean velocity by the cross-section) is always proportional to the applied pressure ΔP. By analogy with electrical Ohm's law, we can define the ratios between Q and ΔP as the microfluidic resistance R and the microfluidic conductance G: R¼
DP Q and G ¼ : Q DP
ð3Þ
For example, the calculations of R and G for the device A, lead to 50.2 mbar h μl− 1 and 19.9 μl h− 1 bar− 1 respectively. Let us note that, R or G have been determined by the least square method applied on four points to be less sensitive to the experimental artifacts described above. 3. Experimental procedure 3.1. Fabrication process The use of an original capillary connection leads to a two steps etching process in silicon. The test-devices were
Fig. 2. SOR algorithm result in the case of a trapezoidal channel with a 21.3 μm depth, 57.5 μm width, 44.9 cm length and 54.7° angle (device A, see Section 5) filled with water at 25 °C with a pressure drop equal to 1 bar.
912
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Fig. 3. Process flow diagram for the test-devices.
fabricated using micromachining technologies as outlined in Fig. 3. Deep cavities (for the capillaries) and trapezoidal microchannels (studied microchannels) were etched in (100)oriented silicon wafers. Long microchannels (up to 44 cm) can be achieved on a 1 cm2 area using a zigzag design. The cavities were used in the final device to connect microchannels and fused silica capillaries. The depth of cavities and microchannels being different, a double etching technique was used. First the silicon wafers were oxidized (800 nm thermal oxide). Then a 200 nm Si3N4 layer was deposited by plasma enhanced chemical vapor deposition. The Si3N4 and oxide layers were patterned using standard photolithography process in order to define microchannels and cavities (Fig. 3a and b). Bulk anisotropic etching in KOH was then used to form 150 μm deep cavities (Fig. 3c). The visible parts of the oxide layer were removed in buffered HF prior to a second etching in KOH used to form 20 μm deep microchannels (Fig. 3d). The remaining Si3N4 and oxide layers were removed and glass substrates (Hoya SD-2) were anodically bonded on the top of silicon wafers to close the channels. Finally the devices were cut with a die-saw and fused silica capillaries were glued inside the cavities (Fig. 3e).
moving out downstream of the channel [12], others use the heat conduction properties of the fluids [13]. In this work, we have measured the fluid flow rate by the displacement of a liquid–air interface in a gauged tube, placed upstream of the microchannel. This method allows adjusting the volume of measurement according to the expected fluid flow, as shown in Fig. 4. A large syringe (60 ml) filled with air is used for applying pressure on the liquid. The total pressure drop ΔPT is measured with 2% accuracy by means of a differential pressure sensor (Honeywell 26PC series) connected to the large syringe. A low volume microsyringe (Hamilton, from 10 μl to 250 μl) is used as gauged tube and connected to the microfluidic system. The flow rate Q is determined through the observation of the meniscus shifting in the syringe for a given time. The
3.2. Microflow measurements One of the main difficulties in the study of microfluidic systems consists in the measurement of low flow rates (typically in the range of few μl h− 1). The size reduction of the microchannel cross-section imposes new experimental techniques of flow measurement. Indeed, friction forces prohibit the use of conventional flow meter for this range of flows. Several original solutions have been developed during last years. Some of them consist in the measurement of the weight of fluid
Fig. 4. Experimental set-up.
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
913
4.1. Study of RC and RS, and experimental set-up validation From Navier–Stokes equations, the microfluidic resistance RC of a circular cross-section capillary is expressed as [14]: Rtube ¼
128 g L p D4
ð4Þ
where L and R are the capillary length and radius respectively. When water is used η is given by [15]: Fig. 5. Schematized view (not to scale) of the pressure drop along the fluid flow.
g¼qm¼
0:1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 120 þ 2:142824 T 8:435 þ 8078:4 þ ðT 8:435Þ2
ð5Þ measurement accuracy ranges from 0.1 μl to 5 μl depending on the microsyringe used. This experimental set-up allows us to measure flow rates from 10 μl h− 1 to 200 μl h− 1 in less than 1 h with accuracy better than 2% (due to the 0.1 μl graduation step of the gauged 10 μl microsyringe, for example). 4. Procedure for the hydrodynamic study of the microfluidic system As shown on the more detailed schematic representation of the experimental set-up (Fig. 5), all the microfluidic elements are in serial configuration (same flow Q), therefore the total microfluidic resistance RT is equal to the sum of each microfluidic resistance. The total pressure drop ΔPT (that is measured) can be split in: 1) ΔPS = RS ⁎ Q, the pressure drop through the gauged tube (Hamilton microsyringe). 2) ΔPC = RC ⁎ Q, the pressure drop through the capillary tube. 3) ΔPJ = RJ ⁎ Q, the pressure drop through the two capillarymicrochannel junctions. 4) ΔPM = RM ⁎ Q, the pressure drop through the microchannel.
where T is the temperature (in Celsius degrees). According to Eq. (4) and the large difference between the syringe and the capillary diameters (typically more than 1.5 mm and 50 μm respectively), RS can be neglected in comparison with the other microfluidic resistances. Moreover, it is obvious that the accuracy on the determination of RC is strongly dependent on that of the diameter. In this study, the inner diameter of the fused silica capillary tube that connects the syringe to the micromachined channel was (50.0 ± 3) μm, as given by the manufacturer (Polymicro Technologies). A diameter value of (52.5 ± 1.0) μm has been measured by means of an optical microscope (Olympus BH-2). In order to validate the experimental set-up, we measured the microfluidic resistance of a (52.5 ± 1.0) μm diameter and (19.9 ± 0.1) cm length capillary at 22.5 °C with pressure range 0– 1 bar. The calculated (from Eq. (4)) and measured microfluidic resistances are (2.8 ± 0.2) mbar h μl− 1 and (2.7 ± 0.1) mbar h μl− 1 respectively; with a 0.9990 correlation coefficient between the applied pressure and the measured flow. This result is in agreement with the literature [14]. 4.2. Study of RJ and RM
Some of these resistances can be calculated (RS and RC) from fluidic theory; others have to be estimated by numerical methods (RJ and RM).
It is obvious that a 2D fluidic model is preferable to a 3D approach in terms of calculation time and convergence limits.
Fig. 6. Geometry of a simulated channel (type C).
914
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Table 1 Geometrical parameters and microfluidic resistances of studied microchannels Type of device Microchannel
Capillary
Device
Temperature (°C) Depth (μm) Width (μm) Length (cm) Estimated RM in mbar h μl− 1 Diameter (μm) Length (cm) Estimated RC in mbar h μl− 1 Estimated RD = RC + RM in mbar h μl − 1 Measured RD in mbar h μl− 1 Correlation coefficient (hydrodynamic flow versus pressure)
However, 2D models are not appropriate to modelize RJ and RM because of the quite complex shapes of the capillary-channel junctions and microchannel bends. In order to work with a 2D approach, the influence of the bends can be taken into account by considering an equivalent straight channel length. This equivalent length obviously depends on the bend shape. Let us note that in this work, the KOH etching ensures a broadening of the channel cross-section in the bends. It is thus reasonable to assume that the equivalent straight channel length is lower than the bend length. Here we can define the bend length as two times the bend size. Fig. 6 shows the geometry (only the main dimensions are shown) of a bend in the case of device C (see Table 1), extracted from microscope views and profilometry measurements. As shown in this figure, the bend length is assumed to be 0.44 mm. A 3D view of the model is shown in Fig. 7a. Two long straight parts (1 mm length) have been added to the bend to define simple boundary conditions (unidirectional velocity and uniform pressure at entrance and exit). The previously proposed assumption has been checked by means of finite element simulation using ANSYS® software. For each device, 3D models have been developed assuming a no-slip condition at the wall. As the number of elements was limited for our software version, only a limited length of the structures was modeled. A scale factor has been applied to the pressure values to get simulated flow values that were close to those measured. During the postprocessing stage, the ANSYS® Parametric Design Language (APDL) has been used to compute the flow and then the microfluidic resistance starting from the velocity on each node of the trapezoidal microchannel exit. For a temperature of 27.6 °C, the 3D simulation gives a microfluidic resistance for a bend equal to 0.0056 mbar h μl− 1. This value has to be compared to the microfluidic resistance of a 0.44 mm length trapezoidal straight channel equal to 0.0114 mbar h μl− 1. This result confirms that the equivalent straight channel length of the bend is lower than the 0.44 mm bend length. Similar analyses have been performed for the other samples (A to E) and gave the same results. In conclusion, the microfluidic resistance for a zigzag channel may be surrounded by two values corresponding to straight channels with extreme lengths. This 2D simplification leads to a dispersion of the total channel length (Table 1).
A
B
C
D
E
25.3 ± 0.3 21.4 ± 0.1 58 ± 0.5 44.4 ± 0.5 48.8 ± 1.4 52.5 ± 1.0 14.2 ± 0.2 1.9 ± 0.2 50.7 ± 1.6 47.5 ± 2.7 0.9964
28.0 ± 0.3 21.4 ± 0.1 94 ± 0.5 44.1 ± 0.7 19.6 ± 0.7 52.5 ± 1.0 9.4 ± 0.2 1.2 ± 0.2 20.8 ± 0.9 20.0 ± 0.3 0.9995
27.6 ± 0.1 21.5 ± 0.1 142 ± 0.5 32.8 ± 0.8 8.3 ± 0.4 52.5 ± 1.0 10.8 ± 0.2 1.4 ± 0.2 9.7 ± 0.6 10.2 ± 0.4 0.9977
27.0 ± 0.2 21.4 ± 0.1 199 ± 0.5 24.9 ± 0.8 4.3 ± 0.2 52.5 ± 1.0 8.9 ± 0.2 1.1 ± 0.2 5.4 ± 0.4 5.3 ± 0.1 0.9990
27.7 ± 0.2 21.3 ± 0.1 267 ± 0.5 20.6 ± 1.0 2.5 ± 0.2 52.5 ± 1.0 9.0 ± 0.2 1.1 ± 0.2 3.6 ± 0.4 3.7 ± 0.1 0.9999
However, let us note that this dispersion is globally negligible because the straight parts of each microchannel are quite longer than the bend length: 8 mm and a maximum of 0.31 mm respectively. Simulations have also been performed to check if the resistance RJ could be neglected according to the short length of the junctions (about 1 mm). The model of a capillary-channel junction is given in Figs. 7b and 8. A particular care has been taken to avoid divergence and oscillation problems during solving, using velocity and pressure capping strategy. Boundary conditions were similar to those previously presented. The simulation shows that a very low pressure drop (around 8 Pa) is obtained near the junction for a 550 Pa total applied pressure. It is clear that the pressure loss at the junction has a negligible influence. These simulation results show that the device microfluidic resistance RD is equal, in a good approximation, to the sum of
Fig. 7. 3D model of a bend of a microchannel (a) and 3D model of a capillarymicrochannel junction (b).
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
915
Fig. 8. Capillary-channel junction model geometry. a) Side view; b) front view.
RC and RM. RC is calculated from Eq. (4) and RM is obtained using the 2D model presented above. 5. Experimental results for trapezoidal microchannels Five different devices have been tested in this work. The geometrical characteristics, numerical flow simulations and flow measurements are shown in Table 1 for all devices. Fig. 9 shows the measured flow rate as a function of the applied pressure for all devices. A good linearity can be pointed out as predicted in case of low Reynolds number value. The linear correlations of the experimental data are always better than 0.995, indicating that the microfluidic systems work in laminar regime. Let us note that the agreement between the theoretical and the experimental results is better than 10%. 6. Conclusion Hydrodynamic flow characteristics in microfluidic devices were studied. Those microfluidic devices were connected to capillaries by an original planar configuration. An experimental set-up dedicated to measurements in the range of few μl h− 1 was realized and validated. A 2D fluidic model based on the Navier–Stokes equations was developed, which is quite simple and less time consuming than an equivalent 3D model. It can be used by a non-specialist audience, and is sufficient for most microfluidic applications.
The obtained experimental results indicate clearly that this 2D model combined to the Navier–Stokes equations and the noslip condition (at the wall) are sufficient to describe and calculate flows in 20 μm depth microchannels etched in silicon with better than a 10% agreement with measurements. According to this result, we propose, in the Appendix, normalized abacuses of the microfluidic resistance for different microchannel's cross-section shapes, as a support to microfluidic design work. These results were used in order to determine the surface charge of capillaries and Lab-On-Chip used in electrophoresis set-up, and published elsewhere [16]. Appendix. Microflow abacus The good agreement between experimental results and the presented 2D fluidic model shows that in the case of microchannels with hydraulic diameters of 40 to 70 μm, the Navier–Stokes equation combined to the Newtonian fluid and no-slip at the wall surface hypothesis are sufficient to describe the hydrodynamic flow. Here, we present some abacuses that can be used to determine microchannel fluidic resistances for the most usual cross-section shapes obtained using microtechnologies, i.e. the rectangular DRIE cross-section (Fig. 10a), the trapezoidal cross-section obtained by KOH wet etching in silicon (Fig. 10b), and the isotropic wet etching (Fig. 10c). Let us note that whatever the cross-section shape, the microchannel fluidic resistance tends
Fig. 9. Measured flow rates versus applied pressure for the five microchannels.
916
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Fig. 10. a) Abacus of the normalized microfluidic resistance RN in mbar h ml− 1 in case of a rectangular w × d cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid. b) Abacus of the normalized microfluidic resistance RN in mbar h μl− 1 in case of a trapezoidal w × d KOH cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid. c) Abacus of the normalized microfluidic resistance RN in mbar h μl− 1 in case of a w × d isotropic wet etching cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid.
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
917
Fig. 10 (continued ).
towards the resistance of an infinite parallel plate microchannel when width is much higher than depth. These abacuses are normalized for one meter microchannel length and a 10− 3 Pa s fluid viscosity (water at 20 °C). The microfluidic resistance RM of a microchannel is calculated from the normalized value RN as following: RM ¼ RN L g
ð6Þ
Where η is expressed in 10− 3 Pa s, and L in meter. For example, in the case of the device D (see Section 5), the cross-section shape of the microchannel is trapezoidal, with the following parameters: w = 200 μm, d = 22 μm, L = 0.25 m. The device was filled with water at 27 °C (η = 0.859 10− 3 Pa s according to Eq. (5)). According to the abacus (Fig. 10b), the normalized resistance is equal to RN = 20 mbar h μl− 1. Then: RM ¼ 2040:2540:859 ¼ 4:3 mbar h Al1 :
ð7Þ
This microchannel is connected to the macroscopic world with a 52.5 μm diameter and 8.9 cm length capillary. Eq. (4) gives: RC ¼ 1:1 mbar h Al1 :
ð8Þ
Then, the microfluidic resistance of the total device is estimated to: RD ¼ RM þ RC ¼ 5:4 mbar h Al1 :
ð9Þ
This estimation is in a good aggreement with the experimental value of (5.3 ± 0.1) mbar h μl− 1 for this device (see Section 5).
References [1] P. Gravesen, J. Branebjerg, O.S. Jensen, Journal of Micromechanics and Microengineering 3 (1993) 168. [2] H.H. Bau, J.N. Pfahler, Experimental observations of liquid flow in micro conduits, in: 39th AIAA Aerospace Sciences Meeting & Exhibit., Reno, NV, January 8–11, 2001. [3] I.P.I. Papautsky, J. Brazzle, T. Ameel, A.B. Frazier, Sensors and Actuators A 73 (1999) 101. [4] J. Pfahler, J. Harley, H. Bau, J. Zemel, Sensors and Actuators A 21–23 (1990) 431. [5] B. Xu, K.T. Ooi, N.T. Wong, W.K. Choi, International Communications in Heat and Mass Transfer 27–8 (2000) 1165. [6] J.T. Cheng, N. Giordano, Physical Review E 65 (2002) 1063. [7] Q. Weilin, G.M. Mala, L. Dongqing, International Journal of Heat and Mass Transfer 43 (2000) 353. [8] Q. Weilin, G.M. Mala, L. Dongqing, International Journal of Heat and Mass Transfer 43 (2000) 3925. [9] H. Herwing, O. Hausner, International Journal of Heat and Mass Transfer 46 (2003) 935. [10] M. Koch, A. Evans, A. Brunnschweiler, Microfluidic technology and applications, Research Studies Press Ltd., R. Pethig, England, 2000. [11] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical Recipes in C, Second Edition, Cambridge University Press, 1999 994 p. [12] G. Schnell, Experimental Thermal and Fluid Science 15 (1997) 406. [13] T.S. Lammerink, N.R. Tas, M. Elwenspoek, J.H. Fluitman, Sensors and Actuators A 37–38 (1993) 45. [14] G.H. Mala, D. Li, International Journal of Heat and Mass Transfer 20 (1999) 142. [15] M. Richter, P. Woias, D. Weis, Sensors and Actuators A 62 (1997) 480. [16] L. Renaud, P. Kleimann, P. Morin, Electrophoresis 25 (2004) 123.
Materials Science and Engineering C 28 (2008) 910 – 917 www.elsevier.com/locate/msec
Theoretical and experimental studies of microflows in silicon microchannels L. Renaud ⁎, C. Malhaire, P. Kleimann, D. Barbier, P. Morin Institut des Nanotechnologies de Lyon, INL, Université de Lyon and INSA de Lyon, Lyon, F-69003, France; Université Lyon 1, Lyon, F-69622, France Available online 18 October 2007
Abstract The determination of fluid flows in silicon microchannels is important for the design of microfluidic systems. In this paper, experimental investigations on the characteristics of low fluid flows (few μl h− 1) in silicon trapezoidal microchannels (21 μm in depth, length and width ranging from 200 to 440 mm and 58 to 267 μm, respectively) are presented. The test-devices have been fabricated using micromachining technologies. A double KOH etching process has been used to achieve microchannels in (100)-oriented silicon wafers as well as deep in-plane cavities used for capillary connections. Silicon has been finally anodically bonded on Pyrex substrates. The experimental set-up, based on the measurement of a differential pressure and a liquid–air interface displacement in a gauged tube, is fully detailed in terms of fluidic connections and measurement principle. The experimental results are in good agreement with the Navier–Stokes theory, solved by a simple iterative method. However, finite element modelling has been used to study complex 3D problems that were found in the devices and the experimental set-up. Finally, we propose abacuses for three different channel cross-sections that may be used to easily compute the flow in a microchannel. © 2007 Elsevier B.V. All rights reserved. Keywords: Microfluidic; Microflow; Navier–Stokes equations; Abacuses; Finite element modelling
1. Introduction The determination of flow velocity profile in microchannels is important for the design of microfluidic systems. In fact, micromachined channels are essential parts in many applications from heat pipes to micro-pumps, μ-TAS and Lab-On-AChip devices [1]. With the device size reduction, many investigations of the liquid flow in microchannels, with different experimental methods, have been presented in the literature to verify the validity of macrofluidic flow theory to the small scale [2]. In fact, the microchannels are characterised by very large surface to volume ratios and surface phenomena that are typically ignored in larger conduits might become important at the microscale. Some significant flow characteristic deviations from theoretical predictions have been experimentally observed but some contradictions among these results appear, according to the work of I.P.I. Papautsky et al. [3]. For example, J. Pfahler et al. [4] have found that experiences carried out with microchannels having hydraulic diameter (Dh) of 76.1 μm were in rough agreement with the prediction from the Navier–Stokes equations. Significant deviations from the theory concerning the ⁎ Corresponding author. Tel.: +33 4 72 17 73 99. E-mail address: [email protected] (L. Renaud). 0928-4931/$ - see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.msec.2007.10.070
friction coefficient (up to 3 times) were only found for smaller channels (Dh = 1.6 and 3.3 μm). B. Xu et al. [5] studied microchannels with Dh ranging from 30 to 344 μm and observed a maximum 20% overestimated relationship between friction factor and Reynolds number. Even for water flow in nanometer depth channels, J.T. Cheng and Giordano validated the Poiseuille flow theory with difference less than 10% between experimental and theoretical flows [6]. On the contrary, Q. Weilin et al. [7,8] observed higher pressure gradient and flow friction factor (up to 40%) in trapezoidal microchannels (hydraulic diameter ranging from 51 to 169 μm) than those given by the conventional laminar flow theory. They proposed a roughness-viscosity model to interpret the experimental results. Those reported deviations from theory could be explained by temperature effects, pressure dependant properties or wall roughness, as noticed by H. Herwig and Hausner [9] or by experimental artifacts: errors due to the fabrication of the samples, fluid leakage, bubbles trapping, and dusts. In this paper, we present experimental investigations on the characteristics of fluid flows in KOH-etched trapezoidal microchannels (Dh ranging from 45.7 to 74.3 μm) for very low Reynolds number value (Re b 1). The experimental set-up is fully detailed and we show some practical aspects in microfluidic experiments like an original capillary connection
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
911
flow is always laminar under low pressure drop (typically a few bar) [7], leading to a unidirectional flow and a uniform absolute pressure in the cross-section. For a permanent flow and a fixed pressure drop ΔP between the inlet and the outlet of the channel, Eq. (1) simplifies in the cross-section to: DP ∂2 w ∂2 w þ 2 þ 2 ¼ 0: gL ∂y ∂z Fig. 1. Geometrical definition of the Navier–Stokes numerical solving.
and microflow measurements. A simple and proven iterative method for solving the Navier–Stokes equations in the case of a 2D flow is presented for a non-specialist audience [10]. Nevertheless, finite element modelling has been used to study the complex three-dimensional problems that were found in the experimental set-up. Finally we give abacuses for three different channel cross-sections in the Appendix that may be used to quickly compute the flow in a microchannel and proposed to support a microfluidic design work. 2. Theory of pressure-driven flows The relationship between the fluid velocity and the absolute pressure for an incompressible viscous liquid is given by the classical fluid dynamics theory and the well-known Navier– Stokes equation: ∂Yv þ ∂t
Y Y
vj
Y
Y
v ¼ j
P Y þ mD v q
ð1Þ
Where Y v stands for the fluid velocity vector with components (u, v, w) each expressed for a set of Euler components (x, y, z, t). P and ρ are the absolute pressure and the voluminal mass, respectively (Fig. 1). The kinematic viscosity ν is the ratio between the dynamic viscosity η and the voluminal mass ρ. In the case of a microfluidic horizontal straight channel (x-direction), the
ð2Þ
The boundary condition at the wall surface is the no-slip condition, assumed for hydrophilic walls (w(y,z) = 0 at the wall). Eq. (2) is a classical elliptic equation that can be solved using a second-order equidistant central difference scheme and a Successive Over Relaxation (SOR) algorithm. The detailed general explanations are given in [11]. Iterations have been performed until convergence of the calculated mean flow. An example of velocity profile in a trapezoidal microchannel crosssection is given in Fig. 2. Since the flow is laminar, the flow rate Q (defined as the product of the mean velocity by the cross-section) is always proportional to the applied pressure ΔP. By analogy with electrical Ohm's law, we can define the ratios between Q and ΔP as the microfluidic resistance R and the microfluidic conductance G: R¼
DP Q and G ¼ : Q DP
ð3Þ
For example, the calculations of R and G for the device A, lead to 50.2 mbar h μl− 1 and 19.9 μl h− 1 bar− 1 respectively. Let us note that, R or G have been determined by the least square method applied on four points to be less sensitive to the experimental artifacts described above. 3. Experimental procedure 3.1. Fabrication process The use of an original capillary connection leads to a two steps etching process in silicon. The test-devices were
Fig. 2. SOR algorithm result in the case of a trapezoidal channel with a 21.3 μm depth, 57.5 μm width, 44.9 cm length and 54.7° angle (device A, see Section 5) filled with water at 25 °C with a pressure drop equal to 1 bar.
912
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Fig. 3. Process flow diagram for the test-devices.
fabricated using micromachining technologies as outlined in Fig. 3. Deep cavities (for the capillaries) and trapezoidal microchannels (studied microchannels) were etched in (100)oriented silicon wafers. Long microchannels (up to 44 cm) can be achieved on a 1 cm2 area using a zigzag design. The cavities were used in the final device to connect microchannels and fused silica capillaries. The depth of cavities and microchannels being different, a double etching technique was used. First the silicon wafers were oxidized (800 nm thermal oxide). Then a 200 nm Si3N4 layer was deposited by plasma enhanced chemical vapor deposition. The Si3N4 and oxide layers were patterned using standard photolithography process in order to define microchannels and cavities (Fig. 3a and b). Bulk anisotropic etching in KOH was then used to form 150 μm deep cavities (Fig. 3c). The visible parts of the oxide layer were removed in buffered HF prior to a second etching in KOH used to form 20 μm deep microchannels (Fig. 3d). The remaining Si3N4 and oxide layers were removed and glass substrates (Hoya SD-2) were anodically bonded on the top of silicon wafers to close the channels. Finally the devices were cut with a die-saw and fused silica capillaries were glued inside the cavities (Fig. 3e).
moving out downstream of the channel [12], others use the heat conduction properties of the fluids [13]. In this work, we have measured the fluid flow rate by the displacement of a liquid–air interface in a gauged tube, placed upstream of the microchannel. This method allows adjusting the volume of measurement according to the expected fluid flow, as shown in Fig. 4. A large syringe (60 ml) filled with air is used for applying pressure on the liquid. The total pressure drop ΔPT is measured with 2% accuracy by means of a differential pressure sensor (Honeywell 26PC series) connected to the large syringe. A low volume microsyringe (Hamilton, from 10 μl to 250 μl) is used as gauged tube and connected to the microfluidic system. The flow rate Q is determined through the observation of the meniscus shifting in the syringe for a given time. The
3.2. Microflow measurements One of the main difficulties in the study of microfluidic systems consists in the measurement of low flow rates (typically in the range of few μl h− 1). The size reduction of the microchannel cross-section imposes new experimental techniques of flow measurement. Indeed, friction forces prohibit the use of conventional flow meter for this range of flows. Several original solutions have been developed during last years. Some of them consist in the measurement of the weight of fluid
Fig. 4. Experimental set-up.
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
913
4.1. Study of RC and RS, and experimental set-up validation From Navier–Stokes equations, the microfluidic resistance RC of a circular cross-section capillary is expressed as [14]: Rtube ¼
128 g L p D4
ð4Þ
where L and R are the capillary length and radius respectively. When water is used η is given by [15]: Fig. 5. Schematized view (not to scale) of the pressure drop along the fluid flow.
g¼qm¼
0:1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 120 þ 2:142824 T 8:435 þ 8078:4 þ ðT 8:435Þ2
ð5Þ measurement accuracy ranges from 0.1 μl to 5 μl depending on the microsyringe used. This experimental set-up allows us to measure flow rates from 10 μl h− 1 to 200 μl h− 1 in less than 1 h with accuracy better than 2% (due to the 0.1 μl graduation step of the gauged 10 μl microsyringe, for example). 4. Procedure for the hydrodynamic study of the microfluidic system As shown on the more detailed schematic representation of the experimental set-up (Fig. 5), all the microfluidic elements are in serial configuration (same flow Q), therefore the total microfluidic resistance RT is equal to the sum of each microfluidic resistance. The total pressure drop ΔPT (that is measured) can be split in: 1) ΔPS = RS ⁎ Q, the pressure drop through the gauged tube (Hamilton microsyringe). 2) ΔPC = RC ⁎ Q, the pressure drop through the capillary tube. 3) ΔPJ = RJ ⁎ Q, the pressure drop through the two capillarymicrochannel junctions. 4) ΔPM = RM ⁎ Q, the pressure drop through the microchannel.
where T is the temperature (in Celsius degrees). According to Eq. (4) and the large difference between the syringe and the capillary diameters (typically more than 1.5 mm and 50 μm respectively), RS can be neglected in comparison with the other microfluidic resistances. Moreover, it is obvious that the accuracy on the determination of RC is strongly dependent on that of the diameter. In this study, the inner diameter of the fused silica capillary tube that connects the syringe to the micromachined channel was (50.0 ± 3) μm, as given by the manufacturer (Polymicro Technologies). A diameter value of (52.5 ± 1.0) μm has been measured by means of an optical microscope (Olympus BH-2). In order to validate the experimental set-up, we measured the microfluidic resistance of a (52.5 ± 1.0) μm diameter and (19.9 ± 0.1) cm length capillary at 22.5 °C with pressure range 0– 1 bar. The calculated (from Eq. (4)) and measured microfluidic resistances are (2.8 ± 0.2) mbar h μl− 1 and (2.7 ± 0.1) mbar h μl− 1 respectively; with a 0.9990 correlation coefficient between the applied pressure and the measured flow. This result is in agreement with the literature [14]. 4.2. Study of RJ and RM
Some of these resistances can be calculated (RS and RC) from fluidic theory; others have to be estimated by numerical methods (RJ and RM).
It is obvious that a 2D fluidic model is preferable to a 3D approach in terms of calculation time and convergence limits.
Fig. 6. Geometry of a simulated channel (type C).
914
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Table 1 Geometrical parameters and microfluidic resistances of studied microchannels Type of device Microchannel
Capillary
Device
Temperature (°C) Depth (μm) Width (μm) Length (cm) Estimated RM in mbar h μl− 1 Diameter (μm) Length (cm) Estimated RC in mbar h μl− 1 Estimated RD = RC + RM in mbar h μl − 1 Measured RD in mbar h μl− 1 Correlation coefficient (hydrodynamic flow versus pressure)
However, 2D models are not appropriate to modelize RJ and RM because of the quite complex shapes of the capillary-channel junctions and microchannel bends. In order to work with a 2D approach, the influence of the bends can be taken into account by considering an equivalent straight channel length. This equivalent length obviously depends on the bend shape. Let us note that in this work, the KOH etching ensures a broadening of the channel cross-section in the bends. It is thus reasonable to assume that the equivalent straight channel length is lower than the bend length. Here we can define the bend length as two times the bend size. Fig. 6 shows the geometry (only the main dimensions are shown) of a bend in the case of device C (see Table 1), extracted from microscope views and profilometry measurements. As shown in this figure, the bend length is assumed to be 0.44 mm. A 3D view of the model is shown in Fig. 7a. Two long straight parts (1 mm length) have been added to the bend to define simple boundary conditions (unidirectional velocity and uniform pressure at entrance and exit). The previously proposed assumption has been checked by means of finite element simulation using ANSYS® software. For each device, 3D models have been developed assuming a no-slip condition at the wall. As the number of elements was limited for our software version, only a limited length of the structures was modeled. A scale factor has been applied to the pressure values to get simulated flow values that were close to those measured. During the postprocessing stage, the ANSYS® Parametric Design Language (APDL) has been used to compute the flow and then the microfluidic resistance starting from the velocity on each node of the trapezoidal microchannel exit. For a temperature of 27.6 °C, the 3D simulation gives a microfluidic resistance for a bend equal to 0.0056 mbar h μl− 1. This value has to be compared to the microfluidic resistance of a 0.44 mm length trapezoidal straight channel equal to 0.0114 mbar h μl− 1. This result confirms that the equivalent straight channel length of the bend is lower than the 0.44 mm bend length. Similar analyses have been performed for the other samples (A to E) and gave the same results. In conclusion, the microfluidic resistance for a zigzag channel may be surrounded by two values corresponding to straight channels with extreme lengths. This 2D simplification leads to a dispersion of the total channel length (Table 1).
A
B
C
D
E
25.3 ± 0.3 21.4 ± 0.1 58 ± 0.5 44.4 ± 0.5 48.8 ± 1.4 52.5 ± 1.0 14.2 ± 0.2 1.9 ± 0.2 50.7 ± 1.6 47.5 ± 2.7 0.9964
28.0 ± 0.3 21.4 ± 0.1 94 ± 0.5 44.1 ± 0.7 19.6 ± 0.7 52.5 ± 1.0 9.4 ± 0.2 1.2 ± 0.2 20.8 ± 0.9 20.0 ± 0.3 0.9995
27.6 ± 0.1 21.5 ± 0.1 142 ± 0.5 32.8 ± 0.8 8.3 ± 0.4 52.5 ± 1.0 10.8 ± 0.2 1.4 ± 0.2 9.7 ± 0.6 10.2 ± 0.4 0.9977
27.0 ± 0.2 21.4 ± 0.1 199 ± 0.5 24.9 ± 0.8 4.3 ± 0.2 52.5 ± 1.0 8.9 ± 0.2 1.1 ± 0.2 5.4 ± 0.4 5.3 ± 0.1 0.9990
27.7 ± 0.2 21.3 ± 0.1 267 ± 0.5 20.6 ± 1.0 2.5 ± 0.2 52.5 ± 1.0 9.0 ± 0.2 1.1 ± 0.2 3.6 ± 0.4 3.7 ± 0.1 0.9999
However, let us note that this dispersion is globally negligible because the straight parts of each microchannel are quite longer than the bend length: 8 mm and a maximum of 0.31 mm respectively. Simulations have also been performed to check if the resistance RJ could be neglected according to the short length of the junctions (about 1 mm). The model of a capillary-channel junction is given in Figs. 7b and 8. A particular care has been taken to avoid divergence and oscillation problems during solving, using velocity and pressure capping strategy. Boundary conditions were similar to those previously presented. The simulation shows that a very low pressure drop (around 8 Pa) is obtained near the junction for a 550 Pa total applied pressure. It is clear that the pressure loss at the junction has a negligible influence. These simulation results show that the device microfluidic resistance RD is equal, in a good approximation, to the sum of
Fig. 7. 3D model of a bend of a microchannel (a) and 3D model of a capillarymicrochannel junction (b).
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
915
Fig. 8. Capillary-channel junction model geometry. a) Side view; b) front view.
RC and RM. RC is calculated from Eq. (4) and RM is obtained using the 2D model presented above. 5. Experimental results for trapezoidal microchannels Five different devices have been tested in this work. The geometrical characteristics, numerical flow simulations and flow measurements are shown in Table 1 for all devices. Fig. 9 shows the measured flow rate as a function of the applied pressure for all devices. A good linearity can be pointed out as predicted in case of low Reynolds number value. The linear correlations of the experimental data are always better than 0.995, indicating that the microfluidic systems work in laminar regime. Let us note that the agreement between the theoretical and the experimental results is better than 10%. 6. Conclusion Hydrodynamic flow characteristics in microfluidic devices were studied. Those microfluidic devices were connected to capillaries by an original planar configuration. An experimental set-up dedicated to measurements in the range of few μl h− 1 was realized and validated. A 2D fluidic model based on the Navier–Stokes equations was developed, which is quite simple and less time consuming than an equivalent 3D model. It can be used by a non-specialist audience, and is sufficient for most microfluidic applications.
The obtained experimental results indicate clearly that this 2D model combined to the Navier–Stokes equations and the noslip condition (at the wall) are sufficient to describe and calculate flows in 20 μm depth microchannels etched in silicon with better than a 10% agreement with measurements. According to this result, we propose, in the Appendix, normalized abacuses of the microfluidic resistance for different microchannel's cross-section shapes, as a support to microfluidic design work. These results were used in order to determine the surface charge of capillaries and Lab-On-Chip used in electrophoresis set-up, and published elsewhere [16]. Appendix. Microflow abacus The good agreement between experimental results and the presented 2D fluidic model shows that in the case of microchannels with hydraulic diameters of 40 to 70 μm, the Navier–Stokes equation combined to the Newtonian fluid and no-slip at the wall surface hypothesis are sufficient to describe the hydrodynamic flow. Here, we present some abacuses that can be used to determine microchannel fluidic resistances for the most usual cross-section shapes obtained using microtechnologies, i.e. the rectangular DRIE cross-section (Fig. 10a), the trapezoidal cross-section obtained by KOH wet etching in silicon (Fig. 10b), and the isotropic wet etching (Fig. 10c). Let us note that whatever the cross-section shape, the microchannel fluidic resistance tends
Fig. 9. Measured flow rates versus applied pressure for the five microchannels.
916
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
Fig. 10. a) Abacus of the normalized microfluidic resistance RN in mbar h ml− 1 in case of a rectangular w × d cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid. b) Abacus of the normalized microfluidic resistance RN in mbar h μl− 1 in case of a trapezoidal w × d KOH cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid. c) Abacus of the normalized microfluidic resistance RN in mbar h μl− 1 in case of a w × d isotropic wet etching cross-section microchannel of 1 m length filled with a η = 10− 3 Pa s Newtonian fluid.
L. Renaud et al. / Materials Science and Engineering C 28 (2008) 910–917
917
Fig. 10 (continued ).
towards the resistance of an infinite parallel plate microchannel when width is much higher than depth. These abacuses are normalized for one meter microchannel length and a 10− 3 Pa s fluid viscosity (water at 20 °C). The microfluidic resistance RM of a microchannel is calculated from the normalized value RN as following: RM ¼ RN L g
ð6Þ
Where η is expressed in 10− 3 Pa s, and L in meter. For example, in the case of the device D (see Section 5), the cross-section shape of the microchannel is trapezoidal, with the following parameters: w = 200 μm, d = 22 μm, L = 0.25 m. The device was filled with water at 27 °C (η = 0.859 10− 3 Pa s according to Eq. (5)). According to the abacus (Fig. 10b), the normalized resistance is equal to RN = 20 mbar h μl− 1. Then: RM ¼ 2040:2540:859 ¼ 4:3 mbar h Al1 :
ð7Þ
This microchannel is connected to the macroscopic world with a 52.5 μm diameter and 8.9 cm length capillary. Eq. (4) gives: RC ¼ 1:1 mbar h Al1 :
ð8Þ
Then, the microfluidic resistance of the total device is estimated to: RD ¼ RM þ RC ¼ 5:4 mbar h Al1 :
ð9Þ
This estimation is in a good aggreement with the experimental value of (5.3 ± 0.1) mbar h μl− 1 for this device (see Section 5).
References [1] P. Gravesen, J. Branebjerg, O.S. Jensen, Journal of Micromechanics and Microengineering 3 (1993) 168. [2] H.H. Bau, J.N. Pfahler, Experimental observations of liquid flow in micro conduits, in: 39th AIAA Aerospace Sciences Meeting & Exhibit., Reno, NV, January 8–11, 2001. [3] I.P.I. Papautsky, J. Brazzle, T. Ameel, A.B. Frazier, Sensors and Actuators A 73 (1999) 101. [4] J. Pfahler, J. Harley, H. Bau, J. Zemel, Sensors and Actuators A 21–23 (1990) 431. [5] B. Xu, K.T. Ooi, N.T. Wong, W.K. Choi, International Communications in Heat and Mass Transfer 27–8 (2000) 1165. [6] J.T. Cheng, N. Giordano, Physical Review E 65 (2002) 1063. [7] Q. Weilin, G.M. Mala, L. Dongqing, International Journal of Heat and Mass Transfer 43 (2000) 353. [8] Q. Weilin, G.M. Mala, L. Dongqing, International Journal of Heat and Mass Transfer 43 (2000) 3925. [9] H. Herwing, O. Hausner, International Journal of Heat and Mass Transfer 46 (2003) 935. [10] M. Koch, A. Evans, A. Brunnschweiler, Microfluidic technology and applications, Research Studies Press Ltd., R. Pethig, England, 2000. [11] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical Recipes in C, Second Edition, Cambridge University Press, 1999 994 p. [12] G. Schnell, Experimental Thermal and Fluid Science 15 (1997) 406. [13] T.S. Lammerink, N.R. Tas, M. Elwenspoek, J.H. Fluitman, Sensors and Actuators A 37–38 (1993) 45. [14] G.H. Mala, D. Li, International Journal of Heat and Mass Transfer 20 (1999) 142. [15] M. Richter, P. Woias, D. Weis, Sensors and Actuators A 62 (1997) 480. [16] L. Renaud, P. Kleimann, P. Morin, Electrophoresis 25 (2004) 123.