第3章单自由度体系5(直接积分法)

……
1. 数值积分概述
一种逐步积分法的优劣，主要由以下四个方面判断：
(1) 收敛性：当Δt→0时，数值解是否收敛于精确解； (2) 计算精度：截断误差与时间步长Δt的关系，若误差
∝ O(Δtn)，则称方法具有n阶精度； (3) 稳定性：随时间步数i的增大，数值解是否变得无穷
大（远离精确解）； (4) 计算效率：所花费的计算时间。
ui
=
ui+1
−
2ui Δt 2
+
ui−1
ui = u(ti ) ui = u(ti ) ui = u(ti )
mu(ti ) + cu(ti ) + ku(ti ) = p(ti )
pi = p(ti )
m ui+1
− 2ui Δt 2
+ ui−1
+ c ui+1 − ui−1 2Δt
⎤ ⎥⎦
ui −1
等效荷载
2、中心差分法
中心差分方法计算中的起步处理方法
ui
=
ui+1 − ui−1 2Δt
ui
=
ui+1
−
2ui Δt 2
+
ui−1
⎡m ⎢⎣Δt 2
+
c 2Δt
⎥⎦⎤ui+1
=
pi
− ⎢⎣⎡k
−
2m Δt 2
⎥⎦⎤ui
−
⎡ ⎢⎣
m Δt 2
−
c 2Δt
⎥⎦⎤ui−1
初始条件为 ：
可以求得t i＋1时刻的位移、速度和加速度。
1. 数值积分概述
一般时域逐步积分法的构造
以上方法也称为平均加速度法，即假设加速度为t i和t i＋1 时间段内的平均值：
a = ui + ui+1 2
也可以假设加速度a为其它形式的变化规律，例如为线
性变化：
ui
=
ui+1 − ui Δt
2Δt
}−1
,
{ui
}
=
{ui
}+1
−
2{ui
Δt 2
}
−
{ui
}−1
(3). 下一步计算用i+1代替i，重复(2)中的计算步骤。
2、中心差分法
中心差分法的特点：
是收敛的；
具有2阶精度，即误差∝O(Δt2) ；
是有条件稳定，稳定条件Δt≤Tn/π； 具有较高的计算效率。
2、中心差分法 优点：
显式积分方法 (explicit integration methods)： ti+1时刻的解是利用ti时刻的动平衡来获得。计 算工作量小，增加的工作量与自由度成线性关 系，如中心差分方法。
1. 数值积分概述
数值积分分类：
按计算稳定性对计算时间间隔Δt 的要求，分
为：
条件稳定（conditionally stable）：
ui+1
=
ui
+
Δtui
+
Δt 2 2
ui
+
Δt 3 6
ui
+"
ui+1
=
ui
+
Δtui
+
Δt 2 2
ui
+"
如果假设在t a，则
i和t
i＋1时刻，即Δt时间段内，体系的加速度为常数
ui = 0
因此， t i＋1时刻体系的速度和位移为：
ui+1
=
ui
+
Δtui
+
Δt 2 2
3、Newmark —β法
将ti+1时刻的位移和速度采用Taylor级数展开，并忽略四 阶以上的高阶小量：
ui+1
=
ui
+
Δtui
+
Δt 2 2
ui
+
βΔt 3ui
ui+1 = ui + Δtui + γΔt 2ui
假定在[ti， ti+1]内加速度为线性变化，即：
ui
(K.J. Bathe, Finite Element Procedures, Prentice-Hall, 1996.)
Two ideas: (1)运动方程并不在任何时间t都得到满足，而仅仅是在
以时间间隔为Δt的离散时间点上得到满足。 (2)在时间间隔Δt内，对位移、速度和加速度的变化作
出某些假定。
将体系ti+1和ti-1时刻的运动采用Taylor级数展开：
ui+1
=
ui
+
Δtui
+
Δt 2 2
ui
+
Δt 3 6
ui
+"
ui−1 = ui − Δtui +
将两式相减，整理得：
Δt 2 2
ui
−
Δt 3 6
ui
+"
ui
=
ui+1 − ui−1 2Δt
+
O(Δt2 )
将两式相加，整理得：
ui+1
=
ui
+
Δtui
+
(
1 2
−
β
)Δt 2ui
+
βΔt 2ui+1
Newmark—β与其它直接积分法的联系：
γ = 1/ 2, β = 1/ 4
ui+1
=
ui
+
Δt 2
(ui
+
ui +1 )
ui+1
=
ui
+
Δtui
+
Δt 2 4
(ui
+
ui +1 )
平均加速 度法
线性加速 度法
γ = 1/ 2, β = 0
=
ui+1 − Δt
ui
于是可得：
ui+1 = ui + (1 − γ )Δtui + γΔtui+1
ui +1
=
ui
+
Δtui
+
(1 2
−
β )Δt 2ui
+
βΔt 2ui+1
3、Newmark —β法 ui+1 = ui + (1− γ )Δtui + γΔtui+1
3、Newmark —β法
Newmark在1959年提出的一种积分方法。 Newmark—β同样将时间离散化，运动方程仅要
求在离散的时间点上满足。假设在ti时刻的运 动均已求得, 然后计算 ti+1时刻的运动。
与中心差分法不同的是，它不是用差分对ti时刻 的运动方程展开，得到外推计算ui+1的公式， 而是通过对加速度的假设，以ti时刻的运动量 为初始值，得到计算ti+1时刻的运动公式。
只有当Δt<tcr时，它的数值解不会无界增大，即
计算结果稳定，如中心差分方法。
无条件稳定 (unconditionally stable)： 计算稳定性对计算时间间隔的大小没有要求。如
平均加速度法。
1. 数值积分概述
数值积分分类： 按计算ti+1时刻的解所需知道前续时刻的个数，分
为：
单步法（single-step integration methods）： 仅需已知ti时刻的解，如Newmark—β法。 多步法 (multi-steps integration methods)： 两步法需要已知ti时刻和ti-1时刻的解。如中心差
一个好的方法首先必须是收敛的、有足够的精度（例如 二阶，满足工程要求）、良好的稳定性、较高的计算 效率。
1. 数值积分概述
数值积分分类： 按计算ti+1时刻的解所需运用平衡方程，可分为：
隐式积分方法（implicit integration methods）： ti+1时刻的解要运用ti+1时刻的动平衡来获得。 计算工作量大，通常增加的工作量与自由度的 平方成正比，例如Newmark—β法、Wilson — θ法。
−
Δtu0
+
Δt 2 2m
( p0
−
cu0
−
ku0 )
(2). 根Ku据i+1i及= ip以i 前时刻的运动，计算i+1时刻的运动
K=
m+
c
Δt2 2Δt
pi
=
pi
−
⎡⎢⎣k
−
2m Δt 2
⎤ ⎥⎦
ui
−
⎡m ⎢⎣ Δt2
−
c 2Δt
⎤ ⎥⎦
ui
−1
如果需要
ui
=
ui+1 − ui−1 2Δt
体系在ti＋1时刻的运动包括：位移、速度和加速度，需 要有三个方程（条件）求这三个量。因此，除体系的 运动方程外，还需补充两个方程（条件） 。
mu(t) + cu(t) + ku(t) = p(t)
1. 数值积分概述
一般时域逐步积分法的构造
两个补充方程可以通过对运动状态的假设得到。
体系在ti及ti以前时刻的运动已知，则ti＋1时刻的运动可采用Taylor 级数展开：
= const,u(τ ) = ui
+
τ
Δt
(ui
+1
−
ui
)
则采用同样的分析步骤可以得到线性加速度法的时域逐 步积分公式。
1. 数值积分概述
一般时域逐步积分法的构造
平均加速度法和线性加速度法的基本假设和补充公式
2、中心差分法（Central Difference Method）
第三章 单自由度体系
直接积分法
主要内容
• 两种直接积分方法 （1）中心差分法 （2）Newmark —β法
• 数值积分的稳定性 • 了解算法阻尼（数值阻尼）现象
1. 数值积分概述（直接积分法，逐步积分法）
(Direct Integration Methods, Step-by-Step Methods)
a
ui+1 = ui + Δta
1. 数值积分概述
一般时域逐步积分法的构造
令t i和t i＋1时间段内的常加速度a=(üi+1+üi)/2，则 得到：
ui+1
=
ui
+
Δtui
+
Δt 2 2
⋅
ui
+ ui+1 2
ui+1
=
ui
+
Δt
⋅
ui
+ ui+1 2
再加上t i＋1时刻的运动方程： mui+1 + cui+1 + kui+1 = pi+1
运动方程： mu+ c(t)u + k(t)u = p(t)
In direct integration the equations of equilibrium are integrated using a numerical step-by-step procedure, the term ‘direct’ meaning that prior to the numerical integration, no transformation of equations into a different form is carried out.
分方法。
1. 数值积分概述
数值积分构造常用手法： （1）基于激励函数插值法 （2）基于有限差分法 （3）基于加速度假定法
1. 数值积分概述
一般时域逐步积分法的构造 所谓的时域逐步积分方法就是构造出根据某一时刻及其
以前时刻的运动，推算下一时刻运动的递推计算公 式。
具体情况可表述为，设体系在ti及ti以前时刻的运动已 知，求ti＋1时刻的运动（ti＝iΔt）。
⎡⎣K
⎤⎦
=
[m]
Δt2
+
[c]
2Δt
（1）由于等效刚度矩阵仅与质量和阻尼有关，当采用
集中质量法进行多自由度计算时，并忽略阻尼的影
响，则
{ui+1}
=
Δt2 diag(mii
)
{
pi
}
从而联立方程解耦。
（2）计算可在单元水平上进行，从而降低计算的存储 要求。
缺点：
由于是条件稳定的积分方法，当体系的最小周期 Tmin → 0 则导致 Δt → 0 ，譬如：当集中质量矩阵中，部分质 量 mii → 0 ，从而使计算工作无法完成。
+ kui
=
pi
⎡m ⎢⎣ Δt 2
+
c 2Δt
⎥⎦⎤u
i
+1
=
pi
− ⎢⎣⎡k −
2m Δt 2
⎥⎦⎤ui
−
⎡ ⎢⎣
m Δt 2
−
c 2Δt
⎥⎦⎤ui−1
等效刚度
Kui +1
=
pi
K=
m
+
c
Δt2 2Δt
pi
=
pi
−
⎡⎢⎣k
−
2m Δt 2
⎤ ⎥⎦
ui
−
⎡m ⎢⎣ Δt2
−
c 2Δt
1. 数值积分概述
常用的数值积分方法： （1）分段解析法； （2）中心差分法； （3）Runge-Kutta法； （4）Houbolt法； （5）平均加速度法； （6）线性加速度法； （7）Newmark —β法； （8）Wilson —θ法； （9）HHT法(Hilber-Hughes-Taylor method)； （10）精细积分法；
ui
=
ui+1 − ui−1 2Δt
ui
=
ui +1
− 2ui Δt 2
+
ui −1
+
O(Δt2 )
ui
=
ui+1
− 2ui Δt 2
+
ui−1
因此，速度和加速度的精度是同一量级的。
2、中心差分法
中心差分方法实际上就是用有限差分代替位移对时间的
求导（即速度和加速度）。
ui
=
ui+1 − ui−1 2Δt
,
ui
=
ui+1
− 2ui Δt 2
− ui−1
(3). 下一步计算用i+1代替i，重复(2)中的计算步骤。
2、中心差分法
对于多自由度，中心差分法计算步骤：
(1). 初始计算
{u0}
和
{u0} 已知
{u−1}
=
{u0}
−
Δt {u0 }
+
Δt 2
2[m]
({
p0 }
−
[c]{u0}
−
[k
u0 = u(0), u0 = u(0)
u0
=
u1 − u−1 2Δt
u0
=
u1
−
2u0 + u−1 Δt 2
u−1
=
u0
−
Δtu0
+
Δt 2 2
u0
u0
=
1 m
( p0
−
cu0
−
ku0 )
2、中心差分法
中心差分法计算步骤：
(1). 初始计算
u0 和 u0 已知
u−1
=
u0ቤተ መጻሕፍቲ ባይዱ
]{u0})
(2). 根据i及i以前时刻的运动，计算i+1时刻的运动
⎡⎣K
⎤⎦{ui+1}
=
{
pi
}
⎡⎣K⎤⎦
{
pi
}
[m]
=+ Δt2
={ pi} −
[c]
2Δt
([k]
−
2Δ[tm2 ] ){ui
}
−
[m]
( Δt2
−
[c]
2Δt
){ui−1}
{ui
}
=
{ui
}+1 − {ui