Matlab第2次作业数值.docx
《Matlab第2次作业数值.docx》由会员分享，可在线阅读，更多相关《Matlab第2次作业数值.docx（23页珍藏版）》请在冰豆网上搜索。
Matlab第2次作业数值
1•设a=（1,2,3）,b=（2,4,3）,分别计算a./b,a.\b,a/b,a\b,分析结果的意义。
解：
>>a=[1,2,3];
>>b=[2,4,3];
>>a./b【意义】矩阵内的元素对应的作-运算
b
ans=
0.50000.50001.0000
>>a.\b【意义】矩阵内的元素对应的作-运算
b
ans=
221
>>a/b【意义】矩阵整体xb=a,求x;
ans=
0.6552
>>a\b
【意义】矩阵整体ax=b,求x;
ans=
000
000
0.66671.33331.0000
2•用矩阵除法解下列线性方程组，并判断解的意义
411%9
1326X22;
153x31
解：
>>a=[41-1;32-6;1-53];
>>b=[9;-2;1];
>>a\b
ans=
2.3830
1.4894
2.0213
【意义】最小二乘法求解ax=b的近似解
Xi
X2
解：
>>a=[41;32;1-5];
>>b=[1;1;1];
>>a\bans=
0.3311
-0.1219
【意义】最小二乘法求解ax=b的近似解
2
11
1
%
X2
1
3
1
21
1
x3
2
1
12
1
3
x4
解:
>>
a=[21-11;1
21
-1;1
121
>>
b=[1;2;3];
>>a\b
ans=
1
0
1
0
【意义】最小二乘法求解欠定方程ax=b的解
3.（人口流动趋势）对城乡人口流动作年度调查，发现有一个稳定的朝向城镇流动的趋势，每年农村居民的5赠居城镇而城镇居民的1%i出，现在总人口的20%u于城镇。
加入城乡总人口保持不变，并且人口流动的这种趋势继续下去，那么
（1）一年以后住在城镇人口所占比例是什么？
两年以后呢？
十年以后呢？
（2）很多年以后呢？
（3）如果现在总人口70%u于城镇，很多年以后城镇人口所占比例是什么？
（4）计算转移矩阵的最大特征值级对应的特征向量，与问题
（2），（3）有何关系？
解:
（1）
>>a=[0.990.05;0.010.95];x0=[0.20.8]';
>>x1=a*x0,x2=aA2*x0,x10=aA10*x0
x1=
0.2380
0.7620x2=
0.2737
0.7263x10=
0.4922
0.5078
（2）
【方法一：
循环的方法】
>>x=x0;fori=1:
1000,x=a*x;
end,x
x=
0.8333
0.1667
【方法二：
累乘的方法】
>>x=aA1000*x0
x=
0.8333
0.1667
【注】若求ax=x，即（a-I）x=0的非零解，得结果如下
>>clear
a=[0.990.05:
0.010.95]:
b=*y^（2,2）:
-0*0606
-0.1M1
错误原因在于，非零解是一组基础解系，不能具体确定x的值。
事实上，所求
X=Kx.（与下（4）所问对应）
（3）
>>x0=[0.70.3]';
>>x=x0;fori=1:
1000,x=a*x;end,x
x=
0.8333
0.1667
（4）>>[v,d]=eig（a）
v=
0.9806-0.7071
0.19610.7071
d=
1.00000
00.9400
>>v（:
1）./x
ans=
1.1767
1.1767
最大特征值1,[0.8333,0.1667]是对应特征值向量之一。
4.（经济预测）在某经济年度内，各经济部门的投入产出表如表3.5（单位：
亿
元）。
消耗部门
最后需求
总产值
工业
农业
第三产业
生产部门
工业
6
2
1
16
25
农业
2.25
1
0.2
1.55
5
第三产
业
3
0.2
1.8
15
20
假设某经济年度工业、农业及第三产业的最后需求均为17亿元，预测该经济年
度工业、农业及第三产业的产出（提示：
对于一个特定的经济系统而言，直接消耗矩阵和Leontief矩阵可视作不变）。
解：
>>B=[621;2.2510.2;30.21.8];x=[25520];
>>C=B/diag（x）
C=
0.24000.40000.0500
0.09000.20000.0100
0.12000.04000.0900
>>A=eye（3,3）-C
0.7600-0.4000-0.0500
-0.09000.8000-0.0100
-0.1200-0.04000.9100>>D=[171717]';x=A\Dx=
37.5696
25.7862
24.7690
所以工业37.5696，农业25.7862，第三产业24.7690.
5.
求下列矩阵的行列式、逆、特征值和特征向量。
5
765
41
1
7
1087
1
32
6
2
6
8109
15
3
5
7910
5
6
1
5
6
4
n阶方阵
1
5■
**°
n分别为5,50和500.
11
+<
6
1
5
解:
（1）
>>a=[41-1;32-6;1-53];det（a）,inv（a）,[v,d]=eig（a）
ans=
-94ans=
0.2553-0.02130.0426
0.1596-0.1383-0.2234
0.1809-0.2234-0.0532v=
0.0185-0.9009-0.3066
-0.7693-0.1240-0.7248
-0.6386-0.41580.6170
d=
-3.052700
03.67600
008.3766
（2）
>>a=[5765;71087;68109;57910];det（a）,inv（a）,[v,d]=eig（a）ans=
1ans=
68.0000-41.0000-17.000010.0000
-41.000025.000010.0000-6.0000
-17.000010.00005.0000-3.0000
d=
0.0102
00
0
0
0.84310
0
0
03.8581
0
0
0030.2887
（3）
【方法一：
循环的方法】
>>n=50;
>>a=zeros（n,n）;
>>fori=1:
n,a（i,i）=5;end
>>fori=1:
（n-1）,a（i,i+1）=6;end
>>fori=1:
（n-1）,a（i+1,i）=1;end
>>a
【方法二：
课上所授稀疏矩阵的方法】通用稀疏矩阵：
>>n=5;>>a1=sparse（1:
n,1:
n,5*ones（1,n）,n,n）;
>>a2=sparse（2:
n,1:
n-1,ones（1,n-1）,n,n）;
>>a3=sparse（1:
n-1,2:
n,6*ones（1,n-1）,n,n）;
>>a=a1+a2+a3;
对角带稀疏矩阵：
>>n=5;
>>a=spdiags（[ones（n,1）,5*ones（n,1）,6*ones（n,1）],[-1,0,1],n,n）>>a
a=
（1,1）5
（2,1）1
（1,2）6
（2,2）5
（3,2）1
（2,3）6
（3,3）5
（4,3）1
（3,4）6
（4,4）5
（5,4）1
（4,5）6
（5,5）5
后同：
>>det（a）
ans=
665
>>inv（a）
ans=
0.3173-0.58651.0286-1.62411.9489
-0.09770.4887-0.85711.3534-1.6241
0.0286-0.14290.5429-0.85711.0286
-0.0075
0.0376
-0.1429
0.4887
-0.5865
0.0015
-0.0075
0.0286
-0.0977
0.3173
>>[v,d]=eig（a）
v=
-0.7843
-0.7843
-0.9237
0.9860
-0.9237
0.5546
-0.5546
-0.3771
-0.0000
0.3771
-0.2614
-0.2614
0.0000
-0.1643
-0.0000
0.0924
-0.0924
0.0628
-0.0000
-0.0628
-0.0218
-0.0218
0.0257
0.0274
0.0257
0.7574
0
0
0
0
0
9.2426
0
0
0
0
0
7.4495
0
0
0
0
05.0000
0
0
0
0
02.5505
将n=50和500带入，类似可求得
6.
（1）判断第5题各小题是否可以相似对角化，如果可以，求出对角矩阵和对应的相似变换矩阵。
解：
若n阶矩阵含有n个不等的特征值，那么他一定可以相似对角化，且此n个特征向量组成的矩阵即相似变换矩阵。
故第
（1）题、第
（2）题、第（3）题的n=5时，可以相似对角化。
对角矩阵就是特征值矩阵，相似变换矩阵就是特征向量矩阵。
（2）判断第5题各小题是否为正定矩阵。
解：
【XX】matlab里如何判断矩阵为正定：
[Dp]=chol（A）,如果A正定,返回的p=0,如果不正定，则返回一个正的p,p—1为A中正定子矩阵的阶次，即D为p—1为D的阶次。
>>a=[41-1;32-6;1-53];
>>[Dp]=chol（a）
D=
2.00000.5000
01.3229
p=
3
P工0，故矩阵不正定。
>>a=[5765;71087;68109;57910];[Dp]=chol（a）
D=
2.23613.13052.68332.2361
00.4472-0.89440
001.41422.1213
0000.7071
0
P=0,故矩阵正定>>n=5;
>>a=spdiags（[ones（n,1）,5*ones（n,1）,6*ones（n,1）],[-1,0,1],n,n）;
>>[Dp]=chol（a）
D=
（1.1）2.2361
（1.2）2.6833
p=
2
Pm0,故矩阵不正定。
7•求下列向量组的秩和它的一个最大线性无关组，并将其余向量用该最大无关组线性表示。
1=7,-6,-7,0.
1=4,-3,1,32=2,-13,5,3=1,-1，-1,1,4=3,-2,3,4,
解：
>>a=[4-313;2-135;1-1-1-1;3-234;7-6-70];
>>rank（a）ans=
3
>>rank（a（1:
3,:
））ans=
>>rank（a（[124],:
））ans=
>>b=a（[124],:
）';c=a（[35],:
）';>>b\cans=
0.50005.0000
-0.50001.0000
-0.0000-5.0000
[a1a2a4]构成了一个极大线性无关组。
由计算结果可以知道,
a3=0.5a-0.5a2,a5=5a1+a2-5a4.
8•假定某天的气温变化记录如下表，试用最小二乘方法找出这一天的气温变化规律.
时刻助）
0
1
J
E
1
5
6
7
s
9
10
11
12
温度
w
15°
14°
uc
14°
U°
16°
23°
25°
2S°
时刻
13
14
15
16
17
18
19
20
21
）2
23
24
31fl
31°
29°
27°
2中
24°
■
2（f
17°
16°
解：
>>a=0:
24;
>>b=[15141414141516182022232528313231292725242220181716];
>>a仁polyfit（a,b,1）;
>>a2=polyfit（a,b,2）;
>>a3=polyfit（a,b,3）;
>>a4=polyfit（a,b,4）;
>>a5=polyfit（a,b,5）;
>>y1=polyval（a1,a）;
>>plot（a,b,'r*',a,y1,'b-'）;
>>y2=polyval（a2,a）;
>>y3=polyval（a3,a）;
>>y4=polyval（a4,a）;
>>y5=polyval（a5,a）;
>>plot（a,b,'r*',a,y2,'b-'）;
>>plot（a,b,'r*',a,y3,'b-'）;
>>plot（a,b,'r*',a,y4,'b-'）;
>>plot（a,b,'r*',a,y5,'b-'）;
>>a5
a5=
0.0001-0.00370.04500.0157-0.830415.0539
拟合曲线如下:
9.一矿脉有13个相邻样本点，人为地设定一原点，现测得各样本点对原点的距离x,与样本点处某种金属含量y的一组数据如下，画出散点图观察二者的关
系,
试建立合适的模型，
如二次曲线
、双曲线、对数曲线等。
X
V
rJ.
4
57
8
10
y
I06J2
109J0
109耶
109.50
11D.00
109,93
110.19
X
JI
14
15
15
18
19
V
i|T
110.59
11060
110.90
110.76
111.00
111.20
二次曲线（y=ax2+bx+c）:
>>x=[23457810111415151819];
>>y=[106.42109.2109.58109.5110109.93110.49110.59110.6110.9110.76111111.2];
>>a=polyfit（x,y,2）;
>>y1=polyval（a,x）;
>>plot（x,y1,'b-',x,y,'r*'）
>>c1=y1-y;
>>dot（c1,c1）ans=
4.1821
112
oo
双曲线（y=-+?
?
：
>>a=[1./x',ones（13,1）]
>>ab=a\y'
ab=
-9.0300
111.4405>>y2=ab
（1）./x+ab
（2）；
>>plot（x,y2,'b-',x,y,'r*'）c2=y2-y;dot（c2,c2）
ans=
对数曲线（y=alogx+b）:
>>a=[log（x'）,ones（13,1）]>>ab=a\y'
ab=
1.5663
106.7113
>>y3=ab
（1）.*log（x）+ab
（2）;
>>plot（x,y3,'b-',x,y,'r*'）c3=y3-y;
dot（c3,c3）
ans=
双曲线拟合效果最佳
10.求下列积分的数值解
2
（1）e2xcos3（x）dx
0
（2）（1xy2）dydx,D为x2y22x
D
解：
（1）
fun.m
functiony=fun（t）y=exp（2*t）.*（（cos（t））.A3）;
>>d=pi/1000;
>>t=0:
d:
2*pi;
>>nt=length（t）;
>>y=fun（t）;
>>sc=cumsum（y）*d;
>>scf=sc（nt）
scf=
9.7505e+04
>>z=trapz（y）*d
z=
9.7054e+04
>>quad（'fun',0,2*pi）
ans=
9.7054e+04
矩形法误差较大，梯形和simpson法结果比较精确
（2）
>>fun=@（x,y）1+x+y42;
>>ymax=@（x）（2.*x-x.A2）.A0.5;
>>ymin=@（x）-（（2.*x-x.A2）.A0.5）;
>>quad2d（fun,0,2,ymin,ymax）
11.用积分法计算下列椭园的周长
解:
所求周长，利用椭圆的参数方程，等价于求积分
2
■■-45cos2d
0
>>fun=@（x）（4+5*（cos（x））.A2）.A0.5;
>>quad（fun,0,2*pi）ans=
15.8654
12.试求下列积分，出现什么问题？
分析原因，设法求出正确的解
1
Ix0.2cos（x）dx
1
解：
>>fun=@（x）（x.A0.2）.*cos（x）;
>>quad（fun,-1,1）
ans=
1.2390+0.4026i
【问题】所得答案为虚数，明显错误。
【错因】在matlab里面,对于区间（-1,0）,xA0.2结果是虚数（数学上是实数才对！
！
具体原因不知啊，求老师课上解答）
【正解】因为x0.2cos（x）为奇函数，积分区间（-1,1）关于原点对称，所以
1=0O
k
13考虑积分1sin（x）dx
0
2k试分别用trapz（取步长h=0.1或n）,quad
和quadl求解1（8）和I（32）
（1）I（8）:
。
发现什么问题？
先求I（8）
>>fun=@（x）abs（sin（x））;
>>d=0.1;
>>t=0:
d:
8*pi;
>>y=fun（t）;
>>z=trapz（y）*d
z=
15.9981
z=trapz（y）*d
z=
1.2311e-14
>>quad（fun,0,8*pi）
ans=
16.0000
d=pi;t=0:
d:
8*pi;
>>quadl（fun,0,8*pi）
ans=
y=fun（t）;16.0000
quad、quadl两种方法所求结果均正确。
在间距为0.1时，trapz所得结果也比较精确。
但间距为pi时，由于所取步长过大，结果错误。
这是因为运用NewtonCotes公式时，误差太大了！
而应取远小于积分长度的间距。
*14.用正交变换化下列二次型为标准形
222f（x1,X2,X3）=X1-4X1X2+4X1x3-2x2+8x2X3-2x3
>>a=[1-22;-2-24;24-2];
>>[v,d]=eig（a）
v=
0.33330.9339-0.1293
0.6667-0.3304-0.6681
-0.66670.1365-0.7327
d=
-7.000000
02.00000
002.0000
>>v'*v
ans=
1.00000.00000.0000
0.00001.00000
0.000001.0000
标准型即为f=-7y1A2+2y2A2+2y3A2
*15.下图是连接三个电压已知终端的电路网，求
由基尔霍夫方程：
1311
?
?
?
?
?
2=10
12…4…3''
11135
?
s+?
?
?
?
=-
3…5…15…3
1471
?
?
?
?
+?
?
=04…60…5''
>>a=[13/12-1/4-1/3;1/31/5-13/15;1/4-47/601/5];
>>b=[10;-5/3;0];
>>a\b
ans=
13.3453
6.4401
8.5420
即a=13.3453V,b=6.4401V,c=8.5420V
123
456
*16.就矩阵A=验证下列性质
780
（i）设1，2,
n为n阶方阵A的特征值，则
i1
aH（A的迹）,
i1
n
i
（1）nA
i1
（ii）设f（x）为A的特征多项式，则f（A）=0解：
（i）
>>A=[123;456;780];
>>trace（A）
ans=
6
>>[v,d]=eig（A）
v=
-0.2998-0.7471-0.2763
-0.70750.6582-0.3884
-0.6400-0.09310.8791
d=
12.122900
0-0.38840
00-5.7345
>>trace（d）
ans=
6.0000>>det（A）
ans=
27.0000>>det（d）ans=
27.0000（ii）
>>root1=[d（1,1）d（2,2）d（3,3）];
>>p=poly（root1）
P=
1.000-6.0000-72.0000-27.0000>>polyvalm（p,A）ans=
1.0e-11*
-0.0373-0.0462-0.0380
-0.0853-0.1133-0.0860
-0.0853-0.1137-0.0625
可以发现，结果约等于0.
升级会员 