Advances in Tribology

Advances in Tribology / 2009 / Article

Review Article | Open Access

Volume 2009 |Article ID 461793 | 10 pages |

An Alternative Simulation Method for Calculation of Microgas Flows under Flying Head Sliders

Academic Editor: Mihai Arghir
Received29 Apr 2009
Revised03 Nov 2009
Accepted07 Dec 2009
Published14 Mar 2010


The precise knowledge of the force and moment generated by the air squeezed under the read-write slider by the rotating disc is an engineering necessity in designing the air bearing surface slider. This paper reviews methods addressing the thin gas film bearings problem. It firstly reviews briefly the relatively well-known two methods of calculations of the microgas flows under flying head sliders, the generalized Reynolds equation, having given a number of useful results of slider design, and the DSMC method, which is precise and appropriate for the flow of complex configurations but is restricted to miniature (~micrometer) size sliders. The main purpose of the paper is to introduce to the reader an alternative method, the information preservation (IP) method, for use in simulation of the flows under air bearing surfaces. Some recent results of IP simulation of slider flows published on conference proceedings are introduced here.

1. Introduction

The problem of thin gas film bearings in the gap between the flying head slider and the magnetic disc now has an increasing interest among scientists concerning computation of actual problems of micro gas flows. To calculate precisely the force and moment generated by the air squeezed by the rotating disc is essential in the design of the head slider. The typical slider length in disk drives is about 1 mm and the width is the same order of the length. The size of the clearances between the slider and the disc is much smaller and is constantly decreasing to increase the recording density. The flying height in the early stage of the disc recording head was of the order of 810 micrometer [1] and 1 micrometer [2] while the drives areal data densities of the order 100 Gbit/in2 and 1 Tbit/in2 in consideration require sliders to have a flying height of 510 nm and less (see e.g., [3]). The flow regime in the modern flying sliders is definitely beyond the slip flow regime, and application of the tools of computation in transitional flow regime is inevitable.

The thin film air flow between the slider air bearing surface and the disc is most appropriately described by the Reynolds equation which is a differential equation relating the pressure 𝑝, density 𝜌, platter velocity 𝑈, and the height of the gap, firstly developed by Reynolds for continuum fluid [4]. The derivation of it can be found in [5]. Burgdorfer [6] introduced the velocity slip correction to the Reynolds equation. Fukui and Kaneko [7, 8] developed a generalization of the equation suitable for the transitional regime. The work in [9] provides a simple and enlightening derivation showing that the equation is the expression of balance of the flow rate of the Poiseuille flow and the flow rate of the Couette flow, and it also shows how the equation can be extended to slip regime and transitional regime. The equation modified to include slip and transitional effects is still called Reynolds equation. It is essentially a mass conservation relation applied to the cross sections of the squeezed air flow and is obtained from the continuity equation by integrating it over the vertical direction with the employment of the momentum equation, and now the generalized Reynolds equation is in routine use to calculate the air bearing parameters in sliders with complicated air bearing surfaces (see e.g., [3, 1013]). The direct simulation Monte-Carlo (DSMC) method [14] is appropriate for simulation of thin film flows in sliders. But the low information-to-noise ratio for low Mach number flows makes the computational process very time-consuming and only simulation results of miniature (micrometer, not authentic millimeter) sizes are available [1517]. The information preservation (IP) method [18, 19] was proposed to overcome the difficulty of low information-to-noise ratio of DSMC method. It was successfully applied to many low-speed microflows [2026] and the IP results for channel flows compared well with DSMC results, experimental data, and exact kinetic theory [27]. The present paper firstly reviews briefly the successful results obtained by the Reynolds equation and the DSMC method, relatively well-known in the literature. But its main purpose is to introduce the IP method and some recent development of the method in calculating the thin film flow under the air bearing surface not familiarized to the readers and show its feasibility as an alternative tool in calculating the micro gas flows under flying head sliders.

2. The Reynolds Equation

The Reynolds equation was first used in the continuum regime, the derivation is enlightening and with its essential assumption clearly revealed, and its framework can be easily used in the slip flow and transitional flow cases; so we give a simplified version of derivation of it. For simplicity the two-dimensional case is considered (see Figure 1), assuming that the head width 𝐷 is larger than the height 𝑜; so the spanwise (z direction) motion can be neglected.

Writing the continuity equation

𝜕𝜌+𝜕𝑡𝜕𝜌𝑢+𝜕𝑥𝜕𝜌𝑣𝜕𝑦=0(1) in the form

𝜕𝜌𝑣𝜕𝑦=𝜕𝜌𝑢+𝜕𝑥𝜕𝜌𝜕𝑡(2) and integrating it over 𝑦 across the whole flow region yields

0𝑑(𝜌𝑣)=0𝜕𝜌𝑢+𝜕𝑥𝜕𝜌𝜕𝑡𝑑𝑦.(3) The left-hand side of (3) vanishes, as there is no fluid flown into or out of the walls. Interchanging the integration and differentiation gives

𝜕𝜕𝑥0𝜕𝜌𝑢𝑑𝑦+𝜕𝑡(𝜌)=0.(4) For thin film flow with the inertial terms neglected the steady momentum Navier-Stokes equation has the following form:

𝜕𝑝=𝜕𝜕𝑥𝜇𝜕𝑦𝜕𝑢𝜕𝑦.(5) Integration across the gap with the nonslip boundary conditions

𝑢𝑦=0=𝑈,𝑢𝑦==0(6) yields the solution of the streamwise velocity component 𝑢:

𝑦𝑢=𝑈1𝜕𝑝/𝜕𝑥𝑦𝑦2𝜇1(7) Substituting (7) into (4) and accomplishing integration over 𝑦, the following equation is attained (we emphasis here that the assumption of 𝜕𝑝/𝜕𝑦=0 is essentially necessary in the derivation):

𝜕𝜕𝑥3𝜌𝜇𝜕𝑝2𝜕𝜕𝑥=6𝜕𝜕𝑡(𝜌)+𝜕𝑥(𝜌𝑈).(8) This is the general form of the Reynolds equation for the two-dimensional case. By introducing 𝑋=𝑥/𝐿,𝐻=/𝑜,𝑃=𝑝/𝑝𝑜, and the bearing number:

Λ=6𝜇𝑈𝐿𝑝𝑜𝑜2.(9) Equation (8) for steady and two-dimensional case can be written in the normalized form [15]:

𝑑𝐻𝑑𝑋3𝑃𝑑𝑃𝑑𝑑𝑋=Λ𝑑𝑋(𝑃𝐻).(10) The first term of (7) is the slip-less solution of the velocity in the Couette flow when the upper plate is stationary and the lower plate moves towards the right with velocity 𝑈 (see [9] compare with its (5.63) for 𝜁=0), and the second term is the slip-less solution of the velocity in the Poiseuille flow when the 𝑥 axis is aligned along the lower plate (see [9], the second term of (7) can be obtained from (5.69) in [9] by a simple translation of the ordinate 𝑦). Equation (10) shows that the flow rate across any cross section is the sum of the flow rate of the Couette flow and the Poiseuille flow and this rate does not change from one cross section to another in steady flow.

In [9] (see (5.73) there) one sees that the flow rate of the Poiseuille flow with slip boundary condition surpasses that of the slip-less case by a factor:

𝑄𝑃,𝑆𝐿𝑄𝑃,𝐶=1+62𝜎𝜎𝐾𝑛.(11) As for the Couette flow the flow rates have a specific feature, and in the slip and the transitional flow cases and the slip-less case they are identical and have the following value independent of the Knudsen number owing to the symmetry of the flow (see Figure 2):

𝑄𝑃,𝐶=𝜌𝑈2.(12) From the flow rate expressions (11) and (12) for Poiseuille and Couette flows in the slip flow case one can conclude that in the slip flow regime the following Reynolds equation is obtained in place of (10):

𝑑𝑑𝑋1+62𝜎𝜎𝐻𝐾𝑛3𝑃𝑑𝑃𝑑𝑑𝑋=Λ𝑑𝑋(𝑃𝐻),(13) where 𝐾𝑛=𝜆/ is local Knudsen number. When the slip boundary conditions

𝑢𝑦=0=𝜍𝑑𝑢𝑑𝑦,𝑢𝑦==𝜁𝑑𝑢𝑑𝑦,𝜁=2𝜎𝜎𝜆(14) instead of the nonslip boundary condition (6) is employed in solving the momentum equation (5), and the resulted velocity profile is substituted into the mass conservation relation (4), one would arrive at the same slip corrected Reynolds equation (13) [6, 15].

Fukui and Kaneko [7] showed that the solution of the linearized Boltzmann equation for the thin film bearing problem can be decomposed into the solutions of the plane Couette flow and the plane Poiseuille flow. On this basis they derived the generalized Reynolds equation for the thin film air bearing problem by employing the flow rates of the fundamental Poiseuille and Couette flows solved by the linearized Boltzmann equation. This generalized Reynolds equation in the isothermal case can be written as [7],

𝑑𝑑𝑋𝑄𝑃,𝑇𝑅(𝐾𝑛)𝐻3𝑃𝑑𝑃𝑑𝑑𝑋=Λ𝑑𝑋(𝑃𝐻),(15) where 𝑄𝑃,𝑇𝑅(𝐾𝑛) is the flow rate in transitional regime (normalized by the slip-less value 𝑄𝑃,𝐶) calculated from the linearized Boltzmann equation for Poiseuille flow. A tabled database of the calculated values of 𝑄𝑃,𝑇𝑅(𝐾𝑛) for 𝜎=1, 𝜎=0.9, 𝜎=0.8, and 𝜎=0.7 is provided in [8], and a fitted formula approximation for diffuse reflection (𝜎=1) by Robert is recorded in [15] (there the second term on the right-hand side is misprinted as 6𝐴𝜋𝐾𝑛):

𝑄𝑃,𝑇𝑅(𝐾𝑛)=1+6𝐴𝐾𝑛+12𝜋𝐾𝑛log(1+𝐵𝐾𝑛),(16) where 𝐴=1.318889 and 𝐵=0.387361. Alexander et al. [15] used the DSMC method to simulate the short head length air bearing problem (𝐿=5𝜇m, 𝑜=50nm=0.05𝜇m, 𝑈=25m/s, 𝜎=1) and found excellent agreement of the DSMC simulation with the generalized Reynolds equations (15) and (16). But at that time they described the latter as continuum hydrodynamic Reynolds equation corrected for slip (in fact only (13) is such an equation); now we have cognized that the generalized Reynolds equation is a global mass conservation equation, but at each cross section its flow is governed by exact kinetic theoretical equation which is appropriate for transitional regime. The comparison made in [15] for the cases (the ratio of the inlet to outlet heights is kept as 2  : 1)

𝐿=1.5𝜇m, 𝑜=15nm=0.015𝜇m, 𝑈=153.9m/𝑠, 𝜎=1.0𝐿=5𝜇m, o=50nm=0.05𝜇m, 𝑈=25m/s, 𝜎=0.7;𝐿=5𝜇m, o=50nm=0.05𝜇m, 𝑈=307.8m/s, 𝜎=1.0

showed good agreement of the generalized Reynolds equation with the results of DSMC simulation; this just confirms that the generalized Reynolds equation is a mass conservation equation in form (although some framework of the N-S equation has been used) but in fact it balances the flow rates of Poiseuille flow and Couette flow calculated from the kinetic theory; so it has the kinetic theoretical merit and can be used to solve the air bearing problem in the entire transitional flow regime.

By analogy with the above derivation of Reynolds equation for the two-dimensional case, it is a simple matter to derive the unsteady and three-dimensional Reynolds equation in the following form (cf., e.g., the 3D Reynolds equation shown in [3]):

𝜎𝜕=𝜕𝜕𝑇(𝑃𝐻)𝜕𝑋𝑄𝑃,𝑇𝑅(𝐾𝑛)𝐻3𝑃𝜕𝑃𝜕𝑋Λ1+𝜕𝑃𝐻𝜕𝑍𝑄𝑃,𝑇𝑅(𝐾𝑛)𝐻3𝑃𝜕𝑃𝜕𝑍Λ3,𝑃𝐻(17) where

Λ𝑇=𝑡𝜔,1=6𝜇𝑈𝐿𝑝𝑜𝑜2,Λ3=6𝜇𝑊𝐿𝑝𝑜𝑜2,𝜎=12𝜇𝜔𝐿2𝑝020.(18)𝑇 is time normalized by 1/𝜔, 𝜔 being angular velocity of the rotating disc, 𝑍, and 𝑊 is the dimensionless coordinate and velocity in the head width direction 𝑧.

Fukui and Kaneko [7] used the generalized Reynolds equation and obtained lubrication characteristic results valid for large Knudsen numbers. By using Fukui-Kaneko’s generalized Reynolds equations some useful results of the slider design were obtained: Hu et al. [10] investigated partial contact air bearing characteristics of tripad sliders. Hu et al. [11] investigated the air bearing dynamics of two configurations of authentic sized sub-ambient pressure sliders and found the way to ensure the reliability of the unloading performance of the types of sliders considered. Tagawa et al. [12] computed the slider’s responses to nanotextured disc surfaces. The computation techniques are developed for solving the slider air bearing problem [3] and new first- and second-order slip models are introduced to the generalized Reynolds equations [13] by Wu and Bogy.

In deriving the Reynolds equation in the form of (17) we have seen that it is essential to use the assumption 𝜕𝑝/𝜕𝑦=0 (see the note in the brackets before (8)). So the Reynolds equation is applicable undoubtedly for cases, where the assumption 𝜕𝑝/𝜕𝑦=0 is valid, for example, for sliders with smooth thin gas layers between the slider and the disc.

In applying the generalized Reynolds equation to the modern authentic sliders having rails on the peripheries (sometimes with thin terraces on them) and on the center part of the slider and fully recessed regions (see, e.g., the NSIC (National Storage Industry Consortium) slider cited in [3]) one should have caution. The depth of the air bearing surface varies drastically at the vertical wall profile regions. This makes the pressure change drastically along the vertical direction and causes perceivable gas flow along the vertical wall direction. The condition 𝜕𝑝/𝜕𝑦=0 exists no more. The extension of the Reynolds equation in the form of (17) to problems of modern configuration air bearing surfaces is questionable and needs further validation. Some comparisons with DSMC or other touch stone results or experimental data are desirable. For comparison with DSMC results even computation of small sized slider with complex configurations is of value.

3. DSMC Method

The direct simulation Monte-Carlo (DSMC) method [14] has been used for simulation thin film flows in sliders. Alexander et al. [15] studied the gas flow under a flat head surface above the moving drive platter and found excellent agreement with Reynolds equation, but their statement that the Fukui Kaneko’s generalized Reynolds equation was a continuum hydrodynamic equation corrected for slip was wrong, and their simulations were restricted to small length sliders (𝐿=5𝜇m). Huang et al. [16] extended the DSMC simulation from two-dimensional to three-dimensional cases and again compared with Fukui–Kaneko’s generalized Reynolds equation (they called it MGL model) and found that, overall, the two solutions agree well with each other and the agreement is better when the spacing decreases to about 5 nm, that is, when the Knudsen number is large. We stress here that the Fukui-Kaneko’s generalized Reynolds equation is a global mass conservation equation, but the detailed flow field is calculated using the exact solution of the linearized Boltzmann equation. So the agreement between the DSMC results and the results of generalized Reynolds equation is not surprising. The size simulated is again very small (𝐿=4𝜇m) in comparison with the authentic size of the modern slider. Wang and Liu [17] by using the DSMC method studied the pressure distribution in head/disc interface at the same position of different disc speeds and at different positions of the same disc speed and the lift force change percentage when air flow velocity changes or slider flying height changes. The size studied was again of the order of 4𝜇m, far smaller than the authentic size of 1 mm. It is the low information-to-noise ratio for low Mach number flows that makes the computational process of DSMC very time-consuming and only simulation results of miniature (not authentic) sizes are available. But the DSMC method is valuable in the microgas flows under flying sliders, even not being able to simulate authentic size, as it has the testing merit, it can be used for small sized but complicated configurations to testify other methods.

4. Information Preservation (IP) Method and Test of It by Channel Flows

The information preservation (IP) method is proposed in [18, 19] to treat the problems encountered by the DSMC method of the huge ratio of the noise to the useful information and the demand of extremely large sample size. This is a method imbedded in the DSMC method in which each simulated molecule is assigned two velocities: thermal velocity 𝐜 and information velocity 𝐮𝑖. The former is just the molecular velocity 𝐜 in the DSMC method and is used to calculate the motion, collision, and the reflection of molecules at the surfaces following the same algorithms and models as the DSMC method. Besides 𝑐 we suppose that each molecule carries the so-called information velocity (IP velocity) 𝑢𝑖 to record the collective velocity of the enormous number of real molecules represented by each simulated molecule. The IP velocities do not produce any influence on the motion of molecules and are used only for summation to obtain the macroscopic velocities; the primitive information is taken from the oncoming flow and the body surface. When the molecules reflect from the surface, collide with each other, experience force action, and enter from boundary, the IP velocities attain new values. The IP method also assigns each cell the IP velocity 𝑈𝑖 and IP density 𝜌. The comprehensive elucidation of the rules governing the renewing of the information velocities 𝑢𝑖 and 𝑈𝑖 and 𝜌 can be found in [9] and [28]. In general 2D and 3D cases the mass and momentum conservation equations for 𝑈𝑖 and 𝜌 are

𝜕𝜌𝜕𝑡𝑑𝑉=𝜌𝑈𝑖𝜌𝑙𝑑𝑆,(19)𝑑𝑈𝑖𝑑𝑡𝑑𝑉=𝑝𝑙𝑑𝑆.(20) The integrals are taken over the whole volume and surface of the cell concerned, respectively. After a time step Δt the cell IP density attains increment according to (21):

Δ𝜌=Δ𝑡Δ𝑉𝜌𝑈𝑖𝑙𝑑𝑆,(21) from where the density and pressure are also renewed: 𝑝=𝑛𝑘𝑇. The increment of the IP velocity 𝐮𝑖 is, according to (22),

Δ𝑢𝑖=Δ𝑡𝜌Δ𝑉𝑝𝑙𝑑𝑆.(22) The IP velocity 𝑈𝑖 of the cell is obtained as the average of 𝐮𝑖 in the cell.

Omitting description of IP simulation of many different flow fields we give here only some detailed simulation of micro channels, for the geometrical forms and flow patterns of two cases of channel and disc flow are similar. Besides there have been two means to test the calculation of the channel flow to show the validity of the IP method. Firstly abundant experimental results of pressure distribution and flow rated through the channel for Knudsen numbers in the transitional regime. And the degenerated Reynolds equation suited for channel flows has been suggested to be derived easily for the generalized Reynolds equation originally derived for the flows under flying sliders and has the merit of kinetic theoretical test stone. So we use both the experimental data and strict theoretical computations to show the validity of the IP method for calculation of the microflows.

The channel flow of fluids has a long history beginning with the famous Poiseuiulle flow, that is, the fluid flow driven by constant pressure gradient. But for gas because of the global mass conservation the pressure distribution along the channel could not be linear and the gas Poiseuille flow with constant pressure gradient is only a hypothetic flow. With the emergence of the MEMS (micro-electromechanical systems) and the technique of manufacture of them, LIGA (abbreviation of German words Lithographie Galvanoformung Abformung) and EDM (electrodischarge machining), fine experiments were carried in new micromachined channels (see, e. g., [2934]). Microchannel/pressure sensor systems with microsized pressure sensors integrated with the flow system were fabricated and exact pressure distributions were obtained. The flow rates were also measured exactly. Nitrogen, argon, and helium have been used as the fluid medium. The Knudsen numbers of 0.5 to 2.5 and larger are obtained and provide database to test the transitional flow. Distinguished nonlinear pressure distributions were also obtained. In [30] the ratio of the pressure difference of the inlet to the outlet is approximately 0.341.7. In [31] for helium pressure difference ratio is approximately 0.591.29. In [34] the pressure difference ratio is even higher; for argon it is 0.942.96. The experimental work provides a solid foundation for testing various methods intended to calculate the microflows.

The generalized Reynolds equation originally is derived for application in the thin film air bearing problem with the lower plate moving with a velocity 𝑈 and the upper plate tilted. Shen suggests to use this Reynolds equation to solve the microchannel flow problem [27] in which the lower plate is stationary and the upper plate is parallel to the lower one. Owing to the steadiness of the lower plate the right-hand side term vanishes; as 𝑈=0 and Λ=0, there is no any contribution of the Couette flow. Owing to the parallelity of the two plates the value 𝐻 is a constant and also can be dropped from the equation. So the generalized Reynolds equation for application to the microchannel problems is degenerated to the form:

𝑑𝑑𝑋𝑄𝑃,𝑇𝑅(𝐾𝑛)𝑃𝑑𝑃𝑑𝑋=0.(23) The values of 𝑃 on the inlet and outlet of the channel are to be specified to make the microchannel problem solvable. For the case of diffuse reflection, the fitted formula approximation of 𝑄𝑃,𝑇𝑅(𝐾𝑛), (16), can be used, and the degenerated Reynolds equation attains the form:

𝑑𝑑𝑋1+6𝐴𝐾𝑛+12𝜋𝑃𝐾𝑛log(1+𝐵𝐾𝑛)𝑑𝑃𝑑𝑋=0.(24) For the ease of integration the local Knudsen number 𝐾𝑛 is most conveniently expressed through 𝑃; for example, for hard sphere model it can be written as

𝜆𝐾𝑛==𝐶𝑃,(25) where

𝜇𝐶=𝑝0𝜋𝑅𝑇02=𝜆0=𝐾𝑛out,(26) for we have for hard sphere

𝜇𝜆=𝑝𝜋𝑅𝑇2(27) (see [9, (2.222)]). 𝑝0 is the pressure at the outlet, 𝑇0 is the temperature of the gas, and 𝜇 is the viscosity of the gas at 𝑇0. The constant 𝐶 has the physical meaning of the Knudsen number 𝐾𝑛out at the outlet of the channel (see (25), at outlet 𝑃=1). Substituting (25) into (24), one arrives at

𝑃+6𝐴𝐶+12𝜋𝐶log1+𝐵𝐶𝑃𝑑𝑃𝑑𝑋=𝐶𝑚,𝑇𝑅,𝐶𝐹𝑅,(28) where 𝐶𝑚,𝑇𝑅,𝐶𝐹𝑅 is an unspecified constant to be determined from the integration and has the physical meaning of the flow rate across the channel normalized by the slip-less flow rate value. Equations (23) and (24) are in fact the generalized Reynolds equation, degenerated for the microchannel flows from where pressure distributions can be calculated. The flow rate across the channel 𝐶𝑚,𝑇𝑅,𝐶𝐹𝑅 has also been integrated in closed [35]

𝐶𝑚,𝑇𝑅,𝐶𝐹𝑅=121𝑃2𝑖+6𝐴𝐾𝑛𝑜1𝑃𝑖+12𝐾𝑛𝑜𝜋×log1+𝐵𝐾𝑛𝑜𝑃𝑖log1+𝐵𝐾𝑛𝑜𝑃𝑖+𝐵𝐾𝑛𝑜log1+𝐵𝐾𝑛𝑜𝑃𝑖+𝐵𝐾𝑛𝑜.(29) And the dimensional corresponding mass flow rate is


We emphasize here, that as the assumption 𝜕𝑝/𝜕𝑦=0 is strictly satisfied here for flat channel, the degenerated Reynolds equation has a strict kinetic theoretical foundation.

The microchannel flow problem was solved by IP method in [26]. The approach, where the pressure 𝑝 is fixed at the ends as the same of the prescribed (experimental) condition and 𝑈𝑖 is allowed to change continuously in the process of computation and finally to reach the steady solution, is used here. In the practice of general IP method the IP process has no reverse influence on the DSMC process. In the solution of channel flow, where the velocities on the boundaries are to be regulated during the simulation, the varying IP velocities on the boundaries are used to continuously adjust the boundary conditions of the DSMC–IP procedure. This enables the DSMC much quickly to have the correct value on the boundaries. To make the calculation process convergent it is essential to use the conservative form of the mass conservation equation. And a superrelaxation technique is employed to speed up the convergence process. For more detailed elucidation of the solution of the channel flow problem by the IP method see [9, 26]. Here only the comparison of some results of pressure distribution given by the IP method and the experimental data and the degenerated Reynolds equation are given (see Figures 3, 4, and 5).

As shown in the above comparisons the IP method gives pressure and mass flow results of microchannels in excellent agreement with the degenerated Reynolds equation and also experimental data; so it can be considered as a reliable tool in dealing with the microscale slow gas flows including the micro flows under flying sliders. The detailed IP method simulation of internal rarefied gas flows can be found in [28] and [9].

5. Microflows under Flying Head Slider Solved by Information Preservation Method

The IP method was applied to simulate the 2D flat slider problem [36]. The conservative scheme of continuum equation and the superrelaxation scheme suggested in [26] are also used here to accelerate the converging process. As the flat air bearing surface is inclined, except rectangular cells, some incised cells appear, and they are combined with the lower adjacent regular cells to avoid extreme small number of simulated molecules in one cell. At the same time a simple weighted averaging method is introduced to smoothen the no-physical oscillation of the cell density:

𝜌𝑖,𝑗=𝑐1𝜌𝑖,𝑗+1.0𝑐1×𝜌𝑖1,𝑗+𝜌𝑖,𝑗1+𝜌𝑖+1,𝑗+𝜌𝑖,𝑗+14,(31) where 𝑐1 is a weighting factor. The pressure distributions for different slider sizes have obtained by the IP method. The 𝐿=5𝜇m result is in good agreement with DSMC method and the Reynolds equation. For 𝐿=25𝜇m good agreement is obtained between the IP and Reynolds equation results. It is remarkable that while DSMC methods only provided results for 𝐿=5𝜇m, the IP simulation result for authentic size head slider 𝐿=1mm is in excellent agreement with Reynolds equation solution (see Figure 6). This is the only comparison of the Reynolds equation with other methods for authentic-sized slider.

In [37] the IP method was used to solve the flows under 3D flat slider and a miniature head slider of authentic complicated configuration. As flows under flat surfaces have been compared between the DSMC method and the Reynolds equation in [16], here only the comparisons of the IP method result with that of the DSMC method are cited from [37]. Figure 7 shows the pressure under the flat slider obtained by the two methods for flying height of 10 nm and the speed of 𝑈=25 m/s. The length of the slider is 4𝜇m and its width is 3.3𝜇m.

Note that the DSMC sample process starts and lasts 0.2 million time steps but for IP sample only 200 time steps are used. The total computation time (the time of convergence process plus the time of sampling process) of the DSMC method is about 100 times longer than the IP method. When the slider surface speed is lower, for example, when 𝑈=2m/s, the DSMC results with 0.6 million sample time steps still have remarkable fluctuations (see Figure 8).

The 3D slider with complex configuration simulated by IP method in [37] is shown in Figure 9. The flow field is divided into 200×160×40 uniform cells, each containing about 20 molecules. The time step used is 7.41×1011s. The size of the slider is chosen as 4𝜇m×3.3𝜇m×0.5𝜇m, so that the DSMC results can be obtained and used for comparison with the present IP results.

The simulated results of the pressure distribution by both IP method and DSMC method are shown in Figure 10. The general agreement is good with minor differences at the corner regions. In a figure (isoheight line, provided through private communication with Mr. Jun Li), one can clearly see that the pressure varies apparently in the vertical direction along the wall, which would violate the basic assumption that the vertical pressure gradient is constant in deriving the Reynolds equation. This indicates that the employ of Reynolds equation for complex configuration should be with caution. Work is under way to calculate the flow field under the complicated configuration slider of authentic size and will soon be published.

6. Conclusions

(1) The precise calculation of the thin gas film bearings characteristics is an engineering necessity in designing the air bearing surface slider. Various calculation methods are desirable to meet the need of design.

(2) The generalized Reynolds equation is now successfully been used to calculate the air bearing parameters in sliders with complicated air bearing surfaces and accomplished a number of design purposes. As the calculation of complex configuration, slider flows by both DSMC and IP method shows the modern complex configuration sliders may violate the condition 𝜕𝑝/𝜕𝑦=0 which is essential in deriving the Reynolds equation. Some comparisons of the Reynolds equation results with DSMC method or experimental data and some modification or improvement are desirable.

(3) The DSMC method is appropriate for simulation of thin film flows under sliders. But the computational process is very time-consuming and only simulation results of sliders of miniature (micrometer, not millimeter, i.e., authentic) sizes are available. The merit of DSMC in the microflows under flying slider is in that its results can be compared with other methods for small size and complicated configurations to check the other methods’ feasibility.

(4) The IP method has been tested in the problem of micro channel flows by DSMC, experimental and kinetic theoretical results. It was used to calculate the 2D and 3D flows under flat sliders. The result of simulation of flow under an authentic complex configuration slider of small size is in good agreement with the DSMC method. It is a promising alternative method in simulating the microflow field under authentic size, authentic configuration air bearing surfaces.


A:Numerical constant
B:Numerical constant
c:Molecular thermal velocity
D:The width of the slider along 𝑧 direction
:Height of the gap between the head slider and the magnetic disc
0:Minimum height of the gap between the head slider and the magnetic disc
𝐻:Dimensionless height /0
𝐿:The length of the slider along 𝑥 direction
𝑝0:Pressure at the front edge of the slider
𝑄:Flow rate
𝑄𝑃:Flow rate of the Poiseuiulle flow
𝑄𝑃,𝐶:Flow rate of the Poiseuiulle flow for continuum flow case
𝑄𝑃,𝑆𝐿:Flow rate of the Poiseuiulle flow for slip flow case
𝑄𝑃,𝑇𝑅:Flow rate of the Poiseuiulle flow for transitional flow case
𝑄𝑃,𝑇𝑅:Dimensionless flow rate of the Poiseuiulle flow for transitional flow case 𝑄𝑃,𝑇𝑅/𝑄𝑃,𝐶
𝑇:Dimensionless time 𝑡𝜔 (sometimes for temperature: 𝑝=𝑛𝑘𝑇)
𝑢:Gas velocity along 𝑥 direction in Reynolds equation
𝑢𝑖:Molecular information velocity in IP method
𝑈:Boundary velocity of the magnetic plate along 𝑥 direction in Reynolds equation
𝑈𝑖:Information velocity of cells in IP method
𝑣:Gas velocity along 𝑦 direction in Reynolds equation
𝑤:Gas velocity along 𝑧 direction (span wise direction) in Reynolds equation
𝑊:Boundary velocity of the magnetic plate along 𝑧 direction (span wise direction) in Reynolds equation
𝑥:Length direction along the moving disc surface
𝑋:Dimensionless coordinate 𝑥/𝐿
𝑦:Direction vertical to the moving disc surface
𝑧:Width direction along the moving disc surface
Z:Dimensionless coordinate 𝑧/𝐷
𝜌:Density of the gas
𝜇:Viscosity of the gas
𝜎:Reflection coefficient
𝜆:Mean free path of the gas
Λ:Bearing number
Λ1:Bearing number calculated by 𝑈𝐵
Λ3:Bearing number calculated by 𝑊𝐵
𝜁:Factor proportional to 𝜆
𝜔:Angular velocity of the rotating disc.


DSMC:Direct simulation Monte Carlo
EDM:Electrodischarge machining
IP:Information preservation
LIGA:Lithographie Galvanoformung Abformung
MEMS:Micro-electromechanical systems
NSIC:National Storage Industry Consortium.


The author would like to express his appreciation for the support of Natural Science Foundation of China through Grants 10621202 and 90205024. He also thanks Mr. Jun Li for his help in formulating this manuscript.


  1. R. K. Brunner, J. M. Harker, K. E. Haughton, and A. G. Osterlund, “A gas film lubrication study—part III: experimental investigation of pivoted slider bearings,” IBM Journal of Research and Development, vol. 8, pp. 260–274, 1959. View at: Google Scholar
  2. R. C. Tseng, “Rarefaction effects of gas-lubricated bearings in a magnetic recording disk file,” ASME Journal of Lubrication Technology, vol. 97, no. 4, pp. 624–629, 1975. View at: Google Scholar
  3. L. Wu and D. B. Bogy, “Numerical simulation of the slider air bearing problem of hard disk drives by two multidimensional upwind residual distribution schemes over unstructured triangular meshes,” Journal of Computational Physics, vol. 172, no. 2, pp. 640–657, 2001. View at: Publisher Site | Google Scholar
  4. O. Reynolds, “On the theory of lubrication and its application to Mr. Beauchamp Tower's experiments, including an experimental determination of the viscosity of the olive oil,” Proceedings of the Royal Society of London, vol. 40, pp. 191–203, 1886. View at: Google Scholar
  5. W. A. Gross, L. A. Matsch, V. Castell, A. Eshel, J. H. Vohr, and M. M. Wildmann, Fluid Film Lubrication, John Wiley & Sons, New York, NY, USA, 1980.
  6. A. Burgdorfer, “The influence of the molecular mean free path on the performance of hydrodynamic gas lubrication bearings,” Journal of Basic Engineering, vol. 81, p. 94, 1959. View at: Google Scholar
  7. S. Fukui and R. Kaneko, “Analysis of ultra-thin gas film lubrication based on linearized Boltzmann equation: first report—derivation of a generalized lubrication equation including thermal creep flow,” Journal of Tribology, vol. 110, no. 2, pp. 253–262, 1988. View at: Google Scholar
  8. S. Fukui and R. Kaneko, “A database for interpolation of Poiseuille flow rates for high Knudsen number lubrication problems,” Journal of Tribology, vol. 112, no. 1, pp. 78–83, 1990. View at: Google Scholar
  9. C. Shen, Rarefied Gas Dynamics, Fundamentals, Simulations and Micro Flows, Springer, New York, NY, USA, 2005.
  10. Y. Hu, P. M. Jones, P. T. Chang, and D. B. Bogy, “Partial contact air bearing characteristics of tripad sliders for proximity recording,” Journal of Tribology, vol. 120, no. 2, pp. 272–279, 1998. View at: Google Scholar
  11. Y. Hu, P. M. Jones, and K. Li, “Air bearing dynamics of sub-ambient pressure sliders during dynamic unload,” Journal of Tribology, vol. 121, no. 3, pp. 553–559, 1999. View at: Google Scholar
  12. N. Tagawa, T. Hayashi, and A. Mori, “Effects of moving three-dimensional nano-textured disk surfaces on thin film gas lubrication characteristics for flying head slider bearings in magnetic disk storage,” Journal of Tribology, vol. 123, no. 1, pp. 151–158, 2001. View at: Publisher Site | Google Scholar
  13. L. Wu and D. B. Bogy, “New first and second order slip models for the compressible Reynolds equation,” Journal of Tribology, vol. 125, no. 3, pp. 558–561, 2003. View at: Publisher Site | Google Scholar
  14. G. A. Bird, Molecular Gas Dynamics and the Direct Simulation of Gas Flow, Clarendon Press, Oxford, UK, 1994.
  15. F. J. Alexander, A. L. Garcia, and B. J. Alder, “Direct simulation Monte Carlo for thin-film bearings,” Physics of Fluids, vol. 6, no. 12, pp. 3854–3860, 1994. View at: Google Scholar
  16. W. Huang, D. B. Bogy, and A. L. Garcia, “Three-dimensional direct simulation Monte Carlo method for slider air bearings,” Physics of Fluids, vol. 9, no. 6, pp. 1764–1769, 1997. View at: Google Scholar
  17. H. C. Wang and T. S. Liu, “Airflow velocity effects on air bearing with grooved disk surface in near-field optical disk drives,” Physics of Fluids, vol. 18, no. 5, Article ID 057103, 2006. View at: Publisher Site | Google Scholar
  18. J. Fan and C. Shen, “Statistical simulation of low-speed unidirectional flows in transitional regime,” in Rarefied Gas Dynamics, R. Brun, R. Campargue, R. Gatignol, and J. C. Lengrand, Eds., vol. 2, pp. 245–252, Cepadues, Toulouse, France, 1999. View at: Google Scholar
  19. J. Fan and C. Shen, “Statistical simulation of low-speed rarefied gas flows,” Journal of Computational Physics, vol. 167, no. 2, pp. 393–412, 2001. View at: Publisher Site | Google Scholar
  20. C. P. Cai, I. D. Boyd, J. Fan, and G. V. Candlir, “Direct simulation method for low-speed micro-channel flows,” Journal of Thermophysics and Heat Transfer, vol. 14, p. 368, 2000. View at: Google Scholar
  21. J. Fan, I. D. Boyd, C.-P. Cai, K. Hennighausen, and G. V. Candler, “Computation of rarefied gas flows around a NACA 0012 airfoil,” AIAA Journal, vol. 39, no. 4, pp. 618–625, 2001. View at: Google Scholar
  22. Q. H. Sun, I. D. Boyd, and G. V. Candler, “Numerical simulation of gas flow over microscale airfoils,” in Proceedings of the 35th Thermophysics Conference, Anaheim, Calif, USA, June 2001, AIAA paper 2001-3071. View at: Google Scholar
  23. H. L. Liu, C. Xie, C. Shen, and J. Fan, “Flow in membrane filter simulated as microchannel flow with diaphragm,” in Rarefied Gas Dynamics, T. J. Bartel and M. A. Gallis, Eds., vol. 585 of AIP Conference Proceedings, pp. 524–530, Springer, New York, NY, USA, 2001. View at: Google Scholar
  24. Q. H. Sun, I. D. Boyd, and J. Fan, “Development of an IP method for subsonic, microscale gas flows,” in Rarefied Gas Dynamics, T. J. Bartel and M. A. Gallis, Eds., vol. 585 of AIP Conference Proceedings, pp. 547–553, Springer, New York, NY, USA, 2001. View at: Google Scholar
  25. J. Z. Jiang, J. Fan, and C. Shen, “Statistical simulation of micro square cavity flows,” in Rarefied Gas Dynamics, A. D. Ketsdever and E. P. Muntz, Eds., vol. 663 of AIP Conference Proceedings, pp. 784–791, Springer, New York, NY, USA, 2002. View at: Google Scholar
  26. C. Shen, J. Fan, and C. Xie, “Statistical simulation of rarefield gas flows in micro-channels,” Journal of Computational Physics, vol. 189, no. 2, pp. 512–526, 2003. View at: Publisher Site | Google Scholar
  27. C. Shen, “Use of the degenerated Reynolds equation in solving the microchannel flow problem,” Physics of Fluids, vol. 17, no. 4, Article ID 046101, 2005. View at: Publisher Site | Google Scholar
  28. C. Shen, “Information preservation (IP) method in simulation of internal rarefied gas flows in MEMS,” in Rarefied Gas Dynamics, M. Capitelli, Ed., vol. 762 of AIP Conference Proceedings, pp. 141–149, Springer, New York, NY, USA, 2004. View at: Google Scholar
  29. J. Q. Liu, Y. C. Tai, K. C. Pong, and C. M. Ho, “Micromachined channel/pressure sensor systems for micro flow studies,” in Proceedings of the 7th International Conference on Solid State Sensors and Actuators (Transducers '93), pp. 995–998, 1993. View at: Google Scholar
  30. K. C. Pong, C. M. Ho, J. Q. Liu, and Y. C. Tai, “Non-linear pressure distribution in uniform microchannels,” in Application of Microfabrication to Fluid Mechanics, vol. 197 of FED, pp. 51–56, ASME, New York, NY, USA, 1994. View at: Google Scholar
  31. J. C. Shin, C. M. Ho, J. Q. Liu, and Y. C. Tai, “Monatomic and polyatomic gas flow through uniform microchannels,” in Proceedings of the ASME International Mechanical Engineering Congress & Exposition, vol. 59, p. 197, 1996. View at: Google Scholar
  32. E. B. Arkilic, M. A. Schmidt, and K. S. Breuer, “Measurement of the TMAC in silicon microchannels,” in Rarefied Gas Dynamics, C. Shen, Ed., pp. 983–988, Peking University Press, Beijing, China, 1997. View at: Google Scholar
  33. E. B. Arkilic, K. S. Breuer, and M. A. Schmidt, “Mass flow and tangential momentum accomodation in silicon micromachined channels,” Journal of Fluid Mechanics, vol. 437, pp. 29–43, 2001. View at: Publisher Site | Google Scholar
  34. Y. Zohar, S. Y. K. Lee, W. Y. Lee, L. Jiang, and P. Tong, “Subsonic gas flow in a straight and uniform microchannel,” Journal of Fluid Mechanics, no. 472, pp. 125–151, 2002. View at: Publisher Site | Google Scholar
  35. C. Shen and J.-Z. Jiang, “Finite length microchannel flow in transitional regime—analytical prediction compared with experiment,” Acta Aerodynamica Sinica, vol. 26, no. 2, pp. 257–262, 2008. View at: Google Scholar
  36. J.-Z. Jiang, C. Shen, and J. Fan, “Statistical simulation of thin-film bearings,” in Rarefied Gas Dynamics, M. Capitelli, Ed., vol. 762 of AIP Conference Proceedings, pp. 180–185, Springer, New York, NY, USA, 2004. View at: Google Scholar
  37. J. Li, C. Shen, and J. Fan, “IP and DSMC simulation of micro gas flows under 3-D head sliders,” in Rarefied Gas Dynamics, T. Abe, Ed., AIP Conference Proceedings, pp. 1003–1008, Springer, New York, NY, USA, 2009. View at: Google Scholar

Copyright © 2009 Ching Shen. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

More related articles

749 Views | 482 Downloads | 1 Citation
 PDF  Download Citation  Citation
 Download other formatsMore
 Order printed copiesOrder

Related articles

We are committed to sharing findings related to COVID-19 as quickly and safely as possible. Any author submitting a COVID-19 paper should notify us at to ensure their research is fast-tracked and made available on a preprint server as soon as possible. We will be providing unlimited waivers of publication charges for accepted articles related to COVID-19. Sign up here as a reviewer to help fast-track new submissions.