数值计算报告讲解

1、数值计算题目名称数值计算学 院机电工程学院专业班级 13微电子制造2班学 号 3113000453姓 名肖铠斌2016-12-30学习数值计算引论，应用MATLA编写Gauss列主元消元法求解线性方程 组、多项式插值、Gauss积分方法、Euler方法求常微分方程初值问题、Newton 迭代法求非线性方程的程序。【关键词】：Gauss列主元消元法求解线性方程组、多项式插值、Gauss 积分方法、Euler方法求常微分方程初值问题、Newt on迭代法求非线性方程目录1 Gauss列主元消元法求解线性方程组错误！未定义书签。1.1问题描述错误！未定义书签。1.2计算程序错误！未定义书签。1.3算

2、例22 Lagrange 多项式插值32.1问题描述32.2计算程序32.3算例33 Gauss积分方法 43.1 .问题描述43.2计算程序43.3算例44 Euler方法求常微分方程初值问题 54.1 .问题描述54.2计算程序44.3算例45Newton迭代法求非线性方程65.1 .冋题描述65.2计算程序65.3算例7参考文献8B的秩，并作比较)生成全零矩阵%返回行向量丫和j ,y向量记录abs(B(p:n,p)的每列的最大值，j向量记录每列最大值的行号把矩阵B中第p行所有列的值全都赋给矩阵C交换A,b 第jk与第k行元素消元计算从n-1循环到1、编写Gauss列主元消元法求解线性方程

3、组的程序，要求附有算例1 、问题描述：计算机中运算的时候常会碰到两个问题。1 )一旦遇到某个主元等于0,消元过程便无法进行下去。2 )在长期使用中还发现，即使消元过程能进行下去，但是当某个 主元的绝对值很小时，求解出的结果与真实结果相差甚远。2、程序：%输岀的量：系数矩阵和增广矩阵的秩RA,RB,方程组中未知量的个数n和有关方程组解及其解X的信息.fun ctionRA,RB, n, X=gaus(A,b)B=A b; n=le ngth(b);RA=ra nk(A);RB=ra nk(B);zhica=RB-RA;%求矩阵 A、if zhica0,disp( RA=RB，此方程组无解.)re

4、turnendif RA=RBif RA=ndisp( RA=RB=n，此方程组有唯一解 .X=zeros( n,1); C=zeros(1, n+1);%for p= 1:n-1Y,j=max(abs(B(p: n,p);C=B(p,:);%B(p,：)= B(j+p-1,:);B(j+p-1,：)=C；%fork=p+1: nm= B(k,p)/B(p,p);%B(k,p: n+1)= B(k,p: n+1)_m* B(p,p: n+1);endendb=B(1: n,n+1);A=B(1: n,1: n);X(n )=b( n)/A( n,n);for q=n-1:-1:1%X(q)=(

5、b(q)-sum(A(q,q+1: n)*X(q+1: n)/A(q,q); %回代计算endelsedisp(RA=RB A=5 2 1;2 8 -3;1 -3 -6;b=8;21;1;RA,RB, n,X =gaus (A,b)RA=RB=n此方程组有唯一解RA =3RB =3n =3X =12-1二 编写多项式插值的程序，要求附有算例1、问题描述：Lagrange差值法满足插值条件的、次数不超过n的多项式是存在而且是唯一的。2、程序:fun cti onyy=lagra nge(x,y,xO) /xle n=len gth(x);/Xyle n=len gth(y);/yif xle n

6、=yle n定义拉格朗日差值函数向量的维数向量的维数war ning(The len gth of x and y is not equa l);/ 提示x.y维数不等endn=len gth(xO);for i=1: n /循环。从1到nz=x0(i);s=0.0;for k=1:xle np=1.0;for j=1:yle nif j=kp=p*(z-x(j)/(x(k)-x(j);endends=p*y(k)+s;endyy(i)=s;end3、算例：已知数据如下表，试用lagrange差值多项式求x分别为0.5626;0.5635;0.5645 时函数的近似值。X0.561600.56

7、2800.564010.56521yi0.827410.826590.825770.81495解：在命令窗口输入执行命令，可得结果: x=0.56160;0.56280;0.56401;0.56521;y=0.82741;0.82659;0.82577;0.81495;x0=0.5626;0.5635;0.5645;yO=lagra nge(x,y,x0)y0 = 0.82650.82680.8231 plot(x,y,o,x0,y0,k)得到图形三编写Gauss积分方法的程序，要求附有算例。1 、问题描述：被积分函数f( E, n , Z )一般是很复杂的，即使能 够得出它的显式，其积分也是

8、很繁的。因此，一般用数值积分来代替函数的定 积分。2 、程序：% f :积分表达式，可以是函数句柄、inline函数、匿名函数、字符串表达式，但是必须可以接受矢量输入%a,b :积分上下限，注意积分区间太大会降低精度，此时建议使用复化求积公式，默认-1 1% n :积分阶数，可以任意正整数，但是不建议设置过大，大不一定能得到更好的精度，默认7% tol :积分精度，默认 1e-6% ql :积分结果% Ak :求积系数% xk :求积节点，满足 ql=sum(Ak.*fun(xk)fun cti onql,Ak,xk=guass(f,a,b ,n ,tol)if nargi n=1%输入参数有

9、效性检验a=-1;b=1; n=7;tol=1e-8;elseif nargin=3n=7;tol=1e-8;elseifnargin=4tol=1e-8;elseifnargin=2|nargin5)；error( The Number of In put Argume nts Is Wrong! endsyms x%计算求积节点p=sym2poly(diff(xA2-1F( n+1), n+1)/(2A n*factorial( n); tk=roots(p);% 求积节点Ak=zeros( n+1,1);%计算求积系数for i=1: n+1xkt=tk;xkt(i)=;pn=poly(

10、xkt);fp=(x)polyval(p n, x)/polyval(p n,tk(i);Ak(i)=quadl(fp,-1,1,tol); end%求积系数xk=(b-a)/2*tk+(b+a)/2;f=fcn chk(f,vectorizefx=f(xk)*(b-a)/2;ql=sum(Ak.*fx);% 积分变量代换，将a,b 变换到-1,1 );%检验积分函数fun有效性%计算变量代换之后积分函数的值%计算积分值3、算例，求1艾0cos(x)的值。编写函数文件:Gaussf.mfun cti onf=gaussf(x)f=cos(x);End在命令窗口输入： ql,Ak,xk=guas

11、s(gaussf, 0,1, 7, 1e-6) ql =0.8415Ak =0.10120.22240.10120.22240.31370.31370.36270.3627 xk =0.01990.10170.98010.89830.23720.76280.40830.5917四编写Euler方法求常微分方程初值问题的程序，要求附有算例1、问题描述：欧拉法的特点：单步，显式，一阶求导精度，截断误差贰阶2、程序：欧拉公式：广yhfQyi)Euler.mfun cti onxx,yy=euler(f,a,b,y0,h)x=a:h:bn=len gth(x);y=zeros(1, n);y(1)=y

12、0;for i=1: n-1;y(i+1)=y(i)+h*feval(f,x(i),y(i);endxx=x;yy=y;plot(xx,yy, rp)xlabel(In depe ndent Variable x),ylabel( Depe ndent Variable y)ti tle(Solution of ODE with modified euler method)3、算例：用欧拉方法解初值问题2 2V =1y 0 乞 X 岂1hx y，取 h=0.1.y(0) =0解：建立函数文件myfu n.mfun cti onf=myfu n(x,y)f=1+xA2+yA2End在窗口输入命令

13、：x y=euler(myfu n,0,1,0,0.1)x =00.10000.20000.30000.40000.50000.60000.70000.80000.90001.0000y =00.10000.20200.31010.42870.56310.71980.90761.13901.43271.8189所得图形:改进后的欧拉公式% 厂 yi hfk y% 厂 Ef(x %)f(x1 %丿fun cti onxx,yy=euler(f,a,b,yO,h) x=a:h:b n=len gth(x);y=zeros(1, n);y(i)=y0;for i=1: n-1;yp=y(i)+h*f

14、eval(f,x(i),y(i); yc=y(i)+h*feval(f,x(i+1),yp); y(i+1)=(yp+yc)/2;end xx=x;yy=y;ey=-1./xplot(xx,yy, rp)严 Q-qEElllvPLJ2cl4A1.5fl-xlabel(In depe ndent Variable x),ylabel(Depe ndent Variable y)title(Solution of ODE with modified euler method)在窗口输入命令：x y=euler(myfu n,0,1,0,0.1)x =00.10000.20000.30000.400

15、00.50000.60000.70000.80000.90001.0000y =00.10100.20610.31960.44690.59430.77100.98991.27191.6529 2.2020U figure 1F-|Ra|文林(E)爲窑旧 查看世)l&AW桌園助 菌口凹鼻W目 0Solution of ODE wrth mrodifed euler method2.5Qi e-0.10 203040 50 60.708091Independent VanaUe x五编写Newt on迭代法求非线性方程的程序，要求附有算例1、问题描述：构造迭代函数的一条重要途径是用近似方程来代替原

16、方程去求根。因此，如果能将非线性方程 f x =0用线性方程去代替，那么， 求近似根问题就很容易解决，而且十分方便。牛顿(Newton)法就是一种将非线性 方程线化的一种方法、程序:Newto n.m%y=n ewt on(a,n,xO,nn, epsi)输入变量：%a n+1元素的一维实数组，输入参数，按升幕存放方程系数。%n整变量，输入参数，方程阶数。%x0实变量，输入参数，初始迭代值。%nn整变量，输入参数，允许的最大迭代次数。%eps1实变量，输入参数，控制根的精度。fun cti ony=n ewt on(a,n,xO,nn, epsi)x(1)=x0;%初始向量b=1;i=1;while (abs(b)eps1*x(i)%abs(b)取模，判断i=i+1;x(i)=x(i-1)-n_f(a, n,x(i-1)/n_df(a, n,x(i-1);%xx-J(xr1)fF(xiJ)b=x(i)-x(i-1);if (inn)error( nn is full);return ;endendy=x(i);i程序中调两个子函数n_f.m和n_df.m文件如下：function y=n_f(a,n,x)%待求根的实数代数方程的函数y=0.0;fori=1:( n+1)y=y+a(i)*xA( n+1-i);endn _df.m:fu