3.3电路一阶微分方程的求解

1、一阶微分方程的求解一阶微分方程的求解3.3电路暂态分析的目的是为了得到电路暂态分析的目的是为了得到电路的时域响应电路的时域响应。 建立动态电路的状态方程，得到一阶微建立动态电路的状态方程，得到一阶微分方程组（或一阶微分方程），再求该分方程组（或一阶微分方程），再求该方程组的解。方程组的解。因此暂态分析的实质就是如何获得并因此暂态分析的实质就是如何获得并且且求解电路的常微分方程求解电路的常微分方程。 一阶微分方程的求解一阶微分方程的求解3.33.3 一阶微分方程的求解 000)(),()( ytytttyfty一阶微分方程的求解可归结为在给定初始条件下，一阶微分方程的求解可归结为在给定初始条件下

2、，求微分方程的初值问题求微分方程的初值问题 基本思想：基本思想： 在初值问题存在唯一解的时间区间内，在若干个时间在初值问题存在唯一解的时间区间内，在若干个时间离散点上，离散点上，用差分方程代替微分方程用差分方程代替微分方程，然后逐点求解差分，然后逐点求解差分方程，得到各时间离散点方程，得到各时间离散点 、 处的函数处的函数 近似近似值值 、 1t2tnt1y2yny)(ty一阶微分方程的求解一阶微分方程的求解3.3 当两相邻离散点之间的间隔较小时，用一阶差商当两相邻离散点之间的间隔较小时，用一阶差商取代一阶导数取代一阶导数 一一. .前向欧拉法前向欧拉法)( )()(11kkkkktyttty

3、ty htytytykkk)( )()(1 hyyykkk1 kktth 1令步长令步长 ，则，则其近似值为：其近似值为： 近似解的误差首先是由差商差商近似代替微商微商引起的，这种近似代替所产生的误差称为截断误差。还有一种误差称为舍入误差，这种误差是由于计算时数值舍入引起的。一阶微分方程的求解一阶微分方程的求解3.3前向欧拉法的几何意义： y(t) y1 y2 y3 y(t3) y(t2) y(t1) y0 y(t0) h h h t0 0 t1 t2 t3 在任一步长内，用一段直在任一步长内，用一段直线代替函数线代替函数 的曲线，的曲线，此直线段的斜率等于该函此直线段的斜率等于该函数在该步长

4、数在该步长起点起点的斜率。的斜率。 )(ty),(1kkkkythfyy 欧拉法的几何意义：过点a0（t0，y0），a1（t1，y1），a n（t n，y n ），斜率分别为f（t0，y0），f（t1，y1），f（tn，y n）所连接的一条折线，所以欧拉法亦称为欧拉折线法欧拉折线法。一阶微分方程的求解一阶微分方程的求解3.3例例1. 1. 应用前向欧拉法解初值问题应用前向欧拉法解初值问题0)1(, 21 ,22 ytetytyt)2(21ntnnnnnetythyy nnhttyyn1 . 010)1(00 取步长取步长h=0.1，并把计算结果与精确解比较，并把计算结果与精确解比较解：据前向欧

5、拉法解：据前向欧拉法又又有：有：271828183. 0)2(0200001 tetythyy684755578. 0)2(1211112 tetythyy【思路思路】 用欧拉法求解常微分方程的初值问题时，首先熟练掌握欧拉公式的用欧拉法求解常微分方程的初值问题时，首先熟练掌握欧拉公式的一般形式，根据具体题目写出找出欧拉公式的迭代式，并根据初始条件和所一般形式，根据具体题目写出找出欧拉公式的迭代式，并根据初始条件和所给步长进行迭代求解。给步长进行迭代求解。一阶微分方程的求解一阶微分方程的求解3.3微分方程微分方程 是一阶线性微分方程，是一阶线性微分方程，可求出其通解：可求出其通解： tetyty

6、22 c)(2tety )-(2eetyt )-()(2eettyntnn 则方程的解为：则方程的解为： 从而有：从而有： 345919876. 0)(1211 eetyt86642536. 0)(2222 eetyt带入初值带入初值 可得可得 0) 1 ( ye c一阶非齐次线性微分方程一阶非齐次线性微分方程一阶微分方程的求解一阶微分方程的求解3.3计算结果列表（计算结果列表（ 为前向欧拉法计算近似值，为前向欧拉法计算近似值， 为精确值）为精确值） n 01.000011.10.2718281830.3459198760.07401969321.20.6847555780.866642536

7、031.2769783441.6072150790.33023673541.42.0935476882.6203595520.52681186451.53.1874451223.9676662950.78022117361.64.6208178465.7209615271.10014368171.76.4663963787.9638734791.497477101ntny)(ntynnyty )(ny)(nty正正一阶微分方程的求解一阶微分方程的求解3.3分析：分析：当步长不是很小时，前向欧拉法的精度不当步长不是很小时，前向欧拉法的精度不是很高。步长取定后，步数越多，

8、误差越是很高。步长取定后，步数越多，误差越大。大。一阶微分方程的求解一阶微分方程的求解3.3二、后向欧拉法)( )()(11 kkktyhtyty欧欧拉拉隐隐式式公公式式 hyyykkk11用一阶差商近似代替用一阶差商近似代替 在一个在一个步长终点步长终点的一阶导的一阶导数，则原微分方程化为：数，则原微分方程化为：),()( ttyfty )(ty00)(yty 对于给定初始条件对于给定初始条件的微分方程的微分方程其近似值：其近似值： 一阶微分方程的求解一阶微分方程的求解3.3在任一步长内，用一段直线在任一步长内，用一段直线代替函数代替函数 的曲线，此直的曲线，此直线段的斜率等于该函数在该线段

9、的斜率等于该函数在该步长步长终点终点的斜率。的斜率。 后向欧拉法的几何意义： y(t) y1 y2 y3 y(t3) y(t2) y(t1) y0 y(t0) h h h t0 0 t1 t2 t3 t y(t) ),(111 kkkkythfyy精确值精确值近似值近似值一阶微分方程的求解一阶微分方程的求解3.3注：后向欧拉法的两种处理方式注：后向欧拉法的两种处理方式 前向前向euler法为显式，后向法为显式，后向euler法法为隐式为隐式须解出须解出yk+1. 可用迭代法可用迭代法yk+1 (n+1) = yk + hf (tk+1，yk+1(n)n = 0，1，2， 解得解得yk+1 ，其

10、中其中yk+1(0) = yk + hf (tk，yk).（结合前（结合前向欧拉法，预报）向欧拉法，预报）一阶微分方程的求解一阶微分方程的求解3.3例例2. 2. 应用后向欧拉法解初值问题应用后向欧拉法解初值问题0)1(,21 ,22 ytetytyt12112111121:)2(11即即ntnnntnnnnnthehtyyetythyynn nnhttyyn1 . 010)1(00 取步长取步长h=0.1，并把计算结果与精确解比较，并把计算结果与精确解比较解：据后向欧拉法解：据后向欧拉法又一阶微分方程的求解一阶微分方程的求解3.3计算结果列表（计算结果列表（ 为后向欧拉法计算近似值，为后向欧

11、拉法计算近似值， 为精确值）为精确值） n 01.000011.10.4442827750.345919876-0.09836289921.21.1068555350.866642536-0.24021299931.32.0409606121.607215079-0.43374553341.43.3084097732.620359552-0.68805022151.54.9809113233.967666295-1.01324502861.67.1415858565.720961527-1.42062432971.79.8866975397.963873479-1.922824060ntny)

12、(ntynnyty )(ny)(nty负负一阶微分方程的求解一阶微分方程的求解3.3三. 梯形法及其预估-矫正法 000)(),()( ytyttytfty)( )( 21)()(11 kkkktytyhtyty)( )( 2)()(11 kkkktytyhtyty用一阶差商近似地代替函数在一个用一阶差商近似地代替函数在一个步长起点和终点的步长起点和终点的一阶导数的平均值一阶导数的平均值 梯形公式梯形公式（欧拉中点公式）（欧拉中点公式）近似值：近似值：)(211 kkkkyyhyy改进欧拉法改进欧拉法一阶微分方程的求解一阶微分方程的求解3.3显然，梯形公式是隐式法，一般求显然，梯形公式是隐式法

13、，一般求 需要解方程，需要解方程，常采用迭代法，常采用迭代法，初值初值由显式的欧拉公式给出：由显式的欧拉公式给出：1 ky),()0(1kkkkythfyy ),(),(2)0(11)1(1 kkkkkkytfytfhyy然后然后将将 替代梯形公式等式右边出现的替代梯形公式等式右边出现的)0(1 ky1 ky),(),(2)n(11)1(1kkkkknkytfytfhyy当步长当步长h足够小，且由前向欧拉法计算的已是较好足够小，且由前向欧拉法计算的已是较好的近似，则的近似，则迭代一、二次即可迭代一、二次即可, 2 , 1 , 0n预报预报校正校正迭代次数迭代次数一阶微分方程的求解一阶微分方程的

14、求解3.3l几何意义几何意义leuler法法折线法折线法l改进改进euler法法平均斜率折线法平均斜率折线法一阶微分方程的求解一阶微分方程的求解3.3例例3. 应用梯形预估应用梯形预估-矫正法解初值问题矫正法解初值问题0)1(,21 ,22 ytetytyt取步长取步长h=0.1，并把计算结果与精确解比较，并把计算结果与精确解比较kktkkkkketythyyy201 . 0)2 . 01(1 ）（解：据前向欧拉法解：据前向欧拉法梯形预估梯形预估- -矫正矫正2205. 0),(),(21121)0(12)0(111）（）（ kkktkktkkkkkkkkkketytetytyytfytfhy

15、y一阶微分方程的求解一阶微分方程的求解3.3计算结果列表（计算结果列表（ 为梯形预估为梯形预估- -矫正法计算矫正法计算近似值，近似值， 为精确值）为精确值）ky)(nty k 01.000011.10.3423777890.3459198760.00354208721.20.8583145370.8666425360.00832799931.31.5927496431.6072150790.01446543641.42.5982982392.6203595520.02206131351.53.9364441143.9676662950.03122218161.65.6789071035.72

16、09615270.04205442471.77.9092092167.9638734790.054664263ktky)(ktykkyty)(一阶微分方程的求解一阶微分方程的求解3.3function t y=trapezia_reckon(odefun,ab,ya,m)% odefun: 微分方程微分方程 a、b：计算区间：计算区间% ya：初值：初值 y(a) m：等分数目：等分数目% t: 离散的时间变量离散的时间变量 y梯形公式的预估校正法解梯形公式的预估校正法解h=(ab(2)-ab(1)/m; 步长步长t=zeros(1,m+1);y=zeros(1,m+1);t=ab(1):h:

17、ab(2); y(1)=ya;for j=1:m k1=feval(odefun,t(j),y(j); k2=feval(odefun,t(j+1),y(j)+h*k1); y(j+1)=y(j)+(h/2)*(k1+k2);end )(2)0(11 jjjjyyhyyfunction y=euler_3_3_2(t,x)y=2/t*x+t2*exp(t ）t y=trapezia_reckon ( euler_3_3_2,1 2,0,10) 一阶微分方程的求解一阶微分方程的求解3.3求解器求解器求解问题求解问题特点特点说明说明ode45非刚性非刚性一步算法；一步算法；4，5阶阶runge-k

18、utta算法算法大部分场合的首选算法大部分场合的首选算法ode23非刚性非刚性一步算法；一步算法；2，3阶阶runge-kutta算法算法使用于精度较低的情形使用于精度较低的情形ode113非刚性非刚性多步法；变阶次的多步法；变阶次的adams -bashforth- moulton 算法算法计算时间比计算时间比ode45短短ode23t刚性刚性采用梯形算法采用梯形算法适合中度刚性问题的求解适合中度刚性问题的求解ode15s刚性刚性多步法；采用了数值差多步法；采用了数值差分算法分算法若若ode45失效时，可尝试失效时，可尝试使用使用;ode23s刚性刚性一步法；一步法；2阶阶rosebrock算法算法当精度较低时，计算时间当精度较低时，计算时间比比ode15s短短ode23tb刚性刚性隐式隐式runge-kutta算法算法当精度较低时，计算时间当精度较低时，计算时间比比ode15s短短不同求解器的特点不同求解器的特点 一阶微分方程的求解一阶微分方程的求解3.3在用常微分方程描述一个电路的暂态过程时，往往又包含着多个变化速度相差十分悬殊的子过程，这样一类过程就认为