数学建模实验第3章

预备：编制计算拉格朗日插值的m文件lagr1.m:
function y=lagr1(x0,y0,x)
n=length(x0);m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=k
p=p*(z-x0(j))/(x0(k)-x0(j));
end
end
s=s+p*y0(k);
end
y(i)=s;
end
预备：对数值函数，用辛普森公式计算定积分的程序simp.m: function s=simp(y,h,m)
s=0;
for k=1:m
s=s+4*y(2*k);
end
for k=1:(m-1)
s=s+2*y(2*k+1);
end
s=(s+y(1)+y(2*m+1))*h/3;
3.1 对函数
2
x
y e-
=, 22
x
-≤≤在n 个节点上（n 不要太大，如5至11）用拉格
朗日、分段线性、三次样条三种插值方法，计算m 个插值点的函数值（m 要适中，如50至100）。

通过数值和图形输出，将三种插值结果与精确值进行比较。

适当增加n ，再作比较，由此作初步分析。

解：
函数插值比较程序如下：
clear;
n=5;%在n个节点上进行插值
m=75;%产生m个插值点，计算函数在插值点处的精确值，将来进行对比
x=-2:4/(m-1):2;
y=exp(-x.^2);
z=0*x;
x0=-2:4/(n-1):2;
y0=exp(-x0.^2);
y1=lagr1(x0,y0,x);% y1为拉格朗日插值
y2=interp1(x0,y0,x);% y2为分段线性插值
y3=spline(x0,y0,x);% y3为三次样条插值
[x' y' y1' y2' y3']
plot(x,z,'k',x,y,'r:',x,y1,'g-.',x,y2,'b',x,y3,'y--')
gtext('Lagr.'), gtext('Pieces. linear'), gtext('Spline'),
gtext('y=exp(-x.^2)')
hold off;%比较插值所得结果与函数在插值点处的精确值
s = ' x y y1 y2 y3'
[x' y' y1' y2' y3']
【MA TLAB 计算结果】
n=5时，得到结果如下：
数值比较如下：
x y y1 y2 y3
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－-2.0000 0.0183 0.0183 0.0183 0.0183
-1.9467 0.0226 -0.0328 0.0370 -0.0082
-1.8933 0.0277 -0.0717 0.0556 -0.0276
-1.8400 0.0339 -0.0990 0.0742 -0.0404
-1.7867 0.0411 -0.1158 0.0929 -0.0468
-1.7333 0.0496 -0.1229 0.1115 -0.0472
-1.6800 0.0595 -0.1211 0.1302 -0.0419
-1.6267 0.0709 -0.1112 0.1488 -0.0314
-1.5733 0.0841 -0.0940 0.1675 -0.0159
-1.5200 0.0992 -0.0702 0.1861 0.0041
-1.4667 0.1164 -0.0406 0.2047 0.0284
-1.4133 0.1357 -0.0058 0.2234 0.0566
-1.3600 0.1573 0.0334 0.2420 0.0883 -1.3067 0.1813 0.0764 0.2607 0.1232 -1.2533 0.2079 0.1226 0.2793 0.1609 -1.2000 0.2369 0.1714 0.2980 0.2011 -1.1467 0.2685 0.2222 0.3166 0.2434 -1.0933 0.3026 0.2745 0.3353 0.2875 -1.0400 0.3391 0.3277 0.3539 0.3330 -0.9867 0.3778 0.3813 0.3763 0.3796 -0.9333 0.4185 0.4349 0.4100 0.4269 -0.8800 0.4610 0.4880 0.4437 0.4746 -0.8267 0.5049 0.5401 0.4774 0.5222 -0.7733 0.5499 0.5910 0.5112 0.5695 -0.7200 0.5955 0.6401 0.5449 0.6162 -0.6667 0.6412 0.6872 0.5786 0.6618 -0.6133 0.6865 0.7320 0.6123 0.7060 -0.5600 0.7308 0.7740 0.6460 0.7484 -0.5067 0.7736 0.8131 0.6797 0.7888 -0.4533 0.8142 0.8490 0.7134 0.8266 -0.4000 0.8521 0.8815 0.7472 0.8617 -0.3467 0.8868 0.9104 0.7809 0.8937 -0.2933 0.9176 0.9355 0.8146 0.9221 -0.2400 0.9440 0.9566 0.8483 0.9467 -0.1867 0.9658 0.9736 0.8820 0.9670 -0.1333 0.9824 0.9865 0.9157 0.9828 -0.0800 0.9936 0.9951 0.9494 0.9937 -0.0267 0.9993 0.9995 0.9831 0.9993 0.0267 0.9993 0.9995 0.9831 0.9993 0.0800 0.9936 0.9951 0.9494 0.9937 0.1333 0.9824 0.9865 0.9157 0.9828 0.1867 0.9658 0.9736 0.8820 0.9670 0.2400 0.9440 0.9566 0.8483 0.9467 0.2933 0.9176 0.9355 0.8146 0.9221
0.3467 0.8868 0.9104 0.7809 0.8937 0.4000 0.8521 0.8815 0.7472 0.8617 0.4533 0.8142 0.8490 0.7134 0.8266 0.5067 0.7736 0.8131 0.6797 0.7888 0.5600 0.7308 0.7740 0.6460 0.7484 0.6133 0.6865 0.7320 0.6123 0.7060 0.6667 0.6412 0.6872 0.5786 0.6618 0.7200 0.5955 0.6401 0.5449 0.6162 0.7733 0.5499 0.5910 0.5112 0.5695 0.8267 0.5049 0.5401 0.4774 0.5222 0.8800 0.4610 0.4880 0.4437 0.4746 0.9333 0.4185 0.4349 0.4100 0.4269
0.9867 0.3778 0.3813 0.3763 0.3796
1.0400 0.3391 0.3277 0.3539 0.3330 1.0933 0.3026 0.2745 0.3353 0.2875 1.1467 0.2685 0.2222 0.3166 0.2434 1.2000 0.2369 0.1714 0.2980 0.2011 1.2533 0.2079 0.1226 0.2793 0.1609 1.3067 0.1813 0.0764 0.2607 0.1232 1.3600 0.1573 0.0334 0.2420 0.0883 1.4133 0.1357 -0.0058 0.2234 0.0566 1.4667 0.1164 -0.0406 0.2047 0.0284 1.5200 0.0992 -0.0702 0.1861 0.0041 1.5733 0.0841 -0.0940 0.1675 -0.0159 1.6267 0.0709 -0.1112 0.1488 -0.0314 1.6800 0.0595 -0.1211 0.1302 -0.0419 1.7333 0.0496 -0.1229 0.1115 -0.0472 1.7867 0.0411 -0.1158 0.0929 -0.0468 1.8400 0.0339 -0.0990 0.0742 -0.0404 1.8933 0.0277 -0.0717 0.0556 -0.0276
1.9467 0.0226 -0.0328 0.0370 -0.0082
2.0000 0.0183 0.0183 0.0183 0.0183
函数图像如下：
n=7时，得到结果如下：
数值比较如下：
x y y1 y2 y3
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－-2.0000 0.0183 0.0183 0.0183 0.0183
-1.9467 0.0226 0.0603 0.0304 0.0115
-1.8933 0.0277 0.0888 0.0424 0.0085
-1.8400 0.0339 0.1071 0.0545 0.0091
-1.7867 0.0411 0.1180 0.0665 0.0133
-1.7333 0.0496 0.1238 0.0786 0.0208
-1.6800 0.0595 0.1267 0.0907 0.0316
-1.6267 0.0709 0.1283 0.1027 0.0454
-1.5733 0.0841 0.1302 0.1148 0.0621
-1.5200 0.0992 0.1336 0.1268 0.0815
-1.4667 0.1164 0.1395 0.1389 0.1036
-1.4133 0.1357 0.1485 0.1509 0.1280
-1.3600 0.1573 0.1612 0.1630 0.1548
-1.3067 0.1813 0.1779 0.1879 0.1837
-1.2533 0.2079 0.1988 0.2257 0.2146 -1.2000 0.2369 0.2240 0.2634 0.2473 -1.1467 0.2685 0.2533 0.3012 0.2817 -1.0933 0.3026 0.2866 0.3390 0.3176 -1.0400 0.3391 0.3234 0.3768 0.3549 -0.9867 0.3778 0.3634 0.4145 0.3934 -0.9333 0.4185 0.4062 0.4523 0.4329 -0.8800 0.4610 0.4511 0.4901 0.4734 -0.8267 0.5049 0.4977 0.5279 0.5146 -0.7733 0.5499 0.5453 0.5656 0.5564 -0.7200 0.5955 0.5934 0.6034 0.5986 -0.6667 0.6412 0.6412 0.6412 0.6412 -0.6133 0.6865 0.6882 0.6699 0.6838 -0.5600 0.7308 0.7337 0.6986 0.7260 -0.5067 0.7736 0.7773 0.7273 0.7671 -0.4533 0.8142 0.8182 0.7560 0.8067 -0.4000 0.8521 0.8560 0.7847 0.8442 -0.3467 0.8868 0.8902 0.8134 0.8790 -0.2933 0.9176 0.9204 0.8421 0.9105 -0.2400 0.9440 0.9461 0.8708 0.9382 -0.1867 0.9658 0.9671 0.8995 0.9614 -0.1333 0.9824 0.9831 0.9282 0.9797 -0.0800 0.9936 0.9939 0.9569 0.9925 -0.0267 0.9993 0.9993 0.9856 0.9991 0.0267 0.9993 0.9993 0.9856 0.9991 0.0800 0.9936 0.9939 0.9569 0.9925 0.1333 0.9824 0.9831 0.9282 0.9797 0.1867 0.9658 0.9671 0.8995 0.9614 0.2400 0.9440 0.9461 0.8708 0.9382 0.2933 0.9176 0.9204 0.8421 0.9105 0.3467 0.8868 0.8902 0.8134 0.8790 0.4000 0.8521 0.8560 0.7847 0.8442
0.4533 0.8142 0.8182 0.7560 0.8067
0.5067 0.7736 0.7773 0.7273 0.7671
0.5600 0.7308 0.7337 0.6986 0.7260
0.6133 0.6865 0.6882 0.6699 0.6838
0.6667 0.6412 0.6412 0.6412 0.6412
0.7200 0.5955 0.5934 0.6034 0.5986
0.7733 0.5499 0.5453 0.5656 0.5564
0.8267 0.5049 0.4977 0.5279 0.5146
0.8800 0.4610 0.4511 0.4901 0.4734
0.9333 0.4185 0.4062 0.4523 0.4329
0.9867 0.3778 0.3634 0.4145 0.3934
1.0400 0.3391 0.3234 0.3768 0.3549
1.0933 0.3026 0.2866 0.3390 0.3176
1.1467 0.2685 0.2533 0.3012 0.2817
1.2000 0.2369 0.2240 0.2634 0.2473
1.2533 0.2079 0.1988 0.2257 0.2146
1.3067 0.1813 0.1779 0.1879 0.1837
1.3600 0.1573 0.1612 0.1630 0.1548
1.4133 0.1357 0.1485 0.1509 0.1280
1.4667 0.1164 0.1395 0.1389 0.1036
1.5200 0.0992 0.1336 0.1268 0.0815
1.5733 0.0841 0.1302 0.1148 0.0621
1.6267 0.0709 0.1283 0.1027 0.0454
1.6800 0.0595 0.1267 0.0907 0.0316
1.7333 0.0496 0.1238 0.0786 0.0208
1.7867 0.0411 0.1180 0.0665 0.0133
1.8400 0.0339 0.1071 0.0545 0.0091
1.8933 0.0277 0.0888 0.0424 0.0085
1.9467 0.0226 0.0603 0.0304 0.0115
2.0000 0.0183 0.0183 0.0183 0.0183 函数图像如下：
n=9时，得到结果如下：
数值比较如下：
x y y1 y2 y3
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－-2.0000 0.0183 0.0183 0.0183 0.0183
-1.9467 0.0226 0.0043 0.0276 0.0209
-1.8933 0.0277 0.0012 0.0369 0.0249
-1.8400 0.0339 0.0057 0.0462 0.0304
-1.7867 0.0411 0.0155 0.0555 0.0375
-1.7333 0.0496 0.0287 0.0648 0.0462
-1.6800 0.0595 0.0441 0.0740 0.0566
-1.6267 0.0709 0.0611 0.0833 0.0689
-1.5733 0.0841 0.0791 0.0926 0.0829
-1.5200 0.0992 0.0980 0.1019 0.0989
-1.4667 0.1164 0.1180 0.1229 0.1168
-1.4133 0.1357 0.1391 0.1509 0.1368
-1.3600 0.1573 0.1616 0.1789 0.1589
-1.3067 0.1813 0.1858 0.2069 0.1832
-1.2533 0.2079 0.2119 0.2349 0.2096
-1.2000 0.2369 0.2402 0.2629 0.2384
-1.1467 0.2685 0.2709 0.2909 0.2695 -1.0933 0.3026 0.3040 0.3189 0.3031 -1.0400 0.3391 0.3396 0.3469 0.3392 -0.9867 0.3778 0.3776 0.3788 0.3778 -0.9333 0.4185 0.4178 0.4227 0.4188 -0.8800 0.4610 0.4599 0.4665 0.4618 -0.8267 0.5049 0.5037 0.5103 0.5063 -0.7733 0.5499 0.5487 0.5542 0.5516 -0.7200 0.5955 0.5945 0.5980 0.5973 -0.6667 0.6412 0.6404 0.6418 0.6429 -0.6133 0.6865 0.6860 0.6857 0.6879 -0.5600 0.7308 0.7306 0.7295 0.7316 -0.5067 0.7736 0.7736 0.7733 0.7737 -0.4533 0.8142 0.8144 0.7994 0.8135 -0.4000 0.8521 0.8524 0.8230 0.8507 -0.3467 0.8868 0.8871 0.8466 0.8848 -0.2933 0.9176 0.9178 0.8702 0.9153 -0.2400 0.9440 0.9443 0.8938 0.9418 -0.1867 0.9658 0.9659 0.9174 0.9639 -0.1333 0.9824 0.9825 0.9410 0.9811 -0.0800 0.9936 0.9937 0.9646 0.9930 -0.0267 0.9993 0.9993 0.9882 0.9992 0.0267 0.9993 0.9993 0.9882 0.9992 0.0800 0.9936 0.9937 0.9646 0.9930 0.1333 0.9824 0.9825 0.9410 0.9811 0.1867 0.9658 0.9659 0.9174 0.9639 0.2400 0.9440 0.9443 0.8938 0.9418 0.2933 0.9176 0.9178 0.8702 0.9153 0.3467 0.8868 0.8871 0.8466 0.8848 0.4000 0.8521 0.8524 0.8230 0.8507 0.4533 0.8142 0.8144 0.7994 0.8135 0.5067 0.7736 0.7736 0.7733 0.7737
0.5600 0.7308 0.7306 0.7295 0.7316
0.6133 0.6865 0.6860 0.6857 0.6879
0.6667 0.6412 0.6404 0.6418 0.6429
0.7200 0.5955 0.5945 0.5980 0.5973
0.7733 0.5499 0.5487 0.5542 0.5516
0.8267 0.5049 0.5037 0.5103 0.5063
0.8800 0.4610 0.4599 0.4665 0.4618
0.9333 0.4185 0.4178 0.4227 0.4188
0.9867 0.3778 0.3776 0.3788 0.3778
1.0400 0.3391 0.3396 0.3469 0.3392
1.0933 0.3026 0.3040 0.3189 0.3031
1.1467 0.2685 0.2709 0.2909 0.2695
1.2000 0.2369 0.2402 0.2629 0.2384
1.2533 0.2079 0.2119 0.2349 0.2096
1.3067 0.1813 0.1858 0.2069 0.1832
1.3600 0.1573 0.1616 0.1789 0.1589
1.4133 0.1357 0.1391 0.1509 0.1368
1.4667 0.1164 0.1180 0.1229 0.1168
1.5200 0.0992 0.0980 0.1019 0.0989
1.5733 0.0841 0.0791 0.0926 0.0829
1.6267 0.0709 0.0611 0.0833 0.0689
1.6800 0.0595 0.0441 0.0740 0.0566
1.7333 0.0496 0.0287 0.0648 0.0462
1.7867 0.0411 0.0155 0.0555 0.0375
1.8400 0.0339 0.0057 0.0462 0.0304
1.8933 0.0277 0.0012 0.0369 0.0249
1.9467 0.0226 0.0043 0.0276 0.0209
2.0000 0.0183 0.0183 0.0183 0.0183 函数图像如下：
n=11时，得到结果如下：
数值比较如下：
x y y1 y2 y3
－－－－－－－－－－－－－－－－－－－－－－－－－－－-2.0000 0.0183 0.0183 0.0183 0.0183 -1.9467 0.0226 0.0296 0.0262 0.0223
-1.8933 0.0277 0.0367 0.0340 0.0273
-1.8400 0.0339 0.0421 0.0419 0.0334
-1.7867 0.0411 0.0474 0.0498 0.0406
-1.7333 0.0496 0.0536 0.0576 0.0492
-1.6800 0.0595 0.0615 0.0655 0.0592
-1.6267 0.0709 0.0715 0.0734 0.0708
-1.5733 0.0841 0.0837 0.0879 0.0842
-1.5200 0.0992 0.0983 0.1092 0.0994
-1.4667 0.1164 0.1152 0.1305 0.1166
-1.4133 0.1357 0.1346 0.1518 0.1359
-1.3600 0.1573 0.1565 0.1731 0.1575
-1.3067 0.1813 0.1808 0.1944 0.1814
-1.2533 0.2079 0.2076 0.2156 0.2079
-1.2000 0.2369 0.2369 0.2369 0.2369
-1.1467 0.2685 0.2687 0.2756 0.2687 -1.0933 0.3026 0.3028 0.3144 0.3030 -1.0400 0.3391 0.3393 0.3531 0.3397 -0.9867 0.3778 0.3780 0.3918 0.3784 -0.9333 0.4185 0.4187 0.4305 0.4190 -0.8800 0.4610 0.4611 0.4692 0.4613 -0.8267 0.5049 0.5049 0.5079 0.5050 -0.7733 0.5499 0.5499 0.5489 0.5499 -0.7200 0.5955 0.5954 0.5923 0.5955 -0.6667 0.6412 0.6411 0.6356 0.6414 -0.6133 0.6865 0.6864 0.6789 0.6868 -0.5600 0.7308 0.7307 0.7222 0.7312 -0.5067 0.7736 0.7736 0.7655 0.7740 -0.4533 0.8142 0.8142 0.8088 0.8145 -0.4000 0.8521 0.8521 0.8521 0.8521 -0.3467 0.8868 0.8868 0.8719 0.8864 -0.2933 0.9176 0.9176 0.8916 0.9168 -0.2400 0.9440 0.9440 0.9113 0.9431 -0.1867 0.9658 0.9658 0.9310 0.9649 -0.1333 0.9824 0.9824 0.9507 0.9817 -0.0800 0.9936 0.9936 0.9704 0.9933 -0.0267 0.9993 0.9993 0.9901 0.9992 0.0267 0.9993 0.9993 0.9901 0.9992 0.0800 0.9936 0.9936 0.9704 0.9933 0.1333 0.9824 0.9824 0.9507 0.9817 0.1867 0.9658 0.9658 0.9310 0.9649 0.2400 0.9440 0.9440 0.9113 0.9431 0.2933 0.9176 0.9176 0.8916 0.9168 0.3467 0.8868 0.8868 0.8719 0.8864 0.4000 0.8521 0.8521 0.8521 0.8521 0.4533 0.8142 0.8142 0.8088 0.8145 0.5067 0.7736 0.7736 0.7655 0.7740
0.5600 0.7308 0.7307 0.7222 0.7312
0.6133 0.6865 0.6864 0.6789 0.6868
0.6667 0.6412 0.6411 0.6356 0.6414
0.7200 0.5955 0.5954 0.5923 0.5955
0.7733 0.5499 0.5499 0.5489 0.5499
0.8267 0.5049 0.5049 0.5079 0.5050
0.8800 0.4610 0.4611 0.4692 0.4613
0.9333 0.4185 0.4187 0.4305 0.4190
0.9867 0.3778 0.3780 0.3918 0.3784
1.0400 0.3391 0.3393 0.3531 0.3397
1.0933 0.3026 0.3028 0.3144 0.3030
1.1467 0.2685 0.2687 0.2756 0.2687
1.2000 0.2369 0.2369 0.2369 0.2369
1.2533 0.2079 0.2076 0.2156 0.2079
1.3067 0.1813 0.1808 0.1944 0.1814
1.3600 0.1573 0.1565 0.1731 0.1575
1.4133 0.1357 0.1346 0.1518 0.1359
1.4667 0.1164 0.1152 0.1305 0.1166
1.5200 0.0992 0.0983 0.1092 0.0994
1.5733 0.0841 0.0837 0.0879 0.0842
1.6267 0.0709 0.0715 0.0734 0.0708
1.6800 0.0595 0.0615 0.0655 0.0592
1.7333 0.0496 0.0536 0.0576 0.0492
1.7867 0.0411 0.0474 0.0498 0.0406
1.8400 0.0339 0.0421 0.0419 0.0334
1.8933 0.0277 0.0367 0.0340 0.0273
1.9467 0.0226 0.0296 0.0262 0.0223
2.0000 0.0183 0.0183 0.0183 0.0183 函数图像如下：
【结果分析】
通过图像及数值的对比都可以看到：随着n 的增大，在区间[-2,2]上，插值函数都越来越逼近于原函数，而且当n 为9、11时，几条插值曲线几乎重合。

下面比较一下三种方法的异同和优劣。

拉格朗日多项式，当n 增大时，并不能保证在所有区间都收敛于原函数——拉格朗日多项式在区间[-2,2] 以外，由于拉格朗日多项式的次数增大，在收敛区间外的点上，高阶导数不为零。

光滑性变差，从而产生了极大的振荡。

也就是说只有已知插值点落在收敛区间以内时，才可采用。

所以影响了这种方法的实用价值。

分段线性插值的曲线不如拉格朗日和三次样条的曲线光滑。

但是当n 趋于无穷时，它总能处处收敛于原函数。

因此，分段线性插值一般应用在需要快速计算而又无特殊要求的情况下。

三次样条插值，当n 趋于无穷时，它也总能处处收敛于原函数。

而且它的曲线更光滑，可以应用于机械加工等领域。

3.2 用1/2
y x =，在0,1,4,9,16x =产生5个节点。

用不同节点构造公式计算5x =处的插值。

与精确值比较并进行分析。

解：
用MATLAB 进行插值：
（1）>> x0=[0 1 4 9 16];y0=x0.^(1/2);x=0:0.1:20;
>> y=spline(x0,y0,5)
y =2.1633
（2）>> x0=[0 1 4 9 ];y0=x0.^(1/2);x=0:0.1:20;
>> y=spline(x0,y0,5)
y = 2
（3）>> x0=[0 1 4 ];y0=x0.^(1/2);x=0:0.1:20;
>> y=spline(x0,y0,5)
y = 1.6667
【结论】随着插值点个数的增加，插值越接近精确值 3.3用梯形和辛普森公式计算由下表数据给出的积分 k 1 2 3 4 5 6 7
x k 0.3 0.5 0.7 0.9 1.1 1.3 1.5
y k 0.3895 0.6598 0.9147 1.1611 1.3971 1.6212 1.8325
已知该表数据为函数y=x+sin(x/3)所产生，将计算值与精确值作比较。

用MA TLAB 作积分：
a.梯形法求积分
》s2=trapz(y)*0.2
s2 =
1.3730
b.用辛普森公式计算积分
》s3=simp(y,0.2,3)
s3 =
1.3743
d.精确值：
3.4选择一些函数用梯形（trapz ）和辛普森（quad ）方法计算积分。

改变步长（对梯形），改变精度要求（对辛普森），进行比较、分析。

如下函数供选择参考：
b. 3sin 2,02;x
y e x x =≤≤
【MA TLAB 源程序】
clear;
a = 0;
b = 2;
fun = inline('exp(3*x) .* sin(2*x)');
format long e
% 以a1 为划分区间数做5次trapz 积分
a1 = [10, 20, 40, 80, 160]; ⎰5.13.0)(dx
x y 1.4323)3cos(321)3sin(|5.13.025
.13.0=-=+⎰x x x x
for i = 1:length(a1)
x1 = linspace(a,b,a1(i)+1);
answer1(i) = trapz(x1,feval(fun,x1));
end
disp(' Trapz :')
disp(' Steplength Answer') disp([(b-a)./a1', answer1'])
% 以a2为精度要求做5次quad积分
a2 = [1e-1, 1e-2, 1e-4, 1e-6, 1e-8];
for i = 1:length(a2)
answer2(i) = quad(fun,a,b,a2(i));
end
disp(' Quad :')
disp(' Tolerance Answer') disp([a2', answer2'])
【运行结果/实验结论】：
Trapz :
Steplength Answer
2.000000000000000e-001 -
3.453120332461184e+001
1.000000000000000e-001 -3.093778888833554e+001
5.000000000000000e-002 -3.003568072046070e+001
2.500000000000000e-002 -2.980992235751527e+001
1.250000000000000e-002 -
2.975346836248791e+001
Quad :
Tolerance Answer
1.000000000000000e-001 -
2.974054419264523e+001
1.000000000000000e-002 -
2.973474126570140e+001
1.000000000000000e-004 -
2.973467945935743e+001
1.000000000000000e-006 -
2.973464947962066e+001
1.000000000000000e-008 -
2.973464908256479e+001
【结果分析】：
计算函数数值积分时有梯形法和辛普森法三种方法。

梯形法和辛普森法比较接近，其精度随着取点数的增加而迅速提高；只是梯形公式的收敛速度明显慢于辛普森公式。

而且梯形公式的步长选取是个问题：既不能太大，否则线性插值与函数真值间的误差较大；也不能太小，否则累加的项太多，累积误差较大。

蒙特卡罗法其结果具有随机性。

即用同一个方法计算时，得到的结果也不完全相同，而是围绕准确值上下波动。

所以蒙特卡罗法的精度较低。

我认为计算一般函数的数值积分时用前两种方法较好。

3.8 求
21x e dx x -+∞⎰的数值积分，使误差在410-以内 解： 先做估计2
81111022x x N N N e dx dx x xe Ne -∞
∞-=<≈⎰⎰,4N = 则有> z=quad('exp(-x.^2)./x',1,4,1e-7)
z =
0.1097
3.9求
20sin x dx x π⎰的数值积分，使误差在210-以内 解：
去掉一个小区间[0,]r ，做如下估计
00sin r r x x dx dx r x
x <=⎰⎰。

所以取210r -=。

z=quad('sin(x)./x',1e-2,2*pi,1e-2)
z =
1.4082
3.10
解： 用MATLAB 计算分段线性插值的方法：
x0=[0,3:2:11,12:15];
y1=[0 1.8 2.2 2.7 3.0 3.1 2.9 2.5 2.0 1.6];
y2=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6];
x=0:0.1:15;
y1_in=interp1(x0,Y1,x);
y2_in=interp1(x0,Y2,x);
[x',y1_in',y1_sp',y2_in',y2_sp']
subplot(2,1,1),plot(x,y1_in,x,y2_in,'b'),title('interp')
trapz(x,y1_in)-trapz(x,y2_in)
利用数值积分可算出机翼面积S=10.7500（平方米）
机翼断面曲线如下机翼断面曲线如下：
用MATLAB计算三次样条插值的方法：
x0=[0,3:2:11,12:15];
y1=[0 1.8 2.2 2.7 3.0 3.1 2.9 2.5 2.0 1.6];
y2=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6];
x=0:0.1:15;
y1_sp=spline(x0,Y1,x);
y2_sp=spline(x0,Y2,x);
[x',y1_in',y1_sp',y2_in',y2_sp']
subplot(2,1,2),plot(x,y1_sp,x,y2_sp,'b'),title('spline') trapz(x,y1_sp)-trapz(x,y2_sp)
利用数值积分可算出机翼面积S=11.3444（平方米）
机翼断面曲线如下：
【实验结果与小结】
比较三种插值算法的结果及所得来的图形。

而由图形可见，三次样条插值出来的曲线要比分段线性插值更光滑。

可见，三次样条较优，分段线性次之。

所以有三次样条得来的机翼断面的面积最精确。

【附表】机翼端面轮廓线数据
3.11
解：
【模型假设】
题中所给数据时真实可用的。

假设国土面积在一定时期没有变化，各种因素对此题无影响。

可以进行实地测量。

【问题的分析与模型的建立】
利用蒙特卡洛算法，在一个面积已知的矩形中随机投入足够多的点，得到点落在该国国土上的概率，乘以矩形面积就是要求的结果。

【模型求解】
程序如下：
x=[7.0 10.5 13.0 17.5 34 40.5 44.5 48 56 61 68.5 76.5 80.5 91 96 101 104 106 111.5 118 123.5 136.5 142 146 150 157 158];
y1=[44 45 47 50 50 38 30 30 34 36 34 41 45 46 43 37 33 28 32 65 55 54 52 50 66 66 68];
y2=[44 59 70 72 93 100 110 110 110 117 118 116 118 118 121 124 121 121 121 122 116 83 81 82 86 85 68];
L=max(x)-min(x);
H=max(y2)-min(y1);
for k=1:10
n=10000;
m=0;
for i=1:n
u(i)=unifrnd(min(x),max(x));
v(i)=unifrnd(min(y1),max(y2));
f1=interp1(x,y1,u(i),'cubic');
f2=interp1(x,y2,u(i),'cubic');
if(v(i)>=f1&v(i)<=f2)
m=m+1;
end
end
area=L*H*m/n/18^2*1600;
end
errorrate=abs(41288-area)/41288;
%format long
%fprintf('area=%.3f\n',area);
fprintf('errorrate=%.3f\n',errorrate);
disp(['area=' num2str(area)])
【计算结果】
area=41834.3822
S= area-41288=546.3822km2
因此的出结论该国国土的近似面积与它的精确值相差546.3822km2
3.12
解：
用MA TLAB进行插值计算
x0=[0 2 4 5 6 7 8 9 10.5 11.5 12.5 14 16 17 18 19 20 21 22 23 24];
y0=[2 2 0 2 5 8 25 12 5 10 12 7 9 28 22 10 9 11 8 9 3];
x=0:1/60:24;
y=spline(x0,y0,x);
plot(x,y),
s=trapz(x,y),
s =211.1414
6.测得活塞中气体压力P和体积V的一组数据如下：
a.求数值微分，在MA TLAB下运行如下命令：
》P=[60 80 100 120 140 160 180]；
》V=[80 69.2 60 52 45 38.6 32.5]；
》pp=spline(V,P)；
》ppd=ppder(pp)；
》dy=ppval(ppd,[60 50]);
dy=
-2.3341 -2.7891
即在V=60, 50(in3)处，V改变1(in3)时P的变化量分别是-2.3341和-2.7891 (lbf/in2)。

b.求数值积分
⎰-7040PdV
》V1=linspace(40,70,1000)；
》P1=spline(V,P,V1)；
》S=trapz(V1,P1)；
S=
3.4144e+003
即V从70减至40(in3)时气体作的功是-3414.4 (lbf in)。

10.炮弹射击的目标为一半径100m的圆形区域，当瞄准目标的中心发射时，在众多因素的影响下，弹着点与目标中心有随机偏差。

可以合理地假设弹着点围绕中心呈二维正态分布，x 方向和y方向的均方差分别为80m和50m，且偏差在x和y之间的相关系数为0.4。

求炮弹命中圆形区域的概率。

用孟特卡罗法计算如下积分
其中Ω是半径为100的圆形区域。

编制如下程序mtc( )
function s=mtc(n)
m=0;ss=0;
for i=1:n
x=(rand-0.5)*200;
y=(rand-0.5)*200;
if x^2+y^2<=100^2
z=1/(2*pi*80*50*sqrt(1-0.4*0.4))*exp(-0.5/(1-0.4*0.4)*(x*x/6400-2 *0.4*x/80*y/50+y*y/2500));
ss=ss+z;
m=m+1;
end
end
s=ss*40000/n;
取随机点数n=1000000，计算mtc (1000000)得结果为0.6986
即炮弹命中的概率是69.9%.。