龙格—库塔法分析Lorenz方程课程设计.doc

龙格—库塔法分析Lorenz方程课程设计.doc《数值分析》课程设计实验报告龙格—库塔法分析Lorenz方程问题叙述考虑著名的Lorenz方程其中s,r,b为变化区域内有一定限制的实参数,该方程形式简单,表面上看并无惊人之处,但由该方程揭示出的许多现象,促使“混沌”成为数学研究的崭新领域,在实际应用中也产生了巨大的影响。问题分析Lorenz方程实际上是一个四元一阶常微分方程,用解析法精确求解是不可能的,只能用数值计算,最主要的有欧拉法、亚当法和龙格-库塔法等。为了得到较高精度的,我们采用经典四阶龙格—库塔方法求解该问题。实验程序及注释(1)算法程序function[T]=Runge_Kutta(f,x0,y0,h,n)%定义算法,其中f为待解方程组,x0是初始自变量,y0是初始函数值,h是步长,n为步数ifnargin<5n=100;%如果输入参数个数小于5,则步数n=100endr=size(y0);r=r(1);%返回初始输出矩阵的行列数,并将值赋给r(1)s=size(x0);s=s(1);%返回初始输入矩阵的行列数,并将值赋给s(1)r=r+s;T=zeros(r,n+1);T(:,1)=[y0;x0];fort=2:n+1%以下是具体的求解过程k1=feval(f,T(1:r-1,t-1));k2=feval(f,[k1*(h/2)+T(1:r-1,t-1);x0+h/2]);k3=feval(f,[k2*(h/2)+T(1:r-1,t-1);x0+h/2]);k4=feval(f,[k3*h+T(1:r-1,t-1);x0+h]);x0=x0+h;T(:,t)=[T(1:r-1,t-1)+(k1+k2*2+k3*2+k4)*(h/6);x0];end(2)主程序functiondy=fun(x)%定义函数s=;%给参数s,r,b赋值r=;b=;dy(1)=s*(x(2)-x(1));%Lorenz方程表达式dy(2)=(r*x(1)-x(3)*x(1)-x(2));dy(3)=x(1)*x(2)-b*x(3);dy=dy';(2)运行程序T=Runge_Kutta('fun',0,[10;10;10],,5000);%调用前面的算法程序plot3(T(1,:),T(2,:),T(3,:));%显示三分量的关系图axis([-2020-5050050])%定义坐标轴长度view(3)%设定观察角度实验数据结果及分析(1)各初始变量相同时的图像分析各初始变量取相同的值[10,10,10],运行上述程序后,得到如下图像:从图中可以看出,各初始变量相同时,曲线总是被吸引回奇怪吸引子附近作来回跳跃。初始变量值取为[-10,-10,-10],[20,20,20]时,依然如此。图像如下:[-10,-10,-10][20,20,20](2)初始值的每个分量变化对图像的影响y分量:[0,2,0][0,5,0][0,15,0][0,20,0]从上面可以看出,随着初始y值的增大,奇怪吸引子中曲线在其附近来回跳跃的两个位置中的一个吸引力变弱,另一个吸引力变强。初始y继续增大到某一特定值,情况又会变回来。这说明在空间存在一些区域,当初始位置位于这些区域外时解将出现奇怪吸引子的性质,而在这些区域以内解将呈现普通吸引子的性质。z分量:[0,0,20]从上图可以看出解的曲线为一直线,