数值分析的MATLAB程序

列主元法 function lianzhuyuan(A,b) n=input('请输入n:') %选择阶数 A=zeros(n,n); %系数矩阵A b=zeros(n,1); %矩阵b X=zeros(n,1); %解X for i=1:n for j=1:n A(i,j)=(1/(i+j-1)); %生成hilbert 矩阵A end b(i,1)=sum(A(i,:)); %生成矩阵b end for i=1:n-1 j=i; top=max(abs(A(i:n,j))); %列主元 k=j; while abs(A(k,j))~=top %列主元所在行 k=k+1; end for z=1:n %交换主元所在行 a1=A(i,z); A(i,z)=A(k,z); A(k,z)=a1; end a2=b(i,1); b(i,1)=b(k,1); b(k,1)=a2; for s=i+1:n %消去算法开始 m=A(s,j)/A(i,j); %化简为上三角矩阵 A(s,j)=0; for p=i+1:n A(s,p)=A(s,p)-m*A(i,p); end b(s,1)=b(s,1)-m*b(i,1); end end X(n,1)=b(n,1)/A(n,n); %回代开始 for i=n-1:-1:1 s=0; %初始化s for j=i+1:n s=s+A(i,j)*X(j,1); end X(i,1)=(b(i,1)-s)/A(i,i); end X 欧拉法 clc clear % 欧拉法 p=10; %贝塔的取值 T=10; %t 取值的上限 y 1=1; %y 1 的初值 r1=1; %y 2 的初值 %输入步长 h 的值 h=inpu t('欧拉法 please inpu t nu mber(h=1 0.5 0.25 0.125 0.0625):h=') ; if h>1 or h<0 break end S1=0:T/h; S2=0:T/h; S3=0:T/h; S4=0:T/h; i=1; % 迭代过程 for t=0:h:T Y=(ex p(-t)); R=(1/(p-1))*ex p(-t)+((p-2)/(p-1))*ex p(-p*t); y =y 1+h*(-y 1); y 1=y ; r=r1+h*(y 1-p*r1); r1=r; S1(i)=Y; S2(i)=R; S3(i)=y ; S4(i)=r; i=i+1; end t=[0:h:T]; % 红线为解析解，'x '为数值解 plot(t,S1,'r',t,S3,'x ') 改进欧拉法 clc clear p=10; %贝塔的取值 T=10; %t 取值的上限 y 1=1; %y 1 的初值 r1=1; %y 2 的初值 %输入步长 h 的值 h=inpu t('改进欧拉法 please inpu t nu mber(h=1 0.5 0.25 0.125 0.0625):h=') ; if h>1 or h<0 break end S1=0:T/h; S2=0:T/h; S3=0:T/h; S4=0:T/h; i=1; % 迭代过程 for t=0:h:T Y=(ex p(-t)); R=(1/(p-1))*ex p(-t)+((p-2)/(p-1))*ex p(-p*t); k1=-y 1; l1=y 1-p*r1; k2=-(y 1+h*k1); l2=y 1+h*k1-p*(r1+h*l1); y =y 1+h*(k1+k2)/2; r=r1+h*(k1+k2)/2; r1=r; y 1=y ; S1(i)=Y; S2(i)=R; S3(i)=y ; S4(i)=r; i=i+1; end t=[0:h:T]; % 红线为解析解，'x '为数值解 plot(t,S1,'r',...