非线性动力学方程的求解方法
1000138生物药剂学与药物动力学_第十一章非线性药物动力学_1002

dt
V Vm
= Km+C
C
因为: X = V C 所以: Vm C dXe = Km+C Vdt
1) C>>Km 时: Vm V CL = C 2) Km>>C 时: CL = Vm V
V Vm C dXe = dt Km+C
Km
3)药物有线性和非线性消除同时存在时:
Vm V
CL =
Km+C
+KV
生物药剂学与药物动力学
Biopharmaceutics and Pharmacokinetics
浙江大学药学院药物制剂研究所
邱利焱 博士, 教授
第十一章
内容概要
非线性药物动力学
概述 非线性药物动力学方程 血药浓度与时间的关系及参数计算
第一节 概述
一、线性动力学的三个基本假设:
1. 吸收速度为零级或一级 2. 药物分布相很快完成 3. 药物在体内消除属一级速度过程
特征
1. 血药浓度与剂量成正比的改变。 2. 血药浓度-时间曲线下面积与剂量成正比。 3. 药物的生物半衰期与剂量无关。 4. 药物的体内过程用一级速度过程或线性过程表示。
二、非线性药物动力学特点
1. 血药浓度与剂量不成正比。 2. 3. 4. 5. 6. AUC与剂量不成正比。 当剂量增加时,消除半衰期延长。 药物的消除不呈一级动力学特征, 即消除过程是非线性的。 其它药物可能与其竞争酶或载体系统影响其动力学过程。 药物代谢物组成比例可能由于剂量变化而变化。
3) ΔC/Δt = Vm- ΔC/Δt Km
C中
ΔC/Δt 对(ΔC/Δt)/ C中作图,斜率-Km,截距Vm
非线性动力系统的数值计算方法及稳定性分析

非线性动力系统的数值计算方法及稳定性分析非线性动力系统是一种非常常见的实际物理系统,例如电路、化学反应、天气系统等,它们的行为通常比线性系统更加复杂。
数值计算非线性动力系统的稳定性与动力学特性是一个非常重要的课题,对于研究和预测实际系统的行为有着非常重要的意义。
在本文中,我们将介绍几种常见的非线性动力系统的数值计算方法及它们的稳定性分析。
一、欧拉法欧拉法是动力系统数值计算中最基本的一种方法。
它的基本思路是将连续的时间离散化,将微分方程转化成差分方程,然后用迭代的方式求解。
欧拉法的迭代公式为:$$y_{n+1}=y_{n}+hf(y_n)$$其中,$h$为步长,$f(y_n)$是微分方程在$y_n$处的导数。
欧拉法是一种比较简单易懂的方法,但是它的稳定性较差,容易产生数值误差。
欧拉法对于初始值的依赖性很强,如果步长$h$选取过大,就会导致解的不稳定。
因此,在使用欧拉法进行数值计算时,我们需要根据实际问题来调整步长,以保证数值解的正确性。
二、龙格-库塔法龙格-库塔法是一种常见的数值积分方法,在动力系统数值计算中也常常被使用。
它的基本思路是利用微分方程的某些性质,选取合适的时间步长和权重,在数值上求得微分方程的积分近似解。
龙格-库塔法通常可以由一些权重系数和步长系数组成,如下:$$Y_{n+1}=Y_n+\frac{1}{6}(k_1+2k_2+2k_3+k_4)$$$$k_1=hf(Y_n)$$$$k_2=hf(Y_n+\frac{1}{2}k_1)$$$$k_3=hf(Y_n+\frac{1}{2}k_2)$$$$k_4=hf(Y_n+k_3)$$其中,$k_1,k_2,k_3,k_4$均为微分方程在相应位置处的导数。
龙格-库塔法比欧拉法更加稳定,适用于多数动力系统的数值计算。
但是,龙格-库塔法在计算一些比较长时间范围内的运动时,需要降低步长以保证解的精度。
同时,权重系数和步长系数也需要根据具体问题调整,才能得到更加准确的数值解。
解非线性方程的牛顿迭代法及其应用

解非线性方程的牛顿迭代法及其应用一、本文概述非线性方程是数学领域中的一个重要研究对象,其在实际应用中广泛存在,如物理学、工程学、经济学等领域。
求解非线性方程是一个具有挑战性的问题,因为这类方程往往没有简单的解析解,需要通过数值方法进行求解。
牛顿迭代法作为一种古老而有效的数值求解方法,对于求解非线性方程具有重要的应用价值。
本文旨在介绍牛顿迭代法的基本原理、实现步骤以及在实际问题中的应用。
我们将详细阐述牛顿迭代法的基本思想,包括其历史背景、数学原理以及收敛性分析。
我们将通过具体实例,展示牛顿迭代法的计算步骤和实际操作过程,以便读者能够更好地理解和掌握该方法。
我们将探讨牛顿迭代法在各个领域中的实际应用,包括其在物理学、工程学、经济学等领域中的典型应用案例,以及在实际应用中可能遇到的问题和解决方法。
通过本文的介绍,读者可以深入了解牛顿迭代法的基本原理和应用技巧,掌握其在求解非线性方程中的实际应用方法,为进一步的研究和应用提供有力支持。
二、牛顿迭代法的基本原理牛顿迭代法,又称为牛顿-拉夫森方法,是一种在实数或复数域上近似求解方程的方法。
其基本原理是利用泰勒级数的前几项来寻找方程的根。
如果函数f(x)在x0点的导数f'(x0)不为零,那么函数f(x)在x0点附近可以用一阶泰勒级数来近似表示,即:这就是牛顿迭代法的基本迭代公式。
给定一个初始值x0,我们可以通过不断迭代这个公式来逼近f(x)的根。
每次迭代,我们都用当前的近似值x0来更新x0,即:这个过程一直持续到满足某个停止条件,例如迭代次数达到预设的上限,或者连续两次迭代的结果之间的差小于某个预设的阈值。
牛顿迭代法的收敛速度通常比线性搜索方法快,因为它利用了函数的导数信息。
然而,这种方法也有其局限性。
它要求函数在其迭代点处可导,且导数不为零。
牛顿迭代法可能不收敛,如果初始点选择不当,或者函数有多个根,或者根是重根。
因此,在使用牛顿迭代法时,需要谨慎选择初始点,并对迭代过程进行适当的监控和调整。
牛顿迭代法解动力学方程不收敛

牛顿迭代法解动力学方程不收敛摘要:1.引言2.牛顿迭代法简介3.动力学方程及其收敛性问题4.牛顿迭代法在解动力学方程中的应用5.牛顿迭代法在解动力学方程中的不收敛问题6.结论正文:1.引言在物理学中,动力学方程是描述物体运动状态的数学模型,广泛应用于各种实际问题中。
然而,在求解动力学方程时,常常会遇到收敛性问题。
牛顿迭代法作为一种求解非线性方程的数值方法,被广泛应用于解动力学方程。
本文将探讨牛顿迭代法在解动力学方程中的不收敛问题。
2.牛顿迭代法简介牛顿迭代法是一种求解非线性方程的数值方法,其基本思想是通过迭代使得函数值逐步逼近零。
对于非线性方程F(x) = 0,牛顿迭代法的迭代公式为:x[n+1] = x[n] - F(x[n])/F"(x[n]),其中F"(x) 表示F(x) 的导数。
牛顿迭代法具有二阶收敛性,即当迭代步长足够小,且初始值足够接近真实解时,可以通过有限次迭代得到精确解。
3.动力学方程及其收敛性问题动力学方程描述了物体在给定力的作用下的运动状态,通常包括质量、速度、加速度等物理量。
求解动力学方程时,通常需要采用数值方法,因为解析解往往难以求得。
然而,在数值求解过程中,可能会遇到收敛性问题。
例如,在迭代过程中,如果迭代步长过大或者初始值与真实解差距过大,可能导致迭代结果发散,无法得到精确解。
4.牛顿迭代法在解动力学方程中的应用由于牛顿迭代法具有二阶收敛性,因此在求解动力学方程时,可以得到较好的数值解。
在实际应用中,可以根据动力学方程的特点,选择合适的牛顿迭代法求解。
例如,对于具有显式解的动力学方程,可以直接使用牛顿迭代法求解;对于具有隐式解的动力学方程,可以通过拟合等方法得到显式解,然后使用牛顿迭代法求解。
5.牛顿迭代法在解动力学方程中的不收敛问题尽管牛顿迭代法具有二阶收敛性,但在求解动力学方程时,仍然可能出现不收敛的情况。
这主要是因为动力学方程的非线性特性和迭代过程中的误差累积。
非线性动力学

t∈R
x∈ Rn
的解,则显然它是不仅是时间的函数,而且也是初值的函数,即解随着初值的改变而改变, 可以将解记为
φ(t, x0 )
当 x0 是 R n 中的某一点时,φ (t, x0 ) 代表了 1 条解轨线,而
{φ(t, x0 ) x0 ∈ D}
则代表了一族轨线。将φ看成是一个映射,即
φ : R× Rn → Rn
运动行为,它在物理上对应了这样的一个观点:在系统的最初阶段,系统由于外界的初始干 扰,将呈现相当复杂的运动形式,但随着时间的延续,运动将进入平稳状态,而这种平稳状 态体现了动态系统的本质结构。
微分方程解的最终形态通常有: (1) 平衡点 (2) 周期解 (3) 拟周期解 (4) 混沌解
6.4.1 平衡点
图 6-7 所示是 2 维线性系统的相轨线,坐标原点是系统的平衡点,图 6-7a、b 中的平衡 点是稳定的,称为稳定结点,图 6-7c 中的平衡点是不稳定的,称为鞍点。
图 6-7 2 维线性系统的相轨线
6.5.2 任意解的稳定性
设 x = ψ (t)是微分方程 x& = F(t, x)
第 6 章 非线性动力学
-0.5
-1
-1.5
0.5
1
1.5
图 6-2 例 1 相图
例2
如图 6-3 所示是微分方程
&y& + 0.2 y& + y = 0
在相平面 (x1, x2 ) ,
x1 = y
x2 = y&
上的轨线图,平衡点为 (0,0),当 t → ∞ 时,解轨线趋于平衡点。
0.6 0.4 0.2
-0.6
-0.4
-0.2 -0.2
使用matlab求解vanderpol方程的研究方法

使用matlab求解vanderpol方程的研究方法研究van der pol方程是一种重要的非线性动力学问题,它描述了一些自振系统的运动行为。
这个方程由物理学家布里尔·范·德·波尔在1920年首次提出,并广泛用于描述电路、生物和化学系统中的自振现象。
van der pol方程的一般形式如下:d^2x/dt^2 - μ(1 - x^2)dx/dt + x = 0其中,x是系统的状态变量,t是时间,而μ是一个非负常数,代表了系统的非线性因素。
为了研究van der pol方程的行为,我们可以使用MATLAB进行求解和分析。
下面将介绍几种常用的方法。
1. 数值解法:通过数值方法求解van der pol方程是最常用的研究方法之一、MATLAB中可以使用ode45函数来进行求解。
该函数基于Adams-Bashforth-Moulton方法,可以自动选择适当的步长,并给出较高的数值精度。
例如,我们可以定义一个匿名函数来表示van der pol方程,并使用ode45来求解:```MATLABmu = 1; % 设置μ的值[t, x] = ode45(f, [0, 10], [0.1; 0]); % 求解方程plot(t, x(:,1)); % 绘制解的图像xlabel('t'); ylabel('x');```这段代码定义了一个匿名函数f,其中t是时间,x是状态变量。
通过ode45函数对该方程进行求解,给定了初始条件[0.1; 0],求解时间范围为[0, 10]。
最后,使用plot函数对解进行可视化。
2. 分岔和稳定性分析:van der pol方程是一个具有非线性耗散的系统,因此它可以表现出各种动力学行为,如周期振荡、混沌和吸引子等。
MATLAB提供了一些函数用于分析方程的分岔和稳定性。
例如,可以使用bifurcation函数分析van der pol方程的分岔图:```MATLABmu = linspace(0, 10, 100); % 设定μ的取值范围x = linspace(-2, 2, 100); % 设定x的取值范围bifplot(bif_x, bif_par); % 绘制分岔图xlabel('μ'); ylabel('x');```这段代码根据设定的μ和x范围,使用bifurcation函数计算系统的分岔图,并使用bifplot函数进行可视化。
第十一章 非线性动力学

可饱和的代谢过程;酶诱导;较高剂量时 的肝中毒;肝血流的变化;代谢物的抑制 作用
二、非线性药物动力学特点与识别
特点:
药物消除为非一级动力学,遵从米氏方程 AUC与剂量不成正比 消除半衰期随剂量增大而延长,剂量增加至一定 程度时,半衰期急剧增大 动力学过程可能会受到合并用药的影响 代谢物的组成比例受剂量的影响
当C0>>Km时, t1/2=C0/(2Vm) 当Km>>C0时, t1/2=0.693Km/Vm
清除率Cl
dX dt Cl C VmC dX dt ( dC dt ) V V Km C Vm V Cl Km C
当C>>Km时, Cl与C成反比:CL=Vm*V/C 当Km>>C时, Cl与C无关: CL=Vm*V/Km
线性动力学
血药浓度与剂量呈正比 ; AUC与剂量呈正比;t1/2、k、 V、Cl与剂量无关
非线性动力学
Dose-dependant PK 动力学参数与剂量有关 存在饱和现象
k
AUC
t1/2
X0
X0
X0
注:图中实线表示非线性,虚线表示线性非线性药代动力学主要见于:
与药物代谢有关的可饱和的酶代谢过程; 与药物吸收、排泄有关的可饱和的载体转 运过程; 与药物分布有关的可饱和的血浆/组织蛋白 结合过程; 酶诱导及代谢产物抑制等其他特殊过程。
五、非线性动力学参数的求算
1. Km及Vm的求算:根据-dC/dt 求算
dC Vm C dt K m C
Lineweaver-Burk方程式: Hanes-Woolf方程式: Eadie-Hofstee方程式:
高速列车轮轨系统非线性动力学模型数值计算

高速列车轮轨系统非线性动力学模型数值计算高速列车轮轨系统非线性动力学模型数值计算是现代运输领域一个重要的研究课题,也是一个复杂的问题。
本文试图从多个角度来探讨这一问题,包括高速列车轮轨系统的研究背景、相关技术的介绍、常用的数值方法和模型分析等方面。
一、研究背景高速列车的轮轨系统是指列车轮子与铁轨之间的接触面,它是列车的运动控制和路面状况诊断的重要部分。
轮轨系统的非线性动力学模型数值计算是指对这一系统进行精确建模,并利用计算机模拟这一系统的运动和响应过程。
在这一领域,有关的研究内容包括车轮和铁轨的动力学特性、轮轨接触力的计算和分析、轮轨系统的振动问题等。
二、相关技术介绍在高速列车轮轨系统的研究中,有几项技术是非常重要的,它们分别是:1. 车轮铁轨滚动接触力的计算和分析技术:轮轨接触力是指列车在高速运行时轮子与铁轨之间所产生的力,这种力对列车的稳定性和安全性有很大的影响。
要精准地计算这种力,需要细致地研究车轮和铁轨的接触面形状、材料参数、载荷等因素。
2. 应力分析和振动分析技术:在高速列车行驶中,车轮和铁轨都会受到很大的压力和振动力,这些力会导致它们产生应力和振动。
要准确地模拟这些过程,需要运用有关的应力分析和振动分析技术,包括有限元分析、模态分析和频域分析等。
3. 数值计算方法:对于高速列车轮轨系统的非线性动力学模型数值计算,需要使用一系列数值计算方法,包括微分方程数值解法、偏微分方程数值解法、常微分方程数值解法和求解线性代数方程组的方法等。
三、常用的数值方法在高速列车轮轨系统的非线性动力学模型数值计算中,有几种常用的数值方法:1. 基于有限元理论的模拟方法:这种方法利用有限元分析的技术,将轮轨系统的各个部分离散化,然后建立数学模型进行模拟。
这种方法具有高效、精确、适用面广等优点,被广泛应用于车轮铁轨的接触力分析和振动分析中。
2. 先进的逆向设计技术:逆向设计技术是指通过反推物体的形状、轮廓、材质和运动特性等信息,来重新设计物体的技术。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
非线性动力学方程的求解方法1、概述在工程实际问题中,我们常常面临这样的选择:我们所遇的问题究竟是静力的还是动力的。
静力问题与动力问题,从力学的角度看就是是否考虑与加速度有关的力,而从数学求解方法看则是一个三维边值问题还是一个四维边值-初值问题。
在这个问题的选择上没有固定的原则,一般取决于我们研究者、分析者对工程问题的判断。
一般认为,实际工程大都是处于动力环境之中,因而属于动力问题。
但是,由于时间、经费等方面的原因的限制,我们不可能把所有的问题都按照动力问题的方法来分析。
对于许多具体的问题,与速度和加速度有关的力足够小,但是又影响结构分析结果的,将采用静力假定来模拟这些力。
线性的动力有限元控制方程如式(1-1)所示。
[]}{}]{[}]{[}{R q K q D qM =++ (1-1) 式中[M ][D ][K ]分别为结构的质量、阻尼和刚度矩阵,{R }为荷载列矢量,}{q、}{q 和}{q 分别是加速度、速度和位移列矢量。
式(1-1)的解法大体上可以分为两类:直接积分法和模态叠加法。
直接积分法在对控制方程进行数值积分之前不对方程做任何形式的变换,直接用数值积分的方法在时域上一步一步地对方程进行积分。
模态叠加法是在求解之前对方程进行某种数学变换,使基底降低,或使矩阵的带宽减小,再进行求解。
这两种方法在形式上不同,但是密切相关。
上述每一类求解方法中又有许多具体的解法,每一种解法又有各自的特点。
因此我们在选择一种方法求解一个问题时,要对该方法的收敛性、稳定性、效率、精度和费用等进行一些分析,讨论它对所求问题的有效性,从而使我们能够针对某一特定的问题,选择合适的方法。
直接积分法基于以下两条:(1)不是在求解时间区间内任意时刻t 都满足式(1-1),而是在相隔△t 上的一些离散时刻满足式(1-1)。
(2)对位移、速度和加速度在每一时间区间△t 内变化的形式进行假设,事实上若把式(1-1)看成一个常系数微分方程组,便可以用任何一种有限差分格式通过位移来近似表示速度和加速度,因此不同的差分格式就得到不同的方法。
从差分格式上看,分显式与隐式两大类方法,所谓显式差分法就是不必对方程进行求解,而是由前一时刻t 的平衡条件直接可求解△t 增量后t+△t 时刻的各参数解。
而隐式差分法则必须对方程进行求解。
显式差分法有中心差分法等。
隐式差分法有Wilson -θ法和Newmark 法。
不同方法解的精度、稳定性、收敛性、效率不同。
2、中心差分法2.1中心差分法简述中心差分法的差分格式为:}]{}{2}[{1}{2t t t t t t q q q tq ∆+∆-+-∆= (2-1) }]{}{[21}{t t t t t q q tq ∆+∆-+-∆= (2-2) 认为}{t q和}{t q 在任意时刻t 上满足平衡方程(1-1),即: []}{}]{[}]{[}{t t t t R q K q D qM =++ (2-3) 将式(2-1)和式(2-2)代入到式(2-3)中,可得:}]){[21][1(}]){[1]([}{}]{[21][1(222t t t t t t q D tM t q M t K R q D t M t ∆-∆+∆-∆-∆--=∆+∆ (2-4)2.2中心差分法的特点(1)}{t t q ∆+的解是基于利用在时刻t 的平衡条件,即}{t t q ∆+是在假定t 时刻式(2-3)成立的条件下来计算的,该积分过程称为显式积分法,因此中心差分法为显式差分法。
但是,用式(2-4)计算}{t t q ∆+时,不仅要用到}{t q 而且还要用到}{t t q ∆-。
所以,计算在时刻t+△t 的解,必须存在一个具体的起始过程,即必须利用初始条件。
这样,已知}{0q 和}{0q,可由式(2-4)求得}{0q,由式(2-1)和式(2-3)可求得}{t q ∆-,如式(2-5)所示。
}{2}{}{}{0200q t q t q q t ∆+∆-=∆- (2-5) 我们称式(2-5)为差分格式(2-4)的初始条件。
(2)在求解时不需对刚度矩阵[K ]进行三角分解。
在应用此法时一般采用集束质量矩阵或称对角质量阵,而阻尼矩阵也通常为对角形式,这样,在运用式(2-4)时,就不需要对等号右端的系数矩阵t D t M ∆+∆2/][/]([2)进行三角分解,从而可以节省计算时间。
如果不考虑系统的阻尼,则0][=D ,这时式(2-4)可简化为:}{}]){[1(2t t t R q M t =∆∆- (2-6) 其中: }]{[1}]){[2]([}{}{22t t t t t q M t q M t K R R ∆-∆-∆--= (2-7) 从式(2-6)和式(2-7)可以看出,若质量阵[M ]为对角矩阵,即在系统中只考虑集束质量的情况,则方程(2-6)实际上是解耦的,即方程组中是相互独立的,因此只需进行矩阵相乘就可得到。
把式(2-6)改写为:}{][}{12t t t R M t q -∆-∆=写成分量形式,即:)(2)(i t ii i R m t q tt ∆=∆+ (2-8) 其中,)(i t t q ∆+和)(i t R 分别为}{t t q ∆+和}{t R 的第i 个分量,而ii m 为质量矩阵的第i 个对角线元素。
对于0=ii m ,我们可以用静凝聚的方法进行处理。
若对总体刚度和质量矩阵不必进行三角分角,就意味着不必形成总体的[K ]和[M ],因此,式(2-7)的计算只需在单元一级进行即可,即:}){2}]({[1})]{([}{}{2t t t e e e e e t t q q M tq K R R -∆--=∆-∑∑ (2-9) 式中[K e ]和[M e ]分别为与待求节点自由度有关的单元刚度矩阵与质量矩阵。
所以,使用式(2-8)和式(2-9)构成的中心差分具有明显的优点。
它不需要计算总体刚度矩阵[K ]和质量矩阵[M ],求解过程在单元一级进行,计算效率较高。
而且,如果所有单元刚度和质量矩阵都相同时,或有大部分的单元都相同时,可以用一个或少数单元矩阵来求解整个系统。
在工程中,常会遇到重复对称式结构,这时用这种方法来求解将具有很高的效率。
(3)加速度和速度的差分格式都具有(△t 2)阶的误差。
(4)可以证明,中心差分法所要求的步长为: πncr T t t ≤∆≤∆ (2-10)其中,T n 为有限无限无集合体的最小周期,n 为单元系统的阶,cr t ∆为临界时间步长。
周期T n 可利用如下任何一种方法来计算:(1) 广义Jacobi 法;(2) 逆迭代法;(3) 正迭代法;(4) Rayleigh 商迭代法;(5) 多项式迭代法;(6) 斯图姆序列对分法;(7) 行列式搜索法;(8) 子空间迭代法。
事实上,T n 为系统最大特征值的倒数。
要求使用的时间步长△t 小于临界时间步长cr t ∆的差分格式,称为条件稳定的。
若使用一个大于cr t ∆的时间步长,则积分将是不稳定的,此时由数值积分在计算机上导致的舍入误差会增大,从而使结果失去意义。
中心差分法是一种条件稳定的方法,这一缺点导致了它在使用上具有很大的局限性。
首选是时间步长选择困难,当矩阵阶数较高时,而且同时原系统的质量分布不均匀,将使得[M ]矩阵的个别对角上的元素很小,从而使系统的高阶特征趋于无限大,T n 趋于零,这样便使积分过程成为不可能。
所以,仅仅因为一个十分小的质量元素,就会导致计算效率的降低。
当这种情况出现时,通常采用静凝聚的方法使矩阵的阶数降低,但这仅仅能在低频段模拟原系统。
因此,中心差分法是条件稳定的非线性动力学分析方法,有其局限性。
3、Wilson -θ法Wilson -θ法是线性加速度方法的扩展,基于以下两条假定:(1)当时间从t 至变化到时间(t t ∆+θ)时,加速度从t q变化为t t q ∆+θ ; (2)在时间区段t ∆θ内,结构的刚度、阻尼、地面运动加速度都无改变。
则由上述两条假定可得:}){}({}{}{t t t t t q q tq q -∆+=∆++θτθτ (3-1) 式中,θ为参数,可根据积分的稳定性和精度进行选取,一般要求大于1。
对式(3-1)两端的τ进行积分,可得:}){}({2}{}{}{2t t t t t t q q tq q q -∆++=∆++θτθττ (3-2) 再对式(3-2)两端的τ进行积分,可得:}){}({6}{21}{}{}{22t t t t t t t q q tq q q q -∆+++=∆++θτθτττ (3-3) 在式(3-2)和式(3-3)中令t ∆=θτ,则可得在(t t ∆+θ)时刻的速度和位移分别如式(3-4)和式(3-5)所示:}){}({2}{}{t t t t t t q qtq q +∆+=∆+∆+θθθ (3-4) }){2}({6}{}{}{22t t t t t t t q qt t q q q +∆+∆+=∆+∆+θθθθ (3-5) 将加速度}{t t q∆+θ 和速度}{t t q ∆+θ 写成由位移}{t t q ∆+θ表达的形式,如式(3-6)和式(3-7)所示:}{2}{6}){}({6}{22t t t t t t t q q t q q t q -∆--∆=∆+∆+θθθθ (3-6) }{2}{2}{}({3}{t t t t t t t q t q q q t q ∆---∆=∆+∆+θθθθ (3-7) Wilson -θ法认为在(t t ∆+θ)时刻,系统满足动力平衡方程(1-1),并由Wilson -θ法的第二条假定,我们把式(3-6)和式(3-7)代入到式(1-1)中,可得到式(3-8)。
[]}{}]{[}]{[}{t t t t t t t t R q K q D qM ∆+∆+∆+∆+=++θθθθ (3-8) 其中,}){}({}{}{t t t t t t R R R R -+=∆+∆+θθ (3-9)求(t t ∆+θ)时刻的位移、速度和加速度是我们所要得到的。
求解的过程如下所示:(1) 将式(3-6)和式(3-7)代入到式(3-8)中,并加以整理可得:}]){[3][)(6]([2t t q D tM t K ∆+∆+∆+θθθ =}){2}{6}{)(6]([}){}({}{2t t t t t t t q q t q t M R R R +∆+∆+-+∆+θθθ +}){2}{2}{3]([t t t q t q q t D ∆++∆θθ (3-10) 由式(3-10)可以求得}{t t q ∆+θ。
(2) 将式(3-6)分别代入式(3-1)、式(3-2)和式(3-3)中,并令t ∆=τ,可得:}){31(}{)(6}){}({)(6}{22tt t t t t t q q t q q t q θθθθ-+∆--∆=∆+∆+ (3-11) }){}({2}{}{t t t t t t q q t q q +∆+=∆+∆+ (3-12) }){2}({6}{}{}{2t t t t t t t q q t q t q q +∆+∆+=∆+∆+ (3-13) Wilson -θ法的特点:(1) 该法是一种隐式积分法,因为刚度矩阵[K ]是未知位移矢量的系数矩阵。