华工数值分析课程设计

1、数值分析课程设计 第一题：1. 设计思路：我打算用选主元法，先算出每一列，然后把买一列加起来就是结果了。2.程序清单：function x=mat(a,b,flag)if nargin<3,flag=0;endn=length(b);a=a,b;for i=1:(n-1) ar,r=max(abs(a(i:n,i);r=r+i-1; if r>i t=a(i,:);a(i,:)=a(r,:);a(r,:)=t; enda(i+1):n,(i+1):(n+1)=a(i+1):n,(i+1):(n+1)-a(i+1):n,i)/a(i,i)*a(i,(i+1):(n+1);a(i+1)

2、:n,i)=zeros(n-i,1);if flag=0,a,endendx=zeros(n,1);x(n)=a(n,n+1)/a(n,n);for i=n-1:-1:1 x(i)=(a(i,n+1)-a(i,(i+1):n)*x(i+1):n)/a(i,i);end3.流程图：我的例子是九阶对角矩阵对角元素为1,2,3,4,5,6,7,8,9,。流程如下图通过比照知道，所编的程序算的结果和inv的结果一致。 我的方法时间约为0.003*9=0.027秒，而inv方法时间约为0.0031秒，明显inv方法快。我第二题：1. 设计思路：我通过limit与求阶方法编程。2. 程序清单：functi

3、on n=qiujie(x0)k=0;syms x;x1=0.99*x-x2;h=(abs(x0-x1)/(abs(x0-x)k;while limit(h,x,x0)=0 k=k+1; h=(abs(x0-x1)/(abs(x0-x)k;endn=kfunction y,m=bdd(x,eps)if nargin=1 eps=1e-100; elseif nargin<1 error returnendx1=0.99*x-x2;n=1;while (abs(x1-x)>=1e-100)&&(n<1e100) x=x1; x1=0.99*x-x2; n=n+1

4、;endy=x1m=n3. 流程图：由所得知，所求阶数为一阶。第三题：1. 设计思路：按照题意用三次样条。从题意分析可知x=4.9与22.9的一阶导数均为0所以带入七组数据即可。2. 程序清单：function s=sp(x,y,dx0,dxn)n=length(x)-1;h=diff(x);d=diff(y)./h;a=h(2:n-1);b=2*(h(1:n-1)+h(2:n);c=h(2:n);u=6*diff(d);b(1)=b(1)-h(1)/2;u(1)=u(1)-3*(d(1)-dx0);b(n-1)=b(n-1)-h(n)/2;u(n-1)=u(n-1)-3*(dxn-d(n);

5、for k=2:n-1 temp=a(k-1)/b(k-1); b(k)=b(k)-temp*c(k-1); u(k)=u(k)-temp*u(k-1);endm(n)=u(n-1)/b(n-1);for k=n-2:-1:1 m(k+1)=(u(k)-c(k)*m(k+2)/b(k);endm(1)=3*(d(1)-dx0)/h(1)-m(2)/2;m(n+1)=3*(dxn-d(n)/h(n)-m(n)/2;for k=0:n-1 s(k+1,1)=(m(k+2)-m(k+1)/(6*h(k+1); s(k+1,2)=m(k+1)/2; s(k+1,3)=d(k+1)-h(k+1)*(2*

6、m(k+1)+m(k+2)/6; s(k+1,4)=y(k+1);end3. 流程图： 可知所求函数为 T=-0.008x3-0.0911x2+22.8, 4.9,9.1 T=0.0291x3-0.1919x2-1,1888x+20.6, 9.1,13.7T=-0.0157x3+0.2093x2-1,1091x+13.9, 13.7,18.3T=0.0039x3-0.0071x2-1,1793x+11.7, 13.7,18.3第五题：1. 设计思路：我设计了精确值与题意的要求，并根据结果进行比拟，只对有相同的数值即对1,1.1，。，2的值进行比照。2. 程序清单：function x,y=eu

7、(fname,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0; for n=1:length(x)-1 y(n+1)=y(n)+h*feval(fname,x(n),y(n);end x=x'y=y'function x,y=eu2(fname,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0;for n=1:length(x)-1 k1=feval(fname,x(n),y(n); y(n+1)=y(n)+h*k1; k2=feval(fname,x(n+1),y(n+1); y(n+1)=y(n)+h*(k1

8、+k2)/2;endx=x'y=y'function x,y=ku(fname,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0;for n=1:length(x)-1 k1=feval(fname,x(n),y(n); k2=feval(fname,x(n)+h/2,y(n)+h/2*k1); k3=feval(fname,x(n)+h/2,y(n)+h/2*k2); k4=feval(fname,x(n+1),y(n)+h*k3); y(n+1)=y(n)+h*(k1+2*k2+2*k3+k4)/6;endx=x'y=y'3.

9、 流程图：精确值如上图。欧拉法，h=0.025欧拉法，h=0.1改良欧拉，h=0.05改良欧拉h=0.1龙格法，h-0.1.比拟表格：X的值11.11.21.31.41.51.61.71.81.92精确值00.34590.86661.60722.62043.96775.72107.963910.793614.323118.6831欧拉h=0.02500.32550.81651.51642.47573.75395.42017.555110.252914.584717.7914欧拉h=0.100.27180.68481.27702.09353.18744.62086.46648.809111.748015.3982改良欧拉h=0.0500.34490.86431.60312.61423.95895.70927.948710.774414.299318.6542改良欧拉h=0.100.34240.85831.59272.59833.93645.67897.909210.724514.237418.5789龙格h=0.100.34590.86661.60722.62033.96765.72097.963810.793514.322918.6829由以