版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领
文档简介
数值积分法仿真第1页,课件共63页,创作于2023年2月Overview数值积分方法的原理是什么?病态系统的特点和仿真算法选取?算法的稳定性分析?第2页,课件共63页,创作于2023年2月第一节数字仿真原理在连续系统的仿真中,数值积分法可分为两大类:单步法:以龙格-库塔法为代表多步法:以Adams法为代表数值积分法的要素:基本特性:稳定性空间特性:精度时间特性:速度第3页,课件共63页,创作于2023年2月数值积分基本原理连续系统的仿真,主要是对一阶微分方程(组)的求解可见仿真关键是对Qm准确,快速的求解第4页,课件共63页,创作于2023年2月步长:将时间t离散t(k)(k=1,2,…n),相邻两点的距离为步长,即h=t(k+1)-t(k)步进法:数值积分法求近似解根据初始值y0,按照离散的时间序列步进求解。 t0t1t2t3…tn y0y1y2y3…tn计算格式:由y(k)计算出y(k+1)(k=0,1,…,n)的递推公式。数值积分基本名词第5页,课件共63页,创作于2023年2月数值积分的基本性能基本性能数值积分算法的性能包含:定性特征:稳定性时间粒度:计算速度空间粒度:计算精度不同的数值积分方法的具有不同的稳定性。同一个模型采用不同的积分算法和不同的积分步长h,稳定性不同。第6页,课件共63页,创作于2023年2月计算速度和计算精度各种数值积分方法的差分方程是对原微分方程的近似逼近,并且因为计算机的字长有限,存在明显的截断误差。这些误差都和计算步距h密切相关,所以计算步距是影响计算精度、速度和稳定性的重要因素。h取得较大,计算时间少,截断误差大;h取得较小,截断误差就会减小,但在给定时间范围内,计算次数必然增加,使误差积累增加。第7页,课件共63页,创作于2023年2月截断误差、累计舍入误差与步长h截断误差、累计舍入误差与步长h关系如图。图中可知,两种误差对步距的要求是矛盾的,但两者之和有一个最小值,步距最好能选在最小值。然而,实际要做到这一点是很困难的。一般只能根据经验确定一个合理步长区,通常将步长h限制在系统的最小时间常数数量级上。第8页,课件共63页,创作于2023年2月引理:泰勒级数:如果f(x)在x0点处任意阶可导,则在该邻域内的n阶泰勒公式为:第二节单步法单步数值积分法的核心就是泰勒级数的近似。第9页,课件共63页,创作于2023年2月2.1一阶欧拉法对于一阶微分方程故一般的一阶欧拉法递推形式为:第10页,课件共63页,创作于2023年2月一阶欧拉法图示第11页,课件共63页,创作于2023年2月2.22阶龙格-库塔对于一阶微分方程第12页,课件共63页,创作于2023年2月2阶龙格-库塔第13页,课件共63页,创作于2023年2月2阶龙格-库塔第14页,课件共63页,创作于2023年2月2阶龙格-库塔故一般的二阶龙格-库塔法递推形式为:2阶龙格-库塔只取到泰勒级数展开式中y的二阶导数项,略去了三阶以上高阶导数项。其截断误差正比于步长h3为纪念提出该方法的德国数学家C.Runge和M.W.Kutta,称这种计算方法为二阶龙格-库塔法。第15页,课件共63页,创作于2023年2月2阶龙格-库塔图示第16页,课件共63页,创作于2023年2月比较第17页,课件共63页,创作于2023年2月高阶龙格-库塔(RK-4)一般在计算精度要求较高的情况下,多使用四阶龙格-库塔法。其计算公式为,其截断误差正比于步长h5第18页,课件共63页,创作于2023年2月高阶龙格-库塔(RK-4)第19页,课件共63页,创作于2023年2月单步法的特点单步法以上介绍的几种数值积分公式,有一个共同的特点,由于采用了泰勒级数展开,在本次计算中,仅仅用到前一步的计算结果,而不需要利用更前面各步的结果。这类计算方法称为单步法。单步法运算有下列优点:(1)需要存储的数据量少,占用的存储空间少。(2)给定初值,就可启动递推公式进行运算(自启动计算能力)(3)容易实现变步长第20页,课件共63页,创作于2023年2月第三节变步长龙格-库塔法步长控制是实现高精度的仿真算法的手段之一。实现步长控制涉及:局部误差估计步长控制策略第21页,课件共63页,创作于2023年2月3.1误差估计通常设法寻找一个低一阶的龙格-库塔公式,两者的结果之差可以设为误差。为减少计算量,Ki通常要求公用。Runge-Kutta-Merson法(RK34)第22页,课件共63页,创作于2023年2月Runge-Kutta-Fehlberg(RK45)计算公式为5阶6级,误差估计低阶公式为4阶五级,具有四阶误差估计和五阶精度,称为RK45法。RK45被公认为对非病态系统仿真的最有效的方法之一。第23页,课件共63页,创作于2023年2月Runge-Kutta-Fehlberg(RK45)iaibijcic*i1016/13525/21621/41/40033/83/329/326656/128251408/2565412/131932/2197-7200/21977296/219728561/564302197/410451439/216-83680/513-847/4104-9/50-1/561/2-8/272-3544/25661859/4104-11/402/550第24页,课件共63页,创作于2023年2月RKF-12第25页,课件共63页,创作于2023年2月RKS-34(1978,Shamping)第26页,课件共63页,创作于2023年2月3.2步长控制步长控制策略一般分为:1)加倍-减半法2)最优步长法第27页,课件共63页,创作于2023年2月步长控制:加倍-减半法加倍-减半法第28页,课件共63页,创作于2023年2月步长控制:最优步长法第29页,课件共63页,创作于2023年2月1.2.2步长控制第30页,课件共63页,创作于2023年2月龙格-库塔方法的一般形式各种龙格-库塔法的公式都由两部分组成,一个是上一步结果,另一个是步长乘以各点导数的加权和。平均斜率第31页,课件共63页,创作于2023年2月第三节线性多步法单步法运算基于泰勒级数展开法,其特点是:(1)需要存储的数据量少,占用的存储空间少。(2)给定初值(t0,y0),就可启动递推公式进行运算(自启动计算能力。(3)容易实现变步长积分。可有效平衡计算速度和精度之间的矛盾。多步法的基本原理是多项式拟合利用一个多项式取匹配变量的若干已知值和各阶导数。第32页,课件共63页,创作于2023年2月线性多步法原理第33页,课件共63页,创作于2023年2月3.1预报公式第34页,课件共63页,创作于2023年2月3.1预报公式第35页,课件共63页,创作于2023年2月3.1预报公式第36页,课件共63页,创作于2023年2月预报举例第37页,课件共63页,创作于2023年2月3.2校正公式第38页,课件共63页,创作于2023年2月3.2校正公式第39页,课件共63页,创作于2023年2月3.2校正公式第40页,课件共63页,创作于2023年2月预报-校正举例第41页,课件共63页,创作于2023年2月预报-校正举例第42页,课件共63页,创作于2023年2月3.3Adams公式根据前面的分析,我们可以将预报和校正公式统一写成:第43页,课件共63页,创作于2023年2月显式Adams系数第44页,课件共63页,创作于2023年2月隐式Adams系数第45页,课件共63页,创作于2023年2月3.4多步法的特点与单步法相比,相同精度下,使用过去多步信息,计算量小。隐式法的精度高,稳定性好,但在计算y(n+k)时需要用到f[y(n+k),t(n+k)],只能采用迭代法计算。缺点之一是不能自启动,需用单步法计算初始值才能启动计算。第46页,课件共63页,创作于2023年2月第四节积分算法的稳定性稳定的系统采用不同的积分算法,其稳定性不同稳定性的测试公式为:当第47页,课件共63页,创作于2023年2月4.1一阶Adams法的稳定性分析第48页,课件共63页,创作于2023年2月一阶Adams法的稳定性分析第49页,课件共63页,创作于2023年2月4.2一般算法的稳定性分析根据上例可得数值积分方法稳定域的一般方法。 设系统测试方程为: 而数值积分公式为: 只有当时,算法稳定。各种数值积分算法的稳定域参见书P96图3.9第50页,课件共63页,创作于2023年2月主要算法的稳定性一阶、二阶Admas法为恒稳算法,其他算法条件稳定。除恒稳法外,其他算法的步长h必须限制在最小时间的数量级对龙格-库塔法,阶次k增大,稳定域略微增大。对Admas法,阶次k增大,稳定域反而缩小。第51页,课件共63页,创作于2023年2月第五节Matlab实现ODE(OrdinaryDifferentaialequation)解法模型描述算法描述算法仿真第52页,课件共63页,创作于2023年2月微分方程模型描述Lorenz曲线filename:mdLorenz.m
functiondx=mdLorenz(t,x)dx=[-8/3*x(1)+x(2)*x(3);-10*x(2)+10*x(3);-x(1)*x(2)+28*x(2)-x(3)];end第53页,课件共63页,创作于2023年2月数值积分算法描述matlab中的数值积分算法函数的格式如下:function[tout,yout]=solver(ModelName,tspan,x0,option)第54页,课件共63页,创作于2023年2月数值积分算法描述%一阶Euler算法,filename:svEulerfunction[tout,yout]=svEuler(odeFcn,tspan,y0)t0=tspan(1);t1=tspan(2);iflength(tspan)<3,h=(t1-t0)/1000;elseh=tspan(3);tout=[t0:h:t1]';N=length(y0);M=length(tout)-1;tout=[t0:h:t1]';yout=[y0';zeros(M,N)];fori=1:Mk1=h*feval(odeFcn,tout(i),y0);y0=y0+k1;yout(i+1,:)=y0';endend第55页,课件共63页,创作于2023年2月数值积分算法描述function[tout,yout]=svRungeKutta4(odeFcn,tspan,y0)t0=tspan(1);t1=tspan(2);iflength(tspan)<3,h=(t1-t0)/1000;elseh=tspan(3);tout=[t0:h:t1]';N=length(y0);M=length(tout)-1;tout=[t0:h:t1]';yout=[y0';zeros(M,N)];fori=1:Mk1=h*feval(odeFcn,tout(i),y0);k2=h*feval(odeFcn,tout(i)+h/2,y0+0.5*k1);k3=h*feval(odeFcn,tout(i)+h/2,y0+0.5*k2);k4=h*feval(odeFcn,tout(i)+h,y0+k3);y0=y0+(k1+2*k2+2*k3+k4)/6;yout(i+1,:)=y0';endend第56页,课件共63页,创作于2023年2月微分方程模型描述在Matlab文件中调用方法为:[t,y]=svEuler(@eqLorenz,[0,100],x0);第57页,课件共63页,创作于2023年2月ode45ODE45Solvenon-stiffdifferentialequations,mediumordermethod.[TOUT,YOUT]=ODE45(ODEFUN,TSPAN,Y0)withTSPAN=[T0TFINAL]integratesthesystemofdifferentialequationsy'=f(t,y)fromtimeT0toTFINALwithinitialconditionsY0.ODEFUNisafunctionhandle.ForascalarTandavectorY,ODEFUN(T,Y)mustreturnacolumnvectorcorrespondingtof(t,y).EachrowinthesolutionarrayYOUTcorrespondstoatimereturnedinthecolumnvectorTOUT.ToobtainsolutionsatspecifictimesT0,T1,...,TFINAL(allincreasingoralldecreasing),useTSPAN=[T0T1...TFINAL].第58页,课件共63页,创作于2023年2月ode45ODE45Solvenon-stiffdifferentialequations,mediumordermethod.
[TOUT,YOUT]=ODE45(ODEFUN,TSPAN,Y0,OPTIONS)solvesasabovewithdefaultintegrationpropertiesreplacedbyvaluesinOPTIONS,anargumentcreatedwiththeODESETfunction.SeeODESETfordetails.Commonlyusedoptionsarescalarrelativeerrortolerance'RelTol'(1e-3bydefault)andvectorofabsoluteerrortolerances'AbsTol'(allcomponents1e-6bydefault).Ifcertaincomponentsofthesolutionmustbenon-negative,useODESETtosetthe'NonNegative'propertytotheindicesofthese第59页,课件共63页,创作于2023年2月ode45ODE45Solvenon-stiffdifferentialequations,mediumordermethod.
ODE45cansolveproblemsM(t,y)*y'=f(t,y)withmassmatrixMthatisnonsingular.UseODESETtosetthe'Mass'propertytoafunctionhandleMASSifMASS(T,Y)returnsthevalueofthemassmatrix.Ifthemassmatrixisconstant,thematrixcanbeusedasthevalueofthe'Mass'option.IfthemassmatrixdoesnotdependonthestatevariableYandthefunctionMASSistobecalledwithoneinputargumentT,set'MStateDependence'to'none'.ODE15SandODE23Tcansolveproblemswithsingularmassmatrices.第60页,课件共63页,创作于2023年2月ode45
ODE45Solvenon-stiffdifferentialequations,mediumordermethod.Example[t,y]=ode45(@vdp1,[020],[20]);plot(t,y(:,1));solvesthesystemy'=vdp1(t,y),usingthedefaultrelativeerrortolerance1e-3andthedefaultabsolutetoleranceof1e-6foreachcomponent,andplotsthefirstcomponentofthesolution.
ClasssupportforinputsTSPAN,Y0,andtheresultofODEFUN(T,Y):float:double,single
温馨提示
- 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
- 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
- 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
- 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
- 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
- 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
- 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。
最新文档
- 质量策划方案20241219
- (案例)标准件项目立项报告
- 2022-2023学年天津市高一(上)期末语文试卷
- 2022-2023学年天津四中高二(上)期末语文试卷
- 重庆2020-2024年中考英语5年真题回-学生版-专题07 阅读理解之说明文
- 2023年碳酸丙烯酯项目筹资方案
- 2023年美容项目筹资方案
- 2023年大功率激光传输石英光纤项目筹资方案
- 复数的加法和减法课件
- 养老院老人心理咨询师职业道德制度
- 2023年上海市中考化学试卷真题(含答案与解析)
- 关于市场化工具及其在中国的运用概述
- 监控维保报告
- 足部健康宣教课件
- 生产准备部进出口报关员关键业绩考核指标(KPI)
- TSZSA 032-2023 SMD 塑料载带技术规范
- 应用文写作《通知》课件 2023-2024学年高教版中职语文基础模块下册
- 车险续保率分析报告
- 基于Android的天气预报系统的设计与实现
- 省级代理合同
- 精神科患者饮食健康宣教知识
评论
0/150
提交评论