数值分析大作业三四五六七

合集下载

北航数值分析大作业3

数值分析第三次作业1.设计方案对Fredholm积分方程，用迭代法进行求解：()'(())u x A u x=，其中11(())()(,)()A u x g x K x y u y dy-=-⋅⎰对于公式中的积分部分用数值积分方法。

复化梯形积分法，取2601个节点，取迭代次数上限为50次。

实际计算迭代次数为18次，最后算得误差为r= 0.97E-10。

复化Simpson积分法，取迭代次数上限为50次，取2*41+1，即83个节点时能满足精度要求。

实际计算迭代次数为17次，最后的误差为r= 0.97E-10。

Guass积分法选择的Gauss—Legendre法，取迭代次数上限为50次，直接选择8个节点，满足精度要求。

实际计算迭代次数为24次，最后算得误差为r= 0.87E-10。

2.全部源程序module integralimplicit nonecontains!//////////复化梯形subroutine trapezoid(m)implicit noneinteger :: i,j,k,mreal*8 :: x(m+1),u(m+1)real*8 :: sum,sum1,g,r,hreal*8 :: e=1.0e-10h=2./mdo i=1,m+1x(i)=-1.+(i-1)*hend dou=0.02do k=1,50do i=1,m+1sum1=0.g=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.)do j=2,msum1=sum1+dexp(x(i)*x(j))*u(j)end dosum=h/2.*(dexp(x(i)*-1.)*u(1)+dexp(x(i)*1.)*u(m+1)+2*sum1)u(i)=g-sumend dor=h/2.*((dexp(x(1)*4)-u(1))**2+(dexp(x(m+1)*4)-u(m+1))**2) do i=2,mr=r+h*(dexp(x(i)*4)-u(i))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(1,file="trapezoid.txt")do i=1,m+1write(1,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(1,'(4x,a2,e9.2)') "r=",rclose(1)returnend subroutine trapezoid!///////////复化simpsonsubroutine simpson(m)implicit noneinteger :: i,j,k,mreal*8 :: x(2*m+1),u(2*m+1)real*8 :: sum,sum1,sum2,g,r,hreal*8 :: e=1.0e-10h=2./(2.*m)do i=1,2*m+1x(i)=-1.+(i-1)*hend dou=0.02do k=1,50do i=1,2*m+1sum1=0.sum2=0.g=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.)do j=1,msum1=sum1+dexp(x(i)*x(2*j))*u(2*j)end dodo j=1,m-1sum2=sum2+dexp(x(i)*x(2*j+1))*u(2*j+1)sum=h/3.*(dexp(x(i)*-1.)*u(1)+dexp(x(i)*1.)*u(2*m+1)+4*sum1+2*sum2) u(i)=g-sumend dor=h/3.*((dexp(x(1)*4)-u(1))**2+(dexp(x(2*m+1)*4)-u(2*m+1))**2)do i=1,mr=r+4.*h/3.*(dexp(x(2*i)*4)-u(2*i))**2end dodo i=1,m-1r=r+2.*h/3.*(dexp(x(2*i+1)*4)-u(2*i+1))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(2,file="simpson.txt")do i=1,2*m+1write(2,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(2,'(4x,a2,e9.2)') "r=",rclose(2)returnend subroutine simpson!///////////Gauss_Legendre法subroutine Gaussimplicit noneinteger,parameter :: m=8integer :: i,j,kreal*8 :: x(m),u(m),a(m)real*8 :: sum,g,rreal*8 :: e=1.0e-10data x /-0.9602898565,-0.7966664774,-0.5255324099,-0.1834346425,&0.1834346425,0.5255324099,0.7966664774,0.9602898565/data a /0.1012285363,0.2223810345,0.3137066459,0.3626837834,&0.3626837834,0.3137066459,0.2223810345,0.1012285363/u=0.02do k=1,50do i=1,mg=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.)do j=1,msum=sum+dexp(x(i)*x(j))*u(j)*a(j)end dou(i)=g-sumend dor=0.do i=1,mr=r+a(i)*(dexp(x(i)*4)-u(i))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(3,file="Gauss.txt")do i=1,mwrite(3,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(3,'(4x,a2,e9.2)') "r=",rclose(3)returnend subroutine Gaussend module!//////////主程序program mainuse integralimplicit noneinteger :: code1=2600integer :: code2=41call trapezoid(code1)call simpson(code2)call Gaussend program3.各种积分方法的节点和数值解（由于数据太多，在打印时用了较计算时少的有效数字）复化Simpson法4.各方法所得曲线（由于所取节点太多，且精度高，所以图中很难看出各曲线的区别。

东南大学数值分析上机作业汇总

东南大学数值分析上机作业汇总-标准化文件发布号：（9456-EUATWK-MWUB-WUNN-INNUL-DDQTY-KII数值分析上机报告院系：学号：姓名：目录作业1、舍入误差与有效数 (1)1、函数文件cxdd.m (1)2、函数文件cddx.m (1)3、两种方法有效位数对比 (1)4、心得 (2)作业2、Newton迭代法 (2)1、通用程序函数文件 (3)2、局部收敛性 (4)（1）最大δ值文件 (4)（2）验证局部收敛性 (4)3、心得 (6)作业3、列主元素Gauss消去法 (7)1、列主元Gauss消去法的通用程序 (7)2、解题中线性方程组 (7)3、心得 (9)作业4、三次样条插值函数 (10)1、第一型三次样条插值函数通用程序： (10)2、数据输入及计算结果 (12)作业1、舍入误差与有效数设∑=-=Nj N j S 2211，其精确值为⎪⎭⎫ ⎝⎛---1112321N N . （1）编制按从小到大的顺序11131121222-⋅⋅⋅+-+-=N S N ，计算N S 的通用程序；（2）编制按从大到小的顺序()12111111222-⋅⋅⋅+--+-=N N S N ，计算N S 的通用程序；（3）按两种顺序分别计算642101010,,S S S ,并指出有效位数；（4）通过本上机你明白了什么？程序：1、函数文件cxdd.mfunction S=cxdd(N) S=0; i=2.0;while (i<=N) S=S+1.0/(i*i-1); i=i+1;endscript 运行结果（省略>>）：S=cxdd(80) S=0.7375772、函数文件cddx.mfunction S=cddx (N) S=0; for i=N:-1:2 S=S+1/(i*i-1); endscript 运行结果（省略>>）： S=cddx(80) S=0.7375773、两种方法有效位数对比精确值函数：function S=jqz(N)S=0.5*(1.5-1.0/N-1.0/(N+1));script运行结果（省略>>）4、心得本题重点体现了数值计算中“大数吃小数”的问题，由于计算机计算的截断特点，从大到小的计算会导致小数的有效数被忽略掉。

数值分析作业答案

. 证明：(1). 假如 A 是对称正定矩阵，则 A 1也是对称正定矩阵(2). 假如 A 是对称正定矩阵，则 A 能够独一地写成 A L T L ，此中 L 是拥有正对角元的下三角矩阵。

证明：(1). 因 A 是对称正定矩阵，故其特点值i皆大于 0 ，所以 A 1的特点值i 1也皆大于 0 。

因此i1也皆大于 0 ，故 A 是可逆的。

又(A 1)T(A T) 1 A 1则 A 1也是对称正定矩阵。

(2).由A是对称正定，故它的全部次序主子阵均不为零，进而有独一的杜利特尔分解~A LU 。

又u111u12u1nu11u11 u221u2 nU DU 0u22u nn1此中 D 为对角矩阵，U0为上三角矩阵，于是~ ~A L U LDU0由 A 的对称性，得~A A T U0T D L T由分解的独一性得U 0T ~ L进而~ ~A LDL T由 A 的对称正定性，假如设D i (i 1,2, ,n) 表示A的各阶次序主子式，则有d1 D1 0 ， d i D i 0 ， i 2,3, , nD i 1故d1 d1 d1d2 d2 d2 Dd n d n d n1 1D2D2所以~11~ ~1~1LL T，A LD2D2 L T LD2(LD2)T~1此中 L L D 2为对角元素为正的下三角矩阵。

. 用列主元消去法解线性方程组12 x1 3x2 3x3 1518 x1 3x2 x3 15x1 x2 x3 6并求出系数矩阵 A 的队列式(即 det A )的值。

解18 3 1 15 ( A b) r1 r2 12 3 3 151 1 1 62m2131m311818 3 1 15 0 1 7 / 375 m32 0 7 /6 17 /18 631/ 618 3 1 150 1 7 / 3 50 0 11/ 3 11所以解为 x3 3 ， x2 2 ， x1 1，det A 66 。

. 用追赶法解三对角方程组Ax b ，此中2 1 0 0 0 1 1 2 1 0 0 0 A012 1 0 ， b 0 。

北航数值分析全部三次大作业

北航数值分析全部三次大作业第一次大作业是关于解线性方程组的数值方法。

我们被要求实现各种常用的线性方程组求解算法，例如高斯消元法、LU分解法和迭代法等。

我首先学习了这些算法的原理和实现方法，并借助Python编程语言编写了这些算法的代码。

在实验中，我们使用了不同规模和条件的线性方程组进行测试，并比较了不同算法的性能和精度。

通过这个作业，我深入了解了线性方程组求解的原理和方法，提高了我的编程和数值计算能力。

第二次大作业是关于数值积分的方法。

数值积分是数值分析中的重要内容，它可以用于计算曲线的长度、函数的面积以及求解微分方程等问题。

在这个作业中，我们需要实现不同的数值积分算法，例如矩形法、梯形法和辛普森法等。

我学习了这些算法的原理和实现方法，并使用Python编写了它们的代码。

在实验中，我们计算了不同函数的积分值，并对比了不同算法的精度和效率。

通过这个作业，我深入了解了数值积分的原理和方法，提高了我的编程和数学建模能力。

第三次大作业是关于常微分方程的数值解法。

常微分方程是数值分析中的核心内容之一，它可以用于描述众多物理、化学和生物现象。

在这个作业中，我们需要实现不同的常微分方程求解算法，例如欧拉法、龙格-库塔法和Adams法等。

我学习了这些算法的原理和实现方法，并使用Python编写了它们的代码。

在实验中，我们解决了一些具体的常微分方程问题，并比较了不同算法的精度和效率。

通过这个作业，我深入了解了常微分方程的原理和方法，提高了我的编程和问题求解能力。

总的来说，北航数值分析课程的三次大作业非常有挑战性，但也非常有意义。

通过这些作业，我在数值计算和编程方面得到了很大的提升，也更加深入地了解了数值分析的理论和方法。

虽然这些作业需要大量的时间和精力，但我相信这些努力将会对我未来的学习和工作产生积极的影响。

BUAA数值分析大作业三

北京航空航天大学2020届研究生《数值分析》实验作业第九题院系：xx学院学号：姓名：2020年11月Q9:方程组A.4一、算法设计方案（一）总体思路1.题目要求∑∑===k i kj s r rsy x cy x p 00),(对f(x, y) 进行拟合，可选用乘积型最小二乘拟合。

),(i i y x 与),(i i y x f 的数表由方程组与表A-1得到。

2.),(**j i y x f 与1使用相同方法求得，),(**j i y x p 由计算得出的p(x,y)直接带入),(**j i y x 求得。

1. ),(i i y x 与),(i i y x f 的数表的获得对区域D ={ （x,y)|1≤x ≤1.24,1.0≤y ≤1.16}上的f (x , y )值可通过xi=1+0.008i ，yj=1+0.008j ，得到),(i i y x 共31×21组。

将每组带入A4方程组，即可获得五个二元函数组，通过简单牛顿迭代法求解这五个二元数组可获得z1~z5有关x,y 的表达式。

再将),(i i y x 分别带入z1~z5表达式即可获得f(x,y)值。

2.乘积型最小二乘曲面拟合2.1使用乘积型最小二乘拟合，根据k 值不用，有基函数矩阵如下：⎪⎪⎪⎭⎫ ⎝⎛=k i i k x x x x B 0000 ， ⎪⎪⎪⎭⎫ ⎝⎛=k j jk y y y y G 0000数表矩阵如下：⎪⎪⎪⎭⎫⎝⎛=),(),(),(),(0000j i i j y x f y x f y x f y x f U记C=[rs c ]，则系数rs c 的表达式矩阵为：11-)(-=G G UG B B B C T TT ）（通过求解如下线性方程，即可得到系数矩阵C 。

UG B G G C B B T T T =)()(2.2计算),(),,(****j i j i y x p y x f (i =1,2,…,31 ; j =1,2,…,21) 的值),(**j i y x f 的计算与),(j i y x f 相同。

北航数值分析大作业3(学硕)

《数值分析》作业三院系：机械学院学号：SY1307145姓名：龙安林2013年11 月24 日1. 算法设计1) 开始；2) 计算数组[][]0.08,0.050.5,0,1,2,,10;0,1,2,,20x i i y j j i j ==+=⋯=⋯（）； 3) 将点[][],0,1,2,,10;0,1,2,,20x i y j i j =⋯=⋯（）,（）带入非线性方程组： 0.5cos 2.670.5sin 1.070.5cos 3.740.5sin 0.79t u v w x t u v w y t u v w x t u v w y +++-=⎧⎪+++-=⎪⎨+++-=⎪⎪+++-=⎩ 得出相应的点,t u （）； 4) 选择拉格朗日插值法，将,t u （）作为中间变量，在题目所给出的二维数表中进行二次代数插值，得到[][],)(z f x i y j =；5) 输出数表：[][][][]()()0,1,2,,10;0,1,2,,20,,,x i y j f x i y j i j =⋯=⋯； 6) 令k=0；7) 以()()(),,,0,1,r r r s x x y y r s ϕψ===…,k 为拟合基函数，将上述数表作为拟合条件，对于给定的k 值，得到矩阵B 、G 、U ；8) 令-1-1(),()T T T A B B B U C AG G G ==，用选主元的LU 分解法分别计算矩阵A 和C 的各列，最后得到系数矩阵C ；9) 以公式：()()()00,k ki j rs r i s j s r p x y C x y ϕψ===∑∑计算每个点的拟合值；10) 利用公式：()()()2102000,,i j i j i j f x y p x y σ===-∑∑计算拟合误差，当σ≤10-7时，循环结束，否则k=k+1，转（6）；11) 令[][]()**0.10.50.2 1,2,81,2,5x i i y j j i j ==+=⋯=⋯；，；，；12) 计算()()()******,,,,,i j i j i jf x y p x y delta x y ，输出数表，观察逼近效果； 13) 结束。

数值分析练习题附答案

1
2-3 对矩阵 A 进行 LDLT 分解和 GGT 分解，求解方程组 Ax=b,其中
16 4 8
1
A=( 4 5 −4) ， b＝(2)
8 −4 22
3
解：（注：课本 P26 P27 根平方法）
设 L=(l i j )，D=diag(di)，对 k=1,2,…,n,
其中��=��－∑��=−11 ��2��
��31=(��31 − ∑0��=1 ��3��1�� )/ ��1=186=12
��32=(��32
−
∑1��=1
��3��2��
�� )/
6.6667
，得 ��3 ＝ 1.78570
−1 209
��4
0
��4
0.47847
(
56
−1
780 (��5) 209)
(200)
(��5) ( 53.718 )
1 −1
4
1 −4
15
��1
25
��2
6.6667再由1源自− 15561
− 56
209

x (k1) 1

1 5
(12

2 x2( k )

x (k) 3
)

2 5
x (k) 2

中科院研究生院信息工程学院课件数值分析数值分析第三次作业及答案

中科院研究⽣院信息⼯程学院课件数值分析数值分析第三次作业及答案6数值分析第三次作业及答案明当h T 0时，它收敛于原初值问题的准确解 y证：梯形公式为 y n ⼗ yn+—[f(X n ，y n )+f(X n^1，y n 』] h由 f (X,y) = —y= y n+ =y n +2( — y n — yn G=L n = 3 Y y n」訓 /乂⼚％l 2+h <12 +h ⼃丫2. (P202(6))写出⽤四阶经典的龙格⼀库塔⽅法求解下列初值问题的计算公式:y n + =y n + — (k 1 +2k 2 +2k 3 +k 4)=0.2214x n +1.2214y n +0.0214 6飞=3y n ⼼+x n )2)” k 2 =3(y n +0.1k 1〃(1+Xn +0.1) )L s =3(yn +0.1k2”(1+Xn +0.1)k 4 =3(yn +0.2k 3”(1+Xn +0.2) 0 2yn + =y n +〒(k 1 +2k 2 +2k 3 +k 4).1. ( P201 (4))⽤梯形⽅法解初值问题〔爲证明其近似解为y n 偌〕:并证⽤ . f 2-h 1 因 yoT " yn F ⼃.⽤上述梯形公式以步长h 经n 步计算到y n ,故有nh ：=x.X◎ T 茹Jf 2—h \n7 l 2+h ⼃1) ]y =x + y, 0 e x ￡1; ly(0) =1;2)l y \3%+x)，O *1； [y(0)=1.解：令h =0.2k 1 = f (X n , y n )= h k2=f (Xn+；；,yn+-k1)=Xn+- + 2 2 2 h k s = f (X n +；, y n +-k 2)=X n +- +y n +-k 2 =1.11(X n + y n )+0.11 2 2 2 2X n +y nh 1)4h h ??yn +；；k i =1.1(Xn +y n )+0.1 2 2 h . . h .................................. .2 ⼋ 2 J 'k 4 = f(X n +h,y n +hk 3)=X n + h + y n +hk 3 =1.222(X n +y n )+0.2223. (P202(7))证明对任意参数t，下列龙格库塔—公式是⼆阶的:r hy n 卄yn+^g+G)；* K i = f (X n, yj；K2 = f (X n +th, y n +thK i);[K3 =f(Xn+(1—t)h,yn+(1—t)hK i).证：由⼀元函数的泰勒展开有2 '''"y(X nG =y(X n) +hy'(X n)⼸[f x(X n,y(X n)) +f y(X n,y(X n))f(X n,y(X n))]中严h'2 3!⼜由⼆元函数的泰勒展开有y n41 =y n +；2(⼼+K3)=y n +；2[(f(X n,y n) + ￡%区『)也+f y(X n, y n)thf (X n, Y n)⼗。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

大作业三1. 给定初值0x 及容许误差，编制牛顿法解方程f (x )=0的通用程序. 解：Matlab 程序如下：函数m 文件：fu.mfunction Fu=fu(x)Fu=x^3/3-x;end函数m 文件：dfu.mfunction Fu=dfu(x)Fu=x^2-1;end用Newton 法求根的通用程序Newton.mclear;x0=input('请输入初值x0:');ep=input('请输入容许误差:');flag=1;while flag==1x1=x0-fu(x0)/dfu(x0);if abs(x1-x0)<epflag=0;endx0=x1;endfprintf('方程的一个近似解为：%f\n',x0);寻找最大δ值的程序：Find.mcleareps=input('请输入搜索精度：');ep=input('请输入容许误差：');flag=1;k=0;x0=0;while flag==1sigma=k*eps;x0=sigma;k=k+1;m=0;flag1=1;while flag1==1 && m<=10^3x1=x0-fu(x0)/dfu(x0);if abs(x1-x0)<ep flag1=0;endm=m+1;x0=x1;endif flag1==1||abs(x0)>=epflag=0;endendfprintf('最大的sigma 值为：%f\n',sigma);2．求下列方程的非零根5130.6651()ln 05130.665114000.0918x x f x x +⎛⎫=-= ⎪-⨯⎝⎭解：Matlab 程序为：（1）主程序clearclcformat longx0=765;N=100;errorlim=10^(-5);x=x0-f(x0)/subs(df(),x0);n=1;while n<Nx=x0-f(x0)/subs(df(),x0);if abs(x-x0)>errorlimn=n+1;elsebreak;endx0=x;enddisp(['迭代次数: n=',num2str(n)])disp(['所求非零根: 正根x1=',num2str(x),' 负根x2=',num2str(-x)])（2）子函数非线性函数ffunction y=f(x)y=log((513+0.6651*x)/(513-0.6651*x))-x/(1400*0.0918);end（3）子函数非线性函数的一阶导数dffunction y=df()syms x1y=log((513+0.6651*x1)/(513-0.6651*x1))-x1/(1400*0.0918);y=diff(y);end运行结果如下：迭代次数: n=5所求非零根: 正根x1=767.3861 负根x2=-767.3861大作业四试编写MATLAB 函数实现Newton 插值，要求能输出插值多项式. 对函数21()14f x x=+在区间[-5，5]上实现10次多项式插值.分析：（1）输出插值多项式。

（2）在区间[-5,5]内均匀插入99个节点，计算这些节点上函数f （x ）的近似值，并在同一张图上画出原函数和插值多项式的图形。

（3）观察龙格现象，计算插值函数在各节点处的误差，并画出误差图。

解：Matlab 程序代码如下：%此函数实现y=1/(1+4*x^2)的n次Newton插值，n由调用函数时指定%函数输出为插值结果的系数向量（行向量）和插值多项式function [t y]=func5(n)x0=linspace(-5,5,n+1)';y0=1./(1.+4.*x0.^2);b=zeros(1,n+1);for i=1:n+1s=0;for j=1:it=1;for k=1:iif k~=jt=(x0(j)-x0(k))*t;end;end;s=s+y0(j)/t;end;b(i)=s;end;t=linspace(0,0,n+1);for i=1:ns=linspace(0,0,n+1);s(n+1-i:n+1)=b(i+1).*poly(x0(1:i)); t=t+s;end;t(n+1)=t(n+1)+b(1);y=poly2sym(t);10次插值运行结果：[b Y]=func5(10)b =Columns 1 through 4-0.0000 0.0000 0.0027 -0.0000Columns 5 through 8-0.0514 -0.0000 0.3920 -0.0000Columns 9 through 11-1.1433 0.0000 1.0000Y =b为插值多项式系数向量，Y为插值多项式。

插值近似值：x1=linspace(-5,5,101);x=x1(2:100);y=polyval(b,x)y =Columns 1 through 122.70033.99944.3515 4.0974 3.4926 2.7237 1.9211 1.1715 0.5274 0.0154 -0.3571 -0.5960Columns 13 through 24-0.7159 -0.7368 -0.6810 -0.5709 -0.4278 -0.2704 -0.1147 0.0270 0.1458 0.2360 0.2949 0.3227Columns 25 through 360.3217 0.2958 0.2504 0.1915 0.1255 0.0588 -0.0027 -0.0537 -0.0900 -0.1082 -0.1062 -0.0830Columns 37 through 48-0.0390 0.0245 0.1052 0.2000 0.3050 0.4158 0.5280 0.6369 0.7379 0.8269 0.9002 0.9549Columns 49 through 600.9886 1.0000 0.9886 0.9549 0.9002 0.8269 0.7379 0.6369 0.5280 0.4158 0.3050 0.2000Columns 61 through 720.1052 0.0245 -0.0390 -0.0830 -0.1062 -0.1082 -0.0900 -0.0537 -0.0027 0.0588 0.1255 0.1915Columns 73 through 840.2504 0.2958 0.3217 0.3227 0.2949 0.2360 0.1458 0.0270 -0.1147 -0.2704 -0.4278 -0.5709Columns 85 through 96-0.6810 -0.7368 -0.7159 -0.5960 -0.3571 0.0154 0.5274 1.1715 1.9211 2.7237 3.4926 4.0974Columns 97 through 994.3515 3.9994 2.7003绘制原函数和拟合多项式的图形代码：plot(x,1./(1+4.*x.^2))hold allplot(x,y,'r')xlabel('X')ylabel('Y')title('Runge现象')gtext('原函数')gtext('十次牛顿插值多项式')绘制结果：误差计数并绘制误差图：hold offey=1./(1+4.*x.^2)-yey =Columns 1 through 12-2.6900 -3.9887 -4.3403 -4.0857 -3.4804 -2.7109 -1.9077 -1.1575 -0.5128 -0.0000 0.3733 0.6130Columns 13 through 240.7339 0.7558 0.7010 0.5921 0.4502 0.2943 0.1401 0.0000 -0.1169 -0.2051 -0.2617 -0.2870Columns 25 through 36-0.2832 -0.2542 -0.2053 -0.1424 -0.0719 -0.0000 0.0674 0.1254 0.1696 0.1971 0.2062 0.1962Columns 37 through 480.1679 0.1234 0.0660 0.0000 -0.0691 -0.1349 -0.1902 -0.2270 -0.2379 -0.2171 -0.1649 -0.0928Columns 49 through 60-0.0271 0 -0.0271 -0.0928 -0.1649 -0.2171 -0.2379 -0.2270 -0.1902 -0.1349 -0.0691 0.0000Columns 61 through 720.0660 0.1234 0.1679 0.1962 0.2062 0.1971 0.1696 0.1254 0.0674 0.0000 -0.0719 -0.1424Columns 73 through 84-0.2053 -0.2542 -0.2832 -0.2870 -0.2617 -0.2051 -0.1169 0.0000 0.1401 0.2943 0.4502 0.5921Columns 85 through 960.7010 0.7558 0.7339 0.6130 0.3733 0.0000 -0.5128 -1.1575 -1.9077 -2.7109 -3.4804 -4.0857Columns 97 through 99-4.3403 -3.9887 -2.6900plot(x,ey)xlabel('X')ylabel('ey')title('Runge现象误差图')输出结果为：大作业五解：Matlab程序为：x = [-520,-280,-156.6,-78,-39.62,-3.1,0,3.1,39.62,78,156.6,280,520]'; y = [0,-30,-36,-35,-28.44,-9.4,0,9.4,28.44,35,36,30,0]';n = 13;%求解Mfor i = 1:1:n-1h(i) = x(i+1)-x(i);endfor i = 2:1:n-1a(i) = h(i-1)/(h(i-1)+h(i));b(i) = 1-a(i);c(i) = 6*((y(i+1)-y(i))/h(i)-(y(i)-y(i-1))/h(i-1))/(h(i-1)+h(i)); enda(n) = h(n-1)/(h(1)+h(n-1));b(n) = h(1)/(h(1)+h(n-1));c(n) = 6/(h(1)+h(n-1))*((y(2)-y(1))/h(1)-(y(n)-y(n-1))/h(n-1)); A(1,1) = 2;A(1,2) = b(2);A(1,n-1) = a(2);A(n-1,n-2) = a(n);A(n-1,n-1) = 2;A(n-1,1) = b(n);for i = 2:1:n-2A(i,i) = 2;A(i,i+1) = b(i+1);A(i,i-1) = a(i+1);endC = c(2:n);C = C';m = A\C;M(1) = m(n-1);M(2:n) = m;xx = -520:10.4:520;for i = 1:51for j = 1:1:n-1if x(j) <= xx(i) && xx(i) < x(j+1)break;endendyy(i) =M(j+1)*(xx(i)-x(j))^3/(6*h(j))-M(j)*(xx(i)-x(j+1))^3/(6*h(j))+(y(j+1)-M( j+1)*h(j)^2/6)*(xx(i)-x(j))/h(j)-(y(j)-M(j)*h(j)^2/6)*(xx(i)-x(j+1))/h(j );end;for i = 52:101yy(i) = -yy(102-i);end;for i = 1:50xx(i) = -xx(i);endplot(xx,yy);hold on;for i = 1:1:n/2x(i) = -x(i);endplot(x,y,'bd');title('机翼外形曲线');输出结果：运行jywx.m文件，得到2. （1）编制求第一型3 次样条插值函数的通用程序；（2）已知汽车门曲线型值点的数据如下：解：（1）Matlab编制求第一型3 次样条插值函数的通用程序：function [Sx] = Threch(X,Y,dy0,dyn)%X为输入变量x的数值%Y为函数值y的数值%dy0为左端一阶导数值%dyn为右端一阶导数值%Sx为输出的函数表达式n = length(X)-1;d = zeros(n+1,1);h = zeros(1,n-1);f1 = zeros(1,n-1);f2 = zeros(1,n-2);for i = 1:n %求函数的一阶差商 h(i) = X(i+1)-X(i);f1(i) = (Y(i+1)-Y(i))/h(i);endfor i = 2:n %求函数的二阶差商 f2(i) = (f1(i)-f1(i-1))/(X(i+1)-X(i-1));d(i) = 6*f2(i);endd(1) = 6*(f1(1)-dy0)/h(1);d(n+1) = 6*(dyn-f1(n-1))/h(n-1); %赋初值A = zeros(n+1,n+1);B = zeros(1,n-1);C = zeros(1,n-1);for i = 1:n-1B(i) = h(i)/(h(i)+h(i+1)); C(i) = 1-B(i);endA(1,2) = 1;A(n+1,n) = 1;for i = 1:n+1A(i,i) = 2;endfor i = 2:nA(i,i-1) = B(i-1);A(i,i+1) = C(i-1);endM = A\d;syms x;for i = 1:nSx(i) =collect(Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(x-X(i))+M(i)/2*(x-X(i))^2+( M(i+1)-M(i))/(6*h(i))*(x-X(i))^3);digits(4);Sx(i) = vpa(Sx(i));endfor i = 1:ndisp('S(x)=');fprintf(' %s (%d,%d)\n',char(Sx(i)),X(i),X(i+1));endS = zeros(1,n);for i = 1:nx = X(i)+0.5;S(i) =Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(x-X(i))+M(i)/2*(x-X(i))^2+(M(i+1)-M (i))/(6*h(i))*(x-X(i))^3;enddisp('S(i+0.5)');disp(' i X(i+0.5) S(i+0.5) ');for i = 1:nfprintf(' %d %.4f %.4f\n',i,X(i)+0.5,S(i));endend输出结果：>> X = [0 1 2 3 4 5 6 7 8 9 10];>> Y = [2.51 3.30 4.04 4.70 5.22 5.54 5.78 5.40 5.57 5.70 5.80]; >> Threch(X,Y,0.8,0.2)S(x)=0.8*x - 0.001486*x^2 - 0.008514*x^3 + 2.51 (0,1)S(x)=0.8122*x - 0.01365*x^2 - 0.004458*x^3 + 2.506 (1,2)S(x)=0.8218*x - 0.01849*x^2 - 0.003652*x^3 + 2.499 (2,3)S(x)=0.317*x^2 - 0.1847*x - 0.04093*x^3 + 3.506 (3,4)S(x)=6.934*x - 1.463*x^2 + 0.1074*x^3 - 5.986 (4,5)S(x)=4.177*x^2 - 21.26*x - 0.2686*x^3 + 41.01 (5,6)S(x)=53.86*x - 8.344*x^2 + 0.427*x^3 - 109.2 (6,7)S(x)=6.282*x^2 - 48.52*x - 0.2694*x^3 + 129.6 (7,8)S(x)=14.88*x - 1.643*x^2 + 0.06076*x^3 - 39.42 (8,9) S(x)=8.966*x - 0.986*x^2 + 0.03641*x^3 - 21.67 (9,10) S(i+0.5)i X(i+0.5) S(i+0.5)1 0.5000 2.90862 1.5000 3.67843 2.5000 4.38154 3.5000 4.98825 4.5000 5.38336 5.5000 5.72377 6.5000 5.59438 7.5000 5.43029 8.5000 5.658510 9.5000 5.7371ans =[ - 0.008514*x^3 - 0.001486*x^2 + 0.8*x + 2.51, - 0.004458*x^3 - 0.01365*x^2 + 0.8122*x + 2.506, - 0.003652*x^3 - 0.01849*x^2 + 0.8218*x + 2.499, - 0.04093*x^3 + 0.317*x^2 - 0.1847*x + 3.506, 0.1074*x^3 - 1.463*x^2 + 6.934*x - 5.986, - 0.2686*x^3 + 4.177*x^2 - 21.26*x + 41.01, 0.427*x^3 - 8.344*x^2 + 53.86*x - 109.2, - 0.2694*x^3 + 6.282*x^2 - 48.52*x + 129.6, 0.06076*x^3 - 1.643*x^2 + 14.88*x - 39.42, 0.03641*x^3 - 0.986*x^2 + 8.966*x - 21.67]大作业六1、炼钢厂出钢时所用的圣刚睡的钢包，在使用过程中由于钢液及炉渣对包衬耐火材料的侵蚀，使其容积不断增大，经试验，钢包的容积与相应的使用次数的数据如下：（使用次数x，容积y）线x b a 1*y 1+=对使用最小二乘法数据进行拟合。