1.6.4 分段三次Hermite 插值
为了利用多项式插值方法而又克服高次插值多项式的缺陷，便引入了分段插值的概念。它的基本思想是把函数整个区间上分成许多段，每段都选用适当的低次插值多项式代替函数，整体上按一定的要求连接起来，构成一个分段的插值函数。
为此，把函数)(x f 的自变量x 在区间],[b a 上用)1(+n 个节点分割成n 段：
b x x x x x a n n =<<<<<=-1210
根据这些节点的取值
i x ，)(x f 在节点上的函数值i i y x f =)(和导数值
i i m x f =')(),,2,1,0(n i =，可以构造一个分段三次插值函数)(x H ，它满足
下述条件：
①i i y x H =)(,i i y x H '=')(),,2,1,0(n i =。
② 在每个小区间],[1+i i x x ),,2,1,0(n i =上，都是一个三次多项式：
3
32210)(x
a x a x a a x H i i i i i +++=
把这样构成的分段三次函数)(x H 称为分段三次Hermite 插值函数，它的
各小段均为三次多项式，而整体上具有一阶连续导数。
由式（1-34）可直接写出分段三次Hermite 插值函数的分段表达式
12
112
1112
1112
111)()(2121)(++++++++++++'???
?
??---+'???? ??---+???
?
??--???? ?
?--++???? ??--???? ??--+=i i i i
i i i i i i i i i i
i i i i i i i i i i
y x
x x x x x y x x x x x x y x
x x x x x x x y x x x x x
x x x x H
也可通过构造基函数给出分段三次Hermite 插值函数的表达式。参照分段线性插值与Hermite 插值基函数公式（1-31）和式（1-32），可得出分段三次
Hermite 插值的基函数为
?????∈∈???
?
??--???? ??--+=]
,[0
],[21)(1102
1010100n x x x x x x x x x x x x x x x h
)1,,1(],(),[0]
,(21],[21)(11012
11112
111-=?????
??????∈∈???? ??--???? ??
--+∈???? ??--???? ??--+=+-++++----n i x x x x x x x x x x x x x x x x x x x x x x x x x x x x h n i i i i i i i i i i i i i i i i i i i （1-38）
?????∈∈?
??? ??--???? ??--+=-----]
,[0
],[21)(1012
111n n n n n
n n n n n x x x x x x x x x x x x x x x h
?????∈∈????
??---=]
,[0
],[)()(1102
10
100n x x x x x x x x x x x x x H
)1,,1(],(),[0],[)(],[)()(11012
1112
11-=?????
?
?????∈∈?
??? ??---∈???? ??---=+-+++---n i x x x x x x x x x x x x x x x x x x x x x x x x H n i i i i i i i i i i i i i i i （1-39）
?????∈∈????
??---=----]
,[0
],[)()(1012
11n n n n n
n n n x x x x x x x x x x x x x H
分段三次Hermite 插值函数为
])()([)(0
∑='+=n
i i i i i x H y x h y x H （1-40）
由余项公式（1-37）可以导出，分段三次Hermite 插值的误差有如下估计
)
(max 384)()()()
4(4
x f h x H x f x R b x a ≤≤≤-= （1-41）
其中
)(max 11
0i i n i x x h -=+-≤≤。
分段三次Hermite 插值函数是插值区间上的光滑函数，它与函数)(x f 在节点处密合程度较好。
1.6.5 分段三次Hermite 插值的MATLAB 实现
【例1-18】给定函数55,11
)(2
≤≤-+=x x
x f ，取插值节点k x k +-=5，10,,1,0 =k ，用分段Hermite 插值计算)(x ?，并画出)(x ?与)(x f 的图形。
解：在编辑窗口输入以下命令： a=-5;b=5;n=10;h=(b-a)/n x=a:h:b;y=1./(1+x.^2); xx=a:0.01:b;yy=1./(1+xx.^2); m=length(xx);z=zeros(1,m);
for i=1:m
z(i)=Hermite_wise(x,y,[],xx(i)); end
plot(x,y,'bo',xx,yy,'r-',xx,z,'-','linewidth',2);grid 命令执行后得到如图1-10所示图形：
分段三次Hermite 插值函数是插值区间上的光滑函数，它与函数)(x f 在
节点处密合程度较好。
图1-10 分段Hermite 插值图形
上例中分段三次Hermite 调用函数如下。 function yi=Hermite_wise(x,y,ydot,xi) % Hermite 分段三次插值的MATLAB 实现 % 分段Hermite 插值公式，其中 % x 为向量，全部的插值节点； % y 为向量，插值节点处的函数值； % ydot 为向量，插值节点处的导数值， % 如果此处值缺省，则用均差代替导数， % 端点用向前、向后均差，中间点用中心均差； % xi 为标量，自变量x ； % yi 为xi 处的函数估计值；
% 如果没有给出y的导数值，则用均差代替导数。
if isempty(ydot)==1
ydot=gradient(y,x);
end
n=length(x);m1=length(y);m2=length(ydot);
% 输入x,y和y的导数的个数必须相同.
if n~=m1|n~=m2|m1~=m2
error('The length of X,Y and Ydot must be equal!');
return;
end
for k=1:n-1
% 插值节点必须互异
if abs(x(k)-x(k+1)) error('The DATA is error!'); return; end if x(k)<=xi&xi<=x(k+1) yi=y(k)*(1-2*(xi-x(k))/(x(k)-x(k+1)))... *(xi-x(k+1))^2/(x(k)-x(k+1))^2+... y(k+1)*(1-2*(xi-x(k+1))/(x(k+1)-x(k)))... *(xi-x(k))^2/(x(k+1)-x(k))^2+... ydot(k)*(xi-x(k))*(xi-x(k+1))^2/(x(k)-x(k+1))^2+... ydot(k+1)*(xi-x(k+1))*(xi-x(k))^2/(x(k+1)-x(k))^2; return; end end 1.7 有关插值的MATLAB 命令 1.7.1 三次插值和三次样条插值的MATLAB 命令 对于三次插值(Hermite 插值)和三次样条插值，MATLAB 中都设有专有命令。 三次插值命令调用格式为： )x y,pchip(x,y k k = 三次样条插值命令调用格式为： )x y,spline(x,y k k = 参数 k k y x y x ,,,的意义及要求与 interp1中的完全一样，插值效果与 interp1中参数"method"分别选用'pchip'和'spline'等价。它允许k x 在区间 max(x)][min(x),外的附近取近似值，因为程序设计中在两个端点处，都使 用了由内部外推的方法。 【例1-19】已知)(x f y =的函数关系中，当x=[-6.2 -5 -3 -1.7 0 1 3 5 6.2 7.5 8.2]时，y=[1 0 -1 0 1 0 -1 0 0.5 -1 0 ]。用“线性插值”、“最近插值”分别求出当xi=-6.2:0.2:8.2时所对应的函数值，并画出曲线图。 解：在编辑窗口输入： x=[-6.2 -5 -3 -1.7 0 1 3 5 6.2 7.5 8.2]; y=[1 0 -1 0 1 0 -1 0 0.5 -1 0 ]; xi=-6.2:0.2:8.2; y1=interp1(x,y,xi,'linear'); y2=interp1(x,y,xi,'nearest'); plot(x,y,'ro',xi,y1,xi,y2,'-', 'linewidth',2); legend('数据节点','线性插值','最近插值'); 执行命令后得到如图1-11所示图形。 图1-11 线性插值与最近插值比较图 【例1-20】已知)(x f y 的函数关系中，当x=[0 2 3 4 5 6]时，y=[-4 1 -2 -4 -4 -5 ]。用“分段三次插值”和“三次样条插值”分别求出当xi=-0.3:0.3:6.8时所对应的函数值，并画出曲线。 解：在编辑窗口输入： x=[0 2 3 4 5 6]; y=[-4 1 -2 -4 -4 -5 ]; xi=-0.3:0.3:6.8; y1=interp1(x,y,xi,'pchip') y2=interp1(x,y,xi,'spline') plot(x,y,'ro',xi,y1,'r',xi,y2,'-','linewidth',2); legend('样本点','分段三次插值','三次样条插值');grid; 图1-12 三次插值与三次样条插值比较图 执行命令得到如图1-12所示的结果。从图形上看，用三次样条插值时曲线的光滑度要好一些。当自变量在限定区间外边界处的曲线上时，三次样条插值和三次插值法的外延趋势不同，差异明显。 xi=-0.3:0.3:6.8时所对应的函数值： y1 = -6.0754 -4.0000 …… -5.4365 -5.7920 y2 = -7.4500 -4.0000 …… -6.0790 -7.7050 其中y1 、y2分别有24个数据点。 y1 、y2两句末尾不加分号，就会输出相应的函数值。 【例7-19】对离散地分布在 x e y x sin 函数曲线上的数据点进行样条插值计算。 解：在编辑窗口输入命令： x=[0 2 4 5 8 12 12.8 17.2 19.9 20]; y=exp(x).*sin(x); xx=0:0.5:20; yy=spline(x,y,xx) plot(x,y,'o',xx,yy) 执行命令得到如图1-12所示的结果 图1-12 三次样条插值 yy = 1.0e+008 * 0 0.0000 …… 1.8156 4.4293 yy 共有41个数据点。 补充知识 1.7.4 三次样条插值的MA TLAB 实现实例 【例1-21】在区间的两个端点处给出条件，称为边界条件。给定两端点 处的导数值n y b S y a S '=''=')(,)(0，然后根据下面的数据点求出其三次样条插值 多项式，并计算当5.3=x 时y 的值。 解：在MATLAB 编辑窗口中输入以下命令： x=1:8; y=[0.84 0.91 0.14 -0.76 -0.96 -0.28 0.66,0.99] [f,f0]=ThrSample1(x,y,0.54,-0.15,3.5) 程序运行后得到如下结果。 f = (-38/25*t+133/25)*(t-5)^2+(-264/25+48/25*t)*(t-4)^2+(-9447/14555*t+3778 8/14555)*(5-t)^2-(25502967261225075/18014398509481984-5100593452245015/ 18014398509481984*t)*(t-4)^2 f0 = -0.9164 程序运行时调用的三次样条函数如下（给定两端点处的一阶导数值）。 function [f,f0] = ThrSample1 (x,y,y_1, y_N,x0) syms t; f = 0.0; f0 = 0.0; if(length(x) == length(y)) n = length(x); else disp('x和y的维数不相等！'); return; end %维数检查 for i=1:n if(x(i)<=x0)&& (x(i+1)>=x0) index = i; break; end end %找到x0所在区间 A = diag(2*ones(1,n)); %求解m的系数矩阵 u = zeros(n-2,1); lamda = zeros(n-1,1); c = zeros(n,1); for i=2:n-1 u(i-1) = (x(i)-x(i-1))/(x(i+1)-x(i-1)); lamda(i) = (x(i+1)-x(i))/(x(i+1)-x(i-1)); c(i) = 3*lamda(i)*(y(i)-y(i-1))/(x(i)-x(i-1))+ ... 3*u(i-1)*(y(i+1)-y(i))/(x(i+1)-x(i)); A(i, i+1) = u(i-1); A(i, i-1) = lamda(i); %形成系数矩阵及向量c end c(1) = 2*y_1; c(n) = 2*y_N; m = followup(A,c); % 用追赶法求解方程组 h = x(index+1) - x(index); % x0所在区间长度 f = y(index)*(2*(t-x(index))+h)*(t-x(index+1))^2/h/h/h + ... y(index+1)*(2*(x(index+1)-t)+h)*(t-x(index))^2/h/h/h + ... m(index)*(t-x(index))*(x(index+1)-t)^2/h/h - ... m(index+1)*(x(index+1)-t)*(t-x(index))^2/h/h; % x0所在区间的插值函数 f0 = subs(f,'t',x0); % x0处的插值 下面函数是ThrSample1 (x,y,y_1, y_N,x0)函数调用运行时调用的函数。function x=followup(A,b) n = rank(A); for(i=1:n) if(A(i,i)==0) disp('Error: 对角有元素为0！'); return; end end; d = ones(n,1); a = ones(n-1,1); c = ones(n-1); for(i=1:n-1) a(i,1)=A(i+1,i); c(i,1)=A(i,i+1); d(i,1)=A(i,i); end d(n,1) = A(n,n); for(i=2:n) d(i,1)=d(i,1) - (a(i-1,1)/d(i-1,1))*c(i-1,1); b(i,1)=b(i,1) - (a(i-1,1)/d(i-1,1))*b(i-1,1); end x(n,1) = b(n,1)/d(n,1); for(i=(n-1):-1:1) x(i,1) = (b(i,1)-c(i,1)*x(i+1,1))/d(i,1); end 【例1-22】在区间的两个端点处给出条件，称为边界条件。给定两端处的 二阶导数值n y b S y a S ''=''''='')(,)(0 ，然后根据下面的数据点求出其三次样条插值多项式，并计算当5.3=x 时 y 的值。 解：在MATLAB 编辑窗口中输入以下命令： x=1:8; y=[0.84 0.91 0.14 -0.76 -0.96 -0.28 0.66 0.99]; [f,f0]=ThrSample2(x,y,-0.84,-0.99,3.5) 程序运行后得到如下结果。 f = (7/25*t-7/10)*(t-4)^2+(-171/25+38/25*t)*(t-3)^2+(-4455276049781723/4503599627370496*t+13365828149345169/4503599627370496)*(4-t)^2-(-7517/2911+7517/11644*t)*(t-3)^2 f0 = -0.3530 程序运行时调用的三次样条函数如下（给定两端点处的二阶导数值）。 function [f,f0] = ThrSample2 (x,y,y2_1, y2_N,x0) syms t; f = 0.0; f0 = 0.0; if(length(x) == length(y)) n = length(x); else disp('x和y的维数不相等！'); return; end %维数检查 for i=1:n if(x(i)<=x0)&& (x(i+1)>=x0) index = i; break; end end %找到x0所在区间for i=1:n if(x(i)<=x0)&& (x(i+1)>=x0) index = i; break; end end %找到x0所在区间A = diag(2*ones(1,n)); %求解m的系数矩阵A(1,2) = 1; A(n,n-1) = 1; u = zeros(n-2,1); lamda = zeros(n-1,1); c = zeros(n,1); for i=2:n-1 u(i-1) = (x(i)-x(i-1))/(x(i+1)-x(i-1)); lamda(i) = (x(i+1)-x(i))/(x(i+1)-x(i-1)); c(i) = 3*lamda(i)*(y(i)-y(i-1))/(x(i)-x(i-1))+ ... 3*u(i-1)*(y(i+1)-y(i))/(x(i+1)-x(i)); A(i, i+1) = u(i-1); A(i, i-1) = lamda(i); %形成系数矩阵及向量c end c(1) = 3*(y(2)-y(1))/(x(2)-x(1))-(x(2)-x(1))*y2_1/2; c(n) = 3*(y(n)-y(n-1))/(x(n)-x(n-1))-(x(n)-x(n-1))*y2_N/2; m = followup(A,c); %用追赶法求解方程组 h = x(index+1) - x(index); %x0所在区间长度 f = y(index)*(2*(t-x(index))+h)*(t-x(index+1))^2/h/h/h + ... y(index+1)*(2*(x(index+1)-t)+h)*(t-x(index))^2/h/h/h + ... m(index)*(t-x(index))*(x(index+1)-t)^2/h/h - ... m(index+1)*(x(index+1)-t)*(t-x(index))^2/h/h; %x0所在区间的插值函数 f0 = subs(f,'t',x0); %x0处的插值 【例1-23】由表中的数据，用程序求三次样条插值函数。 解：用程序进行计算。 在命令窗口输入： >> x=[0,1,2,3];y=[0,0.5,2,1.5];dx0=0.2;dxn=-1; >> s=csfit(x,y,dx0,dxn) 回车得到： s = 0.4800 -0.1800 0.2000 0 -1.0400 1.2600 1.2800 0.5000 0.6800 -1.8600 0.6800 2.0000 即求得： 。 x x x x S x x x x S x x x x S ]3,2[, 268.086.168.0,]2,1[, 5.028.12 6.104.1, ]1,0[,2.018.048.02 3 32322 31∈++-=∈+++-=∈+-= 程序运行时调用的函数如下。 function S=csfit(x,y,dx0,dxn) % x,y 分别为n 个节点的横坐标所组成的向量及纵坐标组成的向量 % dx0,dn 分别为S 的导数在x0,xn 处的值 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)); 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); end m(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); end m(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*m(k+1)+m(k+2))/6; S(k+1,4)=y(k+1); end 特别声明：考试时需带计 算器作辅助计算 1．2015x *=是经四舍五入得到的近似值，则其相对误差* r e ≤-31 104 ?. 2. 01(),(), ,()n l x l x l x 是以01,, ,n x x x 为节点的拉格朗日插值基函数，则 3．设(0)1(1)3(2)4(3)2f =,f =,f =,f =，[0123]f =，，，1 3 - . 4. 利用Simpson 公式求?2 1 2dx x = 7.3 5. 设求积公式1 0()d (),(1)n k k k f x x A f x n ≈≥∑?=是Gauss 型求积公式，则3 n k k k A x == ∑1 .4 6. 数值微分公式(2)(2) ()i i i f x h f x h f x h +≈ --'的截断误差为 2().O h 7. 设1101A ?? = ??? ，则A 的谱半径()A ρ= 1 ，A 的条件数1cond ()A = 4. 8. 用牛顿下山法求解方程3 03 x x -=根的迭代公式是 2 13 3(1),3n n n n x x x x x λ+-=-- 下山条件是 1()().n n f x f x +< 9.对任意初始向量(0)x 及任意向量f ，线性方程组的迭代公式(1)()(0,1,2,)k k k +=+=x Bx f ，迭代序列()k x 收敛于方程组的精确解x *的充分必要条件是()1.ρ 数值分析 第二章 2．当1,1,2x =-时，()0,3,4f x =-,求()f x 的二次插值多项式。 解： 0120121200102021101201220211,1,2, ()0,()3,()4;()()1 ()(1)(2)()()2()()1 ()(1)(2) ()()6 ()()1 ()(1)(1) ()()3 x x x f x f x f x x x x x l x x x x x x x x x x x l x x x x x x x x x x x l x x x x x x x ==-===-=--==-+-----==------= =-+-- 则二次拉格朗日插值多项式为 2 20 ()()k k k L x y l x ==∑ 0223()4() 14 (1)(2)(1)(1)23 537623 l x l x x x x x x x =-+=---+ -+= +- 6．设,0,1,,j x j n =L 为互异节点，求证： （1） 0()n k k j j j x l x x =≡∑ (0,1,,);k n =L （2）0 ()()0n k j j j x x l x =-≡∑ (0,1,,);k n =L 证明 （1） 令()k f x x = 若插值节点为,0,1,,j x j n =L ，则函数()f x 的n 次插值多项式为0 ()()n k n j j j L x x l x == ∑。 插值余项为(1)1() ()()()()(1)! n n n n f R x f x L x x n ξω++=-= + 又,k n ≤Q (1)()0 ()0 n n f R x ξ+∴=∴= 0()n k k j j j x l x x =∴=∑ (0,1,,);k n =L 0 000 (2)()() (())()()(()) n k j j j n n j i k i k j j j i n n i k i i k j j i j x x l x C x x l x C x x l x =-==-==-=-=-∑∑∑∑∑ 0i n ≤≤Q 又 由上题结论可知 ()n k i j j j x l x x ==∑ ()()0 n i k i i k i k C x x x x -=∴=-=-=∑原式 ∴得证。 7设[]2 (),f x C a b ∈且()()0,f a f b ==求证： 21 max ()()max ().8 a x b a x b f x b a f x ≤≤≤≤''≤- 解：令01,x a x b ==，以此为插值节点，则线性插值多项式为 10 101010 ()() ()x x x x L x f x f x x x x x --=+-- =() () x b x a f a f b a b x a --=+-- 1()()0()0 f a f b L x ==∴=Q 又 插值余项为1011 ()()()()()()2 R x f x L x f x x x x x ''=-= -- 011 ()()()()2 f x f x x x x x ''∴= -- 数值计算方法学习心得 ------一个代码的方法是很重要，一个算法的思想也很重要，但 在我看来，更重要的是解决问题的方法，就像爱因斯坦说的内容比 思维本身更重要。 我上去讲的那次其实做了挺充分的准备，程序的运行，pdf文档，算法公式的推导，程序伪代码，不过有一点缺陷的地方，很多细节 没有讲的很清楚吧，下来之后也是更清楚了这个问题。 然后一学期下来，总的来说，看其他同学的分享，我也学习到 许多东西，并非只是代码的方法，更多的是章胜同学的口才，攀忠 的排版，小冯的深入挖掘…都是对我而言比算法更加值得珍惜的东西，又骄傲地回想一下，曾同为一个项目组的我们也更加感到做项 目对自己发展的巨大帮助了。 同时从这些次的实验中我发现以前学到的很多知识都非常有用。 比如说，以前做项目的时候，项目导师一直要求对于要上传的 文件尽量用pdf格式，不管是ppt还是文档，这便算是对产权的一种 保护。 再比如代码分享，最基础的要求便是——其他人拿到你的代码 也能运行出来，其次是代码分享的规范性，像我们可以用轻量级Ubuntu Pastebin，以前做过一小段时间acm，集训队里对于代码的分享都是推荐用这个，像数值计算实验我觉得用这个也差不多了，其 次项目级代码还是推荐github（被微软收购了），它的又是可能更 多在于个人代码平台的搭建，当然像readme文档及必要的一些数据 集放在上面都更方便一些。 然后在实验中，发现debug能力的重要性，对于代码错误点的 正确分析，以及一些与他人交流的“正规”途径，讨论算法可能出 错的地方以及要注意的细节等，比如acm比赛都是以三人为一小组，讨论过后，讲了一遍会发现自己对算法理解更加深刻。 然后学习算法，做项目做算法一般的正常流程是看论文，尽量 看英文文献，一般就是第一手资料，然后根据论文对算法的描述， 就是如同课上的流程一样，对算法进一步理解，然后进行复现，最 后就是尝试自己改进。比如知网查询牛顿法相关论文，会找到大量 可以参考的文献。 最后的最后，想说一下，计算机专业的同学看这个数值分析， 不一定行云流水，但肯定不至于看不懂写不出来，所以我们还是要 提高自己的核心竞争力，就是利用我们的优势，对于这种算法方面 的编程，至少比他们用的更加熟练，至少面对一个问题，我们能思 考出对应问题的最佳算法是哪一个更合适解决问题。 附记: 对课程的一些小建议： 1. debug的能力不容忽视，比如给一个关于代码实现已知错误的代码给同学们，让同学们自己思考一下，然后分享各自的debug方法，一步一步的去修改代码，最后集全班的力量完成代码的debug，这往往更能提升同学们的代码能力。 2. 课堂上的效率其实是有点低的，可能会给学生带来一些负反馈，降低学习热情。 3. 总的来说还是从这门课程中学到许多东西。 数值分析学习心得体会 数值分析习题集 （适合课程《数值方法A 》和《数值方法B 》） 长沙理工大学 第一章 绪 论 1. 设x >0,x 的相对误差为δ,求ln x 的误差. 2. 设x 的相对误差为2％,求n x 的相对误差. 3. 下列各数都是经过四舍五入得到的近似数,即误差限不超过最后一位的半个单位,试指出 它们是几位有效数字: *****123451.1021,0.031,385.6,56.430,7 1.0.x x x x x =====? 4. 利用公式(3.3)求下列各近似值的误差限: ********12412324(),(),()/,i x x x ii x x x iii x x ++其中**** 1234 ,,,x x x x 均为第3题所给的数. 5. 计算球体积要使相对误差限为1％,问度量半径R 时允许的相对误差限是多少? 6. 设028,Y =按递推公式 1n n Y Y -=( n=1,2,…) 计算到100Y .27.982(五位有效数字),试问计算100Y 将有多大误差? 7. 求方程2 5610x x -+=的两个根,使它至少具有四位有效数字27.982). 8. 当N 充分大时,怎样求2 1 1N dx x +∞+?? 9. 正方形的边长大约为100㎝,应怎样测量才能使其面积误差不超过1㎝2 ? 10. 设 212S gt = 假定g 是准确的,而对t 的测量有±0.1秒的误差,证明当t 增加时S 的绝对 误差增加,而相对误差却减小. 11. 序列 {}n y 满足递推关系1101n n y y -=-(n=1,2,…),若0 1.41y =≈(三位有效数字), 计算到 10y 时误差有多大?这个计算过程稳定吗? 12. 计算6 1)f =, 1.4≈,利用下列等式计算,哪一个得到的结果最好? 3 -- 13. ()ln(f x x =,求f (30)的值.若开平方用六位函数表,问求对数时误差有多大?若 第一章典型例题 例3 ln2=0.…，精确到10－3的近似值是多少 解 精确到10－3＝，即绝对误差限是＝, 故至少要保留小数点后三位才可以。ln2 第二章典型例题 例1 用顺序消去法解线性方程组 ??? ??1 -=4+2+4=+2+31 -=4++2321 321321x x x x x x x x x 解 顺序消元 ?? ?? ??????---???→???????????---????→???????????--=-?+-?+-?+1717005.555.00141 25.025.105.555.001412142141231412]b A [)3()2/1()2/3(231312r r r r r r M 于是有同解方程组 ?? ? ??-==--=++17175.555.0142332321x x x x x x 回代得解 x 3=－1, x 2=1,x 1=1,原线性方程组的解为X ＝(1,1,－1)T 例2 取初始向量X (0)=(0,0,0)T ,用雅可比迭代法求解线性方程组 ??? ??5 =+2+23=++1=2-2+321 321321x x x x x x x x x 解 建立迭代格式 ???????+--=+--=++-=+++5223122) (2)(1)1(3 ) (3)(1)1(2 ) (3)(2)1(1k k k k k k k k k x x x x x x x x x (k =1,2,3,…) 第1次迭代,k =0 X (0)＝0，得到X (1)＝(1,3,5)T 第2次迭代，k =1 ???????-=+?-?-=-=+--==+?+?-=3 532123 351515232)2(3) 2(2)2(1x x x X (2)＝(5,－3,－3)T 第3次迭代，k =2 ???????=+-?-?-==+---==+-?+-?-=1 5)3(2521 3)3(511)3(2)3(2)2(3) 3(2)3(1x x x X (3)＝(1,1,1)T 第4次迭代，k =3 ???????=+?-?-==+--==+?+?-=1 512121 311111212)2(3) 2(2)2(1x x x X (4)＝(1,1,1)T 例4 证明例2的线性方程组，雅可比迭代法收敛，而高斯－赛德尔迭代法发散。 证明 例2中线性方程组的系数矩阵为 A ＝?? ?? ? ?????-122111221 于是 D ＝?? ?? ??????100010001 D －1＝D ??????????=022001000L ~ ????? ?????-=000100220U ~ 雅可比迭代矩阵为 （适合课程《数值方法A 》和《数值方法B 》） 第一章 绪 论 1. 设x >0,x 的相对误差为δ,求ln x 的误差. 2. 设x 的相对误差为2％,求n x 的相对误差. 3. 下列各数都是经过四舍五入得到的近似数,即误差限不超过最后一位的半个单位,试指出它们是几位 有效数字: ***** 123451.1021,0.031,385.6,56.430,7 1.0.x x x x x =====? 4. 利用公式(3.3)求下列各近似值的误差限: * * * * * * * * 12412324(),(),()/,i x x x ii x x x iii x x ++其中* * * * 1234,,,x x x x 均为第3题所给的数. 5. 计算球体积要使相对误差限为1％,问度量半径R 时允许的相对误差限是多少? 6. 设028,Y =按递推公式 11783 100 n n Y Y -=- ( n=1,2,…) 计算到100Y .若取783≈27.982(五位有效数字),试问计算100Y 将有多大误差? 7. 求方程2 5610x x -+=的两个根,使它至少具有四位有效数字(783≈27.982). 8. 当N 充分大时,怎样求 2 11N dx x +∞+?? 9. 正方形的边长大约为100㎝,应怎样测量才能使其面积误差不超过1㎝2 ? 10. 设2 12S gt = 假定g 是准确的,而对t 的测量有±0.1秒的误差,证明当t 增加时S 的绝对误差增加, 而相对误差却减小. 11. 序列{}n y 满足递推关系1101 n n y y -=-(n=1,2,…),若02 1.41y =≈(三位有效数字),计算到10 y 时误差有多大?这个计算过程稳定吗? 12. 计算6 (21)f =-,取 2 1.4≈,利用下列等式计算,哪一个得到的结果最好? 3 6 3 11,(322), ,9970 2. (21) (322) --++ 13. 2 ()ln(1)f x x x =- -,求f (30)的值.若开平方用六位函数表,问求对数时误差有多大?若改用另一等 价公式 2 2 ln(1)ln(1)x x x x - -=-+ + 计算,求对数时误差有多大? 14. 试用消元法解方程组{ 10 10 12121010; 2. x x x x +=+=假定只用三位数计算,问结果是否可靠? 15. 已知三角形面积 1sin , 2 s ab c = 其中c 为弧度, 02c π << ,且测量a ,b ,c 的误差分别为,,.a b c ???证 明面积的误差s ?满足 . s a b c s a b c ????≤ ++ 第二章 插值法 1. 根据( 2.2)定义的范德蒙行列式,令 数值分析典型习题 特别声明：考试时需带计 算器作辅助计算 1．2015x *=是经四舍五入得到的近似值，则其相对误差* r e ≤ -31 104 ?. 2. 01(),(),,()n l x l x l x L 是以01,,,n x x x L 为节点的拉格朗日插值基函数，则 3．设(0)1(1)3(2)4(3)2f =,f =,f =,f =，[0123]f =，，，1 3 - . 4. 利用Simpson 公式求?2 1 2dx x = 7.3 5. 设求积公式1 0()d (),(1)n k k k f x x A f x n ≈≥∑?=是Gauss 型求积公式，则3 n k k k A x == ∑1 .4 6. 数值微分公式(2)(2) ()i i i f x h f x h f x h +≈ --'的截断误差为 2().O h 7. 设1101A ?? = ??? ，则A 的谱半径()A ρ= 1 ，A 的条件数1cond ()A = 4. 8. 用牛顿下山法求解方程3 03 x x -=根的迭代公式是 2 13 3(1),3n n n n x x x x x λ+-=-- 下山条件是 1()().n n f x f x +< 9.对任意初始向量(0)x 及任意向量f ，线性方程组的迭代公式(1)()(0,1,2,)k k k +=+=L x Bx f ，迭代序列()k x 收敛于方程组的精确解x *的充分必要条件是()1.ρ 数值分析心得体会 篇一：学习数值分析的经验 数值分析实验的经验、感受、收获、建议班级：计算131 学号：XX014302 姓名：曾欢欢 数值分析实验主要就是学习MATLAB的使用以及对数值分析类容的应用，可以使学生更加理解和记忆数值分析学得类容，也巩固了MATLAB的学习，有利于以后这个软件我们的使用。在做实验中，我们需要具备较好的编程能力、明白MATLAB软件的使用以及掌握数值分析的思想，才能让我们独立自主的完成该作业，如果是上述能力有限的同学，需要借助MATLAB的书以及网络来完成实验。数值分析实验对于我来说还是有一定难度，所以我课下先复习了MATLAB的使用方法以及编写程序的基本类容，借助互联网和同学老师资源完成了数值分析得实验的内容。在实验书写中，我复习了各种知识，所以我认为这门课程是有必要且是有用处的，特别是需要处理大量实验数据的人员，很有必要深入了解学习它，这样在以后的工作学习里面就减少了很多计算问题也提高了实验结果的精确度。 学习数值分析的经验、感受、收获、建议数值分析的内容包括插值与逼近，数值微分与数值积分，非线性方程与线性方程组的数值解法，矩阵的特征值与特征向量计算，常微分方程数值解等。 首先我们必须明白数值分析的用途。通常所学的其他数学类学科都是由公式定理开始，从研究他们的定义，性质再到证明与应用。但实际上，尤其是工程，物理，化学等其它具体的学科。往往我们拿到 手的只是通过实验得到的数据。如果是验证性试验，需要代回到公式 进行分析，验证。但往往更多面对的是研究性或试探性试验，无具体 公式定理可代。那就必须通过插值，拟合等计算方法进行数据处理以得到一个相对可用的一般公式。还有许多计算公式理论上非常复杂，在工程中不实用，所以必须根据实际情况把它转化成多项式近似表 示。学习数值分析，不应盲目记公式，因为公事通常很长且很乏味。其次，应从公式所面临的问题以及用途出发。比如插值方法，就 是就是把实验所得的数据看成是公式的解，由这些解反推出一个近似公式，可以具有局部一般性。再比如说拟合，在插值的基础上考虑实 验误差，通过拟合能将误差尽可能缩小，之后目的也是得到一个具有 一定条件下的一般性的公式。。建议学习本门课程要结合知识与实际，比如在物理实验里面很多 1、解：将)(x V n 按最后一行展开，即知)(x V n 是n 次多项式。 由于 n i i i n n n n n i n x x x x x x x x x x V ...1...1... ......... ...... 1 )(21110 20 0---= ，.1,...,1,0-=n i 故知0)(=i n x V ，即110,...,,-n x x x 是)(x V n 的根。又)(x V n 的最高 次幂 n x 的系数为 )(...1...1... ...... .........1),...,,(101 1 21 11 2 2221 02001101j n i j i n n n n n n n n n n n x x x x x x x x x x x x x x V -== ∏-≤<≤-----------。 故知).)...()()(,...,,()(1101101------=n n n n x x x x x x x x x V x V 6、解：（1）设 .)(k x x f =当n k ,...,1,0=时，有.0)()1(=+x f n 对 )(x f 构造Lagrange 插值多项式， ),()(0 x l x x L j n j k j n ∑== 其 0)()! 1() ()()()(1)1(=+=-=++x w n f x L x F x R n n n n ξ， ξ介于j x 之间，.,...,1,0n j = 故 ),()(x L x f n =即 .,...,1,0,)(0 n k x x l x k j n j k j ==∑= 特别地，当0=k 时， 10) (=∑=n j x j l 。 （2） 0)()1(1) ()1()()(0000=-=??? ? ??-??? ? ??-=--=-===∑∑∑∑k j j i j i k j k i i j i i k j n j k i i j k n j j x x x x i k x l x x i k x l x x ）利用（。 7、证明：以b a ,为节点进行线性插值，得 )()()(1 b f a b a x a f b a b x x P --+--= 因 0)()(==b f a f ，故0)(1=x P 。而 ))()(("2 1 )()(1b x a x f x P x f --= -ξ，b a <<ξ。 故)("max )(8 122)("max )(max 2 2 x f a b a b x f x f b x a b x a b x a ≤≤≤≤≤≤-=??? ??-≤。 14、解：设 ))...()(()(21n n x x x x x x a x f ---=， k x x g =)(，记)() (1 ∏=-=n j j n x x x w ，则 ),()(x w a x f n n =).()(' j n n j x w a x f = 由差商的性质知 [])! 1()(1,..,,1) (' 1 )(')('1 211 11 -== ==-===∑∑∑ n g a x x x g a x w x a x w a x x f x n n n n n j j n k j n n j j n n k j n j j k j ξ， ξ介于n x x ,...,1之间。 当20-≤≤ n k 时，0)()1(=-ξn g ， 当 1-=n k 时，)!1()(1-=-n g n ξ， 故 ???-=-≤≤=-= --=∑1,,20,0)!1()(1) ('1 11 n k a n k n g a x f x n n n n j j k j ξ 16、解：根据差商与微商的关系，有 [] 1! 7! 7!7)(2,...,2,2)7(7 10===ξf f ， [ ] 0! 80 !8)(2,...,2,2)8(8 1 ===ξf f 。 （ 13)(47+++=x x x x f 是7次多项式， 故 ,!7)()7(=x f 0)()8(=x f ）。 25、解：（1） 右边= [][]dx x S x f x S dx x S x f b a b a ??-+-)(")(")("2)(")("2 = [] d x x S x f x S x S x S x f x f b a ?-++-)("2)(")("2)(")(")("2)(" 222 = [] d x x S x f b a ?-)(")(" 22 = [][]dx x S dx x f b a b a 2 2 )(")("??- =左边。 （2）左边= ? -b a dx x S x f x S ))(")(")((" 第一章典型例题 例3…，精确到10－3的近似值是多少？ 解 精确到10－3＝，即绝对误差限是?＝, 故至少要保留小数点后三位才 可以。ln2? 第二章典型例题 例1 用顺序消去法解线性方程组 解 顺序消元 于是有同解方程组 回代得解 x 3=－1, x 2=1,x 1=1,原线性方程组的解为X ＝(1,1,－1)T 例2 取初始向量X (0)=(0,0,0)T ,用雅可比迭代法求解线性方程组 解 建立迭代格式 ??? ????+--=+--=++-=+++5223122)(2)(1)1(3) (3)(1)1(2 )(3)(2)1(1k k k k k k k k k x x x x x x x x x (k =1,2,3,…) 第1次迭代,k =0 X (0)＝0，得到X (1)＝(1,3,5)T 第2次迭代，k =1 X (2)＝(5,－3,－3)T 第3次迭代，k =2 X (3)＝(1,1,1)T 第4次迭代，k =3 X (4)＝(1,1,1)T 例4 证明例2的线性方程组，雅可比迭代法收敛，而高斯－赛德尔迭 代法发散。 证明 例2中线性方程组的系数矩阵为 A ＝?? ?? ? ?????-122111221 于是 D ＝?? ?? ??????100010001 D －1 ＝D ?? ?? ? ?????=022001000L ~ ?? ?? ? ?????-=000100220U ~ 雅可比迭代矩阵为 B 0＝?? ?? ? ?????--=??????????-??????????-=+--022101220022101220100010001)U ~L ~(D 1 得到矩阵B 0的特征根03,2,1=λ，根据迭代基本定理4，雅可比迭代法收敛。 高斯－赛德尔迭代矩阵为 G ＝－U ~ )L ~D (1-+ ＝－?? ?? ??????----=??????????-??????????---=??????????-??????????-2003202200001002201200110010001002201220110011 解得特征根为?1=0，?2,3=2。由迭代基本定理4知，高斯－赛德尔迭代发散。 例5 填空选择题： 1. 用高斯列主元消去法解线性方程组 作第1次消元后的第2，3个方程分别为 。 数值分析 2?当x=1,—1,2时，f(x)=O, 一3,4,求f(x)的二次插值多项式。解： X 0 =1,x j = — 1,x 2 = 2, f(X。)= 0, f (xj = -3, f (x2)= 4; l o(x)=(x-xi^~x2\=-1(x 1)(x-2) (x o -X/X o _x2) 2 (x -x0)(x -x2) 1 l i(x) 0 2(x-1)(x-2) (x i ~x0)(x i ~x2) 6 (x—x0)(x—x,) 1 l2(x) 0 1(x-1)(x 1) (X2 -X°)(X2 - X i) 3 则二次拉格朗日插值多项式为 2 L 2(X)= ' y k 1 k ( x) kz0 = -3l°(x) 4l2(x) 1 4 =(x_1)(x—2) 4 (x-1)(x 1) 2 3 5 2 3 7 x x - 6 2 3 6?设Xj, j =0,1,||(,n 为互异节点，求证： n (1 )7 x：l j(x) =x k(k =0,1川,n); j=0 n (2 )7 (X j -x)k l j(x)三0 (k =0,1川,n); j ￡ 证明 (1)令f(x)=x k n 若插值节点为X j, j =0,1,|l(, n，则函数f (x)的n次插值多项式为L n(x)八x k l j(x)。 j=0 f (n 十)(?) 插值余项为R n(X)二f(X)-L n(X) n1(X) (n +1)! .f(n1)( ^0 R n(X)=O n 二瓦x k l j(x) =x k(k =0,1川,n); j ：o n ⑵、(X j -x)k l j(x) j卫 n n =為(' C?x j(—x)k_L)l j(x) j =0 i =0 n n i k i i =為C k( -x) (、X j l j(x)) i =0 j=0 又70 _i _n 由上题结论可知 n .原式二''C k(-x)k_L x' i=0 =(X -X)k =0 -得证。 7设f (x) c2 la,b 1且f (a) =f (b)二0,求证: max f(x)兰一(b-a) max a $至小一*丘f (x). 解：令x^a,x^b，以此为插值节点，则线性插值多项式为 L i(x^ f(x o) x x f (xj X o —人x -X o X —X o x-b x-a ==f(a) f(b)- a - b x -a 又T f (a) = f (b)二0 L i(x) = 0 1 插值余项为R(x)二f (x) - L,(x) f (x)(x - X Q)(X - xj 1 f(x) = 2 f (x)(x -X g)(X -xj 有限单元法读书报告 摘要：有限单元法以变分原理和加权余量法为基础，其基本求解思想是把计算域划分为有限个互不重叠的单元，在每个单元内，选择一些合适的节点作为求解函数的插值点，将微分方程中的变量改写成由各变量或其导数的节点值与所选用的插值函数组成的线性表达式，借助于变分原理或加权余量法，将微分方程离散求解。采用不同的权函数和插值函数形式，便构成不同的有限元方法。 关键词：有限单元法；插值函数；网格划分；实例分析 1 有限单元法概述 1.1 有限单元法的简介 有限单元法[1]是应用局部的近似解来建立整个定义域的解的一种方法。先把注意力集中在单个单元上，进行上述所谓的单元分析。基本前提是每一单元要尽可能小，以致其边界值在整个边界上的变化也是小的。这样，边界条件就能取某一在结点间插值的光滑函数来近似，在单元内也容易建立简单的近似解。因此，比起经典的近似法，有限元法具有明显的优越性。比如经典的Ritz法，要求选取一个函数来近似描述整个求解区域中的位移，并同时满足边界条件，这是相当困难的。而有限元法采用分块近似，只需对一个单元选择一个近似位移函数，且不必考虑位移边界条件，只须考虑单元之间位移的连续性即可。对于具有复杂几何形状或材料、荷载有突变的实际结构，不仅处理简单，而且合理适宜。 1.2 有限单元法的基本方法简介 有限单元法，是一种有效解决数学问题的解题方法。在有限元方法中，把计算域离散剖分为有限个互不重叠且相互连接的单元，在每个单元内选择基函数，用单元基函数的线形组合来逼近单元中的真解，整个计算域上总体的基函数可以看为由每个单元基函数组成的，则整个计算域内的解可以看作是由所有单元上的近似解构成。在河道数值模拟中[2]，常见的有限元计算方法是由变分法和加权余量法发展而来的里兹法和伽辽金法、最小二乘法等。根据所采用的权函数和插值函数的不同，有限元方法也分为多种计算格式。从权函数的选择来说，有配置法、矩量法、最小二乘法和伽辽金法，从计算单元网格的形状来划分，有三角形网格、四边形网格和多边形网格，从插值函数的精度来划分，又分为线性插值函数和高次插值函数等。不同的组合同样构成不同的有限元计算格式。对于权函数，伽辽金(Galerkin)法是将权函数取为逼近函数中的基函数；最小二乘法是令权函 数值分析典型例题 例1 对下列各数写出具有5位有效数字的近似值。236.478, 0.00234711, 9.000024, 9.0000343 10?. 解：按照定义，以上各数具有5位有效数字的近似值分别为：236.478， 0.0023471， 9.0000， 9.0000310?。 注意: *x =9.000024的5位有效数字是9.0000而不是9，因为9 是1位有效数字。 例2 指出下列各数具有几位有效数字。2.0004， -0.00200， -9000， 9310?， 23 10-?。 解：按照定义，以上各数的有效数字位数分别为5， 3， 4，1，1 例3 已测得某物体行程* s 的近似值s=800m ，所需时间* s 的近似值为t=35s ，若已知m s s s t t 5.0||,05.0||**≤-≤-，试求平均速度v 的绝对误差和相对误差限。 解：因为t s v /=，所以)()(1)()()(2t e t s s e t t e t v s e s v v e -=??+??≈ 从 而 05.00469.035 800 5.0351|)(||||)(|1|)(|22≤≈+?≤+≤t e t s s e t v e 同样v v e v e r )()(≈)()()()(t e s e t e v t t v s e v s s v r r r -=??+??= 所以00205.035 05 .08005.0|)(||)(||)(|≈+≤+≤t e s e v e r r r 因此绝对误差限和相对误差限分别为0.05和0.00205。 例4试建立积分20,,1,05 =+=n dx x x I n n 的递推关系，并研究它的误差 传递。 解：151 --= n n I n I ……………………………………………..…...(1) 5ln 6ln 0-=I ，计算出0I 后可通过（1）依次递推计算出1I ,…,20I 。 但是计算0I 时有误差0e ，由此计算出的1I ,…,20I 也有误差，由（1）可 知近似值之间的递推关系为 151 --= n n I n I ……………………………………………….…..(2) （1）-（2）可得 01)5(5e e e n n n -=-=-，由0I 计算n I 时误差被放大了n 5倍。所以（1）不稳 定。 （1） 可以改写为 n I I n n 51 511+ -=- ……………………………………… （3） 如果能先求出20I ，则依次可以求出19I ,…,0I ，计算20I 时有误差，这样根据（3）计算19I ,…,0I 就有误差，误差传播为 n n n e e ?? ? ??-=-511 ，误差依次减少。 例5 用二分法求解方程012)(23=+--=x x x x f 在区间[0,1]内的1个实根，要求有3为有效数字。 解：因为0)1()0(数值分析典型习题
数值分析课后题答案
数值计算方法学习心得
数值分析习题集及答案[1].(优选)
数值分析典型例题
数值分析习题集及答案
数值分析典型习题资料
数值分析心得体会
数值分析课后答案
数值分析典型例题
数值分析课后题答案
有限单元法读书报告
数值分析典型例题
相关主题
文本预览