1、一阶微分方程的求解一阶微分方程的求解3.3电路暂态分析目标是为了得到电路暂态分析目标是为了得到电路时域响应电路时域响应。建立动态电路状态方程,得到一阶微分建立动态电路状态方程,得到一阶微分方程组(或一阶微分方程),再求该方方程组(或一阶微分方程),再求该方程组解。程组解。所以暂态分析实质就是怎样取得而且所以暂态分析实质就是怎样取得而且求解电路常微分方程求解电路常微分方程。第1页一阶微分方程的求解一阶微分方程的求解3.33.3 一阶微分方程求解一阶微分方程求解可归结为在给定初始条件下,一阶微分方程求解可归结为在给定初始条件下,求微分方程初值问题求微分方程初值问题 基本思想:基本思想:在初值问题存
2、在唯一解时间区间内,在若干个时间离在初值问题存在唯一解时间区间内,在若干个时间离散点上,散点上,用差分方程代替微分方程用差分方程代替微分方程,然后逐点求解差分方,然后逐点求解差分方程,得到各时间离散点程,得到各时间离散点 、处函数处函数 近似值近似值 、第2页一阶微分方程的求解一阶微分方程的求解3.3 当两相邻离散点之间间隔较小时,用一阶差商当两相邻离散点之间间隔较小时,用一阶差商取代一阶导数取代一阶导数 一一.前向欧拉法前向欧拉法令步长令步长 ,则,则其近似值为:其近似值为:近似解误差首先是由差商差商近似代替微商微商引发,这种近似代替所产生误差称为截断误差。还有一个误差称为舍入误差,这种误差
3、是因为计算时数值舍入引发。第3页一阶微分方程的求解一阶微分方程的求解3.3前向欧拉法几何意义:在任一步长内,用一段直在任一步长内,用一段直线代替函数线代替函数 曲线,此曲线,此直线段斜率等于该函数在直线段斜率等于该函数在该步长该步长起点起点斜率。斜率。欧拉法几何意义:过点A0(t0,y0),A1(t1,y1),A n(t n,y n),斜率分别为f(t0,y0),f(t1,y1),f(tn,y n)所连接一条折线,所以欧拉法亦称为欧拉折线法欧拉折线法。第4页一阶微分方程的求解一阶微分方程的求解3.3例例1.1.应用前向欧拉法解初值问题应用前向欧拉法解初值问题取步长取步长h=0.1,并把计算结果
4、与准确解比较,并把计算结果与准确解比较解:据前向欧拉法解:据前向欧拉法又又有:有:【思绪思绪】用欧拉法求解常微分方程初值问题时,首先熟练掌握欧拉公式普通用欧拉法求解常微分方程初值问题时,首先熟练掌握欧拉公式普通形式,依据详细题目写出找出欧拉公式迭代式,并依据初始条件和所给步长形式,依据详细题目写出找出欧拉公式迭代式,并依据初始条件和所给步长进行迭代求解。进行迭代求解。第5页一阶微分方程的求解一阶微分方程的求解3.3微分方程微分方程 是一阶线性微分方程,是一阶线性微分方程,可求出其通解:可求出其通解:则方程解为:则方程解为:从而有:从而有:带入初值带入初值 可得可得 一阶非齐次线性微分方程一阶非
5、齐次线性微分方程第6页一阶微分方程的求解一阶微分方程的求解3.3计算结果列表(计算结果列表(为前向欧拉法计算近似值,为前向欧拉法计算近似值,为准确值)为准确值)n 01.000011.10.2718281830.3459198760.07401969321.20.6847555780.8666425360.18188695831.31.2769783441.6072150790.33023673541.42.0935476882.6203595520.52681186451.53.1874451223.9676662950.78022117361.64.6208178465.720961527
6、1.10014368171.76.4663963787.9638734791.497477101正正第7页一阶微分方程的求解一阶微分方程的求解3.3分析:分析:当步长不是很小时,前向欧拉法精度不当步长不是很小时,前向欧拉法精度不是很高。步长取定后,步数越多,误差越是很高。步长取定后,步数越多,误差越大。大。第8页一阶微分方程的求解一阶微分方程的求解3.3二、后向欧拉法用一阶差商近似代替用一阶差商近似代替 在一个在一个步长终点步长终点一阶导数,一阶导数,则原微分方程化为:则原微分方程化为:对于给定初始条件对于给定初始条件微分方程微分方程 其近似值:其近似值:第9页一阶微分方程的求解一阶微分方程的
7、求解3.3在任一步长内,用一段直线在任一步长内,用一段直线代替函数代替函数 曲线,此直曲线,此直线段斜率等于该函数在该线段斜率等于该函数在该步长步长终点终点斜率。斜率。后向欧拉法几何意义:准确值准确值近似值近似值第10页一阶微分方程的求解一阶微分方程的求解3.3注:后向欧拉法两种处理方式注:后向欧拉法两种处理方式 前向前向Euler法为显式,后向法为显式,后向Euler法法为隐式为隐式须解出须解出yk+1.可用迭代法可用迭代法yk+1(n+1)=yk+hf(tk+1,yk+1(n)n=0,1,2,解得解得yk+1,其中其中yk+1(0)=yk+hf(tk,yk).(结合(结合前向欧拉法,预报)
8、前向欧拉法,预报)第11页一阶微分方程的求解一阶微分方程的求解3.3例例2.2.应用后向欧拉法解初值问题应用后向欧拉法解初值问题取步长取步长h=0.1,并把计算结果与准确解比较,并把计算结果与准确解比较解:据后向欧拉法解:据后向欧拉法又第12页一阶微分方程的求解一阶微分方程的求解3.3计算结果列表(计算结果列表(为后向欧拉法计算近似值,为后向欧拉法计算近似值,为准确值)为准确值)n 01.000011.10.4442827750.345919876-0.09836289921.21.1068555350.866642536-0.24021299931.32.0409606121.6072150
9、79-0.43374553341.43.3084097732.620359552-0.68805022151.54.9809113233.967666295-1.01324502861.67.1415858565.720961527-1.42062432971.79.8866975397.963873479-1.922824060负负第13页一阶微分方程的求解一阶微分方程的求解3.3三.梯形法及其预估-矫正法用一阶差商近似地代替函数在一个用一阶差商近似地代替函数在一个步长起点和终点步长起点和终点一阶导数平均值一阶导数平均值 梯形公式梯形公式(欧拉中点公式)(欧拉中点公式)近似值:近似值:改进欧
10、拉法改进欧拉法第14页一阶微分方程的求解一阶微分方程的求解3.3显然,梯形公式是隐式法,普通求显然,梯形公式是隐式法,普通求 需要解方程,需要解方程,常采取迭代法,常采取迭代法,初值初值由显式欧拉公式给出:由显式欧拉公式给出:然后然后将将 替换梯形公式等式右边出现替换梯形公式等式右边出现当步长当步长h足够小,且由前向欧拉法计算已是很好足够小,且由前向欧拉法计算已是很好近似,则近似,则迭代一、二次即可迭代一、二次即可预报预报校正校正迭代次数迭代次数第15页一阶微分方程的求解一阶微分方程的求解3.3l几何意义几何意义lEuler法法折线法折线法l改进改进Euler法法平均斜率折线法平均斜率折线法第
11、16页一阶微分方程的求解一阶微分方程的求解3.3例例3.应用梯形预估应用梯形预估-矫正法解初值问题矫正法解初值问题取步长取步长h=0.1,并把计算结果与准确解比较,并把计算结果与准确解比较解:据前向欧拉法解:据前向欧拉法梯形预估梯形预估-矫正矫正第17页一阶微分方程的求解一阶微分方程的求解3.3计算结果列表(计算结果列表(为梯形预估为梯形预估-矫正法计算矫正法计算近似值,近似值,为准确值)为准确值)k 01.000011.10.3423777890.3459198760.00354208721.20.8583145370.8666425360.00832799931.31.5927496431
12、.6072150790.01446543641.42.5982982392.6203595520.02206131351.53.9364441143.9676662950.03122218161.65.6789071035.7209615270.04205442471.77.9092092167.9638734790.054664263第18页一阶微分方程的求解一阶微分方程的求解3.3function T Y=Trapezia_reckon(odefun,ab,ya,M)%odefun:微分方程微分方程 a、b:计算区间:计算区间%ya:初值:初值 y(a)M:等分数目:等分数目%T:离散时间
13、变量离散时间变量 Y梯形公式预估校正法解梯形公式预估校正法解h=(ab(2)-ab(1)/M;步长步长T=zeros(1,M+1);Y=zeros(1,M+1);T=ab(1):h: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 Function 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)第19页一阶微分方程
14、的求解一阶微分方程的求解3.3求解器求解器求解求解问题问题特点特点说说明明ode45非非刚刚性性一步算法;一步算法;4,5阶阶Runge-Kutta算法算法大部分场所首选算法ode23非非刚刚性性一步算法;一步算法;2,3阶阶Runge-Kutta算法算法使用于精度较低情形ode113非非刚刚性性多步法;变阶次Adams-Bashforth-Moulton 算法计计算算时间时间比比ode45短短ode23t刚刚性性采取梯形算法适合中度刚性问题求解ode15s刚刚性性多步法;采取了数值差分算法若若ode45失效失效时时,可,可尝试尝试使使用用;ode23s刚刚性性一步法;一步法;2阶阶Roseb
15、rock算法算法当精度当精度较较低低时时,计计算算时间时间比比ode15s短短ode23tb刚刚性性隐隐式式Runge-Kutta算法算法当精度当精度较较低低时时,计计算算时间时间比比ode15s短短不一样求解器特点不一样求解器特点 第20页一阶微分方程的求解一阶微分方程的求解3.3在用常微分方程描述一个电路暂态过程时,往往又包含着多个改变速度相差十分悬殊子过程,这么一类过程就认为含有“刚性(stiff)”,描述这类过程微分方程称为“刚性问题”。比如,电路某一变量以e-t迟缓衰减,而另一变量以e-1000t快速衰减,两变量时间常数相差很大,建立常微分方程就含有“刚性”。刚性问题数值解稳定性通常被最快模式控制,刚性问题解答难度就在于其快变子过程干扰。当我们试图在慢变区间上求解刚性问题时,尽管快变分量值已衰减到微不足道,但这种快速改变干扰仍严重影响数值解稳定性和精度,普通地说,隐型方法比显型方法含有更大稳定性,所以使用隐型方法求解刚性方程组更为适当.在MATLAB中,ode23t、ode15s、ode23s、ode23tb适合求解刚性问题。第21页