2024年4月28日发(作者：卞清昶)
ANALYSIS RESEARCH
虚拟弹簧边界Mindlin板的自由振动和屈曲分析
*
单逸峰 向 维
西南交通大学机械工程学院 成都 610031
摘 要：对于Mindlin板的自由振动和特征屈曲问题，为了用少量结点数求得高精度数值解，并直接求得完整的
模态向量，文中使用虚拟弹簧边界作为边界条件的施加方法，运用微分求积有限单元法推导出Mindlin板的自由
振动和屈曲特征方程。基于Matlab编程计算多种边界、不同厚度的Mindlin板的自振频率和屈曲因子，并讨论
计算方法的收敛性、准确性及计算效率。数值结果表明，虚拟弹簧边界与微分求积有限单元法相结合的数值方法，
计算结果准确，且计算效率有所提高。
关键词：微分求积有限单元法；虚拟弹簧边界；Mindlin板；自由振动；特征屈曲
中图分类号：
TH135
文献标识码：
A
文章编号：1001-0785（2022）1-0054-08
Abstract:
With regard to the free vibration and eigenvalue buckling of Mindlin plate, in order to obtain high-precision
solutions and complete modal vectors with a few nodes, the equations of free vibration and eigenvalue buckling of Mindlin
plate were derived by differential quadrature finite element method with virtual spring boundaries as conditions. According
to Matlab programming, the natural frequencies and buckling factors of Mindlin plates with various boundaries and different
thicknesses were calculated, and the convergence, accuracy and efficiency of the calculation method were discussed. Results
show that the combination of virtual spring boundary and differential quadrature finite element method is accurate and can
improve the calculation efficiency.
Keywords:
differential quadrature finite element method; virtual spring boundaries; Mindlin plate; free vibration; eigenvalue
buckling
0 引言
起重机的箱形梁结构，由上下翼缘板、主副腹板以
及横纵向加劲肋焊接组合而成
。准确求得板的自振频
率和屈曲载荷对起重机的安全工作至关重要。以有限元
法(以下简称FEM)为代表的数值解法易于求得高精度
近似解，然而，FEM需要把研究对象划分成大量微小
单元，使计算过程耗时更多。Xing Y等
[2]
提出的名为
微分求积有限单元法(以下简称DQFEM)的全局数值方
法在计算二维静态问题时表现出高精度和快速收敛性。
在应用数值算法时，Wang Q等
[3，4]
在对壳体的自
由振动分析中多次应用到虚拟弹簧边界。Kim K等
[5]
在
分析功能梯度矩形板的自由振动时，使用虚拟弹簧边界
引入边界条件，所得数值结果准确有效。Du Y等
[6]
在
*基金项目：国家自然科学基金委项目(51705436)
引 用 格 式
[1]
分析球型盖的动态行为时指出，虚拟弹簧边界易于模拟
不同的边界条件。由于采用虚拟弹簧边界不减少整体矩
阵的阶数，故完整的模态向量也可直接求得。Guan X等
[7]
基于经典板理论(以下简称CPT)，将虚拟弹簧边界和
DQFEM共同应用于求解薄板的自由振动问题。但是，
CPT未考虑横向剪切变形，使求得的自然频率偏大
[8]
。
针对板的自由振动和特征屈曲问题，基于考虑了横
向剪切变形的Mindlin板理论，采用DQFEM先后完成
虚拟弹簧边界Mindlin板的自由振动和屈曲分析。数值
结果表明，虚拟弹簧边界与DQFEM相结合的数值方法
所需结点数更少，求得的自振频率精度更高、屈曲因子
偏向于安全保守，且虚拟弹簧边界对DQFEM的计算效
率有所提高。
单逸峰，向维.虚拟弹簧边界Mindlin板的自由振动和屈曲分析[J].起重运输机械，2022（1）：54-61.
54
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
1
22
x
y
2
+
ρ
I
ϕρ
hw
+
ρ
I
ϕ
()
d
x
d
y
(2)
∫∫
A
2
1 虚拟弹簧边界的Mindlin板
Mindlin板理论能够较好地协调计算精度和效率，
因而广泛应用于模拟中厚板的力学行为。其位移场为
u
(
x
,
y
,
z
,
t
)
=−
z
ϕ
x
(
x
,
y
,
t
)
v
(
x
,
y
,
z
,
t
)
=−
z
ϕ
y
(
x
,
y
,
t
)
(1)
w
(
x
,
y
,
z
,
t
)
=
w
(
x
,
y
,
t
)
式中：
w
为沿
z
轴的移动，
φ
x
为
xOz
平面内的转动，
φ
y
为
yOz
平面内的转动。
在板的边界配置3种虚拟弹簧，用以约束边界的每
个广义自由度，如图1所示。弹簧标记
k
t
i
、
k
rx
i
和
k
ry
i
的
含义由表1列出。
=
T
2
1
2
x
y
x
x
2
2
y
x
y
x
D
U
s
d
x
d
y
22
2
A
1
yyy
x
y
2
x
1
x
y
w
2
w
2
22
xy
Gh
x
y
d
x
d
y
A
2
2
w
x
2
w
y
y
x
(3)
Eh
3
D
=
12
(
1
−
υ
2
)
G
=
E
2
(
1
+
υ
)
h
3
I=
12
式中：
D
为弯曲刚度；
G
为剪切模量；
E
为弹性
模量，
E
=10 920 Pa；
υ
为泊松比，
υ
=0.3；
ρ
为密度，
ρ
=1；
A
为板中面的面积；
I
为转动惯量；
к
为剪切修
图1 边界受虚拟弹簧约束的Mindlin板
正系数。
表1 弹簧标记及其含义
弹簧标记
k
ti
(i=1~4)
k
rxi
(i=1,2)
k
ryi
(i=1,2)
k
rx
i
(i
=
3
,
4)
k
ry
i
(i
=
3
,
4)
含义
第
i
条边上的线弹簧
边1和边3上的法向扭转弹簧
边2和边4上的法向扭转弹簧
边2和边4上的切向扭转弹簧
边1和边3上的切向扭转弹簧
z
k
ry2
&k
ry4
y
k
ry2
k
ry2
k
rx1
&k
rx4
k
rx1
k
rx1
&k
rx3
……
(i,N)
……
(M,N)(1,N)(2,N)
…
…
…
…
…
…
…
……
…
……
(M,j)
k
rx2
…
(1, j)
…
…
…
…
……
(M,1)
…
(1,1)(2,1)
…
……
(i,1)
…
k
ry4
k
t3
k
t3
&k
t4
k
ry2
&k
ry4
k
rx2
&k
rx4
k
ry1
k
ry1
&k
ry3
k
t1
&k
t2
k
t2
k
rx3
k
ry1
k
t2
k
rx2
&k
rx3
k
rx3
k
t2
&k
t3
k
ry1
&k
ry4
x
图2 划分Gauss-Lobatto结点的虚拟弹簧边界Mindlin矩形板
Gauss-Lobatto结点已被证明是能有效保证计算精
度的网格划分方法
[9]
，将Mindlin板分别沿
x
轴和
y
轴
划分出
M
个和
N
个Gauss-Lobatto结点。虚拟弹簧边界
离散在各个Gauss-Lobatto结点上，为展示清晰，图2
仅完整表示边2和边3的全部边界弹簧。
边界上，线弹簧的弹性势能
V
t
，扭转弹簧在
xOz
面
内和
yOz
面内的弹性势能
V
rx
、
V
ry
分别为
1
4
l
i
2
Vkwl
d
=
t
2
∑
∫
0
t
iii
i
=
1
1
4
l
i
2
V
rx
=
∑
∫
0
k
rx
i
ϕ
x
i
d
l
i
(4)
2
i
=
1
1
4
l
i
2
V
ry
=
∑
∫
0
k
ry
i
ϕ
y
i
d
l
i
2
i
=
1
式中：
i
为板边编号，
l
i
为第
i
条边的长度。
由式(2)～式(4)可得Mindlin板的Lagrange函数为
1.1 自由振动特征方程
要推导自由振动特征方程，首先要列出系统的能量。
Mindlin板的动能
T
和应变能
U
s
分别为
2022年第1期 /
55
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
L=T−U,(U=U
s
+V
t
+V
rx
+V
ry
)
(5)
式中：
K
为系统的刚度矩阵，
M
为质量矩阵。
w
是由Mindlin板上所有结点的位移所组成的列向
量，即
T
w=
φ
x
T
φ
y
势能
U
中含有充当边界条件的虚拟边界弹簧的弹
性势能，由此推导出的运动方程中便已施加过边界条件，
而无需修改整体矩阵，不同于传统数值方法。
基于DQFEM离散Lagrange函数，由Hamilton变
分原理，令
δ
L=0
，可得Mindlin板的振动特征方程为
2
K
M
w0
(6)
w
T
T
φ
x11
xM1
x1N
xMN
T
x
T
φ
y
(9)
y11
yM1
y1N
yMN
w
T
w
11
w
M1
w
1N
w
MN
KK
s
K
springs
(7)
M=diag
IC
IC
hC
(8)
刚度矩阵
K
中的
K
s
来自于板的应变能，
K
springs
来自
于虚拟弹簧边界的弹性势能。
D
1
A
(1)T(1)(1)T(1)
K
s
＝
D
BCAACB
1
s
A
(1)T
C
(1)T
CB
(1)
1
B
(1)T
CA
(1)
D
2
s
B
(1)T
C
s
CA
(1)
s
CB
(1)
D
3
(10)
4
K
springs
=diag
CK
rx
i
i
1
CK
i
1
4
ry
i
CK
t
i
(11)
i
1
4
式中：
υ
1
=(1-
υ
)/2；
υ
s
=6
κ
(1-
υ
)/
h
2
；
A
(1)
和
B
(1)
分别
为二元函数沿
x
轴和
y
轴方向的一阶权系数矩阵；
C
是
二维Gauss-Lobatto积分系数矩阵；
D
1
、
D
2
、
D
3
由式(12)
求得，即
(1)T(1)
D
1
1
1
s
ACA
D
＝
1
B
(1)T
CB
(1)
(12)
s
2
1
C
D
3
s
s
0
N
M
1
K
rx1
＝
diag(
K
rx1
,
K
rx1
,
,
K
rx1
),
K
rx1
＝
diag(
k
rx1
,0
,
,
0)
N
M
1
diag(
,,,
),diag(0
,,
0
,k
rx2
)
＝＝
KKKKK
rx2rx2rx2rx2rx2
(
N
1)
M
M
(14)
K
r
y1
＝
diag(
k
ry1
,k
ry1
,
,k
ry1
,
0
,
,
0)
(
N
1)
M
M
K
ry2
＝
diag(0
,
,
0
,
k
ry2
,k
ry2
,
,k
ry2
)
K
springs
中的各个分量分别为
N
M
1
K
t1
＝
diag(
K
t1
,
K
t1
,
,
K
t1
),
K
t1
diag(
k
t1
,0
,
,
0)
(
N
1)
M
M
(13)
K
t2
＝
diag(
k
t2
,k
t2
,
,k
t2
,
0
,
,
0)
N
M
1
K
t3
＝
diag(
K
t3
,
K
t3
,
,
K
t3
),
K
t3
diag(0
,
,
0
,k
t3
)
(
N
1)
M
M
K
t4
＝
diag(0
,
,
0,
k
t4
,k
t4
,
,k
t4
)
N
M
1
K
ry3
＝
diag(
K
ry3
,
K
ry3
,
,
K
ry3
),
K
ry3
＝
diag(
k
ry3
,0,
,0)
N
M
1
K
＝＝
diag(,,,),diag(0,,0,
k
ry4
)
KKKK
ry4ry4ry4ry4ry4
(
N
1)
M
M
(15)
K
rx3
＝
diag(
k
rx3
,
k
rx3
,
,
k
rx3
,0,
,0)
(
N
1)
M
M
K
rx4
＝
diag(0,
,0
,
k
rx4
,
k
rx4
,
,
k
rx4
)
为便于同前人文献的数据对比，将式(6)得出的自
由振动频率
ω
代入式(16)，求得无量纲频率参数
ω
为
ω
=
ω
a
ρ
(16)
G
56
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
1.2 屈曲特征方程
虚拟弹簧边界Mindlin矩形板在图3所示的均匀初
始应力作用下，系统的总能量为
T
U
0
U
σ
0
ε
L
d
V
(17)
V
0
T
2 数值方法参数的确定
为使数值方法得到准确的数值结果，需要配置足够
的结点数以使数值解收敛；要使虚拟弹簧边界代替经典
边界条件，需要确定各个虚拟弹簧的刚度取值。
σ
0
x
0
0
xyy
(18)
2.1 结点数
为确定使数值解收敛的最少结点数，针对CCCC和
SSSS边界，厚度比
h
/
b
为0.01和0.1的4种Mindlin方板，
图4为结点数对第1、3、5阶频率参数的
ω
影响情况。
为便于计算，取长和宽方向结点数相同，即
M
=
N
。当模拟CCCC边界时，令所有弹簧刚度为10
10
N/m，
κ
=0.860 1；当模拟SSSS边界时，令
k
x1
=
k
x2
=
k
y1
=
k
y2
=0，
κ
=5/6。图中虚线分别为第1、3、
其他弹簧刚度为10
10
N/m，
5阶自振频率的收敛值，箭头指向前6阶频率收敛时的
最少单边结点数。
由图4可得，使用虚拟弹簧边界代替经典边界条件，
DQFEM仍有良好的收敛性；为保证前6阶无量纲频率
参数收敛，应配置的结点数为12×12。
u
2
v
2
w
2
x
x
x
222
1
u
v
w
L
(19)
ε
=
2
y
y
y
u
u
v
v
w
w
2
x
y
x
y
x
y
式中：
σ
为薄膜力矩阵，
ε
为应变元矩阵
类似地，基于DQFEM离散总能量
U
0
，并令
0
L
δ
U
0
=0
，可得板的屈曲特征方程
K
K
G1
K
G2
K
G3
w
0
(20)
0
x
0
y
0
xy
式中：
K
G1
、
K
G2
、
K
G3
分别为Mindlin板在受
x
轴方
向的单轴压力、受
y
轴方向的单轴压力、受
xOy
面内剪
力时的几何刚度矩阵。
(1)T(1)
K
G1
diag
I
ACA
I
A
(1)T
CA
(1)
I
B
(1)T
CB
(1)
h
A
(1)T
CA
(1)
(21)
h
B
(1)T
CB
(1)
I
B
(1)T
CB
(1)
K
G2
diag
(22)
K
G3
I
A
(1)T
CB
(1)
B
(1)T
CA
(1)
diag
I
A
(1)T
CB
(1)
B
(1)T
CA
(1)
h
A
(1)T
CB
(1)
B
(1)T
CA
(1)
(23)
(a) 模拟CCCC边界，
h
/
b
=0.01 (b) 模拟CCCC边界，
h
/
b
=0.1
为便于数据对比，以
x
轴向受压为例，将临界屈曲
应力
σ
x
代入式(24)，求得屈曲因子
n
*
为
0
hb
2
x
(24)
n
D
π
2
*
0
y
σ
y
0
τ
xy
0
(c) 模拟SSSS边界，
h
/
b
=0.01 (d) 模拟SSSS边界，
h
/
b
=0.1
b
a
σ
x
0
图4 结点数对频率参数的影响
x
O
图3 有面内应力作用的矩形板
2.2 虚拟边界弹簧的刚度取值
为分析不同类型弹簧在模拟经典边界条件时应取的
2022年第1期 /
57
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
刚度值，将弹簧分为3类：法向扭转弹簧
k
rx
i
(
i
=1，2)与
k
ry
i
(
i
=1，2)；切向扭转弹簧
k
rx
i
(
i
=3，4)与
k
ry
i
(
i
=3，4)；
线弹簧
k
t
i
(
i
=1，2，3，4)。针对线弹簧和扭转弹簧，引
入无量纲刚度
k
t
b
3
S
w
D
(25)
S
k
r
b
φ
D
3 自由振动分析数值结果
针对Mindlin板自由振动问题，本小节首先验证数
值结果的准确性，而后讨论虚拟弹簧边界对DQFEM计
算效率的影响。采用12×12的结点划分Mindlin板，边
界弹簧的刚度设置为：固支边(C)上，
S
w
=
S
φ1
=
S
φ2
=10
10
N/m；简支边(S)上，
S
w
=
S
φ2
=10
10
，
S
φ1
=0；自由边(F)上，
S
w
=
S
φ1
=
S
φ2
=0。
用
S
φ1
和
S
φ
2
分别表示法向和切向扭转弹簧的无量
纲刚度。
图5为3类弹簧的刚度对厚度比
h
/
b
为0.01和0.1
的Mindlin方板的第1、3、5阶无量纲频率参数的影响情况。
由图5可得，不论哪种边界弹簧，当其刚度大于某一定值，
频率参数不再随弹簧刚度的增大而变化。这可以解释为
弹簧的刚度太大，以至于板边界被近似严格约束。对于
法向、切向扭转弹簧以及线弹簧，实现这种严格约束的
无量纲刚度分别为10
7
N/m、10
8
N/m和10
10
N/m。
3.1 准确性
为验证本文数值解的准确性，以SSSS、SCSC和
CCCC边界条件的厚度比
h
/
b
为0.01和0.1的Mindlin
板为对象，计算频率参数
ω
。表2~表4将本文数值解
同Hinton
[10]
的精确解、Ferreira
[11]
的有限元解、Dawe D
J等
[12]
的Rayleygh-Ritz解进行对比。Ferreira的角标1
和2分别为15×15和20×20的Q4有限元法。
由表2和表3可知，在SSSS和SCSC边界条件
下得出的前5阶无量纲频率
ω
比FEM求得的数值解
更接近精确解。由表4可知，在CCCC边界条件下求
得的前5阶无量纲频率
ω
比FEM求得的数值解更接近
Rayleygh-Ritz法的数值结果。
由此，使用虚拟弹簧边界替代经典边界条件后，采
用DQFEM求解Mindlin板的自由振动问题的准确性得
(a)
S
φ
2
=
S
w
=102，
h
/
b
=0.01 (b)
S
φ
2
=
S
w
=102，
h
/
b
=0.1
以验证。值得注意的是，DQFEM仅需划分12×12的
结点数，而求得的频率参数比20×20的Q4有限元法
所得结果更加准确。
3.2 计算效率
由于边界条件施加方法的改变，DQFEM的计算效
率可能改变。表5以厚度比
h
/
b
=0.1的Mindlin方板为
(c)
S
φ
1
=
S
w
=102，
h
/
b
=0.01 (d)
S
φ
1
=
S
w
=102，
h
/
b
=0.1
研究对象，分别统计了对整体矩阵消元和使用虚拟弹簧
施加SSSS和CCCC边界条件时DQFEM求得前6阶无
量纲频率及完整的模态向量所花费的时间。
由表5可知，使用虚拟弹簧边界施加边界条件有助
于减少DQFEM的计算时间。对于SSSS边界，消元法
和虚拟弹簧边界的10个算例的平均时间分别为0.422 0
(e)
S
φ
1
=
S
φ
2
==102，
h
/
b
=0.01 (f)
S
φ
1
=
S
φ
2
==102，
h
/
b
=0.1
图5 频率参数随3类弹簧刚度的变化情况
s和0.412 2 s，虚拟弹簧边界使DQFEM的平均计算时
间缩短了2.32%；对于CCCC边界，消元法和虚拟弹簧
58
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
表2 SSSS边界，
κ
=0.833，Mindlin方板前5阶无量纲频率
ω
h/b
计算方法
Present
0.01Hinton
Ferreira
2
Present
0.1Hinton
Ferreira
1
频率阶次
1
0.096 3
0.096 3
0.096 5
0.930 3
0.930
0.934 6
2
0.240 6
0.240 6
0.243 0
2.219 3
2.219
2.254 5
3
0.240 6
0.240 6
0.243 0
2.219 3
2.219
2.254 5
4
0.384 7
0.384 8
0.389 0
3.405 6
3.406
3.459 2
5
0.480 7
0.480 9
0.492 8
4.149 4
4.149
4.303 1
表3 SCSC边界，
κ
=0.822，Mindlin方板前5阶无量纲频
ω
h/b
计算方法
Present
0.01Hinton
Ferreira
1
Present
0.1Hinton
Ferreira
1
频率阶次
1
0.141 1
0.141 1
0.142 4
1.300 1
1.302
1.294 0
2
0.266 8
0.266 8
0.271 0
2.393 9
2.398
2.397 1
3
0.337 6
0.337 7
0.348 4
2.884 5
2.888
2.929 0
4
0.460 4
0.460 8
0.472 2
3.839 1
3.852
3.839 4
5
0.497 7
0.497 9
0.519 1
4.231 3
4.237
4.347 5
表4 CCCC边界，
κ
=0.860 1，Mindlin方板前5阶无量纲频
ω
h/b
计算方法
Present
0.01Ferreira
2
Dawe
Present
0.1Ferreira
2
Dawe
频率阶次
1
0.175 4
0.175 0
0.175 4
1.591 0
1.595 5
1.594 0
2
0.357 4
0.363 5
0.357 6
3.038 9
3.066 2
3.039 0
3
0.357 4
0.363 5
0.357 6
3.038 9
3.066 2
3.039 0
4
0.526 5
0.535 8
0.527 4
4.262 5
4.292 4
4.265 0
5
0.639 9
0.663 4
0.640 2
5.024 7
5.123 2
5.035 0
边界的10个算例的平均时间分别为0.314 9 s和0.305
2 s，虚拟弹簧边界使DQFEM的平均计算时间缩短了
3.08%。
和受双轴向均匀载荷时的临界屈曲因子
n
*
，同时列出
Mizusawa T
[13]
的数值解、Hashemi S等
[14]
的精确解。
剪切修正因子
κ
取
π
2
/12。数值结果与Hashemi S等精确
解的相对误差由式(26)求得，即
4 屈曲分析数值结果
表6和表7列出SSSS、SCSF和SFSF边界厚度比
h
/
b
为0.05和0.1的Mindlin方板，在分别受
x
轴方向
Present-Exact
RelError=
100%
Exact
(26)
2022年第1期 /
59
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
表5 消元法和虚拟弹簧边界DQFEM的计算时间 s
样本
1
2
3
4
5
6
7
8
9
10
平均时间
CCCC
消元法
0.421 8
0.415 6
0.414 3
0.426 1
0.414 0
0.428 1
0.425 1
0.431 3
0.424 7
0.419 2
0.422 0
虚拟弹簧
0.414 3
0.415 0
0.408 9
0.410 5
0.414 2
0.417 7
0.407 5
0.415 6
0.406 4
0.411 4
0.412 2
消元法
0.321 1
0.310 9
0.316 7
0.315 5
0.310 0
0.319 4
0.305 9
0.317 1
0.316 6
0.315 3
0.314 9
SSSS
虚拟弹簧
0.301 9
0.311 2
0.306 2
0.303 2
0.301 0
0.300 8
0.306 2
0.306 7
0.308 3
0.306 5
0.305 2
表7 双轴向受压Mindlin方板的无量纲屈曲因子
n
*
边界条件计算方法
0.05
Present
SSSSHashemi
RelError
Present
SCSFHashemi
RelError
Present
SFSFHashemi
RelError
1.964 0
1.971 8
0.40%
1.108 9
1.111 9
0.27%
0.9187
0.9207
0.22%
h/b
0.1
1.864 5
1.891 9
1.45%
1.052 4
1.064 1
1.10%
0.890 1
0.897 7
0.85%
5 结论
使用虚拟弹簧边界能够很好得模拟板的经典边界条
件而不影响DQFEM的收敛性。结合了虚拟弹簧边界的
DQFEM，与FEM相比，使用更少的结点数即可求得更
高精度的数值解；求得的屈曲因子更加安全保守。同时，
采用虚拟弹簧边界有助于提高DQFEM的计算效率。
对比表6和表7数据可知，将DQFEM和虚拟弹簧
边界共同用于Mindlin板的屈曲分析时，得出的屈曲因
子与Mizusawa通过样条方法得出的数值解几乎完全吻
合，而略小于Hashemi的精确解，但最大相对误差不超
过1.5%。这说明使用虚拟弹簧边界模拟经典边界条件，
并使用DQFEM分析Mindlin板的特征屈曲问题，具有
一定准确性，求得的临界屈曲因子偏向于安全保守。
表6
x
轴方向受压Mindlin方板的无量纲屈曲因子
n
边界条件计算方法
0.05
Present
SSSS
Mizusawa
Hashemi
RelError
Present
SCSF
Mizusawa
Hashemi
RelError
Present
SFSF
Mizusawa
Hashemi
RelError
3.928 0
3.928
3.943 7
0.40%
1.615 3
1.615
1.619 7
0.27%
0.941 2
0.941 2
0.943 1
0.20%
h/b
0.1
3.728 8
3.729
3.783 8
1.45%
1.539 4
1.539
1.555 8
1.05%
0.914 6
0.914 6
0.921 9
0.79%
*
参考文献
[1] 漆静，程文明，濮德璋，等.起重机箱形梁曲形腹板简支
板局部稳定性分析[J].机械设计与制造,2014(4)：
162-164.
[2]
Xing Y
，
Liu B
，
Liu G.A differential quadrature finite
element method
[J].
International Journal of Applied
Mechanics
，2010，2(1)：207-227.
[3]
Wang Q
，
Shi D
，
Qian L
，
et al. Free vibrations of composite
laminated doubly-curved shells and panels of revolution with
general elastic restraints
[J].
Applied Mathematical
Modelling,
2017，46(1)：227-262.
[4]
Wang Q
，
Shao D
，
Qin B.A simple first-order shear
deformation shell theory for vibration analysis of composite
laminated open cylindrical shells with general boundary
conditions
[J].
Composite Structures
，2017，184：
211-232.
60
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
[5]
Kim K
，
Kim K
，
Han C
，
et al. A method for natural
frequency calculation of the functionally graded rectangular
plate with general elastic restraints
[J].
AIP Advances
，2020.
[6]
Du Y
，
Pang F
，
Sun L
，
et al. A unified formulation for
dynamic behavior analysis of spherical cap with uniform and
stepped thickness distribution under different edge constraints
[J].
Thin-Walled Structures
，2019.
[7]
Guan X
，
Tang J
，
Wang Q
，
et al. Application of the differential
quadrature finite element method to free vibration of elastically
restrained plate with irregular geometries
[J].
Engineering
Analysis with Boundary Elements
，2018，90：1-16.
[8] 向维.剪切板自由振动和特征屈曲新解法[D].北京：
北京航空航天大学，2016.
[9] 邢誉峰，李敏.计算固体力学原理与方法[M].北京：
北京航空航天大学出版社，2011.
[10]
Hinton cal methods and software for dynamic
作 者：单逸峰
电子邮箱：****************
收稿日期：2021-06-11
analysis of plates and shells
[J].
Pineridge
，1988.
[11]
Ferreira codes for finite element analysis
[M].
Springer Netherlands
，2009.
[12]
Dawe D J
，
Roufaeil O L. Rayleigh-Ritz vibration analysis
of Mindlin plates
[J].
Journal of Sound & Vibration
，1980，
69(3)：345-359.
[13]
Mizusawa ion of rectangular Mindlin plates with
tapered thickness by the spline strip method
[J].
Computers
and Structures
，1993，30(12)：1 663-1 677.
[14]
Hashemi S
，
Khorshidi K
，
Amabili M. Exact solution for
linear buckling of rectangular Mindlin plates
[J].
Journal
of Sound and Vibration
，2008，315(1，2)：318-342.
三一增程式电动正面吊运机
三一增程式电动正面吊运机广泛应用于各大港口、码头、铁路堆场，以适应性强、不畏严寒酷暑、
科技增效、维修简易、服务及时、数据互联等特点助力客户成为绿色低碳物流的先行者。
1）绿色降本，适应性强：三一增程式电动正面吊运机可以闲时纯电驱动，忙时增程式驱动，一小时
即可满充，作业可自由切换，丝毫不担心电池续航不足出现停机。根据作业数据，三一电动正面吊运机
纯电作业能耗成本较燃油车下降68%，增程器作业能耗成本下降40%，综合作业一年可节省能耗成本超过
25万元，减少二氧化碳排放近200 t/a。
2）巧妙设计，不畏严寒酷暑：在夏季高温时，为了确保电动产品的安全性，该吊运机采用磷酸铁锂
电池系统，热稳定性好，800℃才会分解但不会自燃，可确保高温环境中的人机安全。由于北方冬季寒冷，
针对低温工况增加了隔温保护套，减少外界低温环境对电池的影响，保障冬季电池续航。
3）科技加持，效能提升：针对行驶系统能耗高的反馈，新一代电动正面吊运机由双电动机改为单电
动机直驱，效率优于双电动机。新一代控制系统的物理按键减少，采用新型手柄，设备的动作响应性和
操作实时性随之提高，且具备智能自诊断(IO)，运行周期比上一代提升10倍。第三代控制系统，10寸
大触摸屏，真彩显示，智能控制。
4）维保简易 服务及时：通过设计不断提高设备维保的便捷程度，缩短检修时间。
5）大屏管理 数据互联：为提高产品附加价值，三一开发了C端互联系统，向客户提供电动产品大
数据支持，电量、电耗、总耗电量以及作业量等实时显示，帮助客户更好地安排作业。
2022年第1期 /
61
Copyright©博看网. All Rights Reserved.
2024年4月28日发(作者：卞清昶)
ANALYSIS RESEARCH
虚拟弹簧边界Mindlin板的自由振动和屈曲分析
*
单逸峰 向 维
西南交通大学机械工程学院 成都 610031
摘 要：对于Mindlin板的自由振动和特征屈曲问题，为了用少量结点数求得高精度数值解，并直接求得完整的
模态向量，文中使用虚拟弹簧边界作为边界条件的施加方法，运用微分求积有限单元法推导出Mindlin板的自由
振动和屈曲特征方程。基于Matlab编程计算多种边界、不同厚度的Mindlin板的自振频率和屈曲因子，并讨论
计算方法的收敛性、准确性及计算效率。数值结果表明，虚拟弹簧边界与微分求积有限单元法相结合的数值方法，
计算结果准确，且计算效率有所提高。
关键词：微分求积有限单元法；虚拟弹簧边界；Mindlin板；自由振动；特征屈曲
中图分类号：
TH135
文献标识码：
A
文章编号：1001-0785（2022）1-0054-08
Abstract:
With regard to the free vibration and eigenvalue buckling of Mindlin plate, in order to obtain high-precision
solutions and complete modal vectors with a few nodes, the equations of free vibration and eigenvalue buckling of Mindlin
plate were derived by differential quadrature finite element method with virtual spring boundaries as conditions. According
to Matlab programming, the natural frequencies and buckling factors of Mindlin plates with various boundaries and different
thicknesses were calculated, and the convergence, accuracy and efficiency of the calculation method were discussed. Results
show that the combination of virtual spring boundary and differential quadrature finite element method is accurate and can
improve the calculation efficiency.
Keywords:
differential quadrature finite element method; virtual spring boundaries; Mindlin plate; free vibration; eigenvalue
buckling
0 引言
起重机的箱形梁结构，由上下翼缘板、主副腹板以
及横纵向加劲肋焊接组合而成
。准确求得板的自振频
率和屈曲载荷对起重机的安全工作至关重要。以有限元
法(以下简称FEM)为代表的数值解法易于求得高精度
近似解，然而，FEM需要把研究对象划分成大量微小
单元，使计算过程耗时更多。Xing Y等
[2]
提出的名为
微分求积有限单元法(以下简称DQFEM)的全局数值方
法在计算二维静态问题时表现出高精度和快速收敛性。
在应用数值算法时，Wang Q等
[3，4]
在对壳体的自
由振动分析中多次应用到虚拟弹簧边界。Kim K等
[5]
在
分析功能梯度矩形板的自由振动时，使用虚拟弹簧边界
引入边界条件，所得数值结果准确有效。Du Y等
[6]
在
*基金项目：国家自然科学基金委项目(51705436)
引 用 格 式
[1]
分析球型盖的动态行为时指出，虚拟弹簧边界易于模拟
不同的边界条件。由于采用虚拟弹簧边界不减少整体矩
阵的阶数，故完整的模态向量也可直接求得。Guan X等
[7]
基于经典板理论(以下简称CPT)，将虚拟弹簧边界和
DQFEM共同应用于求解薄板的自由振动问题。但是，
CPT未考虑横向剪切变形，使求得的自然频率偏大
[8]
。
针对板的自由振动和特征屈曲问题，基于考虑了横
向剪切变形的Mindlin板理论，采用DQFEM先后完成
虚拟弹簧边界Mindlin板的自由振动和屈曲分析。数值
结果表明，虚拟弹簧边界与DQFEM相结合的数值方法
所需结点数更少，求得的自振频率精度更高、屈曲因子
偏向于安全保守，且虚拟弹簧边界对DQFEM的计算效
率有所提高。
单逸峰，向维.虚拟弹簧边界Mindlin板的自由振动和屈曲分析[J].起重运输机械，2022（1）：54-61.
54
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
1
22
x
y
2
+
ρ
I
ϕρ
hw
+
ρ
I
ϕ
()
d
x
d
y
(2)
∫∫
A
2
1 虚拟弹簧边界的Mindlin板
Mindlin板理论能够较好地协调计算精度和效率，
因而广泛应用于模拟中厚板的力学行为。其位移场为
u
(
x
,
y
,
z
,
t
)
=−
z
ϕ
x
(
x
,
y
,
t
)
v
(
x
,
y
,
z
,
t
)
=−
z
ϕ
y
(
x
,
y
,
t
)
(1)
w
(
x
,
y
,
z
,
t
)
=
w
(
x
,
y
,
t
)
式中：
w
为沿
z
轴的移动，
φ
x
为
xOz
平面内的转动，
φ
y
为
yOz
平面内的转动。
在板的边界配置3种虚拟弹簧，用以约束边界的每
个广义自由度，如图1所示。弹簧标记
k
t
i
、
k
rx
i
和
k
ry
i
的
含义由表1列出。
=
T
2
1
2
x
y
x
x
2
2
y
x
y
x
D
U
s
d
x
d
y
22
2
A
1
yyy
x
y
2
x
1
x
y
w
2
w
2
22
xy
Gh
x
y
d
x
d
y
A
2
2
w
x
2
w
y
y
x
(3)
Eh
3
D
=
12
(
1
−
υ
2
)
G
=
E
2
(
1
+
υ
)
h
3
I=
12
式中：
D
为弯曲刚度；
G
为剪切模量；
E
为弹性
模量，
E
=10 920 Pa；
υ
为泊松比，
υ
=0.3；
ρ
为密度，
ρ
=1；
A
为板中面的面积；
I
为转动惯量；
к
为剪切修
图1 边界受虚拟弹簧约束的Mindlin板
正系数。
表1 弹簧标记及其含义
弹簧标记
k
ti
(i=1~4)
k
rxi
(i=1,2)
k
ryi
(i=1,2)
k
rx
i
(i
=
3
,
4)
k
ry
i
(i
=
3
,
4)
含义
第
i
条边上的线弹簧
边1和边3上的法向扭转弹簧
边2和边4上的法向扭转弹簧
边2和边4上的切向扭转弹簧
边1和边3上的切向扭转弹簧
z
k
ry2
&k
ry4
y
k
ry2
k
ry2
k
rx1
&k
rx4
k
rx1
k
rx1
&k
rx3
……
(i,N)
……
(M,N)(1,N)(2,N)
…
…
…
…
…
…
…
……
…
……
(M,j)
k
rx2
…
(1, j)
…
…
…
…
……
(M,1)
…
(1,1)(2,1)
…
……
(i,1)
…
k
ry4
k
t3
k
t3
&k
t4
k
ry2
&k
ry4
k
rx2
&k
rx4
k
ry1
k
ry1
&k
ry3
k
t1
&k
t2
k
t2
k
rx3
k
ry1
k
t2
k
rx2
&k
rx3
k
rx3
k
t2
&k
t3
k
ry1
&k
ry4
x
图2 划分Gauss-Lobatto结点的虚拟弹簧边界Mindlin矩形板
Gauss-Lobatto结点已被证明是能有效保证计算精
度的网格划分方法
[9]
，将Mindlin板分别沿
x
轴和
y
轴
划分出
M
个和
N
个Gauss-Lobatto结点。虚拟弹簧边界
离散在各个Gauss-Lobatto结点上，为展示清晰，图2
仅完整表示边2和边3的全部边界弹簧。
边界上，线弹簧的弹性势能
V
t
，扭转弹簧在
xOz
面
内和
yOz
面内的弹性势能
V
rx
、
V
ry
分别为
1
4
l
i
2
Vkwl
d
=
t
2
∑
∫
0
t
iii
i
=
1
1
4
l
i
2
V
rx
=
∑
∫
0
k
rx
i
ϕ
x
i
d
l
i
(4)
2
i
=
1
1
4
l
i
2
V
ry
=
∑
∫
0
k
ry
i
ϕ
y
i
d
l
i
2
i
=
1
式中：
i
为板边编号，
l
i
为第
i
条边的长度。
由式(2)～式(4)可得Mindlin板的Lagrange函数为
1.1 自由振动特征方程
要推导自由振动特征方程，首先要列出系统的能量。
Mindlin板的动能
T
和应变能
U
s
分别为
2022年第1期 /
55
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
L=T−U,(U=U
s
+V
t
+V
rx
+V
ry
)
(5)
式中：
K
为系统的刚度矩阵，
M
为质量矩阵。
w
是由Mindlin板上所有结点的位移所组成的列向
量，即
T
w=
φ
x
T
φ
y
势能
U
中含有充当边界条件的虚拟边界弹簧的弹
性势能，由此推导出的运动方程中便已施加过边界条件，
而无需修改整体矩阵，不同于传统数值方法。
基于DQFEM离散Lagrange函数，由Hamilton变
分原理，令
δ
L=0
，可得Mindlin板的振动特征方程为
2
K
M
w0
(6)
w
T
T
φ
x11
xM1
x1N
xMN
T
x
T
φ
y
(9)
y11
yM1
y1N
yMN
w
T
w
11
w
M1
w
1N
w
MN
KK
s
K
springs
(7)
M=diag
IC
IC
hC
(8)
刚度矩阵
K
中的
K
s
来自于板的应变能，
K
springs
来自
于虚拟弹簧边界的弹性势能。
D
1
A
(1)T(1)(1)T(1)
K
s
＝
D
BCAACB
1
s
A
(1)T
C
(1)T
CB
(1)
1
B
(1)T
CA
(1)
D
2
s
B
(1)T
C
s
CA
(1)
s
CB
(1)
D
3
(10)
4
K
springs
=diag
CK
rx
i
i
1
CK
i
1
4
ry
i
CK
t
i
(11)
i
1
4
式中：
υ
1
=(1-
υ
)/2；
υ
s
=6
κ
(1-
υ
)/
h
2
；
A
(1)
和
B
(1)
分别
为二元函数沿
x
轴和
y
轴方向的一阶权系数矩阵；
C
是
二维Gauss-Lobatto积分系数矩阵；
D
1
、
D
2
、
D
3
由式(12)
求得，即
(1)T(1)
D
1
1
1
s
ACA
D
＝
1
B
(1)T
CB
(1)
(12)
s
2
1
C
D
3
s
s
0
N
M
1
K
rx1
＝
diag(
K
rx1
,
K
rx1
,
,
K
rx1
),
K
rx1
＝
diag(
k
rx1
,0
,
,
0)
N
M
1
diag(
,,,
),diag(0
,,
0
,k
rx2
)
＝＝
KKKKK
rx2rx2rx2rx2rx2
(
N
1)
M
M
(14)
K
r
y1
＝
diag(
k
ry1
,k
ry1
,
,k
ry1
,
0
,
,
0)
(
N
1)
M
M
K
ry2
＝
diag(0
,
,
0
,
k
ry2
,k
ry2
,
,k
ry2
)
K
springs
中的各个分量分别为
N
M
1
K
t1
＝
diag(
K
t1
,
K
t1
,
,
K
t1
),
K
t1
diag(
k
t1
,0
,
,
0)
(
N
1)
M
M
(13)
K
t2
＝
diag(
k
t2
,k
t2
,
,k
t2
,
0
,
,
0)
N
M
1
K
t3
＝
diag(
K
t3
,
K
t3
,
,
K
t3
),
K
t3
diag(0
,
,
0
,k
t3
)
(
N
1)
M
M
K
t4
＝
diag(0
,
,
0,
k
t4
,k
t4
,
,k
t4
)
N
M
1
K
ry3
＝
diag(
K
ry3
,
K
ry3
,
,
K
ry3
),
K
ry3
＝
diag(
k
ry3
,0,
,0)
N
M
1
K
＝＝
diag(,,,),diag(0,,0,
k
ry4
)
KKKK
ry4ry4ry4ry4ry4
(
N
1)
M
M
(15)
K
rx3
＝
diag(
k
rx3
,
k
rx3
,
,
k
rx3
,0,
,0)
(
N
1)
M
M
K
rx4
＝
diag(0,
,0
,
k
rx4
,
k
rx4
,
,
k
rx4
)
为便于同前人文献的数据对比，将式(6)得出的自
由振动频率
ω
代入式(16)，求得无量纲频率参数
ω
为
ω
=
ω
a
ρ
(16)
G
56
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
1.2 屈曲特征方程
虚拟弹簧边界Mindlin矩形板在图3所示的均匀初
始应力作用下，系统的总能量为
T
U
0
U
σ
0
ε
L
d
V
(17)
V
0
T
2 数值方法参数的确定
为使数值方法得到准确的数值结果，需要配置足够
的结点数以使数值解收敛；要使虚拟弹簧边界代替经典
边界条件，需要确定各个虚拟弹簧的刚度取值。
σ
0
x
0
0
xyy
(18)
2.1 结点数
为确定使数值解收敛的最少结点数，针对CCCC和
SSSS边界，厚度比
h
/
b
为0.01和0.1的4种Mindlin方板，
图4为结点数对第1、3、5阶频率参数的
ω
影响情况。
为便于计算，取长和宽方向结点数相同，即
M
=
N
。当模拟CCCC边界时，令所有弹簧刚度为10
10
N/m，
κ
=0.860 1；当模拟SSSS边界时，令
k
x1
=
k
x2
=
k
y1
=
k
y2
=0，
κ
=5/6。图中虚线分别为第1、3、
其他弹簧刚度为10
10
N/m，
5阶自振频率的收敛值，箭头指向前6阶频率收敛时的
最少单边结点数。
由图4可得，使用虚拟弹簧边界代替经典边界条件，
DQFEM仍有良好的收敛性；为保证前6阶无量纲频率
参数收敛，应配置的结点数为12×12。
u
2
v
2
w
2
x
x
x
222
1
u
v
w
L
(19)
ε
=
2
y
y
y
u
u
v
v
w
w
2
x
y
x
y
x
y
式中：
σ
为薄膜力矩阵，
ε
为应变元矩阵
类似地，基于DQFEM离散总能量
U
0
，并令
0
L
δ
U
0
=0
，可得板的屈曲特征方程
K
K
G1
K
G2
K
G3
w
0
(20)
0
x
0
y
0
xy
式中：
K
G1
、
K
G2
、
K
G3
分别为Mindlin板在受
x
轴方
向的单轴压力、受
y
轴方向的单轴压力、受
xOy
面内剪
力时的几何刚度矩阵。
(1)T(1)
K
G1
diag
I
ACA
I
A
(1)T
CA
(1)
I
B
(1)T
CB
(1)
h
A
(1)T
CA
(1)
(21)
h
B
(1)T
CB
(1)
I
B
(1)T
CB
(1)
K
G2
diag
(22)
K
G3
I
A
(1)T
CB
(1)
B
(1)T
CA
(1)
diag
I
A
(1)T
CB
(1)
B
(1)T
CA
(1)
h
A
(1)T
CB
(1)
B
(1)T
CA
(1)
(23)
(a) 模拟CCCC边界，
h
/
b
=0.01 (b) 模拟CCCC边界，
h
/
b
=0.1
为便于数据对比，以
x
轴向受压为例，将临界屈曲
应力
σ
x
代入式(24)，求得屈曲因子
n
*
为
0
hb
2
x
(24)
n
D
π
2
*
0
y
σ
y
0
τ
xy
0
(c) 模拟SSSS边界，
h
/
b
=0.01 (d) 模拟SSSS边界，
h
/
b
=0.1
b
a
σ
x
0
图4 结点数对频率参数的影响
x
O
图3 有面内应力作用的矩形板
2.2 虚拟边界弹簧的刚度取值
为分析不同类型弹簧在模拟经典边界条件时应取的
2022年第1期 /
57
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
刚度值，将弹簧分为3类：法向扭转弹簧
k
rx
i
(
i
=1，2)与
k
ry
i
(
i
=1，2)；切向扭转弹簧
k
rx
i
(
i
=3，4)与
k
ry
i
(
i
=3，4)；
线弹簧
k
t
i
(
i
=1，2，3，4)。针对线弹簧和扭转弹簧，引
入无量纲刚度
k
t
b
3
S
w
D
(25)
S
k
r
b
φ
D
3 自由振动分析数值结果
针对Mindlin板自由振动问题，本小节首先验证数
值结果的准确性，而后讨论虚拟弹簧边界对DQFEM计
算效率的影响。采用12×12的结点划分Mindlin板，边
界弹簧的刚度设置为：固支边(C)上，
S
w
=
S
φ1
=
S
φ2
=10
10
N/m；简支边(S)上，
S
w
=
S
φ2
=10
10
，
S
φ1
=0；自由边(F)上，
S
w
=
S
φ1
=
S
φ2
=0。
用
S
φ1
和
S
φ
2
分别表示法向和切向扭转弹簧的无量
纲刚度。
图5为3类弹簧的刚度对厚度比
h
/
b
为0.01和0.1
的Mindlin方板的第1、3、5阶无量纲频率参数的影响情况。
由图5可得，不论哪种边界弹簧，当其刚度大于某一定值，
频率参数不再随弹簧刚度的增大而变化。这可以解释为
弹簧的刚度太大，以至于板边界被近似严格约束。对于
法向、切向扭转弹簧以及线弹簧，实现这种严格约束的
无量纲刚度分别为10
7
N/m、10
8
N/m和10
10
N/m。
3.1 准确性
为验证本文数值解的准确性，以SSSS、SCSC和
CCCC边界条件的厚度比
h
/
b
为0.01和0.1的Mindlin
板为对象，计算频率参数
ω
。表2~表4将本文数值解
同Hinton
[10]
的精确解、Ferreira
[11]
的有限元解、Dawe D
J等
[12]
的Rayleygh-Ritz解进行对比。Ferreira的角标1
和2分别为15×15和20×20的Q4有限元法。
由表2和表3可知，在SSSS和SCSC边界条件
下得出的前5阶无量纲频率
ω
比FEM求得的数值解
更接近精确解。由表4可知，在CCCC边界条件下求
得的前5阶无量纲频率
ω
比FEM求得的数值解更接近
Rayleygh-Ritz法的数值结果。
由此，使用虚拟弹簧边界替代经典边界条件后，采
用DQFEM求解Mindlin板的自由振动问题的准确性得
(a)
S
φ
2
=
S
w
=102，
h
/
b
=0.01 (b)
S
φ
2
=
S
w
=102，
h
/
b
=0.1
以验证。值得注意的是，DQFEM仅需划分12×12的
结点数，而求得的频率参数比20×20的Q4有限元法
所得结果更加准确。
3.2 计算效率
由于边界条件施加方法的改变，DQFEM的计算效
率可能改变。表5以厚度比
h
/
b
=0.1的Mindlin方板为
(c)
S
φ
1
=
S
w
=102，
h
/
b
=0.01 (d)
S
φ
1
=
S
w
=102，
h
/
b
=0.1
研究对象，分别统计了对整体矩阵消元和使用虚拟弹簧
施加SSSS和CCCC边界条件时DQFEM求得前6阶无
量纲频率及完整的模态向量所花费的时间。
由表5可知，使用虚拟弹簧边界施加边界条件有助
于减少DQFEM的计算时间。对于SSSS边界，消元法
和虚拟弹簧边界的10个算例的平均时间分别为0.422 0
(e)
S
φ
1
=
S
φ
2
==102，
h
/
b
=0.01 (f)
S
φ
1
=
S
φ
2
==102，
h
/
b
=0.1
图5 频率参数随3类弹簧刚度的变化情况
s和0.412 2 s，虚拟弹簧边界使DQFEM的平均计算时
间缩短了2.32%；对于CCCC边界，消元法和虚拟弹簧
58
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
表2 SSSS边界，
κ
=0.833，Mindlin方板前5阶无量纲频率
ω
h/b
计算方法
Present
0.01Hinton
Ferreira
2
Present
0.1Hinton
Ferreira
1
频率阶次
1
0.096 3
0.096 3
0.096 5
0.930 3
0.930
0.934 6
2
0.240 6
0.240 6
0.243 0
2.219 3
2.219
2.254 5
3
0.240 6
0.240 6
0.243 0
2.219 3
2.219
2.254 5
4
0.384 7
0.384 8
0.389 0
3.405 6
3.406
3.459 2
5
0.480 7
0.480 9
0.492 8
4.149 4
4.149
4.303 1
表3 SCSC边界，
κ
=0.822，Mindlin方板前5阶无量纲频
ω
h/b
计算方法
Present
0.01Hinton
Ferreira
1
Present
0.1Hinton
Ferreira
1
频率阶次
1
0.141 1
0.141 1
0.142 4
1.300 1
1.302
1.294 0
2
0.266 8
0.266 8
0.271 0
2.393 9
2.398
2.397 1
3
0.337 6
0.337 7
0.348 4
2.884 5
2.888
2.929 0
4
0.460 4
0.460 8
0.472 2
3.839 1
3.852
3.839 4
5
0.497 7
0.497 9
0.519 1
4.231 3
4.237
4.347 5
表4 CCCC边界，
κ
=0.860 1，Mindlin方板前5阶无量纲频
ω
h/b
计算方法
Present
0.01Ferreira
2
Dawe
Present
0.1Ferreira
2
Dawe
频率阶次
1
0.175 4
0.175 0
0.175 4
1.591 0
1.595 5
1.594 0
2
0.357 4
0.363 5
0.357 6
3.038 9
3.066 2
3.039 0
3
0.357 4
0.363 5
0.357 6
3.038 9
3.066 2
3.039 0
4
0.526 5
0.535 8
0.527 4
4.262 5
4.292 4
4.265 0
5
0.639 9
0.663 4
0.640 2
5.024 7
5.123 2
5.035 0
边界的10个算例的平均时间分别为0.314 9 s和0.305
2 s，虚拟弹簧边界使DQFEM的平均计算时间缩短了
3.08%。
和受双轴向均匀载荷时的临界屈曲因子
n
*
，同时列出
Mizusawa T
[13]
的数值解、Hashemi S等
[14]
的精确解。
剪切修正因子
κ
取
π
2
/12。数值结果与Hashemi S等精确
解的相对误差由式(26)求得，即
4 屈曲分析数值结果
表6和表7列出SSSS、SCSF和SFSF边界厚度比
h
/
b
为0.05和0.1的Mindlin方板，在分别受
x
轴方向
Present-Exact
RelError=
100%
Exact
(26)
2022年第1期 /
59
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
表5 消元法和虚拟弹簧边界DQFEM的计算时间 s
样本
1
2
3
4
5
6
7
8
9
10
平均时间
CCCC
消元法
0.421 8
0.415 6
0.414 3
0.426 1
0.414 0
0.428 1
0.425 1
0.431 3
0.424 7
0.419 2
0.422 0
虚拟弹簧
0.414 3
0.415 0
0.408 9
0.410 5
0.414 2
0.417 7
0.407 5
0.415 6
0.406 4
0.411 4
0.412 2
消元法
0.321 1
0.310 9
0.316 7
0.315 5
0.310 0
0.319 4
0.305 9
0.317 1
0.316 6
0.315 3
0.314 9
SSSS
虚拟弹簧
0.301 9
0.311 2
0.306 2
0.303 2
0.301 0
0.300 8
0.306 2
0.306 7
0.308 3
0.306 5
0.305 2
表7 双轴向受压Mindlin方板的无量纲屈曲因子
n
*
边界条件计算方法
0.05
Present
SSSSHashemi
RelError
Present
SCSFHashemi
RelError
Present
SFSFHashemi
RelError
1.964 0
1.971 8
0.40%
1.108 9
1.111 9
0.27%
0.9187
0.9207
0.22%
h/b
0.1
1.864 5
1.891 9
1.45%
1.052 4
1.064 1
1.10%
0.890 1
0.897 7
0.85%
5 结论
使用虚拟弹簧边界能够很好得模拟板的经典边界条
件而不影响DQFEM的收敛性。结合了虚拟弹簧边界的
DQFEM，与FEM相比，使用更少的结点数即可求得更
高精度的数值解；求得的屈曲因子更加安全保守。同时，
采用虚拟弹簧边界有助于提高DQFEM的计算效率。
对比表6和表7数据可知，将DQFEM和虚拟弹簧
边界共同用于Mindlin板的屈曲分析时，得出的屈曲因
子与Mizusawa通过样条方法得出的数值解几乎完全吻
合，而略小于Hashemi的精确解，但最大相对误差不超
过1.5%。这说明使用虚拟弹簧边界模拟经典边界条件，
并使用DQFEM分析Mindlin板的特征屈曲问题，具有
一定准确性，求得的临界屈曲因子偏向于安全保守。
表6
x
轴方向受压Mindlin方板的无量纲屈曲因子
n
边界条件计算方法
0.05
Present
SSSS
Mizusawa
Hashemi
RelError
Present
SCSF
Mizusawa
Hashemi
RelError
Present
SFSF
Mizusawa
Hashemi
RelError
3.928 0
3.928
3.943 7
0.40%
1.615 3
1.615
1.619 7
0.27%
0.941 2
0.941 2
0.943 1
0.20%
h/b
0.1
3.728 8
3.729
3.783 8
1.45%
1.539 4
1.539
1.555 8
1.05%
0.914 6
0.914 6
0.921 9
0.79%
*
参考文献
[1] 漆静，程文明，濮德璋，等.起重机箱形梁曲形腹板简支
板局部稳定性分析[J].机械设计与制造,2014(4)：
162-164.
[2]
Xing Y
，
Liu B
，
Liu G.A differential quadrature finite
element method
[J].
International Journal of Applied
Mechanics
，2010，2(1)：207-227.
[3]
Wang Q
，
Shi D
，
Qian L
，
et al. Free vibrations of composite
laminated doubly-curved shells and panels of revolution with
general elastic restraints
[J].
Applied Mathematical
Modelling,
2017，46(1)：227-262.
[4]
Wang Q
，
Shao D
，
Qin B.A simple first-order shear
deformation shell theory for vibration analysis of composite
laminated open cylindrical shells with general boundary
conditions
[J].
Composite Structures
，2017，184：
211-232.
60
/ 2022年第1期
Copyright©博看网. All Rights Reserved.
ANALYSIS RESEARCH
[5]
Kim K
，
Kim K
，
Han C
，
et al. A method for natural
frequency calculation of the functionally graded rectangular
plate with general elastic restraints
[J].
AIP Advances
，2020.
[6]
Du Y
，
Pang F
，
Sun L
，
et al. A unified formulation for
dynamic behavior analysis of spherical cap with uniform and
stepped thickness distribution under different edge constraints
[J].
Thin-Walled Structures
，2019.
[7]
Guan X
，
Tang J
，
Wang Q
，
et al. Application of the differential
quadrature finite element method to free vibration of elastically
restrained plate with irregular geometries
[J].
Engineering
Analysis with Boundary Elements
，2018，90：1-16.
[8] 向维.剪切板自由振动和特征屈曲新解法[D].北京：
北京航空航天大学，2016.
[9] 邢誉峰，李敏.计算固体力学原理与方法[M].北京：
北京航空航天大学出版社，2011.
[10]
Hinton cal methods and software for dynamic
作 者：单逸峰
电子邮箱：****************
收稿日期：2021-06-11
analysis of plates and shells
[J].
Pineridge
，1988.
[11]
Ferreira codes for finite element analysis
[M].
Springer Netherlands
，2009.
[12]
Dawe D J
，
Roufaeil O L. Rayleigh-Ritz vibration analysis
of Mindlin plates
[J].
Journal of Sound & Vibration
，1980，
69(3)：345-359.
[13]
Mizusawa ion of rectangular Mindlin plates with
tapered thickness by the spline strip method
[J].
Computers
and Structures
，1993，30(12)：1 663-1 677.
[14]
Hashemi S
，
Khorshidi K
，
Amabili M. Exact solution for
linear buckling of rectangular Mindlin plates
[J].
Journal
of Sound and Vibration
，2008，315(1，2)：318-342.
三一增程式电动正面吊运机
三一增程式电动正面吊运机广泛应用于各大港口、码头、铁路堆场，以适应性强、不畏严寒酷暑、
科技增效、维修简易、服务及时、数据互联等特点助力客户成为绿色低碳物流的先行者。
1）绿色降本，适应性强：三一增程式电动正面吊运机可以闲时纯电驱动，忙时增程式驱动，一小时
即可满充，作业可自由切换，丝毫不担心电池续航不足出现停机。根据作业数据，三一电动正面吊运机
纯电作业能耗成本较燃油车下降68%，增程器作业能耗成本下降40%，综合作业一年可节省能耗成本超过
25万元，减少二氧化碳排放近200 t/a。
2）巧妙设计，不畏严寒酷暑：在夏季高温时，为了确保电动产品的安全性，该吊运机采用磷酸铁锂
电池系统，热稳定性好，800℃才会分解但不会自燃，可确保高温环境中的人机安全。由于北方冬季寒冷，
针对低温工况增加了隔温保护套，减少外界低温环境对电池的影响，保障冬季电池续航。
3）科技加持，效能提升：针对行驶系统能耗高的反馈，新一代电动正面吊运机由双电动机改为单电
动机直驱，效率优于双电动机。新一代控制系统的物理按键减少，采用新型手柄，设备的动作响应性和
操作实时性随之提高，且具备智能自诊断(IO)，运行周期比上一代提升10倍。第三代控制系统，10寸
大触摸屏，真彩显示，智能控制。
4）维保简易 服务及时：通过设计不断提高设备维保的便捷程度，缩短检修时间。
5）大屏管理 数据互联：为提高产品附加价值，三一开发了C端互联系统，向客户提供电动产品大
数据支持，电量、电耗、总耗电量以及作业量等实时显示，帮助客户更好地安排作业。
2022年第1期 /
61
Copyright©博看网. All Rights Reserved.