数值方法计算实习题

数值方法计算实习题要求：1、用Matlab语言或你熟悉的其他算法语言编写程序，使之尽可能具有通用性；2、根据上机计算实践，对所使用的数值方法的特点、性质、有效性、误差和收敛性等方面进行必要的讨论和分析；3、完成计算后写出实验报告，内容包括：课题名称、解决的问题、采用的数值方法、算法程序、数值结果、对实验结果的讨论和分析等；4、特别说明：严禁抄袭，否那么一经发现，所有雷同实验报告最多评为及格。一、下表给出了飞行中鸭子的上部形状的节点数据，试用三次样条插值函数〔自然边界条件〕和20次Lagrange插值多项式对数据进行插值。用图示出给定的数据，以及和。0.91.31.92.12.63.03.94.44.75.06.01.31.51.852.12.62.72.42.152.052.12.257.08.09.210.511.311.61212.613.013.32.32.251.951.40.90.70.60.50.40.25解：>>x=[0.91.31.92.12.63.03.94.44.75.06.07.08.09.210.511.311.61212.613.013.3];>>y=[1.31.51.852.12.62.72.42.152.052.12.252.32.251.951.40.90.70.60.50.40.25];%〔1〕三次样条插值法xi=0.9:0.01:13.3;yi=interp1(x,y,xi,'spline');>>xi=0.9:0.01:13.3;yi=interp1(x,y,xi,'spline');>>title('试验一--三次样条插值图示')>>pp=spline(x,y)pp=form:'pp'breaks:[1x21double]coefs:[20x4double]pieces:20order:4dim:1>>pp.coefsans=0.7735-0.99950.77601.30000.7735-0.07140.34771.5000-2.78941.32091.09741.8500-0.4585-0.35281.29102.10000.4489-1.04050.59442.60000.1738-0.5018-0.02252.70000.0783-0.0325-0.50332.40001.31410.0850-0.47712.1500-1.58121.2676-0.07132.05000.0431-0.15550.26232.1000-0.0047-0.02610.08082.2500-0.0245-0.04010.01462.30000.0175-0.1135-0.13902.2500-0.0128-0.0505-0.33581.9500-0.0201-0.1003-0.53191.40001.2094-0.1485-0.73100.9000-0.82790.9400-0.49350.70000.0122-0.0535-0.13890.6000-0.2960-0.0316-0.19000.5000-0.2960-0.3867-0.35730.4000所以所得方程为%〔2〕用拉格朗日法插值%定义Lagrange程序functionf=Language(x,y,x0)symst;if(length(x)==length(y))n=length(x);elsedisp('xºÍyµÄάÊý²»ÏàµÈ£¡');return;endf=0.0;for(i=1:n)l=y(i);for(j=1:i-1)l=l*(t-x(j))/(x(i)-x(j));end;for(j=i+1:n)l=l*(t-x(j))/(x(i)-x(j));end;f=f+l;simplify(f);if(i==n)if(nargin==3)f=subs(f,'t',x0);elsef=collect(f);f=vpa(f,6);endendend>>Language(x,y)ans=52462.6*t+189995.*t^3-189851.*t^4+136778.*t^5-11.3161*t^12-.277283e-6*t^18+1.18284*t^13-73866.6*t^6+.111076e-4*t^17-.976904e-1*t^14+.427949e-8*t^19-.307453e-10*t^20+30677.6*t^7+2564.20*t^9-9968.98*t^8+.628590e-2*t^15-525.813*t^10-9652.78-.308159e-3*t^16+86.2514*t^11-128683.*t^2二、Wilson矩阵，且向量，那么方程组有准确解。⑴用Matlab内部函数求，的所有特征值和；⑵令，解方程组，并求出向量和，从理论结果和实际计算结果两方面分析方程组解的相对误差与的相对误差的关系；⑶再改变扰动矩阵〔其元素的绝对值不超过0.005〕，重复第2问。解：解：〔1〕A=[10787;7565;86109;75910];b=[32;23;33;31];m=det(A)%A的行列式m=1>>n=30.2887/0.0102%因为A为对称阵，所以A的条件数等于最大特征值除于最小特征值n=2.9695e+003所以A的行列式为1，cond(A)2=2.9695e+003(2)>>B=[1078.17.2;7.085.0465;85.989.899;6.994.9999.98];%令B等于δA+A>>b=[32;23;33;31];>>[rank(B),rank([B,b])]ans=44>>x1=B\b%因为B得秩为4，所以方程组的解唯一，其中x1为的解x1=-81.0000137.0000-34.000022.0000所以δx=x1-x=[-82.0000,136.0000,-35.0000,21.0000]T‖δx‖2=sqrt(82*82+136*136+35*35+21*21)ans=163.9695三、解三对角线性方程组的追赶法及其应用⑴编写解三对角线性方程组的追赶法的通用程序，并应用于方程组，检验程序的正确性；〔解为〕⑵求微分方程边值问题的数值解〔取步长〕，并与精确解比拟〔精确解为〕。说明：离散化微分方程时，解：clearall;a=[2,2,2,2,2];b=[-1,-1,-1,-1];c=[-1,-1,-1,-1];r=[1,0,0,0,0];n=length(a);b=[0,b];u(1)=r(1)/a(1);v(1)=c(1)/a(1);fork=2:n-1u(k)=(r(k)-u(k-1)*b(k))/(a(k)-v(k-1)*b(k));v(k)=c(k)/(a(k)-v(k-1)*b(k));endu(n)=(r(n)-u(n-1)*b(n))/(a(n)-v(n-1)*b(n));x(n)=u(n);fork=n-1:-1:1x(k)=u(k)-v(k)*x(k+1);endfprintf('Èý¶Ô½Ç·½³Ì×éµÄ½âΪ\n')fork=1:nfprintf('x(%1d)=%10.8f\n',k,x(k))end>>zhuiganfa%调用追赶法三对角方程组的解为x(1)=0.83333333x(2)=0.66666667x(3)=0.50000000x(4)=0.33333333x(5)=0.16666667〔2〕四、公元1225年，比萨的数学家Leonardo研究了方程，得到一个根，没有人知道他用什么方法得到这个值。对于这个方程，分别用以下方法：⑴迭代法；⑵迭代法；⑶对⑴的Steffensen加速方法；⑷对⑵的Steffensen加速方法；⑸Newton法。求方程的根〔可取〕，计算到Leonardo所得到的准确度。五、用不同的数值方法计算积functionx=nuwton(fname,dfname,x0,e)ifnargin<4,e=1e-4;endx=x0;x0=x+2*e;whileabs(x0-x)>ex0=x;x=x0-feval(fname,x0)/feval(dfname,x0);end分的近似值，其中⑴取不同的步长，分别用复合梯形公式和复合辛普森公式计算积分，比拟两个公式的计算效果，是否存在一个最小的，使得精度不能再被改善？⑵用龙贝格求积公式，取，并打印出T-表。解：〔1〕①functions=traprl(f,a,b,n)h=(b-a)/n;s=0;fork=1:(n-1)x=a+h*k;s=s+feval('f',x);ends=h*(feval('f',a)+feval('f',b))/2+h*s;>>traprl('f',1,3,1000)ans=-1.4263②functions=simprl(f,a,b,n)h=(b-a)/(2*n);s1=0;s2=0;fork=1:nx=a+h*(1*k-1);s1=s1+feval('f',x);endfork=1:(n-1)x=a+h*2*k;s2=s2+feval('f',x);ends=h*(feval('f',a)+feval('f',b)+4*s1+2*s2)/3;>>simprl('f',1,3,10000)ans=14.4924六、给定矩阵⑴用Matlab的函数“eig”