4-2 五点差分法

§ 4－2 五点差分法求解动压轴承Reynolds 方程
一、差分法基本原理
1、 导数的差商表示：dx
dy
x y x =
∆∆→∆0lim
，用有限的差商代替连续的导数。

有三种表示方法： （1） 前差商：
δ
j
j y y dx dy -=+1 （2） 后差商：δ1
--=j j y y dx dy （3） 中差商：δ
2
12
1-
+-=j j y
y
dx
dy ，
δ
δ22211
11-+-+-=⎪⎪⎭⎫ ⎝⎛+-+=
j j j j j j y y y y y y dx dy ，精度最高。

2、 二阶导数的处理
2
1
11
12
1
2
1222δδ
δ
δ
δ-+-+-++-=
--
-=
⎪
⎭⎫
⎝⎛-⎪⎭⎫ ⎝⎛=⎪⎭
⎫ ⎝⎛=j j j j j j j j j y y y y y y y dx dy dx dy dx dy dx d dx y d
二、有限宽轴承的五点差分近似解
1、 定常工况下（载荷和转速等均不随时间变化），静态等温油膜，方程（7）变为：
ϕλλϕϕ∂∂=⎪⎭⎫ ⎝⎛∂∂∂∂⎪⎭
⎫ ⎝⎛+⎪⎪⎭⎫ ⎝⎛∂∂∂∂H
p H L D p H 332
3 2、 剖分网格
沿ϕ方向的列数用i 编号，沿λ方向的行数用j 编号。

每个节点的位置用()j i ,二维编号表示，如图。

图4 五点差分法网格划分
设在ϕ方向共均匀划分为m 格，m i ~0=，步长m
1
2ϕϕϕ-=∆
在λ方向共均匀划分为n 格，n j ~0=，步长n
2=∆λ 3、 Reynolds 方程差分化
节点()j i ,上之p 为()j i p ,，对于节点()j i ,上的一阶导数
ϕ∂∂p 和λ
∂∂p 可用中差商表示： ϕϕ∆-≈⎪⎪⎭⎫
⎝⎛∂∂-+2,1,1,j i j i j
i p p p ， λλ∆-≈⎪
⎭⎫
⎝⎛∂∂+-21,1,,j i j i j
i p p p 半步长：ϕ
ϕ∆-≈
⎪⎪⎭⎫
⎝⎛∂∂-+j i j i j
i p
p
p ,2
1
,2
1
,
ϕ
ϕϕϕϕ∆⎪⎪⎭⎫ ⎝⎛∂∂-⎪⎪⎭⎫ ⎝⎛∂∂≈⎥⎦⎤
⎢⎣⎡⎪⎪⎭⎫ ⎝
⎛∂∂∂∂-+j
i j i j
i p H p H p H ,2
13,2
13,3
ϕϕ∆-≈⎪⎪⎭⎫ ⎝⎛∂∂+++j i j i j
i j i p p H p H ,,13
,2
1,2
13 ϕϕ∆-≈⎪⎪⎭⎫ ⎝⎛∂∂---j i j i j i j i p p H p H ,1,3
,2
1,2
13 于是： ()
2
,3
,213,21,13,21,13,21,3ϕϕϕ∆⎪⎪⎭⎫ ⎝⎛+-+≈⎥⎦⎤
⎢
⎣⎡⎪⎪⎭⎫ ⎝
⎛∂∂∂∂-+--++j i j i j i j i j i j i j i j i p H H p H p H p H ()
2
,3
21,321,1,3
21,1,3
21,,3λλλ∆⎪⎪⎭⎫ ⎝⎛+-+≈⎥
⎦⎤⎢⎣⎡⎪⎭⎫ ⎝
⎛∂∂∂∂-+--++j
i j i j i j i j i j i j i j i p H H p H p H p H 且ϕ
ϕ∆-≈⎪⎪⎭⎫
⎝⎛∂∂-+j i j i j
i H
H
H ,2
1,2
1,
方程变为：
()()()()()()
()()()()()
j i F j i p j i E j i p j i D j i p j i C j i p j i B j i p j i A ,,,1,,1,,,1,,1,=--+++-++
各系数为：
()()()()()()()()()()⎥
⎦⎤
⎢⎣⎡⎪⎭⎫ ⎝⎛--⎪⎭
⎫ ⎝⎛+∆=+++=⎪
⎭⎫ ⎝
⎛
-⎪⎭⎫ ⎝⎛∆∆=⎪
⎭⎫ ⎝⎛
+⎪⎭⎫ ⎝⎛∆∆=⎪
⎭⎫
⎝⎛-=⎪
⎭⎫
⎝⎛+=j i H j i H j i F j i D j i C j i B j i A j i E j i H L D j i D j i H L D j i C j i H j i B j i H j i A ,21,213,,,,,,21,,21,,,21,,21,3232
33ϕλϕλϕ
适用于全部内节点()1~1,1~1-=-=n j m i
4、 超松弛迭代求解压力()j i p ,
先将边界条件给定的值赋予各边界节点上的压力，并任意估取各个内节点的压力值，逐
次迭代，称为逐点超松弛迭代法（SOR ）
()()()()()()()()()()
()()()
j i p j i p j i E j i F j i p j i D j i p j i C j i p j i B j i p j i A j i p k k k k k k k
,,
,,1,,1,,,1,,1,,1111----+⎥⎦
⎤⎢⎣
⎡---+++-++=β一般取2~1=β 5、 收敛准则
()()
()
δ≤-∑∑∑∑-=-=-=-=-111
1
111
1
1
,,,n j m i k
n j m i k k
j i p j i p j i p ，一般取310-=δ
6、 对称问题的处理
许多轴承都具有λ方向上的左、右对称性，这时宜在λ方向上作偶数分格，即取n 为偶数。

图5 对称问题处理
迭代时只计算一半⎪⎭⎫ ⎝
⎛=-=2~
1,
1~1n j m i 加上一行，令⎪⎭
⎫ ⎝⎛+=⎪⎭⎫ ⎝⎛-12,12,n i p n i p 7、 Reynolds 边界条件的纳入
当轴瓦与轴颈间的间隙既有收敛区，又有开扩区时，油膜可能在开扩区内自然破裂。

常用方法是每行均由起始边向终止边方向逐点计算，如果算出某点压力为负，即取为零，该点之后各点压力均为零。

8、 油膜合力的计算 水平方向ξF ，垂直方向ηF 无量纲值：⎰
⎰ΦΦ-ΦΦ⎪⎭
⎫ ⎝⎛-=Ω=
2
1
sin 112
d d p rl
F F λμψξξ； ⎰⎰ΦΦ-ΦΦ⎪⎭⎫ ⎝⎛-=Ω=
21cos 1
12
d d p rl
F F λμψηη Φ是由轴承上方开始计量的角，θϕ+=Φ
用数值积分法时，先算出各列上p 沿轴向的积分()m i d p G i ~01
1
==⎰
-λ
用Simpson 积分法，设n 为偶数，则
()()()()()[]()m i n i p n i p i p i p i p G i ~0,1,42,21,40,3
1
=∆⋅+-++++≈
λ 再计算Φ方向积分，仍用Simpson 积分法，并设m 为偶数，则
()ϕξ∆⋅Φ+Φ++Φ+Φ+Φ-≈--m m m m G G G G G F sin sin 4sin 2sin 4sin 31
11221100
()ϕ
η∆⋅Φ+Φ++Φ+Φ+Φ-≈--m m m m G G G G G F cos cos 4cos 2cos 4cos 31
11221100 22ηξF F F +=∴
9、 偏位角的修正
对于圆柱轴承（单瓦）来说，ξF 和ηF 即为无量纲油膜合力。

通常载荷方向是给定的，例如垂直向下，则油膜合力方向应垂直向上。

γη
ξ≤F F ，一般取3104-⨯=γ
如果此条件不满足，则需修正偏位角θ，对于圆柱轴承，令⎪⎪⎭
⎫
⎝⎛η
ξθθF F arctg －＝老新
10、 摩擦阻力的计算
每块瓦的油膜对轴颈运动的阻力由压力流阻力2t F 和剪切阻力1t F 叠加而得。

1t F 又由油膜完整区和破裂区的剪切流阻力分别计算叠加得到。

⎰⎰⎰⎰-ΦΦ-ΦΦΦ+Φ=Ω=1
1
21111212121p p H d d H H d d rl F F p t t λλμψ
式中，p Φ表示破裂边的角位置，是λ的函数；p H 是破裂边上的膜厚，是λ的函数；
1Φ是起始角，2Φ是终止角。

()ηξξηe F e F r
F t -=
21
2 ()ηξξηεεμψF F rl F F t t -=Ω=
2
1
22 其中：c
e c
e ηηξξεε=
=
,
在圆柱轴承（单瓦）中：ηθ
F r
e F t 2sin 2=
21t t t F F F +=∴，阻力系数：
F
F f
t
=
ψ
11、 流量的计算
对于每一块瓦：进油一是从供油系统供入的新鲜油，一是从前一块瓦的下游流出的油。

出油一是从轴承两侧漏出，流回油池，一是进入后一块瓦。

（1） 进油量1Q （体积/时间）为：
⎰
-⎥⎦⎤
⎢⎣
⎡⎪⎪⎭⎫ ⎝⎛∂∂-Ω=22
03
101122l l dz p r h rh Q ϕμ 式中：0h 为上游边膜厚；0
⎪⎪⎭⎫
⎝⎛∂∂ϕp 为该处压力导数。

⎥⎦
⎤⎢⎣⎡⎪⎪⎭⎫ ⎝⎛∂∂-⎪⎭⎫ ⎝⎛=Ω⎪⎭⎫ ⎝⎛=⎰⎰--110301102
3113121λϕλψd p H d H d l r d l Q Q 式中：()()()n j i j p j p j p p j
~0,0,2,2,14,03,0==∆-+-≈⎪⎪⎭⎫
⎝⎛∂∂ϕϕ
（2） 下游边出油量
⎥⎦
⎤⎢⎣⎡⎪⎪⎭⎫ ⎝⎛∂∂-⎪⎭⎫ ⎝⎛=Ω⎪⎭⎫ ⎝⎛=⎰⎰--112321122
3223121λϕλψd p H d H d l r d l Q Q 式中：()()()n j m i j m p j m p j m p p j
~0,,2,2,14,3,2==∆-+--≈⎪⎪⎭⎫
⎝⎛∂∂ϕϕ 2H 和2
⎪⎪⎭⎫
⎝⎛∂∂ϕp 是下游边油膜厚度和压力导数。

若在轴瓦末边缘前已发生破裂，则：
λψd H d l r d l Q Q p ⎰-⎪⎭⎫ ⎝⎛=Ω⎪⎭⎫ ⎝⎛=11
2
32221 式中p H 为破裂边膜厚。

（3） 两侧泄漏量
⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝
⎛∂∂-⎪⎭⎫ ⎝⎛∂∂-=Ω⎪
⎭⎫ ⎝⎛=⎰⎰-==λλλλψϕϕλϕϕλd p H d p H r d l Q Q 2121131333361 若在2ϕ前已发生自然破裂，则积分上限改为p ϕ。

式中：()()()0,~0,22,1,40,30
,==∆+-≈⎪
⎭⎫
⎝⎛∂∂j m i i p i p i p p i λλ ()()()n j m i n i p n i p n i p p n
i ==∆---+-≈
⎪⎭⎫
⎝⎛∂∂,~0,22,1,4,3,λλ 若左右对称，则λλψϕϕλd p H r d l Q Q ⎰=⎪⎭⎫ ⎝⎛∂∂-=Ω⎪
⎭⎫
⎝⎛=21133
3331
四、五点差分法计算框图
五、提高计算精度的几种方法
1、 采用不均匀分格：膜厚小，压力变化急剧处分格细；膜厚大处分格粗。

()()()j i H a a a a j i H a a a a j i H a a a a H j i ,2,212,2122
12121122121,--⎪⎭⎫
⎝⎛-+-⎪⎭⎫ ⎝⎛++≈⎪⎪⎭⎫ ⎝⎛∂∂ϕ ()()()()j i p j i p j i p p H j i ,,1,12121,3ααααϕϕ+--++≈⎥⎦⎤
⎢⎣⎡⎪⎪⎭⎫ ⎝
⎛∂∂∂∂ 式中：()()
()
212
2321311,,212a a a j i H a a j i H a +--⎪⎭⎫
⎝⎛+=α ()()
()
212
1321322,,212a a a j i H a a j i H a +-+⎪⎭⎫
⎝⎛-=α ()()()()j i p j i p j i p p H j i ,1,1,2121,3ββββλλ+--++≈⎥⎦⎤
⎢
⎣⎡⎪⎭⎫ ⎝
⎛∂∂∂∂
式中：()()
()
212
2321311,21,2b b b j i H b b j i H b +--⎪⎭⎫ ⎝⎛
+=β ()()
()
212
1321322,21,2b b b j i H b b j i H b +-+⎪⎭⎫ ⎝⎛
-=β 方程
()()()()()()
()()()()()
j i F j i p j i E j i p j i D j i p j i C j i p j i B j i p j i A ,,,1,,1,,,1,,1,=--+++-++中各系数为：
()()()()()()()()()()()
()()()()⎥⎦
⎤
⎢⎣⎡--⎪⎭
⎫
⎝⎛-+-⎪⎭⎫ ⎝⎛++=+++=+⎪⎭⎫
⎝⎛++=⎪⎭
⎫
⎝⎛=⎪⎭⎫
⎝⎛===j i H a a a a j i H a a a a j i H a a a a j i F j i D j i C j i B j i A L D j i E L D j i D L D j i C j i B j i A ,
2,212,2123,,,,,,,,,,212121122121212
212
21
2
2
1ββααββαα
2、 变量代换
令p H 2
3
=∏，则∏的变化比p 平缓，Reynolds 方程变为：
ϕλλϕϕλϕ∂∂=∏⎥⎥⎦
⎤
⎢⎢⎣⎡∂∂⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛∂∂⎪⎭⎫ ⎝⎛+∂∂+⎪⎪⎭⎫ ⎝⎛∂∂-∂∏∂⎪⎭⎫ ⎝⎛+∂∏∂-H H H H L D H H L D H H H H L D 23
22222
22222223212143 边界条件为：
0=p 变为0=∏；
0=∂∂ϕ
p 变为0=∂∏
∂ϕ； in p p =变为in p H 2
3=∏
用差分法解出∏分布后，转化为p 分布，再计算承载力等特性数值。

六、膜厚阶梯变化处的差分方程
某些部位有膜厚的阶梯形变化，使膜厚在此方向上的导数不能取有限值。

可以用取控制体的方法来得到膜厚突变处的差分公式。

1、λ方向有膜厚阶梯变化，如图中A 点。

无量纲迭代式中之各系数为：
()()()()()()()()()()⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛---⎪⎭⎫ ⎝⎛+--⎪⎭⎫ ⎝⎛-++⎪⎭⎫ ⎝⎛++∆=
+++=⎪
⎭⎫ ⎝
⎛
-⎪⎭⎫ ⎝⎛∆∆=⎪
⎭⎫ ⎝⎛
+⎪⎭⎫ ⎝⎛∆∆=⎪⎭⎫ ⎝⎛--+⎪⎭⎫ ⎝⎛+-=
⎪⎭⎫ ⎝⎛-++⎪⎭⎫ ⎝⎛++=
0,210,210,210,2
123,,,,,,21,,21,,20,210,21,2
0,210,21,3232
3333j i H j i H j i H j i H j i F j i D j i C j i B j i A j i E j i H L D j i D j i H L D j i C j i H j i H j i B j i H j i H j i A ϕλϕλϕ
2、Φ方向有膜厚阶梯变化，如图中B 点。

无量纲迭代式中之各系数为：
()()()()()()()()()
()⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛--⎪⎭
⎫ ⎝⎛+∆=
+++=⎪⎭⎫ ⎝⎛--+⎪⎭⎫ ⎝⎛-+⎪
⎭
⎫ ⎝⎛∆∆=⎪
⎭⎫ ⎝⎛
+-+⎪⎭⎫ ⎝⎛++⎪
⎭
⎫ ⎝⎛∆∆=⎪
⎭
⎫
⎝⎛-=⎪
⎭⎫
⎝⎛+=j i H j i H j i F j i D j i C j i B j i A j i E j i H j i H L D j i D j i H j i H L D j i C j i H j i B j i H j i A ,21,2123,,,,,,221,021,0,221,021,0,,21,,21,332
332
33ϕλϕλϕ
3、λ与Φ方向均有膜厚阶梯变化，如图中C 点。

无量纲迭代式中之各系数为：
()()()()()()()()()()⎥⎦
⎤
⎢⎣⎡⎪⎭⎫ ⎝⎛--⎪⎭⎫ ⎝⎛-++⎪⎭⎫ ⎝⎛++∆=
+++=⎪
⎭⎫ ⎝
⎛
-⎪⎭⎫ ⎝⎛∆∆=⎪⎭⎫ ⎝⎛+-+⎪⎭⎫ ⎝⎛++⎪
⎭⎫ ⎝⎛∆∆=⎪
⎭
⎫
⎝⎛-=⎪⎭⎫ ⎝⎛-++⎪⎭⎫ ⎝⎛++=
j i H j i H j i H j i F j i D j i C j i B j i A j i E j i H L D j i D j i H j i H L D j i C j i H j i B j i H j i H j i A ,2120,210,2123,,,,,,21,,221,021,0,,21,2
0,210,21,32
332
333ϕλϕλϕ。