数学实验——线性代数方程组的数值解.docx

数学实验报告实验5 线性代数方程组的数值解法实验5 线性代数方程组的数值解法分1 黄浩 2011011743一、 实验目的1. 学会用MATLAB软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析； 2. 通过实例学习用线性代数方程组解决简化的实际问题。二、 实验内容1. 数学实验第二版（问题1）问题叙述：通过求解线性方程组A1x=b1,A2x=b2，理解条件数的意义和方程组性态对解的影响，其中A1是n阶范德蒙矩阵，即A1=1x0x02x0n-11x1x12x1n-11xn-1xn-12xn-1n-1 , xk=1+0.1k , k=0,1,n-1A2是n阶希尔伯特矩阵，b1，b2分别是A1,A2的行和。（1） 编程构造A1（A2可直接用命令产生）和b1，b2；你能预先知道方程组A1x=b1和A2x=b2的解吗？令n=5，用左除命令求解（用预先知道的解可验证程序）。（2） 令n=5,7,9，计算A1和A2的条件数。为观察他们是否病态，做以下试验：b1，b2不变，A1和A2的元素A1(n,n)，A2(n,n)分别加扰动后求解；A1和A2不变，b1，b2的分量b1(n),b2(n)分别加扰动后求解。分析A与b的微小扰动对解的影响。取10-10,10-8,10-6。（3） 经扰动得到的解记做x，计算误差x-xx，与用条件数估计的误差相比较。模型转换及实验过程：（1）小题.由b1,b2为A1,A2的行和，可知方程组A1x=b1和A2x=b2的精确解均为n行全1的列向量。在n=5的情况下，用matlab编程（程序见四.1），构造A1,A2和b1，b2，使用高斯消去法得到的解x1,x2及其相对误差e1,e2（使用excel计算而得）为：xx1e1x2e211.0000000000000800.0000000000000801.0000000000000200.00000000000002010.999999999999463-0.0000000000005370.999999999999558-0.00000000000044211.0000000000009700.0000000000009701.0000000000020100.00000000000201010.999999999999288-0.0000000000007120.999999999996863-0.00000000000313711.0000000000001800.0000000000001801.0000000000015600.000000000001560由上表可见，当n=5时，所得的解都接近真值，误差在10-12的量级左右。（2）小题分别取n=5,7,9,11,13,15，计算A1和A2的条件数c1和c2，（程序见四.2），结果如下：n579111315c13.574E+058.739E+072.274E+106.518E+122.059E+158.291E+17c24.766E+054.754E+084.932E+115.227E+146.263E+173.675E+17由上表可见，二者的条件数都比较大，可能是病态的。为证实A1和A2是否为病态，先保持b不变，对A(n,n)做扰动，得到该情况下的高斯消元解，（程序见四.3），结果如下：(为使结果清晰简洁，在此仅列出n=5,9,13的情况，n=7,11,15略去)=10-10时：nxx1xx250.9999999284994020.9999999370003031.0000002510851691.0000012599939300.9999996704145830.9999943300273341.0000001916677021.0000088199574730.9999999583331440.999995590021263nxx1xx290.9999975548050940.9999829278035681.0000149244303261.0012291970545280.9999603223369260.9784890657946091.0000600154317051.1577467696177110.9999435079279020.4084498576816771.0000338876411232.2304238808731590.999987348370823-0.4354941236522751.0000026878626961.8788737437097660.9999997511934060.780281609027228nxx1xx2131.0007784033869290.9999998767583270.9937743987274311.0000261185116641.0226872410459610.9987794790026500.9501840256039321.0233044427317511.0734106035353810.7678410435503330.9235081059886532.3669384594613431.057790678168565-4.0983444092691000.96809962810653313.4983690949393491.012769678990370-19.4039025417116410.99638464083603822.9640269119150101.000687237182600-13.9700648285337990.9999212432336596.8528369055604181.0000041151939490.000188901053803=10-8时：nxx1xx250.9999928500302660.9999937027771001.0000251082268421.0001259444579840.9999670418066980.9994332499390881.0000191665851081.0008816112058550.9999958333510870.999559194397072nxx1xx290.9997569218861830.9999249646886321.0014841958519581.0054025423430740.9960527298593050.9054555097425781.0059726416246391.6933262585415040.994376087847693-1.5999734614530311.0033746909787466.4079447885959670.998739698550261-5.3092689110300421.0002678327406924.8628176967655060.9999752006605240.034295576497007nxx1xx2130.9980046593334370.9999998767282631.0164113342717761.0000261243251180.9384398540328910.9987792205553051.1392634494355081.0233092230919600.7883771891288210.7677944673667701.2275767009039882.3672082173739580.822400779656029-4.0993377962873331.10134460129783313.5007796714360190.958029341700176-19.4078052094758501.01230286523370722.9681988664444180.997576920776780-13.9728915924515431.0002879131499596.8539364954676760.9999843910790940.000001889299579=10-6时：nxx1xx250.9992852977929170.9993966095201831.0025097875871361.0120678095963350.9967055393606440.9456948568165061.0019158683868581.0844746671743260.9995835068724450.957762666412834nxx1xx290.9757501420061600.9999223275065541.1480658167117691.0055924194987990.6062136691657740.9021326589030921.5958409910913531.7176938349075140.438948771633566-1.6913518835480071.3366650248743926.5980119249310320.874270161152963-5.5310139213940691.0267194461452354.9985799573668030.9975259772187880.000355009332673nxx1xx2130.7071952366081760.9999998767279653.3904929638964631.000026124383261-7.8992555052967080.99877921797036220.9772809601275481.023309270903186-29.1195663172520300.76779400153262433.1329169063691252.367210915370232-23.874599377714521-4.09934773169287015.07902114649785113.500803780926434-4.782905869975354-19.4078442421848052.68118057389295422.9682405924371480.671624908087923-13.9729198644592591.0386948092493446.8539474930660780.9979195655092280.000000018893025由上表可见：a) 对于希尔伯特阵，随着阶数的增加，微小扰动对解带来的影响越来越大，到了n=9时，已经有了6倍误差的解，到了n=13时，甚至出现了22倍误差的解元素；而随着的增加，解的偏差似乎也有增加的趋势，但仅凭上述表格无法具体判断（在下一小题中具体叙述）。b) 对于范德蒙矩阵，当n10时，微小扰动对范德蒙矩阵的影响急剧增加，更出现了33倍误差的解！此时希尔伯特阵与范德蒙阵的病态程度是类似的。（因此，我认为第一版上的参考答案不够全面）当固定A改变b时，结果与上表类似，微小扰动的影响随阶数的增加而增大，随的增加，解的偏差也有增加的趋势。方程组的解表格在此略去，程序也类似，其具体影响及分析见下一小问。（3）小题在上一小问中，我们解出了当A、b分别作微小扰动时解的变化情况，但是仅凭直观的感觉以及表格横向的比较，并不能很好的找出其内在规律。因此，可以通过计算误差x-xx，来与条件数估计的误差比较，来进一步得出结论。对于这个误差计算公式，我选择的是2-范数，即norm（x）函数进行计算。保持b不变，给A(n,n)以微小的扰动，（程序见四.4），结果如下（其中E1、E2分别为两方程解的误差x-xx）：n=10-10=10-8=10-6条件数E1E2E1E2E1E2Cond(A1)Cond(A2)52.07E-075.12E-062.07E-055.12E-042.07E-034.90E-023.57E+054.77E+0573.19E-062.10E-033.19E-041.89E-013.19E-021.74E+008.74E+074.75E+0893.30E-057.28E-013.30E-033.20E+003.30E-013.31E+002.27E+104.93E+11112.19E-045.94E+002.56E-026.02E+002.55E+006.02E+006.52E+125.23E+14132.51E-024.71E+001.25E-014.71E+001.59E+014.71E+002.06E+156.26E+17152.31E+007.84E+001.79E+007.51E+006.19E+015.96E+008.29E+173.68E+17由上表可得出结论：第一， 随着阶数的增加，误差值有增加的趋势第二， 的变化对误差大小也有影响，但不是绝对的，对于n=5,7,9,11,13，越大，扰动引起的偏差就越大，对于n=15，找不到与E的相关性第三， 实际算得的误差远远小于条件数指出的误差上限第四， 当n10时，微小扰动对范德蒙矩阵的影响急剧增加，到了n=15时，病态程度已与希尔伯特阵类似。再保持A1和A2不变，b1，b2的分量b1(n),b2(n)分别加扰动后，得到该情况下的高斯消元解，计算误差，（程序见四.5）结果如下：n=10-10=10-8=10-6条件数E1E2E1E2E1E2Cond(A1)Cond(A2)52.07E-075.12E-062.07E-055.12E-042.07E-035.12E-023.57E+054.77E+0573.19E-062.10E-033.19E-042.10E-013.19E-022.10E+018.74E+074.75E+0893.33E-059.33E-013.30E-039.33E+013.30E-019.33E+032.27E+104.93E+11111.62E-044.31E+022.57E-024.31E+042.56E+004.31E+066.52E+125.23E+14133.55E-022.49E+041.31E-012.49E+061.59E+012.49E+082.06E+156.26E+17151.06E+013.73E+049.61E-013.73E+067.01E+013.73E+088.29E+173.68E+17由上表可得出结论：第一， 随着阶数的增加，误差值有增加的趋势第二， 的变化对误差大小也有影响，但不是绝对的，对于希尔伯特阵，越大，扰动引起的偏差就越大，对于范德蒙阵，当n=5,7,9,11,13时，越大，扰动引起的偏差就越大，但对于n=15，找不到与E的相关性第三， 实际算得的误差远远小于条件数指出的误差上限第四， 当n10时，微小扰动对范德蒙矩阵的影响急剧增加，但与固定b改变A的情况不同的是，最终范德蒙阵的病态程度也没有赶上希尔伯特阵第五， 将上述两表纵向比较发现，阶数相同偏差相同时，改变b比改变A所造成的结果偏差更大1.2. 数学实验第二版（问题3）问题叙述：已知方程组Ax=b，其中AR(2020)，定义为： 试通过迭代法求解此方程组，认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。实验要求：(1)选取不同的初始向量x(0)和不同的方程组右端向量b，给出迭代误差要求，用雅可比迭代法和高斯-赛德尔迭代法计算，观测得到的迭代向量序列是否均收敛？若收敛，记录迭代次数，分析计算结果并得出你的结论；(2)取定右端向量b和初始向量x(0)，将A的主对角线元素成倍增长若干次，非主对角线元素不变，每次用雅可比迭代法计算，要求迭代误差满足x(k+1)-x(k)10-5，比较收敛速度，分析现象并得出你的结论。模型转换与实验过程：（1）小题 首先，使用matlab设计雅克比迭代法和高斯-赛德尔迭代法的函数，设为JacDD1和GSDD1，（函数见四.6、四.7），通过输入系数矩阵A、常数项向量b、初始解向量x0、步数上限m以及规定误差上限tol，输出结果xJ（或xG）与迭代步数k。首先，保持常数项向量b=(1,1,1,1)T不变，分别设定初始向量x(0)为(0,0,0,0)T、(5,5,5,5)T、(100,100,100,100)T，误差上限tol=10-5，分别使用雅克比迭代法和高斯-赛德尔迭代法计算，（程序见四.8），所得结果如下：（其中xJ，xG分别为使用雅克比迭代法和高斯-赛德尔迭代法得到的解向量）x(0)(0,0,0,0)T(5,5,5,5)T(100,100,100,100)TxJxGxJxGxJxG解 向 量0.4816350110.4816335080.4816371820.4816409430.4816366890.4816398830.5734131580.5734116110.5734168540.5734210570.5734160170.5734194990.6327989970.6327976730.632804280.6328083480.6328030890.6328062750.6520986790.6520977760.6521053730.6521088050.6521038710.6521062630.6609418470.6609414470.6609497970.6609523360.6609480230.6609494090.6643052030.6643053610.664314220.6643156670.664312220.6643125230.6657142460.6657150040.6657241290.6657243440.6657219490.6657211990.6662729030.6662742810.666283440.6662823630.6662811270.6662794310.6664962310.6664982130.6665072060.6665048770.6665048050.6665023220.6665735250.6665760480.6665847190.666581280.6665822740.666579190.6665735250.6665764810.6665847190.6665803960.6665822740.6665787840.6664962310.6664994830.6665072060.6665022740.6665048050.6665010990.6662729030.6662762940.666283440.6662781920.6662811270.6662773790.6657142460.6657176190.6657241290.665718850.6657219490.6657183160.6643052030.6643084110.664314220.664309170.664312220.6643088380.6609418470.6609447570.6609497970.66094520.6609480230.6609450050.6520986790.6521011780.6521053730.6521014210.6521038710.6521013130.6327989970.6328009950.632804280.6328011180.6328030890.6328010640.5734131580.5734145680.5734168540.5734146220.5734160170.5734145980.4816350110.4816358440.4816371820.4816358630.4816366890.481635855迭代次数171119122415由上表可见，随着初值x(0)与真值越来越远离，满足误差上限的解所需的迭代次数越来越多，而且雅克比迭代法比高斯-赛德尔迭代法所需的迭代次数要多，但是二者的结果相近，都在误差上限之内。为了证明两种求解方法的收敛性，我们以x(0)=(5,5,5,5)T时的高斯-赛德尔迭代法为例，将不同迭代次数得到的解列表如下（程序见四.9，其他的情况类似，在此略去）：k=0k=1k=2k=3k=4k=5k=651.5833333330.8097993830.5867626890.5159863770.4929626410.48538750751.8472222220.9783629120.7052563650.6168153850.5877768120.57818087752.0231481481.084426440.7813237890.6819209230.6490964610.63821656752.0744598771.1206981240.8071642940.7035546140.6692025920.65779093752.0976723251.1372138040.8191053420.7135320420.6784443550.66677151852.1058170441.1436399750.8237874540.7173972820.6819885220.67019847552.1091088681.1463317770.8257794470.7190360460.683482730.67163852852.1103362321.1474109130.8265897910.7196988670.6840822740.67221365952.1108151111.1478512890.8269263590.7199736250.6843293130.67244801252.1109972041.1480285630.8270643090.7200860580.684429750.67252176652.111067461.1481001480.8271211030.7201323670.68445080.67243591352.1110943441.1481289050.8271443740.7201513680.6842499120.67203420152.1111046791.148140450.827153910.719918050.6832395860.67110569252.1111086421.1481450730.8271578080.7178314910.6805376480.66952080652.1111101641.1481469220.8242658810.7102614320.675688690.66701021152.1111107471.148147660.805458660.6951998230.6686267110.66268609152.1111109721.1134257320.7566388370.6735497460.6567218090.65311633152.1111110580.9629628860.6907658520.6445060140.6352482410.63333137351.6944444240.7185570680.600404170.5787075670.5745170960.57365217950.7916666590.5333397520.4909645160.4834934290.4820235360.481719644k=7k=8k=9k=10k=11k=120.4828815270.4820500490.4817736970.4816817550.4816511440.4816409430.5749989270.5739417540.573590110.5734730570.5734340650.5734210570.6346027310.6334008590.6330008350.632867610.6328231840.6328083480.6539954920.6527322070.6523115110.6521712730.6521244310.6521088050.6628857660.6615916860.6611604250.6610163940.6609682580.6609523360.6662714780.6649629750.6645261480.6643799480.6643313890.6643156670.6676920040.6663752910.6659342820.6657871530.6657392950.6657243440.6682561060.6669308290.666486750.666341410.666295940.6662823630.6684736390.6671370010.6666960450.6665576910.6665165820.6665048770.6685091180.6671722840.6667505750.666625920.6665908350.666581280.6683888310.6671001890.6667211670.6666159640.666587780.6665803960.6680849140.6669276810.6666119130.6665289910.6665076820.6665022740.6675508740.6666014440.6663582090.666297130.6662819480.6662781920.6666577650.6659471590.6657736310.6657315260.6657213240.665718850.6649462430.6644591990.6643443560.6643171770.6643107140.664309170.6613428570.6610368410.6609663530.6609499580.660946110.66094520.6523289070.6521531220.6521132330.6521040570.6521019220.6521014210.6329187230.6328276480.6328071450.6328024570.6328013720.6328011180.5734671370.5734264380.5734173010.5734152160.5734147340.5734146220.4816544160.4816400440.4816368120.4816360740.4816359030.481635863迭代次数增加由上面的图表可知，随着迭代次数k的增加，所得迭代向量序列之间的偏差越来越小，是逐渐收敛的。然后，保持初始向量x(0)=(2,2,2,2)T不变，分别设定常数项向量b为 (0.01,0.01,0.01,0.01)T、(1,1,1,1)T、(5,5,5,5)T，误差上限tol=10-5，分别使用雅克比迭代法和高斯-赛德尔迭代法计算，所得结果如下（其中xJ，xG分别为使用雅克比迭代法和高斯-赛德尔迭代法得到的解向量）：b(0.01,0.01,0.01,0.01)T(1,1,1,1)T(5,5,5,5)TxJxGxJxGxJxG解 向 量0.0048163670.0048164390.48163670.4816405622.4081724482.4081652010.005734160.0057342440.5734160320.5734205922.8670613492.8670550730.006328030.0063281150.6328031040.6328078783.1639886493.1639849690.0065210380.0065211130.6521038780.6521084383.2604853863.2604853060.006609480.0066095390.6609480160.6609521713.304699743.3047035850.0066431220.0066431640.6643121940.6643157723.321515273.3215231710.0066572190.0066572440.6657219010.665724723.3285594793.3285715020.0066628110.0066628190.6662810590.6662829453.3313520093.3313681320.0066650470.0066650420.6665047220.6665055613.3324681483.3324881660.0066658220.0066658060.6665821830.6665819633.3328543633.3328778060.0066658220.0066657970.6665821830.6665810013.3328543633.3328804820.0066650470.0066650170.6665047220.6665027633.3324681483.3324959760.0066628110.0066627780.6662810590.6662785563.3313520093.3313804540.0066572190.0066571860.6657219010.6657191013.3285594793.3285874180.0066431220.006643090.6643121940.6643093323.321515273.3215416260.006609480.0066094510.6609480160.6609452983.304699743.3047235330.0065210380.0065210140.6521038780.6521014763.2604853863.2605057490.006328030.0063280110.6328031040.6328011463.1639886493.1640049050.005734160.0057341460.5734160320.5734146342.8670613492.8670728070.0048163670.0048163590.48163670.4816358682.4081724482.40817921k251518111510 由上表可见，随着b与（0,0,0,，0）T越来越远离，迭代次数变得越来越少，而且雅克比迭代法比高斯-赛德尔迭代法所需的迭代次数要多，但是二者的结果相近，都在误差上限之内。如果b=（0,0,0,，0）T，那么会出现什么情况呢？根据Ax=b，我们知道，在这种情况下，方程的一个解为x=（0,0,0,，0）T经过试验，迭代向量数列仍是收敛的，不过迭代次数很大，雅克比迭代法需要1042次，而高斯-赛德尔迭代法需要545次。 综上所述，我们可以得到以下几点结论：a) 在迭代收敛的情况下，当初值x(0)远离真值时，满足误差上限的解所需的迭代次数增多。这是因为迭代收敛时，迭代值逐次逼近真值，因而初值与真值越接近，就可以在较少迭代次数的情况下满足误差上限。b) 在迭代收敛的情况下，当b的各个元素远离0时，满足误差上限的解所需的迭代次数减少。这从直观上也很容易理解，b的元素越远离0，解元素的绝对值就可能越大，对于确定的相对误差，其允许的绝对误差就大，迭代次数减少。c) 在其他情况相同的条件下，雅克比迭代法比高斯-赛德尔迭代法所需的迭代次数要多。因矩阵的谱半径不超过（任一种）范数，可求得：BjBj2= 0.4893062, BgBg2= 0.3292286 显然，二者都是收敛的，而且后者的收敛速度快。（2）小题我们固定常数项向量b=初始向量x(0)=(1,1,1,1)T，将A的主对角线元素分别调至3,9,15,21,27,33，即为1,3,5,7,9,11倍，（程序见四.10），得迭代解及迭代次数如下：xii3915212733 解 向 量0.4816367890.1219637270.0703231230.0494345270.0381191020.0310207030.5734161790.1289623430.0727104420.0506266140.0388319270.0314944760.6328033050.1327656850.073965230.0512443130.039198370.0317368180.6521041190.1331791860.0740476840.0512734260.0392118290.0317441330.6609482860.1333090390.0740714350.0512814990.0392154770.0317460970.664312480.1333280570.0740736230.0512820480.0392156720.0317461930.6657221950.1333327940.0740740980.0512821620.0392157110.0317462150.6662813540.133333610.0740741530.0512821740.0392157140.0317462180.6665050150.1333337970.0740741640.0512821760.0392157140.0317462190.6665824740.1333338330.0740741660.0512821770.0392157140.031746220.6665824740.1333338330.0740741660.0512821770.0392157140.031746220.6665050150.1333337970.0740741640.0512821760.0392157140.0317462190.6662813540.133333610.0740741530.0512821740.0392157140.0317462180.6657221950.1333327940.0740740980.0512821620.0392157110.0317462150.664312480.1333280570.0740736230.0512820480.0392156720.0317461930.6609482860.1333090390.0740714350.0512814990.0392154770.0317460970.6521041190.1331791860.0740476840.0512734260.0392118290.0317441330.6328033050.1327656850.073965230.0512443130.039198370.0317368180.5734161790.1289623430.0727104420.0506266140.0388319270.0314944760.4816367890.1219637270.0703231230.0494345270.0381191020.031020703迭代次数1687665 由上表可得出结论：当主对角线的元素增加时，迭代次数逐渐减小，即收敛幅度加快，但是减小的幅度也变得平缓。3. 数学实验第二版（问题5）问题叙述：设国民经济由农业、制造业和服务业三个部门组成，已知某年它们之间的投入产出关系、外部需求、初始投入等如下表所示： 单位：亿元投入 产出农业制造业服务业外部需求总产出农造业301045115200服务业2060070150初始投入3511075总投入1002001501) 如果今年对农业、制造业和服务业的外部需求分别为50，150，100亿元，问这三个部门的总产出分别应为多少？2) 如果三个部门的外部需求分别增加1个单位，问它们的总产出应分别增加多少？模型转换及实验过程：（1）小题.设一定时期内第i个部门的总产出为xi，其中对第j个部门的投入为xij 设消耗系数矩阵A=aij33，其中aij=xijxi则由题中表格，可作出三个部门间投入系数aij的表格：投入 产出农业制造业服务业农业0.150.10.2制造业0.30.050.3服务业0.20.30设产出向量x=x1,x2,x3T，需求向量d=d1,d2,d3T，根据每个部门的总产出等于总投入，可以得到下列关系：x=Ax+d即：I-Ax=d 基于这个模型，我们可以借助matlab求解线性方程组，（程序见四.11）,结果为：d=139.2