第五章 动力反应数值计算
5.2 时间步进法
对于数值算法而言，关注一下三点：
收敛性：随着时间步长的减少，数值解应逼近精确解。
稳定性：由于数值舍如误差的存在，数值解应当稳定。
精 度：数值计算结果应与精确解接近。
本章介绍三种方法：
基于激励函数插值
基于速度和加速度的有限差分
基于假设加速度变化的方法
基于激励函数插值
from math import *k = 10dt = 0.1wn = 6.283x = 0.05wd = wn*sqrt(1-x**2)wdt = wd*dt;wnt = wn*dt;c1 = pow(exp(1),(-x*wnt))c2 = x/sqrt(1-x**2)c3 = sin(wdt)c4 = cos(wdt)c5 = 1 - 2*x*xA = c1*(c2*c3 c4)B = c1*(c3/wd)C = (1/k)*(2*x/wnt c1*((c5/wdt - c2)*c3-(1 2*x/wnt)*c4))D = (1/k)*(1 - 2*x/wnt c1*((2*x*x - 1)/wdt)*c3 (2*x/wdt)*c4)#print("wd is %f\nwdt is %f\nwnt is %f\nc1 is %f\nc2 is %f\nc3 is %f\nc4 is %f\nc5 is %f\n" % (wd, wdt,wnt,c1,c2,c3,c4,c5))print("A is %.6f\nB is %.6f\nC is %.6f\nD is %.6f\n" % (A,B,C,D))A is 0.812939B is 0.090671C is 0.012355D is 0.006766以上内容为求解系数，迭代过程如下用Python实现相当方便。
def pForce(x):
if x > 0.6:
return 0
return 10*sin(pi*x/0.6)a = 0b = 0print("---------------------------------------------------")print(" u v ")for i in range(10):
u = A*a B*b C*pForce(i*dt) D*pForce((i 1)*dt)
v = A1*a B1*b C1*pForce(i*dt) D1*pForce((i 1)*dt)
a = u; b = v
print(" iteration times is %2d % 10.6f % 10.6f" % (i 1,u,v))print("---------------------------------------------------")输出结果大概是这样的：
---------------------------------------------------
u v
iteration times is 1 0.031835 0.935304
iteration times is 2 0.227598 3.067482
iteration times is 3 0.633821 4.854696
iteration times is 4 1.134128 4.730075
iteration times is 5 1.489690 1.931473
iteration times is 6 1.447931 -3.017619
iteration times is 7 0.903468 -7.463887
iteration times is 8 0.057703 -8.876308
iteration times is 9 -0.757919 -6.916763
iteration times is 10 -1.243295 -2.516007
---------------------------------------------------
与书中的答案基本上一致，小数点后地3位开始偶有不同。
小结：
基于激励的插值方法本质上属于直接计算方法，属于显示计算的一种。本文仅给出了数值解（近似解），没有给出解析解（精确解），书中有。
在没有特别说明的情况下，默认计算机语言为Python。
在没有特别说明的情况下，默认结构动力学教材为：
Anil K. Chopra. 谢礼立, 吕大刚等译. 结构动力学理论及其在抗震工程中的应用（第四版）. 北京：高等教育出版社. 2016.09.
首发于 2018-02-28
