Adaptive multilevel method for the air bearing problem in hard disk drives

(1)

www.elsevier.com/locate/triboint

Adaptive multilevel method for the air bearing problem

in hard disk drives

Chung-Jen Lu

, Shean-Son Chiou, Tai-Kuo Wang

Department of Mechanical Engineering, National Taiwan University, No. 1 Roosevelt Rd. Sec. 4, Taipei 10617, Taiwan, ROC Received 10 May 2002; received in revised form 25 December 2003; accepted 7 January 2004

Abstract

An adaptive grid-generating algorithm is constructed and integrated with the multigrid method to form a numerical scheme that suits slider air-bearing simulation of hard disk drives. The relative truncation error, a by-product of the multigrid method, is used in grid adaptation criteria. Finer meshes are constructed over nodes of the current finest grid where the relative truncation error exceeds a predetermined tolerance. The union of these finer meshes forms a new level of grid, which may not cover the entire domain of the coarse grid underneath. The final grid system thus constructed is composed of levels of uniform grids with decreasing mesh sizes. This composite grid structure incorporates with numerical resolution as needed and efficiency of compu-tation. A shaped rail, negative pressure slider is used to demonstrate the effectiveness of this numerical scheme. Compared with the traditional multigrid method, the proposed adaptive multilevel method can significantly reduce the computation work for achieving the same level of accuracy.

Keywords: Gas-lubricated bearings; Numerical analysis; Reynolds equation; Hard disk drives

1. Introduction

Since the introduction of the first magnetic hard disk drive, the head-disk spacing has consistently dimin-ished to meet the demand for higher data storage den-sity. The flying heights of recent sliders are in the range of 20–30 nm or below. Besides the extremely low spa-cing, it is also desirable to have a constant flying height over the entire disk, a low take-off/landing velocity, and a stiff air bearing. In order to meet these strict per-formance requirements, air bearings of significant com-plexity have been proposed [1–3]. New features in air bearing designs include shaped rails, multiple etch depths, and recessed regions. Due to the complicated air-bearing surface topography, the pressure distri-bution under the slider can only be obtained by numerical methods. Therefore, an efficient slider simu-lator is an indispensable tool for the slider design.

The pressure distribution between the slider and the disk is governed by the generalized Reynolds equation. Considerable attention has been drawn to the develop-ment of efficient and stable numerical schemes for solv-ing the generalized Reynolds equation [4–11]. As the modern air-bearing surface (ABS) becomes increasingly complex, the pressure varies dramatically in the area where the ABS height changes suddenly. In order to have a clear picture of the pressure distribution, a fine grid should be used such that the numerical error is below a given border everywhere. This goal can be achieved by using a uniform grid with small mesh size. However, this approach usually leads to a large num-ber of equations and hence significantly deteriorates the efficiency of the numerical simulator. Furthermore, since high resolution is only required over high press-ure-gradient regions, which form a proper subset of the entire domain, it is impractically to employ fine meshes globally. Therefore, it is highly desirable to have a robust adaptive grid generator that concentrates high resolving meshes to the critical parts of the domain.

The adaptive grid generators generally have to fulﬁll two basic goals. First, they have to supply the

_{Corresponding author. Tel.: 23627686; fax:}

+886-2-23631755.

E-mail address: [email protected] (C.-J Lu).

(2)

discretization approach with a domain partitioning method capable of adapting the size of the discretiza-tion cells locally, and second, they should supply the adaptive simulation approach with a method for the evaluation and control of the discrete approximation errors during simulation. The multigrid method is not only an eﬃcient method for solving nonlinear discrete problems but also provides a cheap and reliable way for estimating the discretization error [12–14]. There-fore, it is nature to incorporate the multigrid method with an adaptive grid selection strategy to form an eﬃcient slider simulator.

Lu [10]developed a slider simulator using an adapt-ive multigrid scheme. He employed a geometric pro-gression scheme to generate the initial rectangular grid. After the pressure solution is obtained on the initial grid, the grid lines are redistributed so that the grid is more concentrated in areas where the pressure gradient is high. Their adaptive grid generation scheme has two shortcomings. First, because high resolution is achieved by adjusting the positions of grid lines, local modifi-cation of a discretization cell implies the modifimodifi-cation of all the discretization cells along the grid lines. Hence, unnecessary fine meshes may be introduced to where the pressure is smooth. Second, the number of grid lines of the finest grid is fixed; this may restrict the highest resolution and result in meshes with large aspect ratio. This kind of non-uniform grid may cause lower-order discrete approximation. Wu and Bogy [11]

combined the multigrid method with an unstructured adaptive triangular mesh generation scheme to form an efficient slider simulator. The unstructured triangular mesh makes it easy to describe the complex geometry of shaped sliders and provides local refinement only to the critical regions. However, the use of unstructured irregular meshes requires complicated data structures and huge information overheads for the representation and modification of grid structures. It is also difficult to apply some efficient iteration schemes, e.g. line relax-ation, on such arbitrary grids.

In this study, we present an efficient numerical method, which is suitable for slider air-bearing simulations, based on the multilevel adaptive technique developed by Brandt[12]. A sequence of uniform rectangular grids is used to achieve non-uniform resolution. This capability can be obtained because various grids (levels) need not all extend over the same domain. Finer levels may be confined to increasingly smaller subdomains, so as to provide higher resolution only where desired. Then the multigrid method is used to solve the discrete problem. Because only uniform grids are used, the adaptive grid structures can be defined in terms of simple and man-ageable data structures that naturally fit into the multigrid method. Finally, a tri-rail slider is used to demonstrate the efficiency of this method.

2. Discretization of the governing equation

The pressure distribution of gas-lubricated bearings is governed by the classical Reynolds equation. How-ever, due to the low flying height of the slider in cur-rent hard disk drives, the classical Reynolds equation, which assumes no-slipand continuous flow, is no longer valid. Several modifications have been proposed based on slipboundary conditions and Boltzmann equation. These correction models can be put into the following non-dimensional steady state generalized Reynolds equation @ @X QPH 3@P @X KxPH þ @ @Y QPH 3@P @Y KyPH ¼ 0; ð1Þ

where P¼ p=pa is the dimensionless pressure, H ¼

h=hm the dimensionless bearing height, X ¼ x=L the

dimensionless x-coordinate, and Y¼ y=L the dimen-sionless y-coordinate, in which pa, hm, and L indicate

the ambient pressure, the ﬂying height, and the length of the slider, respectively. Kx¼ 6lUL=pah2m and Ky ¼

6lVL=pah2m are the bearing numbers in the x- and

y-directions, respectively, in which U and V are the x and y velocity components, respectively. Q is the ﬂow fac-tor that assumes diﬀerent forms depending on the type of correction model used[10].

The control volume method is employed to discretize the generalized Reynolds equation. Integration of Eq. (1) over the control volume as shown inFig. 1yields Je Jwþ Jn Js¼ 0;

where J indicates the integrated flux and the lower-case subscript denotes the control surface at which the flux is evaluated. The fluxes crossing the control surfaces can be further represented in terms of the values of clear-ance and pressure at neighboring grid points using the generalized formulation proposed by Patankar[15]. The final discretization form of the equation at grid point P

(3)

can be written as[10]

aPPP¼ aEPEþ aWPW þ aNPNþ aSPSþ fP; ð2Þ

where the capital subscripts refer to the grid points. Since the coeﬃcients a’s in Eq. (2) are functions of the pressure, the discrete equation is nonlinear.

3. Full approximation storage algorithm

The discretized system of equations can be expressed in the matrix form

Lhuh¼ fh; ð3Þ

where Lh is the coeﬃcient matrix formed by the

coeﬃ-cients aP, aE, aW, aN, and aS;uh the vector of pressure on

all grid points, fh the source vector. Note that the

coef-ficient matrix, pressure vector, and source vector all depend on the mesh size h. In this paper, the full approxi-mation storage (FAS) algorithm in Brandt[12], which is well suited for nonlinear equations, is employed for solv-ing Eq. (3). The key idea of the FAS algorithm can be understood by considering the case that only two levels of grid are used. The mesh sizes of the fine and coarse grids are h and H, respectively. A suitable iteration method is used on the fine grid to find an approximate solution ûuhto

Eq. (3). Then the algebraic error in ^uuhor the correction is

^eeh¼ uh ^uuh: ð4Þ

The approximate solution ^uuh will not satisfy Eq. (3)

exactly. The failure is the residual or defect

rh¼ fh Lh^uuh: ð5Þ

Using Eqs. (4) and (5), the fine grid Eq. (3) is replaced by Lhðûuhþ ehÞ Lhûuh¼ rh: ð6Þ

The right-hand side is smooth after a few nonlinear relax-ation sweeps. Then we can transfer the above equrelax-ation to a coarse grid:

LHðuHÞ LH IhH^uuh

¼ I_hHrh; ð7Þ

where IH

h is some restriction operator which transfers

variables or residues from the ﬁne grid to the coarse grid. Let ^uuH denote the approximation solution to Eq. (7).

Then the coarse-grid correction is êeH ¼ ûuH IhHûuh;

Finally, the coarse-grid correction is transferred back to the ﬁne grid to obtain the improved ﬁne grid approxi-mation

^ u

uh ^uuhþ IHh uu^H IhH^uuh;

where Ih

His an interpolation operator which transfers the

correction from the coarse grid to the ﬁne grid. The above procedure can be applied recursively to form a complete multigrid V-cycle[14].

4. Stopping criterion

An essential problem of all iteration processes is when to stop the process. A loose stopping criterion generates useless solutions. On the other hand, a too strict stopping criterion results in an inefficient process. One benefit of the multigrid method is that it provides a natural and optimal stopping criterion. In order to understand this stopping criterion, we analyze the errors introduced by the numerical process first.

Rewrite the generalized Reynolds equation as

LðuÞ ¼ f ; ð8Þ

where L stands for the differential operator, u the exact differential solution, and f the source term. The aim is to have a solution of Eq. (8). But it is only poss-ible to compute an approximate solution ûu_h to the dis-crteized Eq. (3). Let Ih be a projection operator transferring the continuous solution into the discrete function space. The final total error can be expressed as

^

eeh¼ Ihu ^uuh¼ Ihu uhþ uð h ^uuhÞ;

where uh is the exact solution of Eq. (3). The above

equation shows that the total error is composed of the global (discretization) error eh¼ Ihu uh, introduced

by the discretization process, and of the algebraic error uh ^uuh, caused by the numerical solution of the discrete

equation. It is appropriate to require the maximum total error norm of magnitude e:

^ee_h

k k Ihu uh

þ uk h ^uuhk < e: ð9Þ

This condition is satisﬁed if Ih_u_u h

< e=2 and uh ^uuh

k k < e=2 independently. The first condition determines the mesh size of the finest grid. The second condition implies that it is sufficient to solve the discrete problem up to the level of global error. Therefore, a natural stopping criterion for a multigrid iteration is

uh ^uuh

k k Ih_u_u h

:

Since the exact discrete solution uh is unknown in

gen-eral, the residual or defect rh as deﬁned by Eq. (5) is

employed instead to estimate the algebraic error. In this case, the stopping criterion may be expressed

rh

k k c Ih_u_u h

; ð10Þ

where c is a suitable parameter. Because an analytic rep-resentation of the global (discretization) error is imposs-ible in general, another indicator for identifying how well the discrete problem approximates the diﬀerential problem is required. Here, the local error as deﬁned by sh¼ LhðIhuÞ fh

(4)

by

Lhðuhþ ehÞ ¼ fhþ sh:

That means, the discrete solution of the discrete prob-lem, where the right-hand side is modified by adding the local error, coincides with the corresponding projection of the differential solution Ihu. One benefit of the FAS algorithm is that it offers a numerically cheappossibility to estimate the local error with the aid of two consecu-tive grids. Define the relaconsecu-tive truncation error sh

Hby

sh_H ¼ LH IhHuh IhHLhðuhÞ: ð11Þ

This quantity stands for the error that is caused by the substitution of uh restricted to the coarse grid into the

coarse-grid equation. Since LhðuhÞ ¼ fh, Eq. (11) can be

rewritten as

sh_H ¼ LH IhHuh fH: ð12Þ

The local error sH on the coarse grid may be

con-sidered composed of the local error sh on the ﬁne grid

and the relative truncation error on the coarse grid with respect to the ﬁne grid:

sH ¼ IhHshþ shH:

Assuming the existence of asymptotic expansion for the local error, it is easy to show that[16]

sh H ¼

Hp_hp

Hp sH: ð13Þ

Thus, the relative truncation error sh

H equals the local

error on the coarse grid, upto higher order terms and a factor depending on the mesh size ratio and the order of consistency. It follows that the local truncation error may be used as an indicator of the local error. For a standard coarse mesh size, H¼ 2h, Eq. (13) can be rewritten as sH ¼ 2p 2p₁s h H

Since the exact discrete solution uh is unknown, we

cannot exactly compute the local truncation error. But we can have an approximation to it from using ^uuh in

Eq. (12):

s_Hh ^ssh_H¼ LH IhH^uuh fH:

Using this relation, the stopping criterion (10) may be expressed as rh k k c 2 p 2p₁ ^ss h H : ð14Þ

5. Local grid reﬁnement

As discussed in the previous section, the total error is composed of the global (discretization) error and the algebraic error. Eq. (9) is satisfied if both the global and algebraic errors are less than a predetermined tol-erance. Therefore, it is sufficient to solve the discrete problem on the finest grid with the goal that the algebraic error is at the level of the global error. The global error is a function of the mesh size. Let ehðPÞ be

the global error estimator at node P of the grid with mesh size h. Then the mesh size h of the ﬁnest grid is chosen to guarantee that

ehðPÞ < Tol; ð15Þ

where Tol represents the predetermined tolerance. This condition implies that local refinement is required at nodes where the above relation is violated. In this way, finer grids may be confined to increasingly small domain, so as to provide higher resolution only where needed. The procedure for constructing the finer grid is described below.

Let Gl be the currently ﬁnest grid with mesh size

hl¼ h0=2l and P¼ ðxp;ypÞ a node of Gl. If the global

error estimator at node P is larger than the pre-determined tolerance, a grid patch

Elþ1ðPÞ ¼ QjQ ¼ ðx pþ ahlþ1;ypþ ahlþ1Þ; jaj 2

is introduced on the next ﬁner grid Glþ1 with mesh size

hlþ1 ¼ hl=2. In the sequel, this is referred to as the

elementary reﬁnement patch. The complete grid Glþ1 is

obtained by the union of all elementary reﬁnement pat-ches introduced over all Gl nodes, where a better grid

resolution is required. Fig. 2 shows a typical two-grid adaptive structure. The ﬁner grid Glþ1 consists of three

elementary patches introduced over nodes P1, P2, and

P3 of Gl, respectively. Note that the ﬁner grid Glþ1 is

disconnected and only covers part of the entire

(5)

domain. A grid that covers the whole domain is called a global grid. Otherwise, it is local.

The above refinement process can be applied repeat-edly until the requirement on the discretization error is satisfied all over the finest grid. Then the final multi-level adaptive grid structure is a set of uniform grid G¼ Gf 1;G2; . . . ;GMg;

where Gl refers to the uniform grid that covers the

domain Xl on the lth level. The mesh size of Gl is

hl ¼ h1=2l1. On the ﬁrst mg ð1 mg MÞ levels the

grids are global, while the other uniform grids on level l > mg are local ones deﬁned on subdomains for which

Xl Xl1 holds. M denotes the ﬁnest level. By using

properly aligned elementary patches to construct uni-form local grids on diﬀerent levels, the globally required non-uniform resolution of the discrete approximation can be obtained. This leads to important simpliﬁcations both for the discretization and for the handling of the adaptive grid structure.

6. FAS on composite grid

Note that the ﬁnal multilevel grid structure may be composed of grids with diﬀerent domains. In other words, some nodes of Gl1 may not belong to the next

ﬁner grid Gl. Therefore, some modiﬁcations on the

FAS algorithm, which works well on global grids, are required before we can apply it on local grids. The modiﬁed FAS scheme is described brieﬂy on a two-level grid structure. Let Gh and GH denote two

con-secutive grids, and G0

h and GH0 the corresponding sets

of interior points of Gh and GH, respectively (see

Fig. 2). For a reﬁnement local grid Gh, some of the

boundary points, Gh Gh0, are internal boundary

points (gray solid points in Fig. 2). At these points, function values should be interpolated properly from the corresponding points on GH. The discrete equation

together with the boundary conditions can be expres-sed as

Lhuh¼ fh on Gh0

and

uh¼ IHhuH on Gh G0h:

For those interior points of GH which do not belong

to G_h0, GH is the ﬁnest discretization level. A solution

to the discrete problem has to be computed there. On the other hand, the interior points of GH which are

also interior points of Gh are used to calculate

correc-tions to the approximate solution on G0

h. By deﬁning FH¼ LHI_hHuh on G0H\ G0h fH on G0H GH0 \ G0h ;

the discrete problem on the coarse grid can be expres-sed in a compact form as (compared with Eq. (7)) LHuH¼ FH:

This modiﬁed two-grid FAS scheme can be applied recursively as before to form a complete V-cycle multi-grid scheme.

7. Adaptive multilevel method

The grid reﬁnement process and the solution algor-ithm as described above are combined together to form the complete adaptive multilevel method. Fig. 3shows the ﬂow chart of this method. In this method, the glo-bal errors on the entire domain are controlled adap-tively to be under a predetermined tolerance. The

(6)

essential phase of the adaptive control of global errors is their evaluation for the given discrete problem. Here, the s-extrapolation of the relative truncation error on two consecutive grids is employed as the discretization error estimator[12]. Therefore, two global grids G1 and

G2 are required to initialize the adaptive process. The

complete algorithm is described brieﬂy below. Let GM

(M 2) denote the currently ﬁnest grid. The initial guess on GM is obtained by interpolating from the next

coarser grid GM1. Then the V-cycle FAS algorithm is

applied repeatedly until the algebraic error on the cur-rently finest grid is less than the discretization error. After that, new refinement is constructed over the cur-rently finest grid where the discretization error exceeds the predetermined tolerance. These steps are recursively repeated until no further refinement is needed.

8. Results and discussion

The static pressure of a 30% tri-rail slider[17]is calcu-lated using the multilevel numerical scheme developed in this paper. The ﬂying height is 15 nm, the pitch angle is 180 lrad, the sliding velocities in the x- and y-directions are 13.47 and 16.28 m/s, respectively.

Fig. 4 illustrates the air-bearing surface of the 30% tri-rail slider. This slider has a center rail that carries the read–write element. The connected front region enables the eﬃcient generation of the sub-ambient pressure in the central recessed regions. The rail has concave shapes on both sides to minimize the ﬂying height change across the disk. Fig. 5(a) shows the 3-D pressure distribution over the entire slider. Along

the outer shaped rails, the air is pressurized due to the wedge effect provides by the front taper and flat step. Over the central portion of the slider, the flow expands rapidly across the step discontinuity and hence generates sub-ambient pressure. In order to see the sub-ambient pressure clearly, the pressure profile along the line y=L¼ 0:3 is shown in Fig. 5(b). On the front taper and flat step, the pressure rises to the maximum value of about 1.4 and drops right down to 0.36 passing the flat step. In the front recessed region, the pressure increases due to the convergent cross-section. The pressure has a steep rise on the shaped rail and another rapid drop crossing the shaped rail. Then it increases monotonically to the ambient pressure. Note that the pressure gradient in the center recessed region is higher than that on the front taper and flat step.

Figs. 6(a) and (b) show the composite grids at the fourth and fifth level, respectively. From these two figures, it can be seen that the regions with sudden geo-metric change or high pressure gradient have been efficiently covered by fine meshes generated by the adaptation scheme. Almost all the center recessed regions are refined at the fifth level due to the large pressure variations as shown inFig. 5(b). Local refine-ments and the composite grid structure can be seen clearly fromFigs. 6(a) and (b).

In order to see the effects of the mesh size on the numerical results, the normal force obtained at each level of grid is plotted in Fig. 7, where Tol stands for the refinement criterion. The coarsest grid consists of 31 31 nodes. Smaller values of Tol indicate that more elementary refinement patches are introduced at each

(7)

level. As can be seen from the ﬁgure, the normal force converges as the number of levels increases for each value of Tol. The results for Tol¼ 20 and 30 are almost identical. Hence, it can be concluded that the grid struc-ture to the 8th level for Tol¼ 30 is appropriate for describing the pressure distribution under the slider.

Fig. 8 shows the total number of nodes upto the current finest grid for Tol¼ 0 and 30. Tol ¼ 0 indi-cates that uniform global grid is used at each level. In this case, the node number increases exponentially with the level number. As can be seen from the figure, the node number for Tol¼ 30 is much less than that for Tol¼ 0. The total number of nodes upto the 8th level for Tol¼ 30 is about 1/20 of that for Tol ¼ 0. This shows that the proposed adaptive multilevel method is significantly more efficient than the non-adaptive multi-grid method since the computation work increases rap-idly with the node number.

Fig. 5. Pressure distribution of the 30% tri-rail slider: (a) 3-D press-ure proﬁle at the fourth level over the entire slider, and (b) presspress-ure

proﬁle along the line y=L¼ 0:3. Fig. 6. Composite grids at (a) the fourth level and (b) the ﬁfth level.

Fig. 7. Normal force versus level number for various reﬁnement criteria.

(8)

9. Conclusions

Due to the complicated air-bearing surface of slider in hard disk drives, efficient numerical methods are required for analyzing the pressure distribution under the slider. The usual approach is first to discretize the generalized Reynolds equation in some preassigned manner, and then to submit the resulting discrete equa-tions to some numerical solution process. In this study, the discretization and solution processes are combined together, and greatly benefit from each other. An adaptive grid-generating scheme is implanted to dis-cretize the computation domain. The final grid system consists of a sequence of uniform grids (or levels) with decreasing mesh sizes. The domain of any grid may be only a proper part of the domains of the coarse grids. This structure is very flexible and allows an efficient handling of the refinement. The FAS algorithm, which suits well for solving nonlinear equations, is used to obtain solutions on these levels of grids. The grid adap-tation is governed by the relative truncation errors sup-plied by the FAS algorithm, and can naturally be integrated with the FAS algorithm to give high resolution where needed. Calculations on a 30% tri-rail slider show that this multilevel FAS method is much more efficient than the traditional multigrid method using uniform global grids.

Acknowledgements

This research was supported by the National Science Council of the Republic of China under Grant No. NSC-89-2212-E-002-109.

References

[1] White JW. Flying characteristics of the zero-load slider bearing. Journal of Lubrication Technology 1983;105:484–90.

[2] Harde C, Menon A, Crane P, Egbert D. Analysis and perform-ance characteristics of the Seagate advperform-anced air baring slider. IEEE Transactions on Magnetics 1994;30(2):424–32.

[3] White JW. Flying characteristics of the transverse and negative pressure contour (‘‘TNP’’) slider air bearing. ASME Journal of Tribology 1997;119:241–9.

[4] Tang T. Dynamics of air-lubricated slider bearings for non-contact magnetic recording. ASME Journal of Lubrication Technology 1971;93:272–8.

[5] White JW, Nigam A. A factored implicit scheme for the numeri-cal solution of the Reynolds equation at very low spacing. ASME Journal of Lubrication Technology 1980;102:80–4. [6] Ruiz OJ, Bogy DB. A numerical simulation of the head-disk

assembly in magnetic hard disk: 1. Component models. ASME Journal of Tribology 1990;112:593–602.

[7] Wahl MH, Lee PR, Talke FE. An eﬃcient ﬁnite element-based air bearing simulator for pivoted slider bearings using bi-conju-gate gradient algorithms. Tribology Transactions 1996;39(1): 130–138.

[8] Cha E, Bogy DB. A numerical scheme for static and dynamic simulation of subambient pressure shaped rail sliders. ASME Journal of Tribology 1995;117:36–46.

[9] Hu Y. Head-disk-suspension dynamics. Doctoral Dissertation, Department of Mechanical Engineering, University of California, Berkeley, 1996.

[10] Lu S. Numerical simulation of slider air bearings. Doctoral Dis-sertation, Department of Mechanical Engineering, University of California, Berkeley, 1997.

[11] Wu L, Bogy DB. Unstructured adaptive triangular mesh gener-ation techniques and ﬁnite volume schemes for the air bearing problem in hard disk drives. ASME Journal of Tribology 2000;122:761–70.

[12] Brandt A. Multilevel adaptive solutions to boundary-value pro-blems. Mathematics of Computation 1977;31(138):333–90. [13] Hackbusch W, Trottenberg U, editors. Multigrid methods.

Springer-Verlag; 1982.

[14] Briggs WL. A multigrid tutorial. Philadelphia (PA): Society for Industrial and Applied Mathematics; 1987.

[15] Patankar SV. Numerical heat transfer and ﬂuid ﬂow. New York: McGraw-Hill; 1980.

[16] Joppich W, Mijalkovic S. Multigrid methods for process simula-tion. New York: Springer-Verlag; 1993.

[17] Kang TS, Choi DH, Jeong TG. Optimal design of HDD air-lubricated slider bearings for improving dynamic characteristics and operating performance. ASME Journal of Tribology 2001;123:541–7.