基于matlab的数值分析课程设计参考.doc

数值分析试验报告 矩阵的LU分解 题目：求4阶矩阵的LU分解 方法：杜里特尔分解法 程序： function f=LU_decom(A) [m,n]=size(A) L=eye(n); U=zeros(n); flag=ok; for i=1:n U(1,i)=A(1,i); end for r=2:n L(r,1)=A(r,1)/U(1,1); end for i=2:n for j=i:n z=0; for r=1:i-1 z=z+L(i,r)*U(r,j); end U(i,j)=A(i,j)-z; end if abs(U(i,i))eps flag=failure return; end for k=i+1:n m=0; for q=1:i-1 m=m+L(k,q)*U(q,i); end L(k,i)=(A(k,i)-m)/U(i,i); end end L U end 结果 LU_decom(A) m = 4 n = 4 L = 1 0 0 0 2 1 0 0 1 2 1 0 3 3 2 1 U = 2 4 2 6 0 1 2 3 0 0 3 6 0 0 0 1 5.拓展 在编写程序过程中由于角标较多因此在运行过程中出现了不少角标不对的错误 题目；给出函数f(x)=1/(1+25x^2)，求f（x在[-1,1]上取5个和9个等距节点，做最小二乘拟合，得出均方误差方误差。 五个节点时，matlab编码为： 首先建立M文件，并保存 function y=f(x) y=1/(1+25*x^2); end x=[-1 -0.5 0 0.5 1]; for i=1:5 y(i)=f(x(i)); end a=polyfit(x,y,3) syms x f1=a(1)*x^3+a(2)*x^2+a(3)*x+a(4) x=[-1 -0.5 0 0.5 1]; for i=1:5 E(i)=(f(x(i))-(a(1)*x(i)^3+a(2)*x(i)^2+a(3)*x(i)+a(4)))^2; end sum(E) 输出结果为 a = -0.0000 -0.6063 -0.0000 0.5737 f1 = -47781/4*x^3-1600/2639*x^2-30703/*x+1514/2639 （拟合的多项式） ans = 0.3534（均方误差） 九个点的时候，matlab编码为： x=[-1 -0.75 -0.5 -0.25 0 0.25 0.5 0.75 1]; for i=1:9 y(i)=f(x(i)); end a=polyfit(x,y,3) syms x f2=a(1)*x^3+a(2)*x^2+a(3)*x+a(4) x=[-1 -0.75 -0.5 -0.25 0 0.25 0.5 0.75 1]; for i=1:5 E1(i)=(f(x(i))-(a(1)*x(i)^3+a(2)*x(i)^2+a(3)*x(i)+a(4)))^2; end sum(E1) 输出结果为： a = -0.0000 -0.5609 0.0000 0.4855 f2 = -6039/*x^3-12921/85248*x^2+54361/*x+64671/85248 （最小二乘拟合多项式） ans = 0.3350 （均方误差） 用复合梯形公式求积分的值。 function i=combinetraprl(f,a,b,eps) %复化梯形公式求函数f在区间[a，b]上的定积分 %函数名：f %积分下限：a %积分上限：b %积分精度：eps %积分值：i %积分划分的子区间个数：step n=1; h=(b-a)/2; i1=0; i2=(subs(sym(f),findsym(sym(f)),a)+subs(sym(f),findsym(sym(f)),b))/h; while abs(i2-i1)eps n=n+1;