实验6 线性代数方程组的数值解法[实验目的]1. 1. 学会用MATLAB 软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2. 2. 通过实例学习用线性代数方程组解决简化的实际问题。

[实验内容]5-5 输电网络：一种大型输电网络可简化为图5.5（见书）所示电路，其中R 1，R 2，…，R n 表示负载电阻，r 1，r 2，…，r n 表示线路内阻，I 1，I 2，…，I n 表示负载上的电流。

设电源电压为V 。

(1)列出求各负载电阻R 1，R 2，…，R n 的方程；(2)设I 1=I 2=…=I n =I ，r 1=r 2=…=r n =r ，在r=1，I=0.5，V=18，n=10的情况下求R 1，R 2，…，R n 及总电阻R 0。

[问题分析、模型建立及求解](1) 设电源负极为电势为0，电阻R 1上对应节点电压为V 1，对于任意节点，根据KCL 定律列出方程：111++----=k k k k k k k k r V V r V V R V而k kk R V I =，可得： 111111)(++++--++-=k k k k k k k k k k k k R r IR r I r I R r I Ik=2，3，……，n-1；k=1时2221211R r IR r I I +-=，为与上式形式一致，化为22212111111)(R r IR r I r I r V I +--=-k=m （12-≤≤n m ）时 111111)(++++--+--+=m m m n m m m m m m m m R r IR r I r I R r I Ik=n 时n n n n n n n R r IR r I I -=--11设以上方程组的矩阵形式为：b AR =则 []Tn R R R R 21=Tn I I I r V I b ⎥⎦⎤⎢⎣⎡-= 3211⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡---------=------n n nn n n nn n n n n r I r I r I r I r I r I r I r I r I r I r I r I r I r I r I r I r I A 1111114443333233322221222111000000(2)代入参数：5.021=====I I I I n ，121=====r r r r n ，V=18，n=10，⇒⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡-----=5.05.0005.015.005.015.005.015.0005.01A ⇒[]Tb 5.05.05.05.17 -=在命令窗口输入MA TLAB 程序如下：clear all ;n=10; %由题目要求设定A11=sparse(1:n-1,1:n-1,-1,n,n); %定义A 的对角元素，除（n,n) A12=sparse(n,n,-0.5,n,n); %定义(n,n)A1=A11+A12; %对角元素A2=sparse(1:n-1,2:n,0.5,n,n); %输入A 的上次对角元素 A3=sparse(2:n,1:n-1,0.5,n,n); %输入A 的下次对角元素 A=A1+A2+A3;b1=0.5*ones(n,1); %b 的除第一项元素 b2=sparse(1,1,18,n,1); %b 的第一项元素 b=b1-b2; R=A\b输出结果如下：R =26.000017.0000 9.0000 2.0000 -4.0000 -9.0000 -13.0000 -16.0000 -18.0000 -19.0000所以各阻值为(R 1，R 2，…，R 10)=(26，17，9，2，-4，-9，-13，-16，-18，-19)总电阻R 0（即输入等效电阻）为00I V R =，又nII I nk k ==∑=10得到 )(6.35.010180Ω=⨯=R5-6 有5个反应器连接如图5.6（见书），各个Q 表示外部输入、输出及反应器间的流量（m 3/min ），各个c 表示外部输入及反应器内某物质的浓度(mg/m 3)。

假定反应器内的浓度是均匀的，利用质量守恒准则建立模型，求出各反应器内的浓度c 1~c 5，并讨论反应器j 外部输入改变１个单位(mg/min)所引起的反应器i 浓度的变化。

[问题分析、模型建立及求解]当反应容器中反应物浓度稳定时，输入物质质量与输出物质质量平衡，即输入物质质量等于输出物质质量。

由此分析，分别对1~5容器列质量平衡方程，以输入物质为“+”，输出物质为“－”⎪⎪⎪⎩⎪⎪⎪⎨⎧=+-+=+-+-=+-=++--=++-0)(0)(0)()(555542251155544443342240303334312232252423112010133111215c Q Q c Q c Q c Q c Q c Q c Q c Q c Q Q c Q c Q Q Q c Q c Q c Q c Q Q ，代入数据整理为⎪⎪⎪⎩⎪⎪⎪⎨⎧=++=+-+-=+=--=+-0430211816090335065215432322131c c c c c c c c c c c c c由矩阵表示：b Ac =[]Tc c c 51 =，[]Tb 00160050--=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡---=40013211810009100003300106A编写MATLAB 程序如下：A=[-6,0,1,0,0 3,-3,0,0,0 0,1,-9,0,0 0,1,8,-11,2 3,1,0,0,-4]; b=[-50,0,-160,0,0]';c=A\b %计算并输出c pause;dQc01=[1,0,0,0,0]';dQc02=[0,0,1,0,0]';c11=A\(b+dQc01); %01的输入减少1个单位 dc11=c-c2;c12=A\(b-dQc01); %01的输入增加1个单位 dc12=c-c1;c21=A\(b+dQc02); %03的输入减少1个单位 dc21=c-c4;c22=A\(b-dQc02) ; %03的输入增加1个单位 dc22=c-c3;[c,c11,dc11,c12,dc12,c21,dc21,c22,dc22]%外部输入改变引起的 %反应器浓度变化列表比较MATLAB 给出结果如下：c =11.5094 11.5094 19.0566 16.9983 11.5094所以改变前各反应器浓度为(c 1，c 2，c 3，c 4，c 5)=( 11.5094，11.5094，19.0566，16.9983，11.5094) 改变容器01或容器03的外部输入，得到各容器平衡时浓度及其增量如下表：[改进做法]上面这种做法显得有些麻烦，由b Ac =可得b A c 1-=，表明容器平衡浓度c 对外部输入b 是线性的，所以当b 增加1个单位（记作b ∆）时，c 的增量为b A c ∆=∆-1，若外部输入增加一个单位，Tb )0,0,0,0,1(=∆时，c ∆为1-A 的第一列， Tb )0,0,1,0,0(=∆时，c ∆为1-A 的第三列。

外部输入减少1个单位时，取其负值即可。

故改用矩阵求逆的方法来计算：A=[-6,0,1,0,0 3,-3,0,0,0 0,1,-9,0,0 0,1,8,-11,2 3,1,0,0,-4]; b=[-50,0,-160,0,0]';dx=inv(A) %求A 的逆矩阵输出结果为：dx =-0.1698 -0.0063 -0.0189 0 0 -0.1698 -0.3396 -0.0189 0 0 -0.0189 -0.0377 -0.1132 0 0 -0.0600 -0.0746 -0.0875 -0.0909 -0.0455 -0.1698 -0.0896 -0.0189 0 -0.2500红色部分显示为所求，此结果与前面的方法计算的一致，但工作量明显少了许多。

5-8 种群的繁殖与稳定收获：种群的数量因繁殖而增加，因自然死亡而减少，对于人工饲养的种群（比如家畜）而言，为了保证稳定的收获，各个年龄的种群数量应维持不变，种群因雌性个体的繁殖而改变，为方便起见以下种群数量均指其中的雌性。

种群年龄记作k=1，2，…，n ，当年年龄k 的种群数量记作x k ，繁殖率记作b k （每个雌性个体1年繁殖的数量），自然存活率记作s k （s k =1-d k ，d k 为1年的死亡率），收获量记作h k ，则来年年龄k 的种群数量k x ~应为)1,,2,1(~,~11,1-=-==+=∑n k h x s x x b x k k k k nk k k ，要求各个年龄的种群数量每年维持不变就是要使),,1(~n k x x k k ==。

(1)若已知b k ，s k ，给定收获量h k ，建立求个年龄的稳定种群数量x k 的模型（用矩阵、向量表示）。

(2)设n=5，b 1=b 2=b 5=0，b 3=5，b 4=3，s 1=s 4=0.4，s 2=s 3=0.6，如果要求h 1~h 5为500，400，200，200，100，求 (3)要使h 1~h 5均为500， 如何达到？[问题分析、模型建立及求解](1)要使各年龄种群数量每年维持不变即),,1(~n k x x k k ==，依题意得⎪⎩⎪⎨⎧-=-==+=∑)1,,2,1(111n k h x s x x b x k k k k nk kk用矩阵形式表示原方程组为：h Ax =[]Tn x x x x 21=， []Tn h h h h 0121-=nn n n b b b b s s s A ⨯-⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡----=3211211100010001(2)代入题中数据⇒⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡-----=0350114.016.016.014.0A ， []T h 0100200400500=编写MATLAB 程序如下：format bank ;A=[0.4,-1,0,0,0 0,0.6,-1,0,00,0,0.6,-1,0 0,0,0,0.4,-1 -1,0,5,3,0];h=[500,400,200,100,0]'; x=A\hMATLAB 输出结果如下：x =8481.01 2892.41 1335.44 601.27 140.51第5年龄段：x 5=140.5>100=h 5 ，说明收获量h 5可以达到100。

(3)要使h 1~h 5均为500，则h 变为：[]h Ax h ==0500500500500其他参数不变，程序变为：format bank ;A=[0.4,-1,0,0,0 0,0.6,-1,0,0 0,0,0.6,-1,0 0,0,0,0.4,-1 -1,0,5,3,0];h=[500,500,500,500,0]'; x=A\hMATLAB 输出结果如下：x =10981.01 3892.41 1835.44 601.27 -259.49从结果看出，x 5为-259.49，但种群数量不可能为负数，在本题所给条件下，无法使h 1~h 5均为500。