1 1 1 1 1 n=6时: -1.21E-14 6.97E-14 1.18E-13 -6.17E-13 4.59E-13 4.44E-16 0 2.22E-16 0 -1.11E-16 14 13 13 12 13 x 1 1 1 1 1 1 n=7时: Δx -9.28E-13 2.67E-11 -1.82E-10 4.75E-10 -5.26E-10 2.08E-10 rn 0 0 0 0 0 0 有效数字 12 11 10 10 9 10 x 1 1 1 1.00000001 0.99999997 1.00000002 0.99999999 n=8时: Δx -9.26E-12 3.71E-10 -3.59E-09 1.40E-08 -2.57E-08 2.22E-08 -7.30E-09 rn 0 0 2.22E-16 0 0 1.11E-16 0 有效数字 11 10 9 8 8 8 8 x 1 1 0.99999998 1.00000011 0.9999997 1.00000042 0.9999997 1.00000008 n=9时: Δx -2.82E-11 1.53E-09 -2.01E-08 1.09E-07 -2.95E-07 4.19E-07 -2.99E-07 8.44E-08 rn 4.44E-16 0 2.22E-16 0 -2.22E-16 -1.11E-16 1.11E-16 0 有效数字 11 9 8 7 7 7 7 7
x 1 1.00000002 0.99999968 1.00000228 0.99999164 1.00001706 0.99998042 1.00001182 0.99999708 n=10时: Δx -2.75E-10 1.90E-08 -3.20E-07 2.28E-06 -8.36E-06 1.71E-05 -1.96E-05 1.18E-05 -2.92E-06 rn 4.44E-16 0 0 4.44E-16 0 0 -1.11E-16 0 0 有效数字 10 8 7 6 5 5 5 5 6 x 1 1.00000008 0.99999835 1.00001491 0.99992911 1.00019443 0.99968164 1.00030715 0.99983898 1.00003537 n=11时: Δx -9.05E-10 7.76E-08 -1.65E-06 1.49E-05 -7.09E-05 0.000194426 -0.000318365 0.000307149 -0.000161019 3.54E-05 rn 4.44E-16 0 0 4.44E-16 0 0 -1.11E-16 0 0 0 有效数字 9 7 6 5 4 4 4 4 4 5 x 0.999999995 1.000000539 0.999986041 1.000155875 0.999072325 1.003258718 0.992910161 1.009658807 0.991981908 1.003707654 0.999267974 Δx -5.16E-09 5.39E-07 -1.40E-05 0.00015588 -0.0009277 0.00325872 -0.0070898 0.00965881 -0.0080181 0.00370765 -0.000732 rn 0 -4.44E-16 0 0 0 0 0 -2.22E-16 0 0 -2.22E-16 有效数字 8 6 5 4 3 3 2 2 3 3 3
n=12时: x 0.999999965 1.000004409 0.999861744 1.001879711 0.986241983 1.060375974 0.831939173 1.303973363 0.643862006 1.260684018 0.891666924 1.019510753 n=13时: Δx -3.49E-08 4.41E-06 -0.0001383 0.00187971 -0.013758 0.06037597 -0.1680608 0.30397336 -0.356138 0.26068402 -0.1083331 0.01951075 rn 8.88E-16 4.44E-16 0 0 0 2.22E-16 2.22E-16 0 1.11E-16 2.22E-16 1.11E-16 0 有效数字 8 6 4 3 2 2 1 1 1 1 1 2 x 1.000000069 0.999989224 1.000413538 0.993147495 1.061246208 0.669261244 2.149074131 -1.65417119 5.11854808 -3.243049939 3.783004896 -0.051792817 1.174329117 Δx 6.89E-08 -1.08E-05 0.00041354 -0.0068525 0.06124621 -0.3307388 1.14907413 -2.6541712 4.11854808 -4.2430499 2.7830049 -1.0517928 0.17432912 rn 8.88E-16 -4.44E-16 -2.22E-16 2.22E-16 -2.22E-16 2.22E-16 0 0 -1.11E-16 1.11E-16 0 1.11E-16 1.11E-16 有效数字 7 5 4 3 2 1 0 0 0 0 0 0 1 五、分析讨论: 实验的数学原理很容易理解,也容易上手。把运算的结果带入原方程组,可以发现符合的还是比较好。这说明列主元消去法计算这类方程的有效性。当A可逆时,能够将计算进行到底,列主元法就能确保算法的稳定,而且计算量不大。直接三角消去过程,实质上
是将A分解为两个三角矩阵的乘积A=LU,并求解Ly=b的过程。回带过程就是求解上三角方程组Ux=y。所以在实际的运算中,矩阵L和U可以直接计算出,而不需要任何中间步骤,从而在计算过程中将高斯消去法的步骤进行了进一步的简略,大大提高了运算速度,这就是三角分解法。 通过以上的计算比较,2.题方程组具有严重的病态性。当系数矩阵有微小的变化时,wucha =-401.8918 159.5435 124.6330,所得的解与原方程组的解有很大的相对误差。1题方程组中当系数矩阵A和b有微小变化时,wucha =0 0 0 0,所得的解与方程组的解没有相对误差。所以1题方程组是良性的。用MATLAB内部函数inv通过求逆矩阵,然后通过x=inv(A)*b也可以求出方程组的解,但是没有列主元高斯消去法具有良好的稳定性。det函数求方程组系数矩阵的行列式时所得结果和高斯消去法和三角法所得结果相同,具有方便快捷的优点。题四可以看出,条件数越大,有效位数越少,当n=13时,出现有效位数为0的情况。
附:高斯列主消去法源程序代码
function [x,det,index]=Gauss(A,b) % ?ó??D?·?3ì×éμ?áD?÷?aGauss??襷¨£????D % A --- 方程组矩阵 % b --- 方程组右端 % x --- 方程组的解 % det ---方程组行列式
% index --- index=0表示求解失败,index=1表示求解成功|
[n,m]=size(A); nb=length(b); if n~=m
error('The rows and columns of matrix A must be equal!'); return; end if m~=nb
error('The columns of A must be equal the dimension of b!'); return; end
index=1; det=1; x=zeros(n,1); for k=1:n-1 % 选主元 a_max=0; for i=k:n
if abs(A(i,k))>a_max a_max=abs(A(i,k)); r=i; end end
if a_max<1e-20 index=0; return; end % 交换两行 if r>k for j=k:n
z=A(k,j); A(k,j)=A(r,j); A(r,j)=z; end
z=b(k); b(k)=b(r); b(r)=z; det=-det; end
% 消元过程 for i=k+1:n
m=A(i,k)/A(k,k); for j=k+1:n
A(i,j)=A(i,j)-m*A(k,j); end
b(i)=b(i)-m*b(k); end
det=det*A(k,k); end
det=det*A(n,n); % 回代过程
if abs(A(n,n))<1e-20 index=0; return; end
for k=n:-1:1 for j=k+1:n
b(k)=b(k)-A(k,j)*x(j); end
x(k)=b(k)/A(k,k); end