数值分析上机第四次作业

数值分析上机第四次作业实验项目：共轭梯度法求解对称正定的线性方程组实验内容：用共轭梯度法求解下面方程组(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);end clear， 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=pl