用梯形法则算得精确解
x=4:0.1:5;y= 2*pi*exp(x)./x;trapz(x,y) ans =
129.2178
实验七
7.4观察根的变化是否很显著 一、问题描述
7.4一个10次项式
P(x)?x10?a1x9?a2x8?L?a9x?a10
的系数为
[1 a1 a2 …a9 a10]=[1 –55 1320 –18150 157773 –902055 3416930 -8409500 12753576 -10628640 6328800]
试用多项式的求根指令roots求出该10次方程的10个根,然后修改9次项的系数-55为-56,得新的10次方程,求解新的方程,观察根的变化是否很显著。 二、源程序及运行结果
p=[1 -56 1320 -18150 157773 -902055 3416930 -8409500 12753576 -10628640 6328800]; >> x=roots(p) x =
p=[1 -55 1320 -18150 157773 -902055 21.7335 3416930 -8409500 12753576 -10628640 7.3521 + 7.8018i 6328800]; 7.3521 - 7.8018i >> x=roots(p) 3.6433 + 3.5988i x = 3.6433 - 3.5988i
4.9248 + 1.8070i 4.9248 - 1.8070i 1.1808 + 2.0164i 1.1808 - 2.0164i 0.0643
10.6051 + 1.0127i 10.6051 - 1.0127i 8.5850 + 2.7898i 8.5850 - 2.7898i 5.5000 + 3.5058i 5.5000 - 3.5058i 2.4150 + 2.7898i 2.4150 - 2.7898i 0.3949 + 1.0127i 0.3949 - 1.0127i
7.5 将三村短路问题推广为四村短路问题,即已知四个点A,B,C,D的具体位置A(0,
0,),B(0,3),C(8,1),D(10,5),求两个点H1,H2的具体位置,使AH1+BH1+H1H2+H2C+H2D为最短。
解:设所要求的点H1?(x1,y1)构造函数
H2?(x2,y2)
S?x12?y12?x12?(y1?3)2?(x2?x1)2?(y2?y1)2?(x2?8)2?(y2?1)2?(x2?10)2?(y2?5)2
要求使 S 取得最小值的 H1?(x1,y1)H2?(x2,y2)
??S??x?1??S??y?1???S??x2???S??y即必须满足方程组 ?2?2x1?2x1?2(x2?x1)?0?2y1?2(y1?3)?2y(2?y1)?0?2(x2?x1)?2(x2?8)?2x2(?10?)?2(y2?y1)?2(y2?1)?2y(5?)2?
00?3x1?x2?0?3y?y?3?0?12??3x2?x1?18?0?3y2?y1?6?0化简即得 ?
915??2721??x1,y1????,??x2,y2???,?解得
?48??48?
这个解也即使AH1+BH1+H1H2+H2C+H2D为最短的解。
实验八
8.1分别用直接法、雅可比迭代法、赛德尔迭代法求解线性方程组AX=b。
?41?14?A??O?1???(1)
1O41?2???1?????b??O?O????1??1??4??2??10?1 ?10?10,
??5?21??25?21?1?25?21A???OOOOO??1?25?2?1?25(2)
??1?2源代码及结果: (1)
直接法:
e=ones(10,1);
a=spdiags([e,4*e,e],[-1,0,1],10,10); b=[2;1;1;1;1;1;1;1;1;2]; x=transpose(a\\b) x =
0.4793 0.0829 0.1891 0.1608 0.1891 0.0829 0.4793
雅克比:
e=ones(10,1);
a=spdiags([e,4*e,e],[-1,0,1],10,10); b=[2;1;1;1;1;1;1;1;1;2]; x=[0;0;0;0;0;0;0;0;0;0]; y=[0;0;0;0;0;0;0;0;0;0]; for k=1:8 eorr=0; for i=1:10
s=x(i);x(i)=0;
y(i)=(b(i)-a(i,:)*x)/a(i,i); eorr=max(abs(s-x(i)),eorr); x=y; end
x',pause,eorr
??????1???0????1?b??0???2??M???5???0?20?20,?0??20?1 0.1678 0.1678 0.1608 end ans =
0.5000 0.1250 0.2188 0.1953 0.2012 0.1997 0.2001 0.2000 0.2000 0.4500
赛德尔:
e=ones(10,1);
a=spdiags([e,4*e,e],[-1,0,1],10,10); b=[2;1;1;1;1;1;1;1;1;2]; x=[0;0;0;0;0;0;0;0;0;0]; for k=1:8 eorr=0; for i=1:10
s=x(i);x(i)=0;
x(i)=(b(i)-a(i,:)*x)/a(i,i); eorr=max(abs(s-x(i)),eorr); end
x',pause,eorr end ans =
0.4688 0.0781 0.1816 0.1543 0.1615 0.1596 0.0975 0.4756 (2) 直接法:
e=ones(20,1);
a=spdiags([e,-2*e,5*e,-2*e,e],[-2,-1,0,1,2],20,20); b=[1;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; x=transpose(a\\b) x =
Columns 1 through 12
0.2421 0.0944 -0.0218 -0.0314 -0.0069 0.0002 -0.0008 -0.0004 0.0001 0.0001
Columns 13 through 20
0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 雅克比:
e=ones(20,1);
a=spdiags([e,-2*e,5*e,-2*e,e],[-2,-1,0,1,2],20,20);
0.1601 0.0049 0.0000 0.1600 0.0036 -0.0000 b=[1;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; x=[0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; y=[0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; for k=1:8 eorr=0; for i=1:20
s=x(i);x(i)=0;
y(i)=(b(i)-a(i,:)*x)/a(i,i); eorr=max(abs(s-x(i)),eorr); x=y; end
x',pause,eorr end
ans =
Columns 1 through 12
0.2000 0.0800 -0.0080 -0.0192 -0.0061 0.0004 -0.0002 -0.0002 -0.0000 0.0000
Columns 13 through 20
0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 赛德尔:
e=ones(20,1);
a=spdiags([e,-2*e,5*e,-2*e,e],[-2,-1,0,1,2],20,20); b=[1;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; x=[0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0]; for k=1:8 eorr=0; for i=1:20
s=x(i);x(i)=0;
x(i)=(b(i)-a(i,:)*x)/a(i,i); eorr=max(abs(s-x(i)),eorr); end
x',pause,eorr end
ans =
Columns 1 through 12
0.0014 0.0000 0.0018 0.0000