数值分析上机第四次作业.doc

数值分析上机第四次作业实验项目:共轭梯度法求解对称正定得线性方程组实验内容:用共轭梯度法求解下面方程组(1) 迭代20次或满足时停止计算。(2) ,A就是1000阶得Hilbert矩阵或如下得三对角矩阵,Ai,i=4,Ai,i-1=Ai-1,i=-1,i=2,3,、,nb1=3, bn=3, bi=2,i=2,3,n-1迭代10000次或满足时停止计算。(3)*考虑模型问题,方程为 用正方形网格离散化,若取,得到得线性方程组,并用共轭梯度法(CG法)求解,并对解作图。实验要求:迭代初值可以取,计算到停止.本题有精确解,这里表示以为分量得向量,表示在相应点上取值作为分量得向量.实验一:(1)编制函数子程序CGmethod。function x,k=CGmethod(A,b)n=length(A);x=zeros(n,1);r=b-A*x;rho=r*r;k=0;while rho10(-12) & k=10(-7)&k104 k=k+1; if k=1 p=r; else beta=(r1*r1)/(r*r);p=r1+beta*p; end r=r1; w=A*p; alpha=(r*r)/(p*w); x=x+alpha*p; r1=r-alpha*w;end编制主程序shiyan1_2:clear,clcn=1000;A=hilb(n);b=sum(A);x,k=CGmethod_1(A,b)运行结果为:x得值,均接近1,迭代次数k=32实验二 实验目得:用复化Simpson方法、自适应复化梯形方法与Romberg方法求数值积分。实验内容:计算下列定积分(1) (2) (3) 实验要求:(1)分别用复化Simpson公式、自适应复化梯形公式计算要求绝对误差限为,输出每种方法所需得节点数与积分近似值,对于自适应方法,显示实际计算节点上离散函数值得分布图;(2)分析比较计算结果。2、实验目得:高斯数值积分方法用于积分方程求解。实验内容:线性得积分方程得数值求解,可以被转化为线性代数方程组得求解问题。而线性代数方程组所含未知数得个数,与用来离散积分得数值方法得节点个数相同在节点数相同得前提下,高斯数值积分方法有较高得代数精度,用它通常会得到较好得结果。对第二类Fredholm积分方程首先将积分区间a,b等分成n份,在每个子区间上离散方程中得积分就得到线性代数方程组。实验要求:分别使用如下方法,离散积分方程中得积分1、复化梯形方法;2、复化辛甫森方法;3、复化高斯方法。求解如下得积分方程,方程得准确解为,并比较各算法得优劣。实验二1、复化Simpson方法)输入积分区间下限0输入积分区间上限2输入等分份数20输入被积函数(以x为自变量)x6/10-x2+xS =1、1619输入积分区间下限0输入积分区间上限1输入等分份数20输入被积函数(以x为自变量)x*sqrt(x)S =0、4000输入积分区间下限5输入积分区间上限200输入等分份数20输入被积函数(以x为自变量)1/sqrt(x)S = 23、82182、自动变步长Simpson方法函数1:输入积分区间下限0输入积分区间上限2输入为课本得第几个函数(第一个这输入1):1S =1、619(过程省略)i =19函数2:输入积分区间下限0输入积分区间上限1输入为课本得第几个函数(第一个这输入1):2S =0、4(过程省略)i = 17函数3:输入积分区间下限5输入积分区间上限200输入为课本得第几个函数(第一个这输入1):3S=23、8121(过程省略)i = 111编制程序如下:Clear,clcsyms xa=input(输入积分区间下限); b=input(输入积分区间上限);n=input(输入等分份数); ff=input(输入被积函数(以x为自变量);h=(b-a)/n;f=inline(ff,x);sum1=0;sum2=0;for i=0:n-1 sum1=sum1+f(a+i*h+0、5*h);endfor i=1:n-1 sum2=sum2+f(a+i*h);endfor i=0:2*n X(i+1,1)=f(b-a)*i/(n*2)+a);end S=h/6*(f(a)+4*sum1+2*sum2+f(b)function S = zdsps( n )a=0;b=1;h=(b-a)/4;f=inline(x(3/2),x);sum1=0;sum2=0;for i=0:n-1 sum1=sum1+f(a+i*h+0、5*h);endfor i=1:n-1 sum2=sum2+f(a+i*h);endfor i=0:2*n x(i+1,1)=f(b-a)*i/(n*2)+a);end S=h/6*(f(a)+4*sum1+2*sum2+f(b); endfunction S = zpsgs(a,b,n,ff )h=(b-a)/n;sum1=0;sum2=0;sum3=0;sum4=0;if ff=1 f=inline(x6/10-x2+x,x);end if ff=2 f=inline(x(3/2),x); end if ff=3 f=inline(1/sqrt(x),x); endfor i=0:n-1 sum1=sum1+f(a+i*h+0、25*h); sum2=sum2+f(a+i*h+0、75*h); sum4=sum4+f(a+i*h+0、5*h);endfor i=1:n-1 sum3=sum3+f(a+i*h);endfor i=0:4*n x(i+1,1)=f(b-a)*i/(n*4)+a);end S=h/(6*2)*(f(a)+4*sum1+4*sum2+2*(sum3+sum4)+f(b);endclear, clca=input(输入积分区间下限); b=input(输入积分区间上限);ff=input(输入为课本得第几个函数(第一个这输入1):);for i=2:300 S(i)=zpsgs(a,b,(i),ff); S(i+1)=zpsgs(a,b,(i+1),ff); if abs(S(i+1)-S(i)0、5*10(-7) break endendS %所求积分值 i %所分份数 实验三1、对常微分方程初值问题 分别使用Euler显示方法、改进得Euler方法与经典RK法与四阶Adams预测-校正算法,求解常微分方程数值解,并与其精确解进行作图比较。其中多步法需要得初值由经典RK法提供。2、实验目得Lorenz问题与混沌实验内容:考虑著名得Lorenz方程 其中s, r, b为变化区域有一定限制得实参数。该方程形式简单,表面上瞧并无惊人之处,但由该方程揭示出得许多现象,促使“混沌”成为数学研究得崭新领域,在实际应用中也产生了巨大得影响。实验方法:先取定初值Y0=(x, y, z)=(0, 0, 0),参数s=10, r=28, b=8/3,用MATLAB编程数值求解,并与MATLAB函数ods45得计算结果进行对比。实验要求:(1)对目前取定得参数值s, r与b,选取不同得初值Y0进行运算,观察计算得结果有什么特点？解得曲线就是否有界？解得曲线就是不就是周期得或趋于某个固定点？(2)在问题允许得范围内适当改变其中得参数值s, r, b,再选取不同得初始值Y0进行试算,观察并记录计算得结果有什么特点？就是否发现什么不同得现象？3、定义函数子程序为:function z=f(x,y)z=-y+2*cos(x);return主程序为:clear,clcb=pi;a=0;n=100;y(1)=1;h=(b-a)/n;x=a:h:b;for i=1:n y(i+1)=y(i)+h*f(x(i),y(i);endt1=plot(x,y,r-)hold onfor i=1:n K1=f(x(i),y(i); K2=f(x(i+1),y(i)+h*K1); y(i+1)=y(i)+h*(K1+K2)/2;endt2=plot(x,y,b+)for i=1:n K1=f(x(i),y(i); K2=f(x(i)+0、5*h,y(i)+0、5*h*K1); K3=f(x(i)+0、5*h,y(i)+0、5*h*K2); K4=f(x(i),y(i)+h*K3); y(i+1)=y(i)+h*(K1+2*K2+2*K3+K4)/6; endt3=plot(x,y,ko)for i=1:3 K1=f(x(i),y(i); K2=f(x(i)+0、5*h,y(i)+0、5*h*K1); K3=f(x(i)+0、5*h,y(i)+0、5*h*K2); K4=f(x(i),y(i)+h*K3); y(i+1)=y(i)+h*(K1+2*K2+2*K3+K4)/6; endfor i=4:n z=y(i)+h/24*(55*f(x(i),y(i)-59*f(x(i-1),y(i-1)、 +37*f(x(i-2),y(i-2)-9*f(x(i-3),y(i-3); y(i+1)=y(i)+h/24*(9*f(x(i+1),z)+19*f(x(i),y(i)、 -5*f(x(i-1),y(i-1)+f(x(i-2),y(i-2);endt4=plot(x,y,g*)t5=ezplot(sin(x)+cos(x),0,pi)xlabel(x轴,FontWeight,bold)ylabel(y轴, FontWeight,bold)legend(t1,t2,t3,t4,t5,向前Euler法,改进得Euler法,经典四阶Runge-Kutta法,四阶Adams公式,精确解)原图为:局部放大图为:由图可得:四阶Adams公式及改进得欧拉法将为精确。