Gas Velosity and Mass Flowrate Scaling Modeling in Microelectronics ’ Thermal Control Systems

In the present research we investigate pressure driven flow in the transition and free-molecular flow regimes with the objective of developing unified flow models for microchannels. These models are based on a velocity scaling law, which is valid for a wide range of Knudsen number. Simple slip-based descriptions of flowrate in microchannels are corrected for effects in the transition and free-molecular flow regimes with the introduction of a rarefaction factor. The resulting models can predict the velocity distribution, mass flowrate, pressure and shear stress distribution in rectangular microchannels in the entire Knudsen flow regime.


Velocity Scaling
Practical development of micro-devices is creating a need for better understanding of microscale transport phenomena.This is causing a worldwide redirection of transport research and development, from macroscale (conventional) to microscale situations.Significant differences in transport phenomena have been reported at the microscale as compared to the microscale, [1].
From the Direct Simulations Monte Carlo (DSMC) results and solutions of the linearized Boltzmann equation, it is evident that the velocity profiles in pipes and channels remain approximately parabolic for a large range of Knudsen number.This is also consistent with the analysis of the Navier-Stokes and Burnett equations in long channels.Based on this observation, we model the velocity profile as parabolic in the entire Knudsen regime, with a consistent slip condition.We write the dimensional form for velocity distribution in a channel of height h, , where F (dP/dx, μ 0 , h, λ) shows the functional dependence of velocity on the pressure gradient, viscosity, channel height, and local mean free path.Temperature is assumed to be constant, and therefore the dynamic viscosity is also a constant.Here U s is the slip velocity, which satisfies the general slip boundary condition given by where b is the general slip coefficient.Using this boundary condition yields By nondimensionalizing the velocity distribution with the local average velocity, dependence on the local flow conditions F (dP/dx, μ 0 , h, λ) is eliminated.Therefore, the resulting relation is a function of Kn and y only.Assuming diffuse reflection (σ v = 1) for simplicity, we obtain Equation (2) solely describes the shape of the velocity distribution, but it does not properly model the flowrate, which requires additional corrections, as will be shown in the next section.
In Figure 1 we plot the nondimensional velocity variation obtained in a series of DSMC simulations for Kn = 0.1, Kn = 1, Kn = 5, and Kn = 10.We also included the corresponding linearized Boltzmann solutions obtained in [2].It is seen that the DSMC velocity distribution and the linearized Boltzmann solutions agree quite well.We can now use equation ( 2) and compare with the DSMC data by varying the parameter b, which for b = 0 corresponds to Maxwell's first-order and for b = −1 to the second-order boundary condition in the slip regime only.Here we find that for b = −1, equation (2) results in an accurate model of the velocity distribution for a wide range of Knudsen number.From the figure, it is clear that the velocity slip is slightly overestimated with the proposed model for the Kn = 1 case.To obtain a better velocity slip, we varied the value of the parameter b by imposing, for example, b = −1.8 for the Kn = 1 case.Although a better agreement is achieved for the velocity slip, the accuracy of the model in the rest of the channel is destroyed.In Figure 2 we show the nondimensionalized velocity distribution along the centerline and along the wall of the channels for the entire Knudsen number regime considered here, i.e., 0.01 ≤ Kn ≤ 20.We included in the plot data for the velocity slip and centerline velocity from 20 different DSMC runs.It should be noticed that the velocity scaling model is independent of the gas type.The linearized Boltzmann solution of Aoki for a monatomic gas is also shown by triangles.This solution closely matches the DSMC predictions.Maxwell's first-order boundary condition (b=0) (shown by a solid line) erroneously predicts a uniform nondimensional velocity profile for large Knudsen number.The breakdown of slip flow theory based on the first-order slip boundary condition is realized around Kn=0.08 and Kn=0.5 for the wall and the centerline velocity, respectively.This finding is consistent with the commonly accepted limits of the slip flow regime.The prediction using b=-1 is shown by small dashed lines.The corresponding centerline velocity closely follows the DSMC results, while the slip velocity of the model with b=-1 deviates from DSMC in the intermediate range for 0.1 < Kn < 10.One possible reason for this is the effect of the Knudsen layer, a sublayer that is present between the viscous boundary layer and the wall, with a thickness of approximately one mean free path.For small Kn flows the Knudsen layer is thin and does not affect the velocity slip prediction too much.For very large Kn flows, the Knudsen layer covers the entire channel.However, for intermediate Kn values both the fully developed viscous flow (boundary layer) and the Knudsen layer exist in the channel.At this intermediate range, approximating the velocity profile to be parabolic neglects the Knudsen layers.For this reason, the model with b=-1 results in 10% error of the velocity slip at Kn=1.However, the velocity distribution in the rest of the channel is described accurately for the entire flow regime.For a comparison we also included similar predictions by the second-order slip boundary condition of Hsia and Domoto (large dashed line).This boundary condition performs worse than even the first-order Maxwell's boundary condition for large Kn values.Only the general slip boundary condition predicts the scaling of the velocity profiles accurately.

Flowrate Scaling
The volumetric flowrate in a channel is a function of the channel dimensions, fluid properties (μ 0 , λ), and pressure drop, and it can be written as For a channel of height h, using the Navier-Stokes solution and the general slip boundary condition (1) we obtain where Kn = λ /h.
The flowrate for the and free-molecular flows are both linearly dependent on dP/dx, and thus we choose to normalize the flowrate with the pressure gradient.This quantity is computed based on the DSMC simulations and is shown in Figure 3.For comparison we present the Q  / |dP/dx| predictions obtained using Maxwell's first- order slip boundary condition (b = 0, dashed lines) and the general slip boundary condition (b = −1, dashed-dotted lines).In both cases the predictions are erroneous.The general slip boundary condition performs the worst for it is asymptotic to a constant value, while the DSMC data show a considerable increase with Kn.The first-order boundary condition follows the DSMC data, however with a significant error.The model in equation (1) gives good agreement with DSMC data and the linearized Boltzmann solutions for the nondimensional velocity profile, but it does not predict correctly the flowrate.This is expected, since the Navier-Stokes equations are invalid in this regime.In fact, the dynamic viscosity, which defines the diffusion of momentum due to the intermolecular collisions, must be modified to account for the increased rarefaction effects.The kinetic theory description for dynamic viscosity requires where v is the mean thermal speed.
Using mean free path λ in this relation is valid as long as intermolecular collisions are the dominant part of momentum transport in the fluid (i.e., Kn << 1).However, for increased rarefaction, the intermolecular collisions are reduced significantly, and in the free-molecular flow regime, only the collisions of the molecules with the walls should be considered.Therefore, in free-molecular channel flow the diffusion coefficient should be based on characteristic length scale h (channel height) and thus . Since the diffusion coefficient is based on λ in slip or continuum flow regimes and h in the free-molecular flow regime, we propose to model the variation of diffusion coefficient with the following hybrid formula: where μ 0 is the dynamic viscosity of the gas at a specified temperature and μ is the generalized diffusion coefficient.
The variable diffusion coefficient model presented above is based on a simple analysis.
In general, the increased rarefaction effects in our flowrate model can be taken into account by introducing a correction expressed as rarefaction factor C r (Kn), which is a function of the Knudsen number.The flowrate is then obtained as where C r (Kn) is a general function of Knudsen number.A possible model for C r is suggested by equation ( 4) as where α is a parameter.If we assume that α is constant in the entire Knudsen regime, the flowrate in the slip flow regime will be erroneously enhanced, resulting in , where c M  corresponds to continuum mass flowrate.This model becomes inaccurate for a nonzero value of α in the slip flow regime.Moreover, in the free-molecular flow regime, for very long channels (L >> λ >> h) there are no physical values for α, since the flowrate increases logarithmically with Kn. finite-length channels the flowrate is asymptotic to a constant value proportional to log e (L/h).Therefore, for finite-length two-dimensional channel flows, the coefficient α should smoothly vary from zero in the slip flow regime to an appropriate constant value in the free-molecular flow regime.The physical meaning of this behavior is that the dynamic viscosity remains the standard diffusion coefficient in the early slip flow regime.
The value of  increases slowly with Kn in the slip flow regime, and therefore the effect of change of the diffusion coefficient is second-order in Kn.For this reason the experimental slip flow results are accurately predicted by the slip flow theory, which does not require change of the diffusion coefficient length scale from  to channel height h.Variation of  as a function of Kn can be represented accurately with the following relation: where  o is determined to result in the desired freemolecular flowrate.Note that the values for  1 and  are the only two undetermined parameters of the model.

Conclusions
We developed a unified flow model that can accurately predict the volumetric flowrate, velocity profile, and pressure distribution in the entire Knudsen regime for rectangular microchannels in microelectronics thermal control systems.The new model is based on the hypothesis that the velocity distribution remains parabolic in the transition flow regime, which is supported by the asymptotic analysis of the Burnett equations.The general velocity slip boundary condition and the rarefaction correction factor are the two primary components of this unified model.
The general slip boundary condition gives the correct nondimensional velocity profile, where the normalization is obtained using the local channel averaged velocity.This eliminates the flowrate dependence in modeling the velocity profile.For channel flows, we obtain b=-1 in the slip flow regime.Evidence based on comparisons of the model with the DSMC and Boltzmann solutions shows that b=-1 is valid in the entire Knudsen regime.
In order to model the flowrate variations with respect to the Knudsen number Kn, we introduced the rarefaction correction factor as C r =1+ Kn.This form of the correction factor was justified using two independent arguments: first, the apparent diffusion coefficient; and second, the ratio of intermolecular collisions to the total molecular collisions.We must note that  cannot be a constant.Physical considerations to match the slip flowrate require → 0 for Kn ≤ 0.1, while →  o in the free molecular flow regime.The variation of  between zero and a known   value is approximated using equation (7) which introduced two empirical parameters   and  to the new model.
Therefore, the unified model employs two empirical parameters (  and ) and two known parameters b=-1 and  o .Although this empiricism is not desired, the  value in C r varies from zero in the slip flow regime to an order-one value of  o as Kn → ∞.
of velocity distribution, the average velocity in the channel (

Figure 1 -
Figure 1 -Velocity profile comparisons of the model (equation (2)) with DSMC and linearized Boltzmann solutions [2].Maxwell's first-order boundary condition is shown with dashed lines (b = 0), and the general slip boundary condition (b = −1) shown with solid lines.

Figure 2 -
Figure 2 -Velocity scaling at wall and centerline of the channels for slip and transition flows.The linearized Boltzmann solution of Aoki is shown by triangles, and the DSMC simulations are shown by points.Theoretical predictions of velocity scaling for different values of b, and Hsia and Domoto's second-order slip boundary condition are also shown.

Figure 3 -
Figure 3 -Volumetric flowrate (per channel width) per absolute value of the pressure gradient in [m 3 /(sPa)] as a function of Kn.The solid line represents the proposed model.