- Email: [email protected]
Study on dynamic characteristics for high speed water-lubricated spiral groove thrust bearing considering cavitating effect
Journal Pre-proof Study on dynamic characteristics for high speed water-lubricated spiral groove thrust bearing considering cavitating effect Xiaohui Lin, Ruiqi Wang, Shaowen Zhang, Shuyun Jiang PII:
S0301-679X(19)30539-0
DOI:
https://doi.org/10.1016/j.triboint.2019.106022
Reference:
JTRI 106022
To appear in:
Tribology International
Received Date: 15 August 2019 Revised Date:
14 October 2019
Accepted Date: 15 October 2019
Please cite this article as: Lin X, Wang R, Zhang S, Jiang S, Study on dynamic characteristics for high speed water-lubricated spiral groove thrust bearing considering cavitating effect, Tribology International (2019), doi: https://doi.org/10.1016/j.triboint.2019.106022. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.
Study on Dynamic Characteristics for High Speed Water-Lubricated Spiral Groove Thrust Bearing Considering Cavitating Effect Xiaohui Lin, Ruiqi Wang, Shaowen Zhang, Shuyun Jiang1 School of Mechanical Engineering Southeast University, Nanjing 211189, China
Abstract: The purpose of this paper was to study dynamic characteristics for high-speed water-lubricated spiral groove thrust bearing considering cavitating effect. A lubrication model was established with multiple gas-liquid interface effects, the evolution of bubble volume caused by breakage and coalescence of bubbles was quantitatively
described
using
the
population
balance
equation.
Dynamic
characteristics of the spiral groove thrust bearing (SGTB) were predicted. An experimental rig was developed to measure the dynamic coefficients of the SGTB. The result shows that, when spiral angle is in small range, direct stiffness coefficients of the SGTB with the cavitating effect are larger than those without cavitating effect. The cavitating effect has less influence on the damping coefficients than the stiffness coefficients.
Keywords: spiral groove thrust bearing; water lubrication; cavitating effect; dynamic characteristics
1. Introduction Currently, water has been used as a lubricant of fluid bearing. Comparing with oil-lubricated bearing, the water-lubricated bearing has a lower temperature rise, a higher limiting speed as well as no pollution. Comparing with conventional gas-lubricated bearing, the water-lubricated bearing has a larger stiffness, a larger load-carrying capacity. Spiral groove thrust bearings (SGTB) have been widely used 1
Corresponding author. Tel.: +86 25 52090533; fax: +86 25 52090504.
E-mail address: [email protected] (S. Jiang).
1
in high speed mechanical spindle support systems due to its high load carrying capacity[1-3]. Therefore, it is of great practical significance and engineering value to develop a water-lubricated SGTB. The cavitation will occur when the water-lubricated bearing runs at a high speed. In addition, the stiffness and damping coefficients are of main dynamic parameters for the bearing. Thus, it is necessary to study the dynamic characteristics of high-speed water-lubricated SGTB with consideration of cavitating effect. In the past few decades, there has been an increasing amount of literature on the dynamic characteristics of the gas-lubricated SGTB, for example, Malanoski and Pan [4] firstly obtained the axial stiffness and damping coefficients for gas-lubricated SGTB by using the narrow groove theory. Lund [5] firstly put forward the classic perturbation method to determine the dynamic coefficients of gas bearing. Zhou [6] obtained the axial stiffness and damping coefficients for spiral grooved self-acting gas bearings. Ren Liu etal. [7] studied the effects of gas-rarefaction on dynamic characteristics of gas lubricated SGTB. Meanwhile, the existing literature on the water-lubricated SGTB is extensive and focuses particularly on its static characteristics, for example: Yoshimoto et al. [8,9] investigated the static and dynamic characteristics of water-lubricated conical SGTB. Cabrera et al.[10] investigated the static characteristics of water-lubricated rubber journal bearings. Makoto Gohara et al. [11] investigated the static characteristics of a water-lubricated hydrostatic thrust bearing with a membrane restrictor. Liu et al.[12] and Durazo-Cardenas et al. [13] investigated water-lubricated hydrostatic journal bearings with a porous restrictor. It should be mentioned that the water-lubricated models proposed in the above literature are mostly based on single-phase flow. However, gas-liquid interface effect exists in the cavitating flow of water-lubricated bearings under high speed, and the interfacial effects of cavitating flow cannot be quantitatively described by single-phase flow model. In addition, the thermal and inertial effects can't be ignored when at high speed. In this paper, a complete water-lubricated model for water-lubricated SGTB considering cavitation effects was developed based on two-phase flow theory. 2
Multiple interface effects were considered and the influence of the gas nuclei in the water liquid on the size distribution of the cavitation bubble was also considered in the model [14,15]. The parameter weighted average method and Newton-SOR iterative method are used to solve the population balance equation of bubble. The stiffness and damping coefficients of the water-lubricated SGTB with three degrees of freedom are calculated using the small perturbation method. An experimental setup is developed to measure the dynamic characteristics of the SGTB. The influence of cavitation on the dynamic characteristics of the water-lubricated SGTB are analyzed.
2. Model governing equation 2.1 Pressure control equation of cavitating flow
Fig. 1
Schematic of SGTB
Fig.1 shows a schematic of the inward pumping spiral groove thrust bearing, its spiral groove line is a logarithmic spiral, which is expressed as
r =rbeθ cot β
(1)
For the water-lubricated SGTB in this study, Reynolds number is lower than 5000, the water flow is in laminar flow state, and the bubbles will generate under high speed due to high inception cavitation number of water, so the lubricant state is two-phase flow, the momentum transfer occurs on the interface of two phases. In addition, the inertial effect is even more highlighted under the high speed condition. Thus, the following assumptions are given for cavitating flow lubrication: 3
(1) The cavitating flow is regarded as two-phase flow, and the bubble is as discrete phase, and the shape of bubble is approximate to sphere. (2) Since the circumferential velocity of bearing is greater than the radial velocity of bearing, the momentum transfer between interface of two phases in radial direction is ignored. （3）The inertial effect in the r direction is only considered. (4) The water is regarded as laminar fluid. Based on two-phase flow theory [16] and the above assumptions, the N-S equations of the water lubricant in the polar coordinate system is expressed as
Cw
∂ ∂u 1 ∂ ( Cw p ) − M (θ , z ) µ = ∂z ∂z r ∂θ
（1）
Cw
∂ ∂v ∂ Cw ρ w u 2 = C p − µ ( w ) ∂z ∂z ∂r r
（2）
where, the interface momentum transfer effect is composed of three parts:
M (θ , z)=M1 (θ , z )+M 2 (θ , z )+M 3 (θ )
（3）
where, M 1 (θ , z ) represents the item of mass transfer，The expression of M 1 (θ , z ) is ∞
M 1 (θ , z )= − Cπ ρ wu (u − ub ) ∫ (V )
23
0
f eq (θ ,V )dV
（4 ）
M 2 (θ , z ) represents the item of viscous drag, The expression of M 2 (θ , z ) is ∞
M 2 (θ , z )=Cπ CDs ρ w (u − ub )2 ∫ (V ) 0
23
f eq (θ ,V )dV
（5 ）
M 3 (θ ) represents the item of surface tension, The expression of M 3 (θ ) is M 3 (θ )=
σ
∞
(V ) 2∫
13
0
f eq (θ ,V )dV
（6 ）
where
Cπ =
4π 3 3 4π
23
4
Cw =1 − C g
Cg =
∫
∞
0
Vf eq (θ , V )dV ∞
1 + ∫ Vf eq (θ , V )dV
（7 ）
0
The continuous equation of water lubricant in the polar coordinate system is expressed as
∫
h
0
h ∂ (C u ) ∂( rCw v ) w dz + ∫ dz = 0 0 ∂r ∂θ
（8 ）
The generalized Reynolds equation for cavitation flow lubrication can be derived from equations (1); (2); (8) in polar coordinates. (The derivation process is given in Appendix A).
∂ ∂ h3 ∂ ( Cw p ) + r ∂r 12µ ∂r ∂θ
h3 1 ∂ ( Cw p ) = 12µ r ∂θ
hU ∂ ∂ h 3 ρ wu 2 ∂ h 3 M (θ ) ( Cw ) + rCw + 2 ∂θ ∂r 12µ r ∂θ 12µ
(9)
The spiral groove is converted into a regular fan shape by the following coordinate transformation(as shown in Fig. 2).
ζ =r η = θ − ψ (ζ )
ψ (ζ )=
1 ζ ln( ) cot β ζ b
5
Fig. 2 Coordinates transformation of spiral curve
The generalized Reynolds equation (9) in the (ζ − η ) coordinates system is as follows
h 3 ∂Cw p ∂ h3 ∂Cw p ∂ 2 2 + ζ (ζ ψ ′ (ζ ) + 1) ∂ζ µ ∂ζ ∂η µ ζ∂η ∂ h 3 ∂Cw p ∂ h 3 ∂Cw p ′ ′ ζ ψ ζ ζψ ζ − ( ) − ( ) ∂ζ µ ∂η µ ∂ζ ∂η ∂ h3 ρ wu 2 ∂ ∂ h3 = ζψ ′(ζ ) ( 6CwUh ) + M (η ) − Cw µ ζ ∂η ∂η µ ∂η +
∂ ∂ζ
(10)
h3 2 Cw ρ w u µ
where
ψ ′(ζ )=
1 ζ cot β
The circumferential velocity of the fluid is
u=
Uz 1 1 ∂ + ( Cw p ) − M (η ) ( z 2 − hz ) h 2 µCw ζ ∂η
(11)
Boundary conditions for Eq. (10) are
2π p( η ) = p η + N p ζ =r = 0
(12)
out
p ζ = r = pin in
Because film thickness of the groove–ridge boundary is discontinuous, the 6
generalized Reynolds equation (10) is not satisfied at the groove–ridge boundary, the pressure of the groove–ridge boundary should be calculated using finite volume flow balance condition.
Fig.3 Finite volume separated by boundary
The finite volume including the groove-ridge boundary is established, as shown in Fig. 3. Finite volume flow balance condition is represented by the equation (13)
QAB + QBF + QFC + QCD + QDE + QEA = 0
(13)
where Qη = ∫
ζ2
ζ1
γ qη d ζ
η2
Qζ = ∫ qζ ζ dη
(η =AB,CD )
(14)
(ζ =BF ,FC ,DE,EA )
(15)
η1
γ = 1+
1 cot 2 β
Uη h h3 1 ∂ − ( Cw p ) − M (η ) 12 µ C w ζ ∂η 1 2 qη = 3 2 γ + 1 h 1 ∂Cw p − ζψ ′(ζ ) ∂C w p − ρ wu cot β 12 µ C ∂ζ ζ ∂η ζ w
qζ = −
h3 1 ∂Cw p ∂C p ρ u 2 − ζψ ′(ζ ) w − w 12 µ Cw ∂ζ ζ∂η ζ 7
(16)
(17)
2.2 The temperature control equation of cavitating flows The temperature field of cavitating flows was calculated using energy equation, the following assumptions were given :
（1）Considering convection heat transfer was much greater than heat conduction under the high speed condition, the heat conduction is neglected in the energy equation.
（2）Considering the water film is extremely thin, the variation of temperature or specific heat at constant volume of water in the direction of film thickness is neglected.
（3）Considering the radial velocity is much smaller than the circumferential velocity, the convective heat transfer in the radial direction is neglected.
（4）The viscous dissipation of the film thickness is only considered in the energy equation. According the above assumption (2), a average temperature of the water in the direction of film is defined as h
Tm = ∫ c v Tvθ dz / c v qη 0
Hence, the energy equation of cavitating flows in the (ζ − η ) coordinates system is simplified as h h ∂ ∂u ∂v ρ w c v qη ( CwTm ) = ∫0 Cw µ dz + ∫0 Cw µ dz ζ∂η ∂z ∂z 2
2
(18)
+ Ev (η ) (u − ub ) uh +Eσ (η )σ uh 2
where, the left-hand side of equation (18) represents the convective heat transfer term, the first and second items on the right-hand side of equation (18) represent the heat dissipation term, the third and fourth items on the right-hand side of equation (18) indicate the energy flowing into the water due to the gas-liquid interface effect. where ∂u U 1 1 ∂ ( Cw p ) = + − M (η ) ( 2 z − h ) ∂z h 2 µ C w ζ ∂η 8
(19)
∂v 1 ∂ ( Cw p ) ∂ ( Cw p ) ρ wu 2 ′ = − − ( ) ψ ζ (2z − h) ∂z 2 µ C w ∂ζ C w ∂η r
(20)
23 ∞ 3 23 Ev (η ) = 4π CDs ρ w ∫ V f ( η , V ) dV ( ) eq 0 4π 13 ∞ 3 13 Eσ (η ) = 4π ∫ V ) f eq (η ,V )dV ( 0 4π
The boundary conditions for Eq. (18) are
Tm ζ = r = T0 in
∂Tm ∂ζ
=0 ζ = rout
（21）
∂Tm ∂T = m =0 ∂η η =0 ∂η η =2π The viscosity–temperature relationship of fluid is described by an exponential function (22)
µ = µ0e−α
T
(Tm −T0 )
(22)
2.3 The bubble velocity control equation The velocity of the bubble in the cavitation flow is described by the bubble motion equation [17]
d ( ρgVub ) = ∑ Fi dt i
(23)
The forces acting on the bubble mainly include : (1) The force of relative acceleration between bubble and liquid
F1 = − V
∂p ∂x
(2) The viscous drag between the bubble and the liquid due to relative motion
1 F1 = ρ wπ R 2CD u − ub (u − ub ) 2 The inertial force can be ignored because the mass of the bubble is small, so the equation (23) is simplified as
9
1 ∂p ρ wπ R 2CD u − ub (u − ub ) − V =0 2 r∂ θ
(24)
The bubble velocity can be obtained by solving equation (24)
2.5 The population balance equation of bubble. The bubble volume of cavitating flow changes randomly due to breakup and coalescence effect. The breakup effect causes large volume bubbles to break down into small volume bubbles, and the coalescence effect causes small volume bubbles to be grown up into large volume bubbles. The random state of bubble volume can be described by using a bubble number density function f (η ,V ) , the bubble expansion due to change of gas density and mass transfer is ignored, the evolution of the bubble number density function f (η ,V ) is described by the following population balance equation of bubble: ∞ 6 ∂f ∂f 1 + ub =∫ h (ν ,V ) b (ν ) f (ν )dν ∂t ζ∂η V π 3ν 2 3 b 1 V − b (V ) f (V ) + ∫ c(ν ,V − ν ) f (ν ) f (V − ν )dν 2 0 13
(25)
∞
− f (V ) ∫ c(ν ,V ) f (ν )dν + S 0
where, the first term on the right-hand side of Eq. (25) represents breakup source term of bubbles, the second term on the right-hand side of Eq. (25) represents breakup sink term of bubbles, the third term on the right-hand side of Eq. (25) represents coalescence source term of bubbles, the fourth term on the right-hand side of Eq. (25) represents coalescence sink term of bubble, the fifth term on the right-hand side of Eq. (13) represents gas nucleus (microbubble) source term. The gas nucleus source term represents bubbles with volume V formed due to expansion of gas nucleus in the water under low pressure state. The expression of the source term S was
S = Cρ
∂ ( χ (V ) ) ∂V
( pv − p )
where
3 2 3 2 1 2 Cρ = 4π 4π 3ρ w 10
χ (V ) = ∫ (ϕ( ν ,V ) f i (ν )ν 2 3 )dν V
vl
where the ϕ (ν ,V ) was the redistribution function of bubble volume V ; the redistribution function of bubble volume V was given as [18] 2 V −ν ϕ( ν ,V ) = 2 exp −c2 V ν max
c1
The gas nuclei distribution is given as [19]
f i (ν ) = N i
κ α −1 −(κ ν ) −α e ν Γ (α − 1)
where, c1 =0.1 , c2 =0.5 , κ =6 , α =4 The breakage frequency is given as [20] b( ν )=
κ 1ε 1 3 κ 2σ exp − 29 5 9 ρ ε 2 3 ( 6ν π ) ( 6ν π ) w
The daughter size redistribution function is given as[20]
hb ( v,V ) =
2.4π ( 6V π )
23
ν
(2V − ν )2 exp −4.5 ν2
where the kernel function of coalescence is defined as respectively [21]
c(ν ,V ) = 0.0906
2 3
ν +V 2
2 3
ub − u
1 0.6 − C g
The initial condition for Eq. (25) is
f
t =0
=0
(26)
The boundary condition for Eq. (25)
f f
V =0
=0
V =∞
=0
are
(27)
2.6 Perturbation of Generalized Reynolds Equation of Cavitating Flows In practice, the SGTB has an axial displacement and two angular displacements, so the bearing should be modeled under three degrees of freedom. Supposing the origin of Cartesian coordinate is fixed in the steady equilibrium position of the bearing 11
surface, as shown in Fig.1, the independent motions of the bearing are the axial movement along z-axis, the angular movement around x-axis and the angular movement around y-axis respectively, where ( ∆z, ∆ϕ x , ∆ϕ y ) represent perturbations quantity of the thrust plate under three degrees of freedom. Assuming the initial position of the thrust plate is ( hR , ϕ x 0 , ϕ y 0 ) , the film thickness and pressure in quasi-equilibrium can be expressed as
h = h0 + ∆z + r cos θ∆ϕ y − r sin θ∆ϕx
(
where, h0 = hR + r −ϕ x 0 sinθ + ϕ y 0 cos θ
(28)
)
p = p0 + ∑ pq ∆q + ∑ pq& ∆q& q&
q
( q = z, ϕ , ϕ ) x
y
(29)
Substitute equations (28) and (29) into equation (10), perturbed generalized Reynolds equations of cavitating flows can be obtained
Rey ( pq ) = Fq (η , ζ ) where, the operators Rey (
(q = z, ϕ x , ϕ y , z&, ϕ&x , ϕ& y )
(30)
) was defined as
ζ h03 ∂ ( ) 12 µ ∂ζ ∂ h03 ∂ ( ) ∂ ζ h03 ∂ ( ) ′ ′ −ζψ (ζ ) − ψ (ζ ) ∂ζ 12 µ ∂η ∂η 12 µ ∂ζ Rey (
)=
∂( ) ∂ ∂ h03 ζ 2ψ ′2 (ζ ) + 1) ( + ∂η 12 µζ ∂η ∂ζ
The expressions of Fq (η , ζ ) were defined as follows ∂ CwU η ∂ h02 M (η ) + ∂η 2 ∂η 4µ ∂ ∂ h02 h02 2 + ′ C u − ( ) C ρ ψ ζ ρ wuη2 w w w η ∂η 4µ ∂ζ 4µ Fq (η , ζ ) = 2 2 − ∂ h0 ζ 2ψ ′2 (ζ ) + 1 ∂Cw p0 + ∂ h0 ζ ∂Cw p0 ) ∂η ∂ζ 4µ ∂ζ ∂η 4µζ ( 2 2 ∂ h0 ∂Cw p0 ∂ h0 ∂Cw p0 +ζψ ′(ζ ) + ψ ′(ζ ) ζ ∂ 4 ∂ ∂ ∂ζ ζ µ η η 4µ
12
q=z
∂ h02 ∂C p sin (η ′ ) (ζ 2ψ ′2 (ζ ) + 1) w 0 ∂η ∂η 4µ ∂ h2 ∂Cw p0 2 0 + ζ sin (η ′) ∂ζ ∂ζ 4µ 2 − 1 ∂ h0 ζ sin (η ′) ∂Cw p0 cot β ∂ζ 4µ ∂η 1 ∂ h02 2 ∂Cw p0 − ζ cot β ∂η 4µ ζ sin (η ′ ) ∂ζ Fq (η ,ζ ) = CU − ∂ w η ζ sin (η ′ ) ∂η 2 2 − ∂ h0 M (η ′ ) ζ sin η ′ ( ) ∂η 4µ 2 2 ∂ h0 ρ wuη − ζ sin (η ′) Cw 4µ ∂ζ 2 2 ∂ h0 ρ wuη + 1 ′) C sin ζ η ( w ζ cot β ∂η 4µ
∂ CwUη ζ cos (η ′ ) ∂η 2 ∂ h 2 M (η ′ ) 0 + ζ cos (η ′ ) ∂η 4µ 2 2 ∂ h0 ρ wuη ′ ζ η + C cos ( ) w ∂ζ 4µ 2 2 1 ∂ h0 ρ wuη − ζ cos (η ′) Cw 4µ ζ cot β ∂η Fq (η ,ζ ) = 2 ∂ ζ h0 1 ∂Cw p0 cos (η ′) + cot β ∂η ∂ζ 4µ 1 ∂ ζ h02 ∂C p + cos (η ′) w 0 ∂ζ cot β ∂η 4µ 2 − ∂ ζ h0 cos (η ′) (ζ 2ψ ′2 (ζ ) + 1) ∂Cw p0 ∂η 4µ ∂η ∂Cw p0 ∂ ζ 2 h02 − ∂ζ 4µ cos (η ′ ) ∂ζ
13
q = ϕx
q = ϕy
Fq (η , ζ ) =Cwζ
q = z&
Fq (η , ζ ) = − Cwζ 2 sin (η ′ )
q = ϕ& x
Fq (η , ζ ) =Cwζ 2 cos (η ′)
q = ϕ& y
where
η ′=η +ψ (ζ ) The boundary conditions for Eq. (30) are
pq pq pq
=0
ζ =rin
ζ =rout
η =0
(q = z, ϕ x , ϕ y , z&, ϕ&x , ϕ& y )
=0
= pq
(31)
η =2π
2.7 Calculation of stiffness coefficients and damping coefficients The dynamic characteristic parameters of the water-lubricated SGTB includes the stiffness coefficients and damping coefficients, which are calculated using the perturbation method. Once the perturbation pressure are determined by solving equation (30), the stiffness coefficient and damping coefficient can be obtained by integrating the perturbed pressure. For the SGTB with three degrees of freedom, totally nine stiffness coefficients and nine damping coefficients should be calculated using equations (32) and (33). K zz K zϕ K zϕ x y Kϕ x z Kϕ xϕ x Kϕ xϕ y Kϕ y z Kϕ yϕ x Kϕ yϕ y
−1 2π rout = ∫ ∫ ζ sin (η + f (ζ ) ) pz pϕ x pϕ y ζ d ζ d η 0 rin −ζ cos (η + f (ζ ) )
(32)
C zz Czϕ C zϕ x y Cϕ x z Cϕ xϕ x Cϕ xϕ y Cϕ y z Cϕ yϕ x Cϕ yϕ y
−1 2π rout = ∫ ∫ ζ sin (η + f (ζ ) ) pz& pϕ&x pϕ& y ζ d ζ d η 0 rin −ζ cos (η + f ( ζ ) )
(33)
{
{
}
}
The dimensionless stiffness coefficients and dimensionless damping coefficients are defined as follows
14
K zz∗ K z∗ϕ K z∗ϕ x y Kϕ∗ z Kϕ∗ ϕ Kϕ∗ ϕ x x x y x ∗ ∗ ∗ Kϕ z Kϕ ϕ Kϕ ϕ y x y y y
Czz∗ Cz∗ϕ Cz∗ϕ x y ∗ ∗ ∗ Cϕ z Cϕ ϕ Cϕ ϕ x x x y x ∗ ∗ ∗ Cϕ z Cϕ ϕ Cϕ ϕ y x y y y
K zz h3 K = 3 g ϕx z routUθ µ0 rout Kϕ y z r out Czz h3 C = 4 g ϕx z rout µ0 rout Cϕ y z r out
K zϕx
K zϕ y
rout
rout
Kϕxϕx Kϕxϕ y 2 rout
2 rout
Kϕ yϕx Kϕ yϕ y 2 rout
2 rout
Czϕx
Czϕ y
rout
rout
Cϕ xϕx Cϕxϕ y 2 rout
2 rout
Cϕ yϕx Cϕ yϕ y 2 rout
2 rout
3. Numerical calculation of the model The equations (10), (18) and equation (30) were discretized using the finite difference method, and the corresponding linear equations were obtained. The pressure distribution p and pq , the fluid velocity u and bubble velocity ub were calculated using the SOR iterative method, and the temperature field Tm were calculated using the stepping approach. The population balance equation (25) is an integral-differential equation of parabolic type. Since the coalescence source and sink terms in the population balance equation are nonlinear, the discretization equation of population balance equation is a system of nonlinear equations. Defining the solution domain node index
q = i + ( j − 2) Ni with 1 ≤ i ≤ Ni and 2 ≤ j ≤ N j − 1 . At time t , a finite difference
(
)
discrete operator L% ( f )q at q node is defined as
L%
(( f ) ) = ( K1)
N j −1
t
q
t
q
f + ( K 2 )q f − ∑ ( K 3)q f qt + ( K 4 ) q f qt t q2
t q1
N j −1
j −1
−∑ ( K 5 ) q f f + f e =1
e
t qe
t ql
t q
∑ ( K 6) e =1
qe
e
e= j
f − Sq = 0 t qe
e
( q = 1,...N )
(35)
q
Applying time-differencing through the λ weighted average method, the discretization equation for the equation (25) at q node can be written as 15
( f )q
t +1
−( f
)q +∆t ( (1 − λ ) L% ( ( f )q ) +λ L% ( ( f )q t
t +1
t
)) =0
(36)
where 0 ≤ λ ≤ 1 , the equation (36) can be rewritten as Fq
(( f ) ) = ( f ) t +1
t +1
q
q
+∆t λ L%
(( f ) ) +∆t (1 − λ ) L% (( f ) ) − ( f ) t +1
t
t
q
q
q
=0
( q = 1,...N ) q
(37) Equation (37) is a system of nonlinear equations for
( f )q
t +1
, which is solved
numerically by the Newton-SOR iterative method. In order to improve convergence, the downhill condition is included in the iterative procedure, the iterative procedure is written as
(
t +1( k )
Fq′ f q
t +1( k +1)
fq
(
) ∆f = − F ( f q
t +1( k )
= fq
)
q
t +1( k ) q
)
+ ωk ∆f
(
Fq f qt +1( k +1) ≤ F f qt +1( k )
(38)
)
1] . where the downhill factor ωk ∈ [ 0， The nine stiffness coefficients and nine damping coefficients were calculated by integrating (32) and (33) numerically. Fig. 4 shows the flow chart of the whole algorithm.
16
Start
Assign calculation parameters, the space coordinate
t =1 Calculate the interface items M (η )
Solve the generalized Reynolds equation ( 4 ) with SOR
Judge convergence of pressure at t vη = α vη ( k ) + (1 − α )vη (k −1)
vb = β vb ( k ) + (1 − β )vb ( k −1)
Solve equation ( 5) to obtain new fluid velocity
Judge convergence of fluid velocity at t
t = t +1
Solve the energy equation ( 9 ) with step method
Judge convergence of temperature at t Solve equation (12 ) to obtain new bubble velocity
Judge convergence of bubble velocity at t Solve the transport equation of bubble volume (13 ) with N-SOR
Judge convergence of bubble disturbution at t Output pressure, velocity, viscosity p, vη , µ and other parameters and bubble balance distribution f eq (η , V )
Fig.4 Flow diagram of solution algorithm
17
Small disturbance initial parameters ( ∆z, ∆ϕx , ∆ϕ y )
4. Experimental setup and procedure 4.1 Experimental setup An experimental setup is developed to measure the direct stiffness coefficient of the high-speed water-lubricated SGTB. Fig. 5 presents the schematic view and real scene image of the experiment setup. The grooved thrust disc is installed on the end of the high-speed motorized spindle, the rotational speed is controlled with an inverter; the thrust pad is fixed in the bearing housing; the external load is applied to the bearing housing by using a loading unit, and the external load is applied by adding a dead weight. The thrust pad displacement due to water film formation is detected by a non-contact displacement probe using eddy current; signals collected from the sensors are fed into a data acquisition system and interfaced with a computer. The water lubricant with ordinary pressure is supplied to the test bearing by a water circulation unit.
(a) Schematic view
18
(b) Real scene image Fig. 5 Experimental setup to measure the stiffness coefficient of SGTB 5.2 Specimens Fig. 6 shows the specimens of inwardly pumping SGTB. The material of the test thrust discs is stainless steel, the disc surfaces are machined using a grinding processing, and then the spiral grooves are etched on the surfaces using laser processing. The material of the test thrust pad is plexiglass, and a hole is fabricated in the center of pad, which is acted as outlet of water lubricant for the SGTB.
Fig. 6 Experimental specimens of SGTB 5.3 Measurement of stiffness coefficients of SGTB (1) Film thickness The motorized spindle is started to drive the disc up to a given rotational speed; an 19
axial load is applied on the test bearing. The film thickness is obtained by detecting the axial displacement of pad with the eddy current probe. Here, the film thickness refers to the film thickness in ridge region, namely, the nominal film thickness. The data is measured three times for each condition. (2) Axial stiffness coefficient The motorized spindle is started to drive the dis,c up to a given rotational speed; an initial axial load is applied on the test bearing, and a small load (about 5 percent of the applied load) is added on the bearing, then the variation of film thickness due to this small increment of load is measured. Thus, the axial stiffness coefficient of the bearing at this condition is obtained with dividing the small load by the variation of film thickness.
5. Results and Discussion 5.1 Validation of the present model A comparison of the predicting stiffness coefficients and the experimental ones is performed to verify the proposed model. Considering that it is difficult to test the angular stiffness coefficients of the SGTB, only the axial stiffness coefficient Kzz is measured in this study. The specifications of the spiral groove thrust bearing in this test are listed in Table 1. Table 1 The parameters of the spiral groove thrust bearing
Item
Value
Inner radius rin (mm)
7.5
Outer radius rout (mm)
20
Base circle radius rb (mm)
12
Water viscosity µ ( Pa.s )
0.001
Groove number N
12
Groove depth hg ( µm)
40
Spiral angles β ( o )
20
20
Pressure of supply water pin ( MPa)
0.1
Specific heat at constant volume of water cv ( J kg .Co )
4200
Density of water ρw ( kg m3 )
1000
Bearing capacity W ( N )
20
Surface tension of bubble σ ( N m )
0.3
Groove-to-land ratio λ b
0.5
Groove-to-dam ratio λl
0.6
Fig.7 gives a comparison between the predicting axial stiffness coefficient Kzz∗ and the experimental ones under three rotational speeds. The results show the predicting values generally agree with the experimental values, and the results with the cavitating flow are closer to the experimental values than those with non-cavitating flow, which indicates that the proposed model can be used to analyze the dynamic coefficients for high speed water-lubricated spiral groove thrust bearing considering cavitating effect.
(a) 9000 r/min
21
(b) 12000 r/min
(c) 15000 r/min
Fig.7 The predicting stiffness coefficients vs. the experimental ones It can be seen that the predicting values of axial stiffness coefficient are always lower than the experimental values. A possible explanation for this might be that the centrifugal effect on the water lubricant causes an error for the experimental values. In this study, the water lubricant is supplied by a pump and injected to the bearing clearence from the outside of bearing through three nozzles. The inwardly pumping bearing is employed, so the water is pumped into the bearing from the outside, and then flows through the bearing, finally reach the inside of the bearing. The water in the region near the inlet is disturbed by the centrifugal effect, and part of the water is thrown out from the thrust pad rim. Therefore, the actual flowrates of the bearing are less than the predicting values, which means that the water film thickness is less than the theoratical one. Therefore, it should be resonable that the the predicting values of 22
axial stiffness coefficient are always lower than the experimental values. In addition, the centrifugal effect on the water will increase with the increasing of rotary speed, as a result, the difference value between the predicting value and the experimental one gets larger with the inceasing of rotary speed.
5.2 The effect of cavitation on the dynamic characteristics of SGTB 1） ）Spiral angle Fig.8 shows the effect of cavitation on the stiffness coefficients of the SGTB with different spiral angles. It can be seen that the cavitation affects the stiffness coefficient Kzz∗ obviously; and the direct stiffness coefficients Kzz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ with x x
y y
cavitating effect are larger than those without cavitating effect when the spiral angles is less than 45 degrees. This result may be explained by the fact that the hydrodynamic effect is enhanced due to the interfacial effect of the cavitating flow, and the load-carrying capacity of the SGTB with cavitating flow is higher than that with non-cavitating flow, which means that a larger force is required to produce a unit displacement for the SGTB with the cavitation flow, so the stiffness coefficients Kzz∗ , Kϕ∗xϕ x , Kϕ∗yϕ y increase. The deviation of the direct stiffness coefficients between the
cavitating flow and the non-cavitating flow decreases with the increasing of spiral angles, which means the cavitation effect becomes weak as the spiral angle increases. In addition, no matter whether of the cavitating flow or the non-cavitating flow, the direct stiffness coefficients Kzz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the cross-coupled stiffness coefficients x x
y y
Kϕ∗yϕx , Kϕ∗xϕ y decrease with the increasing of the spiral angles. The cross-coupled stiffness
coefficients Kϕ∗ ϕ , Kϕ∗ ϕ are of opposite sign ( Kϕ∗ ϕ = − Kϕ∗ ϕ ). The x y
y x
x
y
y x
cross-coupled
stiffness coefficients Kϕ∗ z , K z∗ϕ , K z∗ϕ , Kϕ∗ z are approximately equal to zero for x
x
y
y
cavitating flow and non-cavitating flow. The cavitation effect hardly affect the cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ . x y
y x
23
2.0
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
cavitating flow non-cavitating flow
2.5
K*zz
1.5 1.0 0.5
K*zϕϕ
K zϕ y *
x
0.0 0
10
20
30
40
50
60
70
cavitating flow non-cavitating flow
0.7 0.5
K*ϕxϕx
0.3 0.1 -0.1
K* ϕ x z
-0.3 -0.5
K*ϕxϕy
-0.7 -0.9 -1.1 0
10
Spiral angle β (°)
30
40
50
60
70
(b)
Dimensionless stiffness coefficients
(a)
20
Spiral angle β (°)
cavitating flow non-cavitating flow
1.2 1.0
K*ϕyϕx 0.8 0.6 0.4
K*ϕyϕy
0.2
K*ϕyz
0.0 0
10
20
30
40
50
Spiral angle β (°)
60
70
(c) Fig.8 Effect of cavitation on the stiffness coefficient of SGTB with different spiral angles
( ω = 15000rpm , hg = 40 µ m , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 ) Fig. 9 shows the effect of cavitation on the damping coefficients of the SGTB with different spiral angles. It can be seen that cavitation effect hardly affect the nine damping coefficients, this result may be explained by the fact that the cavitation effect has a great influence on the pressure p0 of equilibrium position, while the perturbed pressure components ( p z& , pϕ& , pϕ& ) associating with damping coefficient are only x
y
related to the volume fraction Cg of bubbles per unit volume of water rather than the pressure p0 . In other word, the damping coefficients are independent of the pressure p0 at equilibrium position, furthermore, the volume fraction Cg of bubbles per unit volume of water is usually small. Therefore, the effect of cavitation on the 24
damping coefficient is very weak. No matter whether the cavitating flow or the non-cavitating flow, the direct damping coefficients C zz∗ , Cϕ∗ ϕ , Cϕ∗ ϕ decrease gradually with the increasing of the x x
y
y
spiral angles, the cross-coupled damping coefficients C z∗ϕ , Cz∗ϕ , Cϕ∗ z Cϕ∗ ϕ Cϕ∗ ϕ Cϕ∗ z are x
y
x
x
y
y x
y
cavitating flow non-cavitating flow
0.12
Dimensionless damping coefficients
Dimensionless damping coefficients
approximately equal to zero.
0.08
C*zz 0.04
C*zϕϕx
C*zϕϕy
0.00 0
10
20
30
40
50
60
Spiral angle β (°)
70
cavitating flow non-cavitating flow
0.030 0.025 0.020
C*ϕxϕx
0.015 0.010 0.005
C*ϕxϕy
0.000
C*ϕxz
-0.005 0
10
30
40
50
60
70
(b)
Dimensionless damping coefficients
(a)
20
Spiral angle β (°)
cavitating flow non-cavitating flow
0.030 0.025 0.020
C*ϕyϕy
0.015 0.010 0.005
C*ϕyϕx
0.000
C*ϕyz
-0.005 0
10
20
30
40
50
Spiral angle β (°)
60
70
(c) Fig.9 Effect of cavitation on the damping coefficients of SGTB with different spiral angles
( ω = 15000rpm , hg = 40 µ m , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
2) Groove depth Fig. 10 shows the effect of cavitation on the stiffness coefficients of the SGTB with different groove depths. It can be seen that the direct stiffness coefficients Kzz∗ , Kϕ∗xϕx , Kϕ∗yϕ y with cavitating effect are larger than those without cavitating effect, the 25
deviation of the direct stiffness coefficients between the cavitating flow and non-cavitating flow increases with the increasing of groove depth. A possible explanation for this is that the hydrodynamic pressure decreases when the groove depth increases, leading to an increase in the number of bubbles per unit volume of water lubricant. The increasing of the number of bubbles will cause the hydrodynamic effect of interface to be enhanced. Furthermore, the hydrodynamic effect of interface due to the bubble number is more dominant than the hydrodynamic effect due to the groove depth. Therefore, as the groove depth increases, the direct stiffness coefficients with cavitating effect are larger than those without cavitating effect. The effect of cavitation on the cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ is very weak for x y
y x
different groove depths. In addition, no matter whether the cavitating flow or the non-cavitating flow, the direct stiffness coefficients K zz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the x x
y y
cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ decrease gradually with the increasing y x
y x
of groove depth, the cross-coupled stiffness coefficients K z∗ϕ , K z∗ϕ , Kϕ∗ z , Kϕ∗ z are x
y
x
y
cavitating flow non-cavitating flow
12
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
approximately equal to zero.
10 8
K*zz
6 4 2
K*zϕϕx
K*zϕϕy
0 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(a)
cavitating flow non-cavitating flow
4
K*ϕxϕx
2
K*ϕxz 0
K*ϕxϕy
-2 0
5
10
15
20
(b)
26
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless stiffness coefficients
cavitating flow non-cavitating flow
3.5 3.0
K*ϕyϕy
2.5 2.0 1.5
K*ϕyϕx
1.0 0.5
K*ϕyz
0.0 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(c) Fig.10 Effect of cavitation on stiffness coefficients of SGTB with different groove depths
( ω = 15000rpm , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
Fig. 11 shows the effect of cavitation on the damping coefficients of the SGTB with different groove depths. It can be also seen that the effect of cavitation hardly affect the damping coefficients, the cause for this result is identical with the section above. For the cavitating flow and the non-cavitating flow, the direct damping coefficients C zz∗ , Cϕ∗xϕ x , Cϕ∗yϕ y decrease nonlinearly with the increasing of groove depth, the
cross-coupled damping coefficients C z∗ϕ , Cz∗ϕ , Cϕ∗ z Cϕ∗ ϕ Cϕ∗ ϕ Cϕ∗ z are approximately x
y
x
x
y
y x
y
0.20
cavitating flow non-cavitating flow
0.15
C*zz 0.10
0.05
C*zϕϕx
C*zϕϕy
0.00 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless damping coefficients
Dimensionless damping coefficients
equal to zero.
(a)
cavitating flow non-cavitating flow
0.06 0.05 0.04
C* ϕ x ϕ x 0.03 0.02 0.01
C*ϕxϕy
C*ϕxz
0.00 0
5
10
15
20
(b)
27
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless damping coefficients
cavitating flow non-cavitating flow
0.06 0.05 0.04
C*ϕyϕy
0.03 0.02
C* ϕ y ϕ x
0.01
C* ϕ y z
0.00 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(c) Fig.11 Influence of cavitation effect on damping coefficients for different groove depth
( ω = 15000rpm , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
3) Rotational speed Fig. 12 shows the effect of cavitation on the stiffness coefficients of the SGTB with different rotational speeds. It can be seen that the direct stiffness coefficient K zz∗ with cavitating effect are larger than those without cavitating effect, and the direct stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ with cavitating effect are slightly larger than those without x x
y y
cavitating effect in the whole speed range. A possible explanation for this is that the the number of bubbles per unit volume of water lubricant increases with the increasing of rotary speed, so the hydrodynamic effect of interface is enhanced. Consequently, the load-carrying capacity of the SGTB with cavitating flow is larger than that with non-cavitating flow. Furthermore, no matter whether the cavitating flow or the non-cavitating flow, the direct stiffness coefficients K zz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the x x
y y
cross-coupled stiffness coefficient Kϕ∗ ϕ increase approximately linearly with the y x
increasing
of
rotational
speed,
and
the
cross-coupled
stiffness
coefficient Kϕ∗ ϕ decreases approximately linearly with the increasing of rotational x y
speed, the cross-coupled stiffness coefficients K z∗ϕ , K z∗ϕ , Kϕ∗ z , Kϕ∗ z are approximately x
equal to zero. 28
y
x
y
4.0
cavitating flow non-cavitating flow
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
4.5
3.5
K zz *
3.0 2.5 2.0 1.5 1.0 0.5
K*zϕϕy
K*zϕϕx
0.0 8000 10000 12000 14000 16000 18000 20000 22000
cavitating flow non-cavitating flow
1.5 1.0 0.5
K*ϕxϕx 0.0
K*ϕxz
-0.5 -1.0 -1.5
K*ϕxϕy -2.0 8000 10000 12000 14000 16000 18000 20000 22000
rpm
rpm
(b)
Dimensionless stiffness coefficients
(a)
2.0
cavitating flow non-cavitating flow K*ϕyϕx
1.5
1.0
0.5
K*ϕyϕy
K*ϕyz 0.0
8000 10000 12000 14000 16000 18000 20000 22000
rpm
(c) Fig.12 Effect of cavitation on stiffness coefficients of SGTB with different rotational speeds
( hg = 40 µ m , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 ) Fig. 13 shows the effect of cavitation on the damping coefficients of the SGTB with different rotational speeds. It can be also seen that the cavitation effect hardly affect the damping coefficients, the cause for this result is identical with the section above, which is not repeated here. For the cavitating flow and the non-cavitating flow, the direct damping coefficients Czz∗ , Cϕ∗ ϕ , Cϕ∗ ϕ and cross-coupled damping coefficients x x
Cϕ∗xϕ y , Cϕ∗yϕ x
y y
are approximately constant with the increasing of rotational speed, and
cross-coupled damping coefficients Cϕ∗ ϕ , Cϕ∗ ϕ are of opposite sign ( Cϕ∗ ϕ = −Cϕ∗ ϕ ). The x
y
y
x
x
y
y
x
cross-coupled damping coefficients C z∗ϕ , C z∗ϕ , Cϕ∗ z , Cϕ∗ z are approximately equal to zero. x
y
29
x
y
C zz *
0.08 0.07
cavitating flow non- cavitating flow
0.06 0.05
Dimensionless damping coefficients
Dimensionless damping coefficients
0.09
0.0004 0.0003 0.0002 0.0001 0.0000
-0.0001
0.04
C* ϕx ϕx
C*ϕyϕy
C*ϕxz
-0.0002
0.03
cavitating flow non-cavitating flow C*ϕyz
C*ϕyϕx
C*ϕxϕy
-0.0003
0.02 8000 10000 12000 14000 16000 18000 20000 22000
-0.0004 8000 10000 12000 14000 16000 18000 20000 22000
rpm
rpm
(b) Dimensionless damping coefficients
(a)
0.02
cavitating flow non-cavitating flow
0.01
0.00
C*zϕϕy
C*zϕϕx -0.01
-0.02 8000 10000 12000 14000 16000 18000 20000 22000
rpm
(c) Fig.13 Effect of cavitation on damping coefficients of SGTB with different rotational speeds
( hg = 40 µ m , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
6. Conclusions A theoretical model of cavitation flow lubrication considering various gas-liquid interface effects has been established. The dynamic coefficients of the high-speed water-lubricated hydrodynamic SGTB is calculated using the model. An experimental setup has been developed to verify the proposed model. Through numerical simulation and discussion, the following conclusions can be drawn: (1) The direct stiffness coefficients of SGTB with cavitating flow are larger than those of SGTB with non-cavitating flow when the spiral angles is less than 45 degrees, or the groove depth reaches a greater value. The effect of cavitation on the direct 30
stiffness coefficients is gradually weakened with the increasing of spiral angle; however, the effect of cavitation on the direct stiffness coefficients is gradually enhanced with the increasing of groove depth; (2) The influence of cavitation on the direct stiffness coefficients K zz∗ depends on the rotational speed. The direct stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ with cavitating effect are x x
y y
larger than those without cavitating effect slightly in the whole range of speed. (3) The cavitation hardly affect the cross-coupled stiffness coefficients and the damping coefficients of SGTB.
Acknowledgment The authors gratefully acknowledge the support of the National Science Foundation of China through Grant Nos. 51635004 and 11472078.
31
Nomenclature b ( v ) The bubble breakage frequency function （ ） 1 s
Cij (i = z, ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine damping coefficients（
N .s N .s N .m.s , , N .s, ） m rad rad
Cij∗ (i = z , ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine dimensionless damping coefficients
Cw the water volume fraction (dimensionless)
Cg the bubble volume fraction (dimensionless) Cπ Constant coefficient (dimensionless)
C Ds drag coefficient (dimensionless) f eq (θ , V ) the equilibrium distribution function of bubble volume（
1 m6
）
feq (η ,V ) the equilibrium distribution function of bubble volume at coordinates system（
1 m6
(ζ ,η )
）
f eq* the dimensionless equilibrium distribution function of bubble volume (dimensionless)
f (V ) the distribution of bubble volume (
1 ) m6
h the film thickness ( m )
hg the groove depth ( m )
hR steady film thickness at ridge area ( m ) h0 Static film thickness ( m ) hb ( v,V ) bubble volume redistribution function (
1 ) m
Ki , j (i = z, ϕ x , ϕ y , j = z, ϕ x , ϕ y ) nine stiffness coefficients（
N N N.m N.m ; ; ; ） m rad m rad
K ij∗ (i = z , ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine dimensionless stiffness coefficients
(dimensionless)
32
M (θ ) the interface item ( N3 ) m
M (η ) the interface item at (ζ ,η ) coordinates system ( N3 ) m
N the number of spiral groove (dimensionless)
Ni the number density of nuclei ( p the water film pressure (
p0 Static pressure (
1 ) m3
N ) m2
N ) m2
pq , pq& perturbed pressure components (
N N N .s N .s , , , ) m3 m2 m3 m2
pv the vaporization pressure of water (
N ) m2
pout the pressure at outer radius ( pin the pressure at inter radius (
N ) m2
N ) m2 3
QAB ,QBF ,QFC ,QCD ,QDE ,QEA the volume flow ( m ) s
3
Qζ the flow rate in the ζ direction ( m ) s
3
Qη the flow rate in the η direction ( m ) s
2
qζ the ζ direction flows of water in spiral coordinates ( m ) s
2
qη the η direction flows of water in spiral coordinates ( m ) s
r the radial coordinate in polar coordinates ( m )
rb the base circle radius of spiral line ( m ) rin thrust disc inner radius ( m )
rout
thrust disc outer radius ( m )
R the statistical average radius of bubbles ( m )
c (ν ,V ) the bubble coalescence kernel function ( m ) 3
s
S the source terms
t the time ( s )
33
T Water film temperature ( K )
Tm the average temperature of water film ( K ) T0 the ambient temperature ( C o ) U ( = routω ) the circumferential velocity of thrust disc ( ) m s
u the circumferential velocity of water ( m ) s
m s
u the circumferential average velocity of water ( )
v the circumferential velocity of water at (ζ ,η ) coordinates system ( ) m s
m s
ub the bubble velocity ( ) V the bubble volume ( m3 )
Vmax maximum volume of bubbles ( m3 ) z the dimensionless coordinate of the film thickness (dimensionless) Greek Letters
α T the temperature coefficient of viscosity (
1 ) K
cv the specific heat capacity of water ( J / ( kg ⋅ K ) )
σ the surface tension of bubble (
N ) m
µ the dynamic viscosity of water (
Ns ) m2
κ , κ 2 the adjustable parameters (dimensionless) 1
ρ g the density of bubble ( kg3 ) m
ρw the density of water ( kg ) m3 2
ε the eddy diffusivity ( m ) s
λb the groove-to-land ratio (dimensionless) λl the groove-to-dam ratio (dimensionless) θ the circumferential coordinate in polar coordinates (deg) 34
ν the bubble volume ( m3 )
β the spiral angle (deg) η the spiral coordinates (deg)
ζ the spiral coordinates ( m ) ∆z, ∆ϕx , ∆ϕ y perturbations quantity ( m ) ϕ x 0 , ϕ y 0 the initial position of the thrust plate ( m )
Appendix A. The derivation of the generalized Reynolds equation
Cw
∂ ∂u 1 ∂ ( Cw p ) − M (θ , z ) µ = ∂z ∂z r ∂θ
Cw
∂ ∂v ∂ Cw ρ w u 2 µ = C p − ( w ) ∂z ∂z ∂r r
（A1）
（A2）
Equations (A1) and (A2) can be integrated,
∂u 1 ∂ ( Cw p ) z − F1( z )+c1 Cw µ = ∂z r ∂θ
Cw µ
(A3)
∂vr ∂ C ρ = ( Cw p ) z − w w F2 ( z ) + d1 ∂z ∂r r
(A4)
Where z
F1 ( z ) = ∫ M (θ , z′)dz′ 0
z
F2 ( z )= ∫ u 2 dz ′ 0
Further integration of equations (A3) and (A4) yields
( Cw µ u ) =
z 1 ∂ z2 ( Cw p ) − ∫0 F1( z′)dz′+c1z + c2 (A5) r ∂θ 2
∂ z 2 Cw ρ w C w µ v = ( Cw p ) − ∂r 2 r
∫
z
0
F2 ( z )dz ′ + d1 z + d 2
Boundary conditions are
35
(A6)
u z =h = U , u z =0 = 0 v z = h = 0, v z =0 = 0
(A7)
Substitute boundary conditions (A7) into equations (A5) and (A6) , the flow field distribution of the cavitating flow are obtained
( Cw u ) = Cw v =
h CwU 1 ∂ 1 z z + ( Cw p ) ( z 2 − hz ) − ( ∫0 F1 ( z′)dz′ − ∫0 F1 ( z′)dz′ ) (A8) h 2r µ ∂θ µ h
z z 2 hz 1 C ρ 1 ∂ z h ′ ′ ′ ′ F（ ( Cw p ) − − w w ∫0 F（ 2 z ）dz − 2 z ）dz ∫ h 0 µ ∂r 2 2 µ r
(A9)
Since the water film is very thin, the velocity u in the interface term and the inertia term can be approximated by the average velocity u of the film thickness. The u is defined as 1 h u = ∫ udz h 0 and F1 ( z ) = M (θ ) z
∫
z
∫
h
∫
z
0
0
0
F1 ( z ′)dz ′=M (θ )
z2 2
F1 ( z ′)dz ′=M (θ )
h2 2
′ ′ 2 F（ 2 z ）dz =u
z2 2
Equations (A8) and (A9) are simplified to
( Cw u ) = Cw v =
CwU 1 ∂ 1 + ( Cw p ) ( z 2 − hz ) − M (θ ) ( z 2 − hz ) （A10） h 2r µ ∂θ 2µ
z 2 hz 1 Cw ρ wu 2 2 1 ∂ C p ( w ) − − ( z − hz ) (A11) µ ∂r 2 2 µ 2r
Where
36
M (θ ) = (Cπ CDs ρ w (u − ub )2 − Cπ ρ wu (u − ub ) ) ∫ (V ) ∞
f eq (θ ,V )dV
23
0
+
σ
∞
(V ) 2∫ 0
13
f eq (θ ,V )dV
The continuity equation of the two-phase flow aqueous phase in polar coordinates is:
∫
h
0
h ∂ (C u ) ∂( rCw v ) w dz + ∫ dz = 0 0 ∂r ∂θ
（A12）
Substitute equations (A10) and (A11) into equations(A12) , The generalized Reynolds equation (A13) of the cavitating flow is given ∂ ∂r
h3 ∂ ∂ h3 1 r C p ( ) + w 12 µ ∂ r ∂ θ 12 µ r hU θ ∂ h 3 ρ wu 2 ∂ C + rC ( w) w 2 ∂θ ∂r 12 µ r
∂ (C w p ) = ∂θ
（A13）
∂ h3 M (θ ) + ∂ θ 12 µ
References [1] Murata S. Exact Two-Dimensional Analysis of Circular Disk Spiral-Groove Bearing. ASME Journal of Tribology 1979;101(4):424-430. [2] Elrod HG. Some Refinements of the Theory of the viscous Screw Pump. ASME Journal of Tribology 1973;95(1):82-93. [3] Muijderman EA. Spiral Groove Bearings. New York: Springer-Verlag; 1966. [4] Manaloski SB, Pan CHT. The Static and Dynamic Characteristic of the Spiral-Grooved
Thrust
Bearing.
ASME
Journal
of
Basic
Engineering
1965;87:547–558. [5] Lund L. Calculation of Stiffness and Damping Properties of Gas Bearings. ASME Journal of Lubrication Technology 1968;90:793–803. [6] Zhou H. A High Order, Isoparametric, Finite Element Method for the Calculation of Static and Perturbed Frequency Characteristics of Spiral-Grooved Self Acting Gas Bearings. Proc. 8th International Gas Bearing Symposium, BHRA Fluid Engineering, Cranfield, England; 1981 pp. 229–236. [7] Liu R, Wang XL, Zhang XQ. Effects of gas rarefaction on dynamic characteristics of
micro
spiral-grooved
thrust
bearing. 37
Journal
of
Tribology
2012;
134(2):222011-222017. [8]Yoshimoto S, Anno Y, Tamura M, Kakiuchi Y, Kimura K. Axial load carrying capacity of water-lubricated hydrostatic conical bearings with spiral grooves (on the case of rigid surface bearings). Trans ASME 1965; 118(4):893–899. [9] Yoshimoto S, Oshima S, Danbara S, Shitara T. Stability of water-lubricated, hydrostatic, conical bearings with spiral grooves for high-speed spindles. Trans ASME 2002;124(2):398–405. [10] Cabrera DL, Woolley NH, Allanson DR, Tridimas YD. Film pressure distribution in water-lubricated rubber journal bearings. Proc. Institution Mech. Eng., Part J: J. Eng. Tribol. 2005;219(2):125–132. [11] Gohara M, Somaya K, Miyatake M, Yoshimoto S. Static characteristics of a water-lubricated hydrostatic thrust bearingusing a membrane restrictor. Tribology International 2014;75:111–116. [12] Feng L, Bin L, Xiaofeng Z. Numerical design method for water-lubricated hybrid sliding bearings. International Journal of Precision Engineering & Manufacturing 2008;9(1):47-50. [13] Durazo-Cardenas IS, Corbett J, Stephenson DJ. The performance of a porous ceramic hydrostatic journal bearing. Proceedings of the Institution of Mechanical Engineers, Part J: Journal of Engineering Tribology 2010;224(1):81-89. [14] Lin XH, Jiang SY, Zhang CB, Liu X. Thermohydrodynamic analysis of high speed water-lubricated spiral groove thrust bearing considering effects of cavitation, inertia and turbulence. Tribology International 2018;119:645-658. [15] Lin XH, Jiang SY, Zhang CB, Liu X. Thermohydrodynamic Analysis of High-Speed Water-Lubricated Spiral Groove Thrust Bearing Using Cavitating Flow Model. Journal of Tribology Trans ASME 2018;140.051703:1-12. [16] Ni J, Wang G, Zhang H. The Basic Theory of Solid Liquid Two-Phase Flow and Its Latest Applications. Science Press, Beijing, China; 1991 pp.35-52. [17] Guo L. Two-Phase and Multiphase Flow Dynamics. Xi’an Jiaotong University Press, Xi’an, China; 2002 pp. 430. [18] Bork Z, Jakobsen HA. Evaluation of breakage and coalescence kernels for 38
vertical bubbly flows using a combined multifluid-population balance model solved by least squares method. Procedia Engineering 2012;42:682-694. [19] Billet ML. Cavitation Nuelei Measurements,Proe.of ZndIntorn Symp. on Cavitation Ineeption, Dee, 1984. p.33-42. [20] Coulaloglou CA, Tavlarides LL. Description of interaction processes in agitated liquid-liquid dispersions. Chemical Engineering Science 1977;32(11):1289-1297. [21] Millies M, Mewes D. Interfacial Area Density in Bubbly Flow. Chemical Engineering and Processing: Process Intensification 1999;38(4-6): 307–319.
39
Highlights
1) The dynamic characteristics of high-speed water-lubricated SGTB considering cavitating effect is studied. 2) The hydrodynamic lubrication model of cavitating flow for the water-lubricated SGTB was established. 3) The cavitation effect on dynamic characteristics of water-lubricated SGTB is analyzed.
Declaration of Interest Statement
This is to declare that there is no interest conflict about this paper.
S0301-679X(19)30539-0
DOI:
https://doi.org/10.1016/j.triboint.2019.106022
Reference:
JTRI 106022
To appear in:
Tribology International
Received Date: 15 August 2019 Revised Date:
14 October 2019
Accepted Date: 15 October 2019
Please cite this article as: Lin X, Wang R, Zhang S, Jiang S, Study on dynamic characteristics for high speed water-lubricated spiral groove thrust bearing considering cavitating effect, Tribology International (2019), doi: https://doi.org/10.1016/j.triboint.2019.106022. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.
Study on Dynamic Characteristics for High Speed Water-Lubricated Spiral Groove Thrust Bearing Considering Cavitating Effect Xiaohui Lin, Ruiqi Wang, Shaowen Zhang, Shuyun Jiang1 School of Mechanical Engineering Southeast University, Nanjing 211189, China
Abstract: The purpose of this paper was to study dynamic characteristics for high-speed water-lubricated spiral groove thrust bearing considering cavitating effect. A lubrication model was established with multiple gas-liquid interface effects, the evolution of bubble volume caused by breakage and coalescence of bubbles was quantitatively
described
using
the
population
balance
equation.
Dynamic
characteristics of the spiral groove thrust bearing (SGTB) were predicted. An experimental rig was developed to measure the dynamic coefficients of the SGTB. The result shows that, when spiral angle is in small range, direct stiffness coefficients of the SGTB with the cavitating effect are larger than those without cavitating effect. The cavitating effect has less influence on the damping coefficients than the stiffness coefficients.
Keywords: spiral groove thrust bearing; water lubrication; cavitating effect; dynamic characteristics
1. Introduction Currently, water has been used as a lubricant of fluid bearing. Comparing with oil-lubricated bearing, the water-lubricated bearing has a lower temperature rise, a higher limiting speed as well as no pollution. Comparing with conventional gas-lubricated bearing, the water-lubricated bearing has a larger stiffness, a larger load-carrying capacity. Spiral groove thrust bearings (SGTB) have been widely used 1
Corresponding author. Tel.: +86 25 52090533; fax: +86 25 52090504.
E-mail address: [email protected] (S. Jiang).
1
in high speed mechanical spindle support systems due to its high load carrying capacity[1-3]. Therefore, it is of great practical significance and engineering value to develop a water-lubricated SGTB. The cavitation will occur when the water-lubricated bearing runs at a high speed. In addition, the stiffness and damping coefficients are of main dynamic parameters for the bearing. Thus, it is necessary to study the dynamic characteristics of high-speed water-lubricated SGTB with consideration of cavitating effect. In the past few decades, there has been an increasing amount of literature on the dynamic characteristics of the gas-lubricated SGTB, for example, Malanoski and Pan [4] firstly obtained the axial stiffness and damping coefficients for gas-lubricated SGTB by using the narrow groove theory. Lund [5] firstly put forward the classic perturbation method to determine the dynamic coefficients of gas bearing. Zhou [6] obtained the axial stiffness and damping coefficients for spiral grooved self-acting gas bearings. Ren Liu etal. [7] studied the effects of gas-rarefaction on dynamic characteristics of gas lubricated SGTB. Meanwhile, the existing literature on the water-lubricated SGTB is extensive and focuses particularly on its static characteristics, for example: Yoshimoto et al. [8,9] investigated the static and dynamic characteristics of water-lubricated conical SGTB. Cabrera et al.[10] investigated the static characteristics of water-lubricated rubber journal bearings. Makoto Gohara et al. [11] investigated the static characteristics of a water-lubricated hydrostatic thrust bearing with a membrane restrictor. Liu et al.[12] and Durazo-Cardenas et al. [13] investigated water-lubricated hydrostatic journal bearings with a porous restrictor. It should be mentioned that the water-lubricated models proposed in the above literature are mostly based on single-phase flow. However, gas-liquid interface effect exists in the cavitating flow of water-lubricated bearings under high speed, and the interfacial effects of cavitating flow cannot be quantitatively described by single-phase flow model. In addition, the thermal and inertial effects can't be ignored when at high speed. In this paper, a complete water-lubricated model for water-lubricated SGTB considering cavitation effects was developed based on two-phase flow theory. 2
Multiple interface effects were considered and the influence of the gas nuclei in the water liquid on the size distribution of the cavitation bubble was also considered in the model [14,15]. The parameter weighted average method and Newton-SOR iterative method are used to solve the population balance equation of bubble. The stiffness and damping coefficients of the water-lubricated SGTB with three degrees of freedom are calculated using the small perturbation method. An experimental setup is developed to measure the dynamic characteristics of the SGTB. The influence of cavitation on the dynamic characteristics of the water-lubricated SGTB are analyzed.
2. Model governing equation 2.1 Pressure control equation of cavitating flow
Fig. 1
Schematic of SGTB
Fig.1 shows a schematic of the inward pumping spiral groove thrust bearing, its spiral groove line is a logarithmic spiral, which is expressed as
r =rbeθ cot β
(1)
For the water-lubricated SGTB in this study, Reynolds number is lower than 5000, the water flow is in laminar flow state, and the bubbles will generate under high speed due to high inception cavitation number of water, so the lubricant state is two-phase flow, the momentum transfer occurs on the interface of two phases. In addition, the inertial effect is even more highlighted under the high speed condition. Thus, the following assumptions are given for cavitating flow lubrication: 3
(1) The cavitating flow is regarded as two-phase flow, and the bubble is as discrete phase, and the shape of bubble is approximate to sphere. (2) Since the circumferential velocity of bearing is greater than the radial velocity of bearing, the momentum transfer between interface of two phases in radial direction is ignored. （3）The inertial effect in the r direction is only considered. (4) The water is regarded as laminar fluid. Based on two-phase flow theory [16] and the above assumptions, the N-S equations of the water lubricant in the polar coordinate system is expressed as
Cw
∂ ∂u 1 ∂ ( Cw p ) − M (θ , z ) µ = ∂z ∂z r ∂θ
（1）
Cw
∂ ∂v ∂ Cw ρ w u 2 = C p − µ ( w ) ∂z ∂z ∂r r
（2）
where, the interface momentum transfer effect is composed of three parts:
M (θ , z)=M1 (θ , z )+M 2 (θ , z )+M 3 (θ )
（3）
where, M 1 (θ , z ) represents the item of mass transfer，The expression of M 1 (θ , z ) is ∞
M 1 (θ , z )= − Cπ ρ wu (u − ub ) ∫ (V )
23
0
f eq (θ ,V )dV
（4 ）
M 2 (θ , z ) represents the item of viscous drag, The expression of M 2 (θ , z ) is ∞
M 2 (θ , z )=Cπ CDs ρ w (u − ub )2 ∫ (V ) 0
23
f eq (θ ,V )dV
（5 ）
M 3 (θ ) represents the item of surface tension, The expression of M 3 (θ ) is M 3 (θ )=
σ
∞
(V ) 2∫
13
0
f eq (θ ,V )dV
（6 ）
where
Cπ =
4π 3 3 4π
23
4
Cw =1 − C g
Cg =
∫
∞
0
Vf eq (θ , V )dV ∞
1 + ∫ Vf eq (θ , V )dV
（7 ）
0
The continuous equation of water lubricant in the polar coordinate system is expressed as
∫
h
0
h ∂ (C u ) ∂( rCw v ) w dz + ∫ dz = 0 0 ∂r ∂θ
（8 ）
The generalized Reynolds equation for cavitation flow lubrication can be derived from equations (1); (2); (8) in polar coordinates. (The derivation process is given in Appendix A).
∂ ∂ h3 ∂ ( Cw p ) + r ∂r 12µ ∂r ∂θ
h3 1 ∂ ( Cw p ) = 12µ r ∂θ
hU ∂ ∂ h 3 ρ wu 2 ∂ h 3 M (θ ) ( Cw ) + rCw + 2 ∂θ ∂r 12µ r ∂θ 12µ
(9)
The spiral groove is converted into a regular fan shape by the following coordinate transformation(as shown in Fig. 2).
ζ =r η = θ − ψ (ζ )
ψ (ζ )=
1 ζ ln( ) cot β ζ b
5
Fig. 2 Coordinates transformation of spiral curve
The generalized Reynolds equation (9) in the (ζ − η ) coordinates system is as follows
h 3 ∂Cw p ∂ h3 ∂Cw p ∂ 2 2 + ζ (ζ ψ ′ (ζ ) + 1) ∂ζ µ ∂ζ ∂η µ ζ∂η ∂ h 3 ∂Cw p ∂ h 3 ∂Cw p ′ ′ ζ ψ ζ ζψ ζ − ( ) − ( ) ∂ζ µ ∂η µ ∂ζ ∂η ∂ h3 ρ wu 2 ∂ ∂ h3 = ζψ ′(ζ ) ( 6CwUh ) + M (η ) − Cw µ ζ ∂η ∂η µ ∂η +
∂ ∂ζ
(10)
h3 2 Cw ρ w u µ
where
ψ ′(ζ )=
1 ζ cot β
The circumferential velocity of the fluid is
u=
Uz 1 1 ∂ + ( Cw p ) − M (η ) ( z 2 − hz ) h 2 µCw ζ ∂η
(11)
Boundary conditions for Eq. (10) are
2π p( η ) = p η + N p ζ =r = 0
(12)
out
p ζ = r = pin in
Because film thickness of the groove–ridge boundary is discontinuous, the 6
generalized Reynolds equation (10) is not satisfied at the groove–ridge boundary, the pressure of the groove–ridge boundary should be calculated using finite volume flow balance condition.
Fig.3 Finite volume separated by boundary
The finite volume including the groove-ridge boundary is established, as shown in Fig. 3. Finite volume flow balance condition is represented by the equation (13)
QAB + QBF + QFC + QCD + QDE + QEA = 0
(13)
where Qη = ∫
ζ2
ζ1
γ qη d ζ
η2
Qζ = ∫ qζ ζ dη
(η =AB,CD )
(14)
(ζ =BF ,FC ,DE,EA )
(15)
η1
γ = 1+
1 cot 2 β
Uη h h3 1 ∂ − ( Cw p ) − M (η ) 12 µ C w ζ ∂η 1 2 qη = 3 2 γ + 1 h 1 ∂Cw p − ζψ ′(ζ ) ∂C w p − ρ wu cot β 12 µ C ∂ζ ζ ∂η ζ w
qζ = −
h3 1 ∂Cw p ∂C p ρ u 2 − ζψ ′(ζ ) w − w 12 µ Cw ∂ζ ζ∂η ζ 7
(16)
(17)
2.2 The temperature control equation of cavitating flows The temperature field of cavitating flows was calculated using energy equation, the following assumptions were given :
（1）Considering convection heat transfer was much greater than heat conduction under the high speed condition, the heat conduction is neglected in the energy equation.
（2）Considering the water film is extremely thin, the variation of temperature or specific heat at constant volume of water in the direction of film thickness is neglected.
（3）Considering the radial velocity is much smaller than the circumferential velocity, the convective heat transfer in the radial direction is neglected.
（4）The viscous dissipation of the film thickness is only considered in the energy equation. According the above assumption (2), a average temperature of the water in the direction of film is defined as h
Tm = ∫ c v Tvθ dz / c v qη 0
Hence, the energy equation of cavitating flows in the (ζ − η ) coordinates system is simplified as h h ∂ ∂u ∂v ρ w c v qη ( CwTm ) = ∫0 Cw µ dz + ∫0 Cw µ dz ζ∂η ∂z ∂z 2
2
(18)
+ Ev (η ) (u − ub ) uh +Eσ (η )σ uh 2
where, the left-hand side of equation (18) represents the convective heat transfer term, the first and second items on the right-hand side of equation (18) represent the heat dissipation term, the third and fourth items on the right-hand side of equation (18) indicate the energy flowing into the water due to the gas-liquid interface effect. where ∂u U 1 1 ∂ ( Cw p ) = + − M (η ) ( 2 z − h ) ∂z h 2 µ C w ζ ∂η 8
(19)
∂v 1 ∂ ( Cw p ) ∂ ( Cw p ) ρ wu 2 ′ = − − ( ) ψ ζ (2z − h) ∂z 2 µ C w ∂ζ C w ∂η r
(20)
23 ∞ 3 23 Ev (η ) = 4π CDs ρ w ∫ V f ( η , V ) dV ( ) eq 0 4π 13 ∞ 3 13 Eσ (η ) = 4π ∫ V ) f eq (η ,V )dV ( 0 4π
The boundary conditions for Eq. (18) are
Tm ζ = r = T0 in
∂Tm ∂ζ
=0 ζ = rout
（21）
∂Tm ∂T = m =0 ∂η η =0 ∂η η =2π The viscosity–temperature relationship of fluid is described by an exponential function (22)
µ = µ0e−α
T
(Tm −T0 )
(22)
2.3 The bubble velocity control equation The velocity of the bubble in the cavitation flow is described by the bubble motion equation [17]
d ( ρgVub ) = ∑ Fi dt i
(23)
The forces acting on the bubble mainly include : (1) The force of relative acceleration between bubble and liquid
F1 = − V
∂p ∂x
(2) The viscous drag between the bubble and the liquid due to relative motion
1 F1 = ρ wπ R 2CD u − ub (u − ub ) 2 The inertial force can be ignored because the mass of the bubble is small, so the equation (23) is simplified as
9
1 ∂p ρ wπ R 2CD u − ub (u − ub ) − V =0 2 r∂ θ
(24)
The bubble velocity can be obtained by solving equation (24)
2.5 The population balance equation of bubble. The bubble volume of cavitating flow changes randomly due to breakup and coalescence effect. The breakup effect causes large volume bubbles to break down into small volume bubbles, and the coalescence effect causes small volume bubbles to be grown up into large volume bubbles. The random state of bubble volume can be described by using a bubble number density function f (η ,V ) , the bubble expansion due to change of gas density and mass transfer is ignored, the evolution of the bubble number density function f (η ,V ) is described by the following population balance equation of bubble: ∞ 6 ∂f ∂f 1 + ub =∫ h (ν ,V ) b (ν ) f (ν )dν ∂t ζ∂η V π 3ν 2 3 b 1 V − b (V ) f (V ) + ∫ c(ν ,V − ν ) f (ν ) f (V − ν )dν 2 0 13
(25)
∞
− f (V ) ∫ c(ν ,V ) f (ν )dν + S 0
where, the first term on the right-hand side of Eq. (25) represents breakup source term of bubbles, the second term on the right-hand side of Eq. (25) represents breakup sink term of bubbles, the third term on the right-hand side of Eq. (25) represents coalescence source term of bubbles, the fourth term on the right-hand side of Eq. (25) represents coalescence sink term of bubble, the fifth term on the right-hand side of Eq. (13) represents gas nucleus (microbubble) source term. The gas nucleus source term represents bubbles with volume V formed due to expansion of gas nucleus in the water under low pressure state. The expression of the source term S was
S = Cρ
∂ ( χ (V ) ) ∂V
( pv − p )
where
3 2 3 2 1 2 Cρ = 4π 4π 3ρ w 10
χ (V ) = ∫ (ϕ( ν ,V ) f i (ν )ν 2 3 )dν V
vl
where the ϕ (ν ,V ) was the redistribution function of bubble volume V ; the redistribution function of bubble volume V was given as [18] 2 V −ν ϕ( ν ,V ) = 2 exp −c2 V ν max
c1
The gas nuclei distribution is given as [19]
f i (ν ) = N i
κ α −1 −(κ ν ) −α e ν Γ (α − 1)
where, c1 =0.1 , c2 =0.5 , κ =6 , α =4 The breakage frequency is given as [20] b( ν )=
κ 1ε 1 3 κ 2σ exp − 29 5 9 ρ ε 2 3 ( 6ν π ) ( 6ν π ) w
The daughter size redistribution function is given as[20]
hb ( v,V ) =
2.4π ( 6V π )
23
ν
(2V − ν )2 exp −4.5 ν2
where the kernel function of coalescence is defined as respectively [21]
c(ν ,V ) = 0.0906
2 3
ν +V 2
2 3
ub − u
1 0.6 − C g
The initial condition for Eq. (25) is
f
t =0
=0
(26)
The boundary condition for Eq. (25)
f f
V =0
=0
V =∞
=0
are
(27)
2.6 Perturbation of Generalized Reynolds Equation of Cavitating Flows In practice, the SGTB has an axial displacement and two angular displacements, so the bearing should be modeled under three degrees of freedom. Supposing the origin of Cartesian coordinate is fixed in the steady equilibrium position of the bearing 11
surface, as shown in Fig.1, the independent motions of the bearing are the axial movement along z-axis, the angular movement around x-axis and the angular movement around y-axis respectively, where ( ∆z, ∆ϕ x , ∆ϕ y ) represent perturbations quantity of the thrust plate under three degrees of freedom. Assuming the initial position of the thrust plate is ( hR , ϕ x 0 , ϕ y 0 ) , the film thickness and pressure in quasi-equilibrium can be expressed as
h = h0 + ∆z + r cos θ∆ϕ y − r sin θ∆ϕx
(
where, h0 = hR + r −ϕ x 0 sinθ + ϕ y 0 cos θ
(28)
)
p = p0 + ∑ pq ∆q + ∑ pq& ∆q& q&
q
( q = z, ϕ , ϕ ) x
y
(29)
Substitute equations (28) and (29) into equation (10), perturbed generalized Reynolds equations of cavitating flows can be obtained
Rey ( pq ) = Fq (η , ζ ) where, the operators Rey (
(q = z, ϕ x , ϕ y , z&, ϕ&x , ϕ& y )
(30)
) was defined as
ζ h03 ∂ ( ) 12 µ ∂ζ ∂ h03 ∂ ( ) ∂ ζ h03 ∂ ( ) ′ ′ −ζψ (ζ ) − ψ (ζ ) ∂ζ 12 µ ∂η ∂η 12 µ ∂ζ Rey (
)=
∂( ) ∂ ∂ h03 ζ 2ψ ′2 (ζ ) + 1) ( + ∂η 12 µζ ∂η ∂ζ
The expressions of Fq (η , ζ ) were defined as follows ∂ CwU η ∂ h02 M (η ) + ∂η 2 ∂η 4µ ∂ ∂ h02 h02 2 + ′ C u − ( ) C ρ ψ ζ ρ wuη2 w w w η ∂η 4µ ∂ζ 4µ Fq (η , ζ ) = 2 2 − ∂ h0 ζ 2ψ ′2 (ζ ) + 1 ∂Cw p0 + ∂ h0 ζ ∂Cw p0 ) ∂η ∂ζ 4µ ∂ζ ∂η 4µζ ( 2 2 ∂ h0 ∂Cw p0 ∂ h0 ∂Cw p0 +ζψ ′(ζ ) + ψ ′(ζ ) ζ ∂ 4 ∂ ∂ ∂ζ ζ µ η η 4µ
12
q=z
∂ h02 ∂C p sin (η ′ ) (ζ 2ψ ′2 (ζ ) + 1) w 0 ∂η ∂η 4µ ∂ h2 ∂Cw p0 2 0 + ζ sin (η ′) ∂ζ ∂ζ 4µ 2 − 1 ∂ h0 ζ sin (η ′) ∂Cw p0 cot β ∂ζ 4µ ∂η 1 ∂ h02 2 ∂Cw p0 − ζ cot β ∂η 4µ ζ sin (η ′ ) ∂ζ Fq (η ,ζ ) = CU − ∂ w η ζ sin (η ′ ) ∂η 2 2 − ∂ h0 M (η ′ ) ζ sin η ′ ( ) ∂η 4µ 2 2 ∂ h0 ρ wuη − ζ sin (η ′) Cw 4µ ∂ζ 2 2 ∂ h0 ρ wuη + 1 ′) C sin ζ η ( w ζ cot β ∂η 4µ
∂ CwUη ζ cos (η ′ ) ∂η 2 ∂ h 2 M (η ′ ) 0 + ζ cos (η ′ ) ∂η 4µ 2 2 ∂ h0 ρ wuη ′ ζ η + C cos ( ) w ∂ζ 4µ 2 2 1 ∂ h0 ρ wuη − ζ cos (η ′) Cw 4µ ζ cot β ∂η Fq (η ,ζ ) = 2 ∂ ζ h0 1 ∂Cw p0 cos (η ′) + cot β ∂η ∂ζ 4µ 1 ∂ ζ h02 ∂C p + cos (η ′) w 0 ∂ζ cot β ∂η 4µ 2 − ∂ ζ h0 cos (η ′) (ζ 2ψ ′2 (ζ ) + 1) ∂Cw p0 ∂η 4µ ∂η ∂Cw p0 ∂ ζ 2 h02 − ∂ζ 4µ cos (η ′ ) ∂ζ
13
q = ϕx
q = ϕy
Fq (η , ζ ) =Cwζ
q = z&
Fq (η , ζ ) = − Cwζ 2 sin (η ′ )
q = ϕ& x
Fq (η , ζ ) =Cwζ 2 cos (η ′)
q = ϕ& y
where
η ′=η +ψ (ζ ) The boundary conditions for Eq. (30) are
pq pq pq
=0
ζ =rin
ζ =rout
η =0
(q = z, ϕ x , ϕ y , z&, ϕ&x , ϕ& y )
=0
= pq
(31)
η =2π
2.7 Calculation of stiffness coefficients and damping coefficients The dynamic characteristic parameters of the water-lubricated SGTB includes the stiffness coefficients and damping coefficients, which are calculated using the perturbation method. Once the perturbation pressure are determined by solving equation (30), the stiffness coefficient and damping coefficient can be obtained by integrating the perturbed pressure. For the SGTB with three degrees of freedom, totally nine stiffness coefficients and nine damping coefficients should be calculated using equations (32) and (33). K zz K zϕ K zϕ x y Kϕ x z Kϕ xϕ x Kϕ xϕ y Kϕ y z Kϕ yϕ x Kϕ yϕ y
−1 2π rout = ∫ ∫ ζ sin (η + f (ζ ) ) pz pϕ x pϕ y ζ d ζ d η 0 rin −ζ cos (η + f (ζ ) )
(32)
C zz Czϕ C zϕ x y Cϕ x z Cϕ xϕ x Cϕ xϕ y Cϕ y z Cϕ yϕ x Cϕ yϕ y
−1 2π rout = ∫ ∫ ζ sin (η + f (ζ ) ) pz& pϕ&x pϕ& y ζ d ζ d η 0 rin −ζ cos (η + f ( ζ ) )
(33)
{
{
}
}
The dimensionless stiffness coefficients and dimensionless damping coefficients are defined as follows
14
K zz∗ K z∗ϕ K z∗ϕ x y Kϕ∗ z Kϕ∗ ϕ Kϕ∗ ϕ x x x y x ∗ ∗ ∗ Kϕ z Kϕ ϕ Kϕ ϕ y x y y y
Czz∗ Cz∗ϕ Cz∗ϕ x y ∗ ∗ ∗ Cϕ z Cϕ ϕ Cϕ ϕ x x x y x ∗ ∗ ∗ Cϕ z Cϕ ϕ Cϕ ϕ y x y y y
K zz h3 K = 3 g ϕx z routUθ µ0 rout Kϕ y z r out Czz h3 C = 4 g ϕx z rout µ0 rout Cϕ y z r out
K zϕx
K zϕ y
rout
rout
Kϕxϕx Kϕxϕ y 2 rout
2 rout
Kϕ yϕx Kϕ yϕ y 2 rout
2 rout
Czϕx
Czϕ y
rout
rout
Cϕ xϕx Cϕxϕ y 2 rout
2 rout
Cϕ yϕx Cϕ yϕ y 2 rout
2 rout
3. Numerical calculation of the model The equations (10), (18) and equation (30) were discretized using the finite difference method, and the corresponding linear equations were obtained. The pressure distribution p and pq , the fluid velocity u and bubble velocity ub were calculated using the SOR iterative method, and the temperature field Tm were calculated using the stepping approach. The population balance equation (25) is an integral-differential equation of parabolic type. Since the coalescence source and sink terms in the population balance equation are nonlinear, the discretization equation of population balance equation is a system of nonlinear equations. Defining the solution domain node index
q = i + ( j − 2) Ni with 1 ≤ i ≤ Ni and 2 ≤ j ≤ N j − 1 . At time t , a finite difference
(
)
discrete operator L% ( f )q at q node is defined as
L%
(( f ) ) = ( K1)
N j −1
t
q
t
q
f + ( K 2 )q f − ∑ ( K 3)q f qt + ( K 4 ) q f qt t q2
t q1
N j −1
j −1
−∑ ( K 5 ) q f f + f e =1
e
t qe
t ql
t q
∑ ( K 6) e =1
qe
e
e= j
f − Sq = 0 t qe
e
( q = 1,...N )
(35)
q
Applying time-differencing through the λ weighted average method, the discretization equation for the equation (25) at q node can be written as 15
( f )q
t +1
−( f
)q +∆t ( (1 − λ ) L% ( ( f )q ) +λ L% ( ( f )q t
t +1
t
)) =0
(36)
where 0 ≤ λ ≤ 1 , the equation (36) can be rewritten as Fq
(( f ) ) = ( f ) t +1
t +1
q
q
+∆t λ L%
(( f ) ) +∆t (1 − λ ) L% (( f ) ) − ( f ) t +1
t
t
q
q
q
=0
( q = 1,...N ) q
(37) Equation (37) is a system of nonlinear equations for
( f )q
t +1
, which is solved
numerically by the Newton-SOR iterative method. In order to improve convergence, the downhill condition is included in the iterative procedure, the iterative procedure is written as
(
t +1( k )
Fq′ f q
t +1( k +1)
fq
(
) ∆f = − F ( f q
t +1( k )
= fq
)
q
t +1( k ) q
)
+ ωk ∆f
(
Fq f qt +1( k +1) ≤ F f qt +1( k )
(38)
)
1] . where the downhill factor ωk ∈ [ 0， The nine stiffness coefficients and nine damping coefficients were calculated by integrating (32) and (33) numerically. Fig. 4 shows the flow chart of the whole algorithm.
16
Start
Assign calculation parameters, the space coordinate
t =1 Calculate the interface items M (η )
Solve the generalized Reynolds equation ( 4 ) with SOR
Judge convergence of pressure at t vη = α vη ( k ) + (1 − α )vη (k −1)
vb = β vb ( k ) + (1 − β )vb ( k −1)
Solve equation ( 5) to obtain new fluid velocity
Judge convergence of fluid velocity at t
t = t +1
Solve the energy equation ( 9 ) with step method
Judge convergence of temperature at t Solve equation (12 ) to obtain new bubble velocity
Judge convergence of bubble velocity at t Solve the transport equation of bubble volume (13 ) with N-SOR
Judge convergence of bubble disturbution at t Output pressure, velocity, viscosity p, vη , µ and other parameters and bubble balance distribution f eq (η , V )
Fig.4 Flow diagram of solution algorithm
17
Small disturbance initial parameters ( ∆z, ∆ϕx , ∆ϕ y )
4. Experimental setup and procedure 4.1 Experimental setup An experimental setup is developed to measure the direct stiffness coefficient of the high-speed water-lubricated SGTB. Fig. 5 presents the schematic view and real scene image of the experiment setup. The grooved thrust disc is installed on the end of the high-speed motorized spindle, the rotational speed is controlled with an inverter; the thrust pad is fixed in the bearing housing; the external load is applied to the bearing housing by using a loading unit, and the external load is applied by adding a dead weight. The thrust pad displacement due to water film formation is detected by a non-contact displacement probe using eddy current; signals collected from the sensors are fed into a data acquisition system and interfaced with a computer. The water lubricant with ordinary pressure is supplied to the test bearing by a water circulation unit.
(a) Schematic view
18
(b) Real scene image Fig. 5 Experimental setup to measure the stiffness coefficient of SGTB 5.2 Specimens Fig. 6 shows the specimens of inwardly pumping SGTB. The material of the test thrust discs is stainless steel, the disc surfaces are machined using a grinding processing, and then the spiral grooves are etched on the surfaces using laser processing. The material of the test thrust pad is plexiglass, and a hole is fabricated in the center of pad, which is acted as outlet of water lubricant for the SGTB.
Fig. 6 Experimental specimens of SGTB 5.3 Measurement of stiffness coefficients of SGTB (1) Film thickness The motorized spindle is started to drive the disc up to a given rotational speed; an 19
axial load is applied on the test bearing. The film thickness is obtained by detecting the axial displacement of pad with the eddy current probe. Here, the film thickness refers to the film thickness in ridge region, namely, the nominal film thickness. The data is measured three times for each condition. (2) Axial stiffness coefficient The motorized spindle is started to drive the dis,c up to a given rotational speed; an initial axial load is applied on the test bearing, and a small load (about 5 percent of the applied load) is added on the bearing, then the variation of film thickness due to this small increment of load is measured. Thus, the axial stiffness coefficient of the bearing at this condition is obtained with dividing the small load by the variation of film thickness.
5. Results and Discussion 5.1 Validation of the present model A comparison of the predicting stiffness coefficients and the experimental ones is performed to verify the proposed model. Considering that it is difficult to test the angular stiffness coefficients of the SGTB, only the axial stiffness coefficient Kzz is measured in this study. The specifications of the spiral groove thrust bearing in this test are listed in Table 1. Table 1 The parameters of the spiral groove thrust bearing
Item
Value
Inner radius rin (mm)
7.5
Outer radius rout (mm)
20
Base circle radius rb (mm)
12
Water viscosity µ ( Pa.s )
0.001
Groove number N
12
Groove depth hg ( µm)
40
Spiral angles β ( o )
20
20
Pressure of supply water pin ( MPa)
0.1
Specific heat at constant volume of water cv ( J kg .Co )
4200
Density of water ρw ( kg m3 )
1000
Bearing capacity W ( N )
20
Surface tension of bubble σ ( N m )
0.3
Groove-to-land ratio λ b
0.5
Groove-to-dam ratio λl
0.6
Fig.7 gives a comparison between the predicting axial stiffness coefficient Kzz∗ and the experimental ones under three rotational speeds. The results show the predicting values generally agree with the experimental values, and the results with the cavitating flow are closer to the experimental values than those with non-cavitating flow, which indicates that the proposed model can be used to analyze the dynamic coefficients for high speed water-lubricated spiral groove thrust bearing considering cavitating effect.
(a) 9000 r/min
21
(b) 12000 r/min
(c) 15000 r/min
Fig.7 The predicting stiffness coefficients vs. the experimental ones It can be seen that the predicting values of axial stiffness coefficient are always lower than the experimental values. A possible explanation for this might be that the centrifugal effect on the water lubricant causes an error for the experimental values. In this study, the water lubricant is supplied by a pump and injected to the bearing clearence from the outside of bearing through three nozzles. The inwardly pumping bearing is employed, so the water is pumped into the bearing from the outside, and then flows through the bearing, finally reach the inside of the bearing. The water in the region near the inlet is disturbed by the centrifugal effect, and part of the water is thrown out from the thrust pad rim. Therefore, the actual flowrates of the bearing are less than the predicting values, which means that the water film thickness is less than the theoratical one. Therefore, it should be resonable that the the predicting values of 22
axial stiffness coefficient are always lower than the experimental values. In addition, the centrifugal effect on the water will increase with the increasing of rotary speed, as a result, the difference value between the predicting value and the experimental one gets larger with the inceasing of rotary speed.
5.2 The effect of cavitation on the dynamic characteristics of SGTB 1） ）Spiral angle Fig.8 shows the effect of cavitation on the stiffness coefficients of the SGTB with different spiral angles. It can be seen that the cavitation affects the stiffness coefficient Kzz∗ obviously; and the direct stiffness coefficients Kzz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ with x x
y y
cavitating effect are larger than those without cavitating effect when the spiral angles is less than 45 degrees. This result may be explained by the fact that the hydrodynamic effect is enhanced due to the interfacial effect of the cavitating flow, and the load-carrying capacity of the SGTB with cavitating flow is higher than that with non-cavitating flow, which means that a larger force is required to produce a unit displacement for the SGTB with the cavitation flow, so the stiffness coefficients Kzz∗ , Kϕ∗xϕ x , Kϕ∗yϕ y increase. The deviation of the direct stiffness coefficients between the
cavitating flow and the non-cavitating flow decreases with the increasing of spiral angles, which means the cavitation effect becomes weak as the spiral angle increases. In addition, no matter whether of the cavitating flow or the non-cavitating flow, the direct stiffness coefficients Kzz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the cross-coupled stiffness coefficients x x
y y
Kϕ∗yϕx , Kϕ∗xϕ y decrease with the increasing of the spiral angles. The cross-coupled stiffness
coefficients Kϕ∗ ϕ , Kϕ∗ ϕ are of opposite sign ( Kϕ∗ ϕ = − Kϕ∗ ϕ ). The x y
y x
x
y
y x
cross-coupled
stiffness coefficients Kϕ∗ z , K z∗ϕ , K z∗ϕ , Kϕ∗ z are approximately equal to zero for x
x
y
y
cavitating flow and non-cavitating flow. The cavitation effect hardly affect the cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ . x y
y x
23
2.0
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
cavitating flow non-cavitating flow
2.5
K*zz
1.5 1.0 0.5
K*zϕϕ
K zϕ y *
x
0.0 0
10
20
30
40
50
60
70
cavitating flow non-cavitating flow
0.7 0.5
K*ϕxϕx
0.3 0.1 -0.1
K* ϕ x z
-0.3 -0.5
K*ϕxϕy
-0.7 -0.9 -1.1 0
10
Spiral angle β (°)
30
40
50
60
70
(b)
Dimensionless stiffness coefficients
(a)
20
Spiral angle β (°)
cavitating flow non-cavitating flow
1.2 1.0
K*ϕyϕx 0.8 0.6 0.4
K*ϕyϕy
0.2
K*ϕyz
0.0 0
10
20
30
40
50
Spiral angle β (°)
60
70
(c) Fig.8 Effect of cavitation on the stiffness coefficient of SGTB with different spiral angles
( ω = 15000rpm , hg = 40 µ m , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 ) Fig. 9 shows the effect of cavitation on the damping coefficients of the SGTB with different spiral angles. It can be seen that cavitation effect hardly affect the nine damping coefficients, this result may be explained by the fact that the cavitation effect has a great influence on the pressure p0 of equilibrium position, while the perturbed pressure components ( p z& , pϕ& , pϕ& ) associating with damping coefficient are only x
y
related to the volume fraction Cg of bubbles per unit volume of water rather than the pressure p0 . In other word, the damping coefficients are independent of the pressure p0 at equilibrium position, furthermore, the volume fraction Cg of bubbles per unit volume of water is usually small. Therefore, the effect of cavitation on the 24
damping coefficient is very weak. No matter whether the cavitating flow or the non-cavitating flow, the direct damping coefficients C zz∗ , Cϕ∗ ϕ , Cϕ∗ ϕ decrease gradually with the increasing of the x x
y
y
spiral angles, the cross-coupled damping coefficients C z∗ϕ , Cz∗ϕ , Cϕ∗ z Cϕ∗ ϕ Cϕ∗ ϕ Cϕ∗ z are x
y
x
x
y
y x
y
cavitating flow non-cavitating flow
0.12
Dimensionless damping coefficients
Dimensionless damping coefficients
approximately equal to zero.
0.08
C*zz 0.04
C*zϕϕx
C*zϕϕy
0.00 0
10
20
30
40
50
60
Spiral angle β (°)
70
cavitating flow non-cavitating flow
0.030 0.025 0.020
C*ϕxϕx
0.015 0.010 0.005
C*ϕxϕy
0.000
C*ϕxz
-0.005 0
10
30
40
50
60
70
(b)
Dimensionless damping coefficients
(a)
20
Spiral angle β (°)
cavitating flow non-cavitating flow
0.030 0.025 0.020
C*ϕyϕy
0.015 0.010 0.005
C*ϕyϕx
0.000
C*ϕyz
-0.005 0
10
20
30
40
50
Spiral angle β (°)
60
70
(c) Fig.9 Effect of cavitation on the damping coefficients of SGTB with different spiral angles
( ω = 15000rpm , hg = 40 µ m , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
2) Groove depth Fig. 10 shows the effect of cavitation on the stiffness coefficients of the SGTB with different groove depths. It can be seen that the direct stiffness coefficients Kzz∗ , Kϕ∗xϕx , Kϕ∗yϕ y with cavitating effect are larger than those without cavitating effect, the 25
deviation of the direct stiffness coefficients between the cavitating flow and non-cavitating flow increases with the increasing of groove depth. A possible explanation for this is that the hydrodynamic pressure decreases when the groove depth increases, leading to an increase in the number of bubbles per unit volume of water lubricant. The increasing of the number of bubbles will cause the hydrodynamic effect of interface to be enhanced. Furthermore, the hydrodynamic effect of interface due to the bubble number is more dominant than the hydrodynamic effect due to the groove depth. Therefore, as the groove depth increases, the direct stiffness coefficients with cavitating effect are larger than those without cavitating effect. The effect of cavitation on the cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ is very weak for x y
y x
different groove depths. In addition, no matter whether the cavitating flow or the non-cavitating flow, the direct stiffness coefficients K zz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the x x
y y
cross-coupled stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ decrease gradually with the increasing y x
y x
of groove depth, the cross-coupled stiffness coefficients K z∗ϕ , K z∗ϕ , Kϕ∗ z , Kϕ∗ z are x
y
x
y
cavitating flow non-cavitating flow
12
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
approximately equal to zero.
10 8
K*zz
6 4 2
K*zϕϕx
K*zϕϕy
0 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(a)
cavitating flow non-cavitating flow
4
K*ϕxϕx
2
K*ϕxz 0
K*ϕxϕy
-2 0
5
10
15
20
(b)
26
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless stiffness coefficients
cavitating flow non-cavitating flow
3.5 3.0
K*ϕyϕy
2.5 2.0 1.5
K*ϕyϕx
1.0 0.5
K*ϕyz
0.0 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(c) Fig.10 Effect of cavitation on stiffness coefficients of SGTB with different groove depths
( ω = 15000rpm , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
Fig. 11 shows the effect of cavitation on the damping coefficients of the SGTB with different groove depths. It can be also seen that the effect of cavitation hardly affect the damping coefficients, the cause for this result is identical with the section above. For the cavitating flow and the non-cavitating flow, the direct damping coefficients C zz∗ , Cϕ∗xϕ x , Cϕ∗yϕ y decrease nonlinearly with the increasing of groove depth, the
cross-coupled damping coefficients C z∗ϕ , Cz∗ϕ , Cϕ∗ z Cϕ∗ ϕ Cϕ∗ ϕ Cϕ∗ z are approximately x
y
x
x
y
y x
y
0.20
cavitating flow non-cavitating flow
0.15
C*zz 0.10
0.05
C*zϕϕx
C*zϕϕy
0.00 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless damping coefficients
Dimensionless damping coefficients
equal to zero.
(a)
cavitating flow non-cavitating flow
0.06 0.05 0.04
C* ϕ x ϕ x 0.03 0.02 0.01
C*ϕxϕy
C*ϕxz
0.00 0
5
10
15
20
(b)
27
25
30
35
Groove Depth hg (µm)
40
45
Dimensionless damping coefficients
cavitating flow non-cavitating flow
0.06 0.05 0.04
C*ϕyϕy
0.03 0.02
C* ϕ y ϕ x
0.01
C* ϕ y z
0.00 0
5
10
15
20
25
30
35
Groove Depth hg (µm)
40
45
(c) Fig.11 Influence of cavitation effect on damping coefficients for different groove depth
( ω = 15000rpm , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
3) Rotational speed Fig. 12 shows the effect of cavitation on the stiffness coefficients of the SGTB with different rotational speeds. It can be seen that the direct stiffness coefficient K zz∗ with cavitating effect are larger than those without cavitating effect, and the direct stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ with cavitating effect are slightly larger than those without x x
y y
cavitating effect in the whole speed range. A possible explanation for this is that the the number of bubbles per unit volume of water lubricant increases with the increasing of rotary speed, so the hydrodynamic effect of interface is enhanced. Consequently, the load-carrying capacity of the SGTB with cavitating flow is larger than that with non-cavitating flow. Furthermore, no matter whether the cavitating flow or the non-cavitating flow, the direct stiffness coefficients K zz∗ , Kϕ∗ ϕ , Kϕ∗ ϕ and the x x
y y
cross-coupled stiffness coefficient Kϕ∗ ϕ increase approximately linearly with the y x
increasing
of
rotational
speed,
and
the
cross-coupled
stiffness
coefficient Kϕ∗ ϕ decreases approximately linearly with the increasing of rotational x y
speed, the cross-coupled stiffness coefficients K z∗ϕ , K z∗ϕ , Kϕ∗ z , Kϕ∗ z are approximately x
equal to zero. 28
y
x
y
4.0
cavitating flow non-cavitating flow
Dimensionless stiffness coefficients
Dimensionless stiffness coefficients
4.5
3.5
K zz *
3.0 2.5 2.0 1.5 1.0 0.5
K*zϕϕy
K*zϕϕx
0.0 8000 10000 12000 14000 16000 18000 20000 22000
cavitating flow non-cavitating flow
1.5 1.0 0.5
K*ϕxϕx 0.0
K*ϕxz
-0.5 -1.0 -1.5
K*ϕxϕy -2.0 8000 10000 12000 14000 16000 18000 20000 22000
rpm
rpm
(b)
Dimensionless stiffness coefficients
(a)
2.0
cavitating flow non-cavitating flow K*ϕyϕx
1.5
1.0
0.5
K*ϕyϕy
K*ϕyz 0.0
8000 10000 12000 14000 16000 18000 20000 22000
rpm
(c) Fig.12 Effect of cavitation on stiffness coefficients of SGTB with different rotational speeds
( hg = 40 µ m , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 ) Fig. 13 shows the effect of cavitation on the damping coefficients of the SGTB with different rotational speeds. It can be also seen that the cavitation effect hardly affect the damping coefficients, the cause for this result is identical with the section above, which is not repeated here. For the cavitating flow and the non-cavitating flow, the direct damping coefficients Czz∗ , Cϕ∗ ϕ , Cϕ∗ ϕ and cross-coupled damping coefficients x x
Cϕ∗xϕ y , Cϕ∗yϕ x
y y
are approximately constant with the increasing of rotational speed, and
cross-coupled damping coefficients Cϕ∗ ϕ , Cϕ∗ ϕ are of opposite sign ( Cϕ∗ ϕ = −Cϕ∗ ϕ ). The x
y
y
x
x
y
y
x
cross-coupled damping coefficients C z∗ϕ , C z∗ϕ , Cϕ∗ z , Cϕ∗ z are approximately equal to zero. x
y
29
x
y
C zz *
0.08 0.07
cavitating flow non- cavitating flow
0.06 0.05
Dimensionless damping coefficients
Dimensionless damping coefficients
0.09
0.0004 0.0003 0.0002 0.0001 0.0000
-0.0001
0.04
C* ϕx ϕx
C*ϕyϕy
C*ϕxz
-0.0002
0.03
cavitating flow non-cavitating flow C*ϕyz
C*ϕyϕx
C*ϕxϕy
-0.0003
0.02 8000 10000 12000 14000 16000 18000 20000 22000
-0.0004 8000 10000 12000 14000 16000 18000 20000 22000
rpm
rpm
(b) Dimensionless damping coefficients
(a)
0.02
cavitating flow non-cavitating flow
0.01
0.00
C*zϕϕy
C*zϕϕx -0.01
-0.02 8000 10000 12000 14000 16000 18000 20000 22000
rpm
(c) Fig.13 Effect of cavitation on damping coefficients of SGTB with different rotational speeds
( hg = 40 µ m , β = 20o , h0 = 15µ m , ϕ x 0 = ϕ y 0 = 0 )
6. Conclusions A theoretical model of cavitation flow lubrication considering various gas-liquid interface effects has been established. The dynamic coefficients of the high-speed water-lubricated hydrodynamic SGTB is calculated using the model. An experimental setup has been developed to verify the proposed model. Through numerical simulation and discussion, the following conclusions can be drawn: (1) The direct stiffness coefficients of SGTB with cavitating flow are larger than those of SGTB with non-cavitating flow when the spiral angles is less than 45 degrees, or the groove depth reaches a greater value. The effect of cavitation on the direct 30
stiffness coefficients is gradually weakened with the increasing of spiral angle; however, the effect of cavitation on the direct stiffness coefficients is gradually enhanced with the increasing of groove depth; (2) The influence of cavitation on the direct stiffness coefficients K zz∗ depends on the rotational speed. The direct stiffness coefficients Kϕ∗ ϕ , Kϕ∗ ϕ with cavitating effect are x x
y y
larger than those without cavitating effect slightly in the whole range of speed. (3) The cavitation hardly affect the cross-coupled stiffness coefficients and the damping coefficients of SGTB.
Acknowledgment The authors gratefully acknowledge the support of the National Science Foundation of China through Grant Nos. 51635004 and 11472078.
31
Nomenclature b ( v ) The bubble breakage frequency function （ ） 1 s
Cij (i = z, ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine damping coefficients（
N .s N .s N .m.s , , N .s, ） m rad rad
Cij∗ (i = z , ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine dimensionless damping coefficients
Cw the water volume fraction (dimensionless)
Cg the bubble volume fraction (dimensionless) Cπ Constant coefficient (dimensionless)
C Ds drag coefficient (dimensionless) f eq (θ , V ) the equilibrium distribution function of bubble volume（
1 m6
）
feq (η ,V ) the equilibrium distribution function of bubble volume at coordinates system（
1 m6
(ζ ,η )
）
f eq* the dimensionless equilibrium distribution function of bubble volume (dimensionless)
f (V ) the distribution of bubble volume (
1 ) m6
h the film thickness ( m )
hg the groove depth ( m )
hR steady film thickness at ridge area ( m ) h0 Static film thickness ( m ) hb ( v,V ) bubble volume redistribution function (
1 ) m
Ki , j (i = z, ϕ x , ϕ y , j = z, ϕ x , ϕ y ) nine stiffness coefficients（
N N N.m N.m ; ; ; ） m rad m rad
K ij∗ (i = z , ϕ x , ϕ y , j = z , ϕ x , ϕ y ) nine dimensionless stiffness coefficients
(dimensionless)
32
M (θ ) the interface item ( N3 ) m
M (η ) the interface item at (ζ ,η ) coordinates system ( N3 ) m
N the number of spiral groove (dimensionless)
Ni the number density of nuclei ( p the water film pressure (
p0 Static pressure (
1 ) m3
N ) m2
N ) m2
pq , pq& perturbed pressure components (
N N N .s N .s , , , ) m3 m2 m3 m2
pv the vaporization pressure of water (
N ) m2
pout the pressure at outer radius ( pin the pressure at inter radius (
N ) m2
N ) m2 3
QAB ,QBF ,QFC ,QCD ,QDE ,QEA the volume flow ( m ) s
3
Qζ the flow rate in the ζ direction ( m ) s
3
Qη the flow rate in the η direction ( m ) s
2
qζ the ζ direction flows of water in spiral coordinates ( m ) s
2
qη the η direction flows of water in spiral coordinates ( m ) s
r the radial coordinate in polar coordinates ( m )
rb the base circle radius of spiral line ( m ) rin thrust disc inner radius ( m )
rout
thrust disc outer radius ( m )
R the statistical average radius of bubbles ( m )
c (ν ,V ) the bubble coalescence kernel function ( m ) 3
s
S the source terms
t the time ( s )
33
T Water film temperature ( K )
Tm the average temperature of water film ( K ) T0 the ambient temperature ( C o ) U ( = routω ) the circumferential velocity of thrust disc ( ) m s
u the circumferential velocity of water ( m ) s
m s
u the circumferential average velocity of water ( )
v the circumferential velocity of water at (ζ ,η ) coordinates system ( ) m s
m s
ub the bubble velocity ( ) V the bubble volume ( m3 )
Vmax maximum volume of bubbles ( m3 ) z the dimensionless coordinate of the film thickness (dimensionless) Greek Letters
α T the temperature coefficient of viscosity (
1 ) K
cv the specific heat capacity of water ( J / ( kg ⋅ K ) )
σ the surface tension of bubble (
N ) m
µ the dynamic viscosity of water (
Ns ) m2
κ , κ 2 the adjustable parameters (dimensionless) 1
ρ g the density of bubble ( kg3 ) m
ρw the density of water ( kg ) m3 2
ε the eddy diffusivity ( m ) s
λb the groove-to-land ratio (dimensionless) λl the groove-to-dam ratio (dimensionless) θ the circumferential coordinate in polar coordinates (deg) 34
ν the bubble volume ( m3 )
β the spiral angle (deg) η the spiral coordinates (deg)
ζ the spiral coordinates ( m ) ∆z, ∆ϕx , ∆ϕ y perturbations quantity ( m ) ϕ x 0 , ϕ y 0 the initial position of the thrust plate ( m )
Appendix A. The derivation of the generalized Reynolds equation
Cw
∂ ∂u 1 ∂ ( Cw p ) − M (θ , z ) µ = ∂z ∂z r ∂θ
Cw
∂ ∂v ∂ Cw ρ w u 2 µ = C p − ( w ) ∂z ∂z ∂r r
（A1）
（A2）
Equations (A1) and (A2) can be integrated,
∂u 1 ∂ ( Cw p ) z − F1( z )+c1 Cw µ = ∂z r ∂θ
Cw µ
(A3)
∂vr ∂ C ρ = ( Cw p ) z − w w F2 ( z ) + d1 ∂z ∂r r
(A4)
Where z
F1 ( z ) = ∫ M (θ , z′)dz′ 0
z
F2 ( z )= ∫ u 2 dz ′ 0
Further integration of equations (A3) and (A4) yields
( Cw µ u ) =
z 1 ∂ z2 ( Cw p ) − ∫0 F1( z′)dz′+c1z + c2 (A5) r ∂θ 2
∂ z 2 Cw ρ w C w µ v = ( Cw p ) − ∂r 2 r
∫
z
0
F2 ( z )dz ′ + d1 z + d 2
Boundary conditions are
35
(A6)
u z =h = U , u z =0 = 0 v z = h = 0, v z =0 = 0
(A7)
Substitute boundary conditions (A7) into equations (A5) and (A6) , the flow field distribution of the cavitating flow are obtained
( Cw u ) = Cw v =
h CwU 1 ∂ 1 z z + ( Cw p ) ( z 2 − hz ) − ( ∫0 F1 ( z′)dz′ − ∫0 F1 ( z′)dz′ ) (A8) h 2r µ ∂θ µ h
z z 2 hz 1 C ρ 1 ∂ z h ′ ′ ′ ′ F（ ( Cw p ) − − w w ∫0 F（ 2 z ）dz − 2 z ）dz ∫ h 0 µ ∂r 2 2 µ r
(A9)
Since the water film is very thin, the velocity u in the interface term and the inertia term can be approximated by the average velocity u of the film thickness. The u is defined as 1 h u = ∫ udz h 0 and F1 ( z ) = M (θ ) z
∫
z
∫
h
∫
z
0
0
0
F1 ( z ′)dz ′=M (θ )
z2 2
F1 ( z ′)dz ′=M (θ )
h2 2
′ ′ 2 F（ 2 z ）dz =u
z2 2
Equations (A8) and (A9) are simplified to
( Cw u ) = Cw v =
CwU 1 ∂ 1 + ( Cw p ) ( z 2 − hz ) − M (θ ) ( z 2 − hz ) （A10） h 2r µ ∂θ 2µ
z 2 hz 1 Cw ρ wu 2 2 1 ∂ C p ( w ) − − ( z − hz ) (A11) µ ∂r 2 2 µ 2r
Where
36
M (θ ) = (Cπ CDs ρ w (u − ub )2 − Cπ ρ wu (u − ub ) ) ∫ (V ) ∞
f eq (θ ,V )dV
23
0
+
σ
∞
(V ) 2∫ 0
13
f eq (θ ,V )dV
The continuity equation of the two-phase flow aqueous phase in polar coordinates is:
∫
h
0
h ∂ (C u ) ∂( rCw v ) w dz + ∫ dz = 0 0 ∂r ∂θ
（A12）
Substitute equations (A10) and (A11) into equations(A12) , The generalized Reynolds equation (A13) of the cavitating flow is given ∂ ∂r
h3 ∂ ∂ h3 1 r C p ( ) + w 12 µ ∂ r ∂ θ 12 µ r hU θ ∂ h 3 ρ wu 2 ∂ C + rC ( w) w 2 ∂θ ∂r 12 µ r
∂ (C w p ) = ∂θ
（A13）
∂ h3 M (θ ) + ∂ θ 12 µ
References [1] Murata S. Exact Two-Dimensional Analysis of Circular Disk Spiral-Groove Bearing. ASME Journal of Tribology 1979;101(4):424-430. [2] Elrod HG. Some Refinements of the Theory of the viscous Screw Pump. ASME Journal of Tribology 1973;95(1):82-93. [3] Muijderman EA. Spiral Groove Bearings. New York: Springer-Verlag; 1966. [4] Manaloski SB, Pan CHT. The Static and Dynamic Characteristic of the Spiral-Grooved
Thrust
Bearing.
ASME
Journal
of
Basic
Engineering
1965;87:547–558. [5] Lund L. Calculation of Stiffness and Damping Properties of Gas Bearings. ASME Journal of Lubrication Technology 1968;90:793–803. [6] Zhou H. A High Order, Isoparametric, Finite Element Method for the Calculation of Static and Perturbed Frequency Characteristics of Spiral-Grooved Self Acting Gas Bearings. Proc. 8th International Gas Bearing Symposium, BHRA Fluid Engineering, Cranfield, England; 1981 pp. 229–236. [7] Liu R, Wang XL, Zhang XQ. Effects of gas rarefaction on dynamic characteristics of
micro
spiral-grooved
thrust
bearing. 37
Journal
of
Tribology
2012;
134(2):222011-222017. [8]Yoshimoto S, Anno Y, Tamura M, Kakiuchi Y, Kimura K. Axial load carrying capacity of water-lubricated hydrostatic conical bearings with spiral grooves (on the case of rigid surface bearings). Trans ASME 1965; 118(4):893–899. [9] Yoshimoto S, Oshima S, Danbara S, Shitara T. Stability of water-lubricated, hydrostatic, conical bearings with spiral grooves for high-speed spindles. Trans ASME 2002;124(2):398–405. [10] Cabrera DL, Woolley NH, Allanson DR, Tridimas YD. Film pressure distribution in water-lubricated rubber journal bearings. Proc. Institution Mech. Eng., Part J: J. Eng. Tribol. 2005;219(2):125–132. [11] Gohara M, Somaya K, Miyatake M, Yoshimoto S. Static characteristics of a water-lubricated hydrostatic thrust bearingusing a membrane restrictor. Tribology International 2014;75:111–116. [12] Feng L, Bin L, Xiaofeng Z. Numerical design method for water-lubricated hybrid sliding bearings. International Journal of Precision Engineering & Manufacturing 2008;9(1):47-50. [13] Durazo-Cardenas IS, Corbett J, Stephenson DJ. The performance of a porous ceramic hydrostatic journal bearing. Proceedings of the Institution of Mechanical Engineers, Part J: Journal of Engineering Tribology 2010;224(1):81-89. [14] Lin XH, Jiang SY, Zhang CB, Liu X. Thermohydrodynamic analysis of high speed water-lubricated spiral groove thrust bearing considering effects of cavitation, inertia and turbulence. Tribology International 2018;119:645-658. [15] Lin XH, Jiang SY, Zhang CB, Liu X. Thermohydrodynamic Analysis of High-Speed Water-Lubricated Spiral Groove Thrust Bearing Using Cavitating Flow Model. Journal of Tribology Trans ASME 2018;140.051703:1-12. [16] Ni J, Wang G, Zhang H. The Basic Theory of Solid Liquid Two-Phase Flow and Its Latest Applications. Science Press, Beijing, China; 1991 pp.35-52. [17] Guo L. Two-Phase and Multiphase Flow Dynamics. Xi’an Jiaotong University Press, Xi’an, China; 2002 pp. 430. [18] Bork Z, Jakobsen HA. Evaluation of breakage and coalescence kernels for 38
vertical bubbly flows using a combined multifluid-population balance model solved by least squares method. Procedia Engineering 2012;42:682-694. [19] Billet ML. Cavitation Nuelei Measurements,Proe.of ZndIntorn Symp. on Cavitation Ineeption, Dee, 1984. p.33-42. [20] Coulaloglou CA, Tavlarides LL. Description of interaction processes in agitated liquid-liquid dispersions. Chemical Engineering Science 1977;32(11):1289-1297. [21] Millies M, Mewes D. Interfacial Area Density in Bubbly Flow. Chemical Engineering and Processing: Process Intensification 1999;38(4-6): 307–319.
39
Highlights
1) The dynamic characteristics of high-speed water-lubricated SGTB considering cavitating effect is studied. 2) The hydrodynamic lubrication model of cavitating flow for the water-lubricated SGTB was established. 3) The cavitation effect on dynamic characteristics of water-lubricated SGTB is analyzed.
Declaration of Interest Statement
This is to declare that there is no interest conflict about this paper.