病态矩阵例子
Hilbert矩阵病态性分析
Hilbert矩阵病态性分析200820303013 王刚一、问题叙述Hilbert矩阵是著名的病态矩阵，用迭代法解矩阵方程时，如果系数矩阵是Hilbert矩阵，则求解结果误差较大。
本文就研究Hilbert矩阵的病态性，和Hilbert 矩阵的阶数与迭代求解误差大小的关系。
二、问题分析MATLAB中有专门的Hilbert矩阵及其准确逆矩阵的生成函数，从而可以通过MATLAB求解出系数矩阵是Hilbert矩阵的矩阵方程比较准确的解，再根据条件数估计求解误差的大小。
三、实验程序及注释m=input('input m:='); %输入矩阵的阶数N=[m];for k=1:length(N)n=N(k); %矩阵的阶H=hilb(n); %产生n阶Hilbert矩阵disp(H) %输出n阶Hilbert矩阵Hi=invhilb(n); %产生完全准确的n阶逆Hilbert矩阵b=ones(n,1); %生成n阶全1向量x_approx=H\b; %利用左除H求近似解x_exact=Hi*b; %利用准确逆Hilbert矩阵求准确解ndb=norm(H*x_approx-b);nb=norm(b);ndx=norm(x_approx - x_exact);nx=norm(x_approx);er_actual(k)=ndx/nx; %实际相对误差K=cond(H); %计算Hilbert矩阵的条件数er_approx(k)=K*eps; %最大可能的近似相对误差er_max(k)=K*ndb/nb; %最大可能的相对误差enddisp('Hilbert矩阵阶数'),disp(N)format short edisp('实际误差 er_actual'),disp(er_actual),disp('')disp('近似的最大可能差 er_approx'),disp(er_approx),disp('') disp('最大可能误差 er_max'),disp(er_max),disp('')四、实验数据结果及分析程序运行后，输入矩阵阶数4时的输出为：input m:=41.0000 0.5000 0.3333 0.25000.5000 0.3333 0.2500 0.20000.3333 0.2500 0.2000 0.16670.2500 0.2000 0.1667 0.1429Hilbert矩阵阶数4实际误差er_actual1.0284e-013近似的最大可能误差er_approx3.4447e-012最大可能误差er_max4.7732e-011当输入不同的矩阵阶数时，误差大小如表1所示。
病态矩阵
病态矩阵
病态矩阵的最根本的特征：病态矩阵的条件数（||A||*||A^(-1)||定义为A 的条件数，||A||为A的范数，MATLAB中范数的求法c = cond(A,p)，p为几范数，默认情况下为二范数）很大（可以用机器精度的一半来衡量,比如双精度浮点数的精度大约是10^(-16)，那么条件数大于10^8就算比较坏了.），即最大特征值与最小特征值的比值很大。
病态矩阵对扰动特别敏感。
病态矩阵矩阵不能求逆。
例如，有一个变量是其他三个变量之和，这个变量也存在于模型中，这个矩阵就是病态矩阵。
解决方法：
1.利用正则化方法可以解决这一问题
最经典的正则化方法就是Tikhonov正则化，另外一种方法就是Total Variation
2.奇异值分解
将病态矩阵分解为A=SVD
3.模拟退火算法
/view/3a75ad82bceb19e8b8f6baad.html#opennewwindow
4./link?url=ripw0P7tMkcEl2x5zmkOJ-_iDFdfOLQ0Ox9ao cKY3sAJMFXEBQa__Yq6OVb3R41__nvgiWEQhA0bxXWzZzpQiFH7_IlNWlvqyrkKjLqvQ qC。
病态矩阵的例子
>> norm(x'-ones(8,1))
ans =
9.3866e-007
>> norm(A*x'-b)
ans =
7.0278e-008
1.0000
1.0000
就是说，对于系数矩阵为 8 阶的 Hilbert 矩阵，用 Matlab 中 的”inv”语句，大致能达到 6-7 位有效数字，这时条件数为 10 的 10 次方量级。
Results may be inaccurate. RCOND = 2.632091e-017. >> rank(A) ans =
11
>> A=hilb(100); >> rank(A) ans =
18 >> A=hilb(1000); >> rank(A) ans =
24
>> rank(hilb(2000)) ans =
0.0000
0.0000
0.0000
0.0262 0.2981 1.6959
>> cond(A)
ans =
1.5258e+010
0.0001
0.0015
>> b=A*ones(8,1);
>> x=(inv(A)*b)'
x=
1.0000
1.0000
1.0000
1.0000 1.0000 1.0000
n=10; >> A=hilb(10); >> cond(A) ans =
1.6025e+013
> b=A*ones(10,1); >> solution=inv(A)*b; >> solution' ans =
测量数据处理中病态矩阵和部分有偏估计方法
文章编号:0494 0911(2010)09 0009 03中图分类号:P207 文献标识码:B测量数据处理中病态矩阵和部分有偏估计方法吴 杰,苗恒严(宿迁学院,江苏宿迁223800)Ill posed M atrix and Parti al B iased E stimati on M ethod i nM easure m ent Data Processi ngW U J i e ,M IAO H engyan摘要:从扰动角度、空间几何角度、超椭球的角度对病态方程进行研究,加深对病态方程的认识,并给出部分有偏估计方法。
数值试验结果表明,该方法效果较好,可以作为相关计算的参照。
关键词:测量平差;病态性;条件数;空间几何分析;复共线性;部分有偏估计方法收稿日期:2010 02 02基金项目:宿迁学院科研基金资助项目(2010KY14)作者简介:吴 杰(1975 ),男,江苏宿迁人,硕士,讲师,现主要从事工程测量及数据处理方面的研究。
一、引 言在测量数据处理中,最小二乘法是被广泛采用的估计方法,但是当模型存在病态时,观测数据的微小扰动常常会造成估值与其真值存在很大偏差,估值表现不稳定。
病态问题危害非常严重,且存在于变形观测数据处理、大地测量、GPS 数据处理等许多领域。
目前处理病态方程效果较好的方法有有偏估计中的岭估计和截断奇异值分解方法,但都存在不同程度的缺点。
本文拟从不同角度去认识病态方程的本质,同时探讨用部分有偏估计方法对病态方程进行解算,并给出符合相应精度要求的解。
二、病态方程的本质1.最小二乘原理[1]Guass M ar kov 模型为L n 1=A m n X n 1+ n 1E ( )=0,D ( )= 20P-1(1)用误差方程表示V =A X ^-L(2)按最小二乘原理V TPV =m i n ,求得法方程N X ^=W(3)其中,N =A T PA,W =A TPL 。
最小二乘解为X ^=(A T PA )-1A T PL =N -1W (4)2.从扰动的角度认识病态性在式(3)中,当N 、W 有微小扰动时,参数估值X ^也有扰动(N + N )(X ^+ X ^)=W + W (5)化简得X ^=-N -1 N (X ^+ X ^)+N -1W (6)将上式两端取范数,并应用向量范数的三角不等式及矩阵和向量范数相容条件有! X ^!∀!N -1!! N !(!X ^!+! X ^!)+!N-1!! W !(7)整理得1-K! N !!N !! X ^!!X ^!∀K ! N !!N !+! W !!W !(8)其中,K =!N !!N -1!,当N 的扰动 N 很小时,使!N-1!! N !<1,有! X ^!!X ^!∀K 1-K! N !!N !! N !!N !+! W !!W !(9)条件数K 度量了法方程的病态程度,一般认为条件数小于100时为良态,条件数在100~1000时为中等程度的病态;条件数超过1000时存在严重病态。
矩阵病态线性代数方程组的求解
实验一病态线性代数方程组的求解1.估计Hilbert矩阵2-条件数与阶数的关系运行tiaojianshu.m输入m=10 可以得到如下表的结果2.选择不同维数，分别用Guass消去（LU分解），Jacobi迭代，GS 迭代，SOR迭代求解，比较结果。
说明：Hx=b，H矩阵可以由matlab直接给出，为了设定参考解，我们先设x为分量全1的向量，求出b,然后将H和b作为已知量，求x，与设定的参考解对比。
对于Jacobi迭代，GS迭代，SOR迭代，取迭代初值x0为0向量，迭代精度eps=1.0e-6，迭代次数<100000, SOR迭代中w=1.2和0.8分别计算。
取不同的n值，得到如下结果:对于Guass法，可以看出来，随着n的增大，求解结果误差变大，这是因为随着n增大，系数矩阵的条件数变大，微小的扰动就容易造成很大的误差。
最后得不到精确解。
对于Jacobi迭代，计算结果为Inf，说明是发散的。
对于GS迭代和SOR迭代，结果是收敛的，但是可以看出迭代次数比较多，并且对于不同维数GS和SOR收敛速度不一样，有时候GS快，有时SOR快。
对SOR取不同的w迭代速度也不一样，存在一个最优的松弛因子w。
并且可以知道，迭代次数多少跟初值x0也有关系。
3.讨论病态问题求解的算法。
通过上面的实验分析，可以看出，求解病态矩阵的时候要小心，否则可能得不到所要求的精确度。
可以采用高精度运算，用双倍多倍字长，使得由于误差放大而损失若干有效数字位之后，还能保留一些有效位。
另外可以通过对原方程作某些预处理，降低系数矩阵的条件数，因为cond(aA)=cond(A),所以不能通过将每一个方程乘上相同的常数来达到这个目标，可考虑将矩阵的每一行和每一列分别乘上不同的常数，亦即找到可逆的对角阵D1和D2将方程组化为D1AD2y=D1b，x=D2y这称为矩阵的平衡问题，但是这样计算量比原问题本身要多。
或者通过变分原理将求解线性方程组的问题转化为等价的求解无约束函数最优化问题的极小值等等，可以参考[1]郑洲顺，黄光辉，杨晓辉求解病态线性方程组的混合算法实验一所编程序如下：1.求条件数tiaojianshu.mm=input ('input m:=') ;N=[1:m];for i=1:mn=N(i);H=hilb(n);k=cond(H);disp('矩阵的阶数')disp(n)disp('矩阵')disp(H)disp('矩阵的条件数')disp(k)end2.Guass法①Guass.mfunction guass(n)n=input('请输入系数矩阵的维数n=');H=hilb(n);x_exact=ones(n,1);b=H*x_exact;x=Doolittle(H,b)a=input('是否继续,继续请按1,结束请按0:');if a==1guass(n);else end②Doolittle.mfunction x=Doolittle(A,b)% LU分解求Ax=b 的解N=size(A);n=N(1);L=eye(n,n);%L的对角元素为1U=zeros(n,n); %U为零矩阵U(1,1:n)=A(1,1:n)%U第一行L(1:n,1)=A(1:n,1)/U(1,1)%L第一列for k=2:nfor i=k:nU(k,i)=A(k,i)-L(k,1:(k-1))*U(1:(k-1),i)endfor j=(k+1):nL(j,k)=(A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy=solveDownTriangle(L,b);%调用下三角系数矩阵求解函数x=solveUpTriangle(U,y);%调用上三角系数矩阵求解函数③solveDownTriangle.mfunction x=solveDownTriangle(A,b)%求下三角系数矩阵的线性方程组Ax=bN=size(A);n=N(1);for i=1:nif (i>1)s=A(i,1:(i-1))*x(1:(i-1),1);elses=0;endx(i,1)=(b(i)-s)/A(i,i);end④solveUpTriangle.mfunction x=solveUpTriangle(A,b)% 求上三角系数矩阵的线性方程组Ax=bN=size(A);n=N(1);for i=n:-1:1if (i<n)s=A(i,(i+1):n)*x((i+1):n,1);elses=0;endx(i,1)=(b(i)-s)/A(i,i);end3.Jacobi法function [x,n]=jacobi(a,x0)a=input('请输入系数矩阵的维数n=');x0=input('请输入迭代初始向量x0');eps=1.0e-6;%解的精度控制m=10000;%迭代步数控制H=hilb(a);%生成h矩阵x_exact=ones(a,1);%求出x精确值b=H*x_exact;D=diag(diag(H));L=-tril(H,-1);U=-triu(H,1);B=D\(L+U);f=D\b;x=x0;n=0;tol=1;while tol>=epsx=B*x0+f;n=n+1;tol=norm(x-x0);x0=x;if (n>=m)disp('迭代次数太多，可能不收敛');return;endenddisp(x)disp(n)4.GS法function [x,n]=gauseidel(a,x0)a=input('请输入系数矩阵的维数n=');x0=input('请输入迭代初始向量x0');eps=1.0e-6;%解的精度控制m=100000;%迭代步数控制H=hilb(a);%生成h矩阵x_exact=ones(a,1);%求出x精确值b=H*x_exact;D=diag(diag(H));L=-tril(H,-1);U=-triu(H,1);G=(D-L)\U;f=(D-L)\b;x=x0;n=0;tol=1;while tol>=epsx=G*x0+f;n=n+1;tol=norm(x-x0);x0=x;if (n>=m)disp('迭代次数太多，可能不收敛');return;endenddisp(x)disp(n)5.SOR法function [x,n]=SOR(a,x0,w)a=input('请输入系数矩阵的维数n=');x0=input('请输入迭代初始向量x0=');w=input('请输入w=');eps=1.0e-6;%解的精度控制m=100000;%迭代步数控制H=hilb(a);%生成h矩阵x_exact=ones(a,1);%求出x精确值b=H*x_exact;if (w<=0||w>=2)error;return;endD=diag(diag(H));L=-tril(H,-1);U=-triu(H,1);B=inv(D-L*w)*((1-w)*D+w*U);f=w*inv((D-L*w))*b;x=x0;n=0;tol=1;while tol>=epsx=B*x0+f;n=n+1;tol=norm(x-x0);x0=x;if (n>=m)disp('迭代次数太多，可能不收敛');return;endenddisp(x)disp(n)。
病态矩阵 预处理 matlab
病态矩阵在数学和计算机科学中是一个重要的概念，它在多种领域都有着广泛的应用。
本文将结合 matlab 中的预处理方法，对病态矩阵进行深入探讨和分析，希望能够为相关领域的研究和实践提供一定的参考价值。
一、病态矩阵的定义和特性1. 病态矩阵的定义病态矩阵是指条件数(矩阵的最大奇异值与最小奇异值的比值)非常大的矩阵，这意味着它对输入数据的微小扰动非常敏感。
具体来说，当一个矩阵的条件数远大于1时，就可以认为它是病态矩阵。
2. 病态矩阵的特性病态矩阵在数学和计算机科学中有着重要的地位，它的特性主要包括：（1）对输入数据敏感：病态矩阵对输入数据的微小扰动非常敏感，可能导致计算结果的误差被放大。
（2）计算困难：病态矩阵在求逆、求解线性方程组等计算过程中可能会出现数值不稳定的情况，导致计算困难和误差的积累。
二、病态矩阵的预处理方法1. 规范化规范化是一种常用的病态矩阵预处理方法，它通过将矩阵的每一行除以该行的范数，使得矩阵的每一行都具有单位范数。
这样可以有效地减小矩阵的条件数，提高计算的稳定性和准确性。
2. 正交化正交化是另一种常用的病态矩阵预处理方法，它通过将矩阵进行正交变换，使得矩阵的列向量之间相互正交。
这样可以提高矩阵的条件数，改善计算的稳定性和精度。
3. 奇异值分解奇异值分解是一种更为复杂但有效的病态矩阵预处理方法，它通过将矩阵分解为奇异值矩阵、左奇异向量和右奇异向量，从而减小矩阵的条件数，提高计算的准确性。
三、matlab 中的病态矩阵预处理实现matlab 是一种强大的科学计算软件，它提供了丰富的病态矩阵预处理工具和函数，可以方便地进行病态矩阵的处理和分析。
下面我们主要介绍 matlab 中常用的病态矩阵预处理函数和实现方法。
1. cond 函数matlab 中的 cond 函数可以用来计算矩阵的条件数，从而帮助用户评估矩阵的病态程度。
通过 cond 函数，用户可以及时发现矩阵的病态特性，为后续的预处理和计算提供参考。
希尔伯特矩阵的病态性问题
希尔伯特矩阵的病态性问题（学号：200921100102 姓 名：陈丹丹 ）一、问题叙述希尔伯特矩阵是著名的病态矩阵。
n 阶Hilbert 矩阵为⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡−++=)12/(1)1/(1/1)1/(13/12/1/12/11n n n n n H n L L L L L L L 其条件数 Cond(H n ) 如下表所示n2 3 4 5 6 Cond(H n )1.9281e+001 5.2406e+002 1.5514e+004 4.7661e+005 1.4951e+007随着n 的增大，矩阵条件数迅速增加。
猜测：希尔伯特矩阵条件数以指数规律增长。
即，设矩阵阶数为n ，有Cond(H n ) ≈ exp( a n + b )对表中的条件数做对数变换，问题转化为线性拟合问题ln[Cond(H n )] = a n + b ，( n = 2，3，4，5，6)线性函数的两个系数分别为a = 3.3935，b = – 3.8811故指数拟合函数为：Cond(H n ) ≈ exp(3.3935 n –3.8811 )拟合函数的残差向量为r 2 r 3 r 4 r 5 r 69.9860e-001 -2.0231e+001 -6.8994e+002 -5.7825e+003 5.9012e+005 我们讨论以下三个问题：（1）由于线性拟合所得残差向量较大，下面用实验确定：选择用几次多项式拟合较好（2）设D 是由H n 的对角线元素开方构成的对角矩阵，令11ˆ−−=D H D H nn ，不难看出该矩阵仍然是对称正定矩阵，而且其对角元素全为1。
将H n 变换为n H ˆ的技术称为预处理（Preconditioning ）。
试分析函数)](/)ˆ(ln[nn H Cond H Cond 随n 变化的规律，绘制该函数的曲线。
（3）如果只用左预处理，即D 取希尔伯特矩阵对角元构成对角矩阵，用D 的逆阵左乘H 后条件数如何。
几何的角度解释矩阵出现病态
矩阵的范数范数的定义若X是数域K上的线性空间，泛函║·║: X->R 满足：1. 正定性：║x║≥0，且║x║=0 <=> x=0；2. 正齐次性：║AB║│≤│A║B║；3. 次可加性(三角不等式)：║A+B║≤║A║+║B║ 。
那么║·║称为X上的一个范数。
（注意到║x+y║≤║x║+║y║中如令y=-x，再利用║-x║=║x║可以得到║x║≥0，即║x║≥0在定义中不是必要的。
）在AX=B式中，当A、B有微小扰动时，参数估值Ｘ也有扰动：（A＋ΔA）（Ｘ＋ΔＸ）＝B＋ΔB （1）AΔＸ＝ΔA（Ｘ＋ΔＸ）+ΔBΔＸ＝－A-1ΔA（Ｘ＋ΔＸ）＋A－１ΔB将上式两端取范数，并应用向量范数的三角不等式║AB║│≤│A║B║；║A+B║≤║A║+║B║ 。
||ΔＸ||≤||A-1||||ΔA||（||Ｘ|| ＋||ΔＸ||）＋||A-1| ||ΔB|| 把含有ΔＸ的项移到式子的左边有：(1-||A-1||||ΔA||)||ΔＸ||≤||A-1||||ΔA|| ||Ｘ||+||A-1| ||ΔB||由于有1/||Ｘ||≤||A||/||B||将上式两端同时乘以1/||Ｘ||得：(1-||A-1||||ΔA||)||ΔＸ||/||X||≤||A-1||||ΔA|| +||A-1||||ΔB||A||/||B||;设K=||A||||A-||将上式整理的：（1-K||ΔA||/||A||）||ΔＸ||/||X||≤K(||ΔA||/||A||+||ΔB||/||B||);即有：||ΔＸ||/||X||≤k/（1-K||ΔA||/||A||）(||ΔA||/||A||+||ΔB||/||B||);问题与实验1:试从几何的角度解释矩阵出现病态的原因，并用‘有说服力’的例子来支持你的观点；线性方程组解的敏感性的几何解释(2x2矩阵)线性方程组求解：两直线求交点下面两图分别反映了良态问题和病态问题两种情况。
题目：主要研究希尔伯特矩阵的病态性
题目：主要研究希尔伯特矩阵的病态性Hilbert 矩阵的定义为： H n =⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-++)1n 2/(1)1n /(1/1）1n /(13/12/1n /12/11 n 该矩阵是一个对称正定矩阵，其条件数随n 的增加迅速增加。
画出ln （cond （H n ）2）~n 之间的曲线，观察他们的关系。
其中cond （H n ）2表示矩阵在2范数下的条件数，它可按如下公式计算： （1） cond （H n ）2=nλλ1，1λ和n λ分别是H n 的最大特征值和最小特征值。
（2） 设D 是H n 对角元素开方构成的对角矩阵，令nH∧=D 1- H n D 1-，不难看出nH∧依然还是对称矩阵，，并且对角元素都是1，画出ln （cond （H n ）2）~之间的曲线，观察他们之间的关系。
（3） 对于4≤n ≤12，给定不同的右端项b ，分别用x 1=H 1-n ；x 2= D 1-（nH ∧1-D 1-b ）；x 3= H n \b ，求解H n x=b ，比较结果x 1 、x 2、 x 3。
程序如下：首先，建立eiggg-M 文件，求解矩阵的最大特征值和最小特征值。
function [max,min]=eiggg(s) ss=size(s) if ss==1; max=s(1); min=s(1); elsemax=s(1); min=s(1); for i=1:1:ss; if s(i)>=max max=s(i); endif s(i)<=min min=s(i); end endend再建立Hilbert.m文件，给出Hilbert矩阵的定义：function H=Hilbert(n)H=zeros(n,n);for i=1:1:n;for j=1:1:n;H(i,j)=1/(i+j-1);endend然后求出D的逆矩阵function H=Hilbert2(n)H=zeros(n,n);for j=1:1:n;H(j,j)=1/(2*j-1);endD=sqrt(H);D1=inv(D);H=D1*H*H;编译主程序：clear;clc;N=4;mv2=zeros(N,1);for i=1:1:N;H=Hilbert2(i);s=eig(H);[max,min]=eiggg(s);mv1=max/min;mv2(i)=log(mv1);endn=1:1:N;figure(1),subplot(221),plot(n,mv2(n));title('ln(cond(Hn)2)与n之间的关系');H1=Hilbert(N);H2=Hilbert2(N);b=[1;2;3;4;];;x1=inv(H1)*b;x2=H2*b;x3=H1\b;sss=size(x1)i=1:1:sss(1);subplot(222),plot(i,x1);subplot(223),plot(i,x2);subplot(224),plot(i,x3);运行之后的结果：ss =1 1ss =2 1ss =3 1ss =4 1sss =4 1得出图像：由此可见，其条件数随n的增加而增加，则该矩阵不收敛，故则证得希尔伯特矩阵式的病态性。
计算方法Hilbert矩阵病态问题
Solution （1）猜想：()()ln cond Hn n 成线性比例关系。
验证：猜想部分正确。
根据不同n 建立对应的Hilbert 矩阵，计算该矩阵的2-范数的条件数，并绘制曲线。
编程实现见附录1中solution(1)。
实际编程中，分别取N=5，=10N ，=20N ，N=30，可以发现()()ln cond Hn n 曲线先是按线性规律变化；当n 达到大约12以后，()()ln cond Hn 的取值在40~50之间产生振荡。
各结果如下：N=5时 N=10时ln(cond(Hn))-nnl n (c o n d (H n ))ln(cond(Hn))-nnl n (c o n d (H n ))ln(cond(Hn))-nnl n (c o n d (H n ))ln(cond(Hn))-nnl n (c o n d (H n ))N=20时N=30时Solution （2）猜想：通过对Hilbert 矩阵的预处理，能改善Hilbert 的病态性质。
验证：根据不同n 建立对应的Hilbert 矩阵，计算该矩阵的2-范数的条件数，其次计算经过预处理后的H^矩阵的2-范数的条件数，并绘制()()()ln ^/cond Hn cond Hn n 曲线。
编程实现见附录1中的solution(2)。
实际编程中，分别取N=5，=12N ，=50N ，N=100，可以发现所得到的变化曲线()()()ln ^/cond Hn cond Hn n 在12n ≤时，()()()ln ^/cond Hn cond Hn 的取值较为平缓的减小，且均小于0。
而随着Hilbert 矩阵阶数的增大，()()()ln ^/cond Hn cond Hn 的取值在[]-7,3区间中振荡，主要集中在[]-31，，且多数取值()()()ln ^/0cond Hn cond Hn ≤。
则对多数的Hilbert 矩阵，其条件数在经过预处理后，有()()^cond Hn cond Hn ≤，即预处理在一定程度上改善了原Hilbert 矩阵的病态性质。
几何角度解释矩阵出现病态
矩阵的范数范数的定义若X是数域K上的线性空间，泛函║·║: X->R 满足：1. 正定性：║x║≥0，且║x║=0 <=> x=0；2. 正齐次性：║AB║│≤│A║B║；3. 次可加性(三角不等式)：║A+B║≤║A║+║B║ 。
那么║·║称为X上的一个范数。
（注意到║x+y║≤║x║+║y║中如令y=-x，再利用║-x║=║x║可以得到║x║≥0，即║x║≥0在定义中不是必要的。
）在AX=B式中，当A、B有微小扰动时，参数估值Ｘ也有扰动：（A＋ΔA）（Ｘ＋ΔＸ）＝B＋ΔB （1）AΔＸ＝ΔA（Ｘ＋ΔＸ）+ΔBΔＸ＝－A-1ΔA（Ｘ＋ΔＸ）＋A－１ΔB将上式两端取范数，并应用向量范数的三角不等式║AB║│≤│A║B║；║A+B║≤║A║+║B║ 。
||ΔＸ||≤||A-1||||ΔA||（||Ｘ|| ＋||ΔＸ||）＋||A-1| ||ΔB|| 把含有ΔＸ的项移到式子的左边有：(1-||A-1||||ΔA||)||ΔＸ||≤||A-1||||ΔA|| ||Ｘ||+||A-1| ||ΔB||由于有1/||Ｘ||≤||A||/||B||将上式两端同时乘以1/||Ｘ||得：(1-||A-1||||ΔA||)||ΔＸ||/||X||≤||A-1||||ΔA|| +||A-1||||ΔB||A||/||B||;设K=||A||||A-||将上式整理的：（1-K||ΔA||/||A||）||ΔＸ||/||X||≤K(||ΔA||/||A||+||ΔB||/||B||);即有：||ΔＸ||/||X||≤k/（1-K||ΔA||/||A||）(||ΔA||/||A||+||ΔB||/||B||); 问题与实验1:试从几何的角度解释矩阵出现病态的原因，并用‘有说服力’的例子来支持你的观点；线性方程组解的敏感性的几何解释(2x2矩阵)线性方程组求解：两直线求交点下面两图分别反映了良态问题和病态问题两种情况。
病态矩阵的特征值
病态矩阵的特征值嘿，朋友们！今天咱来聊聊病态矩阵的特征值这个神奇的玩意儿。
你说这病态矩阵啊，就好像是一个有点“脾气古怪”的家伙。
咱正常的矩阵就像是个乖巧听话的孩子，啥都好商量，可病态矩阵呢，就有那么些让人捉摸不透的地方。
咱想想啊，特征值这东西就像是矩阵的“性格特点”。
对于一般的矩阵，它的特征值都挺正常，能让咱清楚地了解它。
可病态矩阵的特征值啊，就有点像那调皮捣蛋的小孩子，时不时就给你闹出点幺蛾子来。
比如说吧，正常矩阵的特征值计算起来顺顺利利的，结果也很明确。
但病态矩阵的特征值呢，可能就会变得特别奇怪，让你算半天都摸不着头脑。
这就好比你要找个东西，在正常地方一下子就找到了，可在病态矩阵这里啊，就好像那东西藏得特别深，让你一顿好找。
你说这病态矩阵的特征值为啥就这么特殊呢？这就好像有些人的性格就是很特别一样，没办法呀！它就是这么与众不同。
有时候你想用常规方法去理解它，嘿，它还就不按套路出牌。
咱再打个比方，正常矩阵就像是一条笔直的大路，你顺着走就行。
可病态矩阵呢，就像是一条弯弯曲曲的小道，你得小心翼翼地走，不然一不小心就会迷路。
那遇到病态矩阵的特征值咱咋办呢？这可得有点耐心和技巧了。
不能着急，得慢慢琢磨，就像解开一个复杂的谜题一样。
而且啊，还得多尝试几种方法，可不能在一棵树上吊死。
你想想看，如果一直被病态矩阵的特征值给难住了，那多憋屈啊！咱可不能轻易被它打败，得想办法搞定它。
这就跟咱生活中遇到困难一样，不能退缩，得迎头而上。
反正啊，病态矩阵的特征值虽然有点麻烦，但咱也不能怕它。
咱要勇敢地去面对，去探索，去找到解决问题的办法。
咱可不能被它这点小古怪给吓住了，对吧？总之，病态矩阵的特征值虽然有点让人头疼，但也是个很有意思的东西。
它让我们知道，数学的世界里可不是只有那些规规矩矩的东西，还有很多奇怪又有趣的现象等着我们去发现呢！所以，大胆地去和病态矩阵的特征值过过招吧！。
奇异矩阵的例子
奇异矩阵的例子
1. 嘿，你知道吗？就像一个拼图少了关键一块就拼不完整一样，矩阵里也有这样奇怪的存在呢！比如说一个矩阵，它的行列式值为 0，这可就成了奇异矩阵啦！就好像走路突然没了方向一样让人不知所措。
2. 哇塞，来看看这个例子。
当一个矩阵没办法找到唯一解的时候，哎呀，这不就是奇异矩阵嘛！就好比你找东西怎么找都找不到，急得团团转呀！
3. 嘿呀，想想看，如果一个矩阵在解方程的时候老是出问题，那它大概率就是奇异矩阵咯！这不就像一辆车总是熄火，怎么都跑不起来嘛！
4. 哎呀呀，比如有个矩阵，不管你怎么努力去求解，都没办法得到确定的结果，那它妥妥的是奇异矩阵呀！这多像在迷雾中找不到出路一样让人郁闷呢！
5. 诶，你瞧这个，一个矩阵它自己有些“毛病”，导致很多操作都不能正常进行，这不是奇异矩阵还能是什么呀！这就和一个生病的人没办法正常工作一样让人无奈呀！
6. 哇哦，有一种矩阵，它就像一个不听话的小孩，总是不按常理出牌，那就是奇异矩阵啦！真是让人又爱又恨呢！
7. 嘿，你想想，如果一个矩阵总是给你制造难题，让你头疼不已，那毫无疑问是奇异矩阵嘛！就好像一只调皮的小猫总在捣乱一样。
8. 哈哈，还有那种怎么看都很奇怪的矩阵，一深究，嘿，就是奇异矩阵呀！简直就像一个行为古怪的人一样让人捉摸不透呀！
我的观点结论就是：奇异矩阵真的很特别，让人又好奇又苦恼呀！。
病态矩阵的例子
>> norm(x'-ones(8,1))
ans =
9.3866e-007
>> norm(A*x'-b)
ans =
7.0278e-008
1.0000
1.0000
就是说，对于系数矩阵为 8 阶的 Hilbert 矩阵，用 Matlab 中 的”inv”语句，大致能达到 6-7 位有效数字，这时条件数为 10 的 10 次方量级。
Results may be inaccurate. RCOND = 2.632091e-017. >> rank(A) ans =
11
>> A=hilb(100); >> rank(A) ans =
18 >> A=hilb(1000); >> rank(A) ans =
24
>> rank(hilb(2000)) ans =
病态矩阵的例子
>> Eigenvalue=eig(A) Eigenvalue =
0.0000 0.0003 0.0114 0.2085 1.5671 >> cond(A) ans = 4.7661t;> x=inv(A)*b; >> x solution =x
0.0000
0.0000
0.0000
0.0262 0.2981 1.6959
>> cond(A)
ans =
1.5258e+010
0.0001
0.0015
>> b=A*ones(8,1);
>> x=(inv(A)*b)'
矿区坐标转换中病态矩阵与坐标粗差的处理算法
矿区坐标转换中病态矩阵与坐标粗差的处理算法矿区坐标转换是指将采矿现场的实际坐标转换成地理坐标的过程。
在矿区坐标转换过程中，由于各种原因，可能会出现病态矩阵和坐标粗差的情况。
本文将介绍病态矩阵和坐标粗差的处理算法。
1. 病态矩阵的处理算法病态矩阵是指矩阵的条件数非常大，即矩阵的行列式很接近于零的情况。
在矿区坐标转换过程中，病态矩阵的存在会导致计算结果非常不稳定，甚至可能无法得到准确的解。
病态矩阵的处理算法一般有以下几种方法：（1）数据平差法：利用最小二乘法对数据进行平差处理，通过求解最优解来降低矩阵的条件数。
（2）奇异值分解法：将矩阵分解为奇异值矩阵的乘积形式，通过对奇异值进行处理来减小矩阵的条件数。
（3）正则化方法：在目标函数中引入正则化项，通过调整正则化参数来平衡模型的复杂度和数据的拟合程度，从而降低矩阵的条件数。
2. 坐标粗差的处理算法坐标粗差是指由于观测误差或其他原因导致的坐标值明显偏离真实值的情况。
在矿区坐标转换过程中，坐标粗差的存在会对计算结果产生较大的影响，因此需要采取适当的处理方法。
（1）3σ原则：根据正态分布的性质，坐标观测值的误差在三个标准差范围内的概率约为99.7%，所以可以将超出三个标准差范围的观测值认为是坐标粗差，并进行剔除或修正。
（2）最小二乘法：通过采用最小二乘法进行平差处理，可以将坐标粗差的影响最小化，并得到尽可能准确的结果。
（3）加权最小二乘法：根据坐标观测值的精度估计，为每个观测值赋予不同的权重，通过加权最小二乘法进行计算，可以有效降低坐标粗差的影响。
病态矩阵和坐标粗差是矿区坐标转换中常见的问题，需要通过合适的处理算法来解决。
根据具体情况选择适用的方法，可以提高坐标转换的准确性和稳定性。
几何的角度解释矩阵出现病态
矩阵的范数范数的定义若X是数域K上的线性空间，泛函║·║: X->R 满足：1. 正定性：║x║≥0，且║x║=0 <=> x=0；2. 正齐次性：║AB║│≤│A║B║；3. 次可加性(三角不等式)：║A+B║≤║A║+║B║ 。
那么║·║称为X上的一个范数。
（注意到║x+y║≤║x║+║y║中如令y=-x，再利用║-x║=║x║可以得到║x║≥0，即║x║≥0在定义中不是必要的。
）在AX=B式中，当A、B有微小扰动时，参数估值Ｘ也有扰动：（A＋ΔA）（Ｘ＋ΔＸ）＝B＋ΔB （1）AΔＸ＝ΔA（Ｘ＋ΔＸ）+ΔBΔＸ＝－A-1ΔA（Ｘ＋ΔＸ）＋A－１ΔB将上式两端取范数，并应用向量范数的三角不等式║AB║│≤│A║B║；║A+B║≤║A║+║B║ 。
||ΔＸ||≤||A-1||||ΔA||（||Ｘ|| ＋||ΔＸ||）＋||A-1| ||ΔB|| 把含有ΔＸ的项移到式子的左边有：(1-||A-1||||ΔA||)||ΔＸ||≤||A-1||||ΔA|| ||Ｘ||+||A-1| ||ΔB||由于有1/||Ｘ||≤||A||/||B||将上式两端同时乘以1/||Ｘ||得：(1-||A-1||||ΔA||)||ΔＸ||/||X||≤||A-1||||ΔA|| +||A-1||||ΔB||A||/||B||;设K=||A||||A-||将上式整理的：（1-K||ΔA||/||A||）||ΔＸ||/||X||≤K(||ΔA||/||A||+||ΔB||/||B||);即有：||ΔＸ||/||X||≤k/（1-K||ΔA||/||A||）(||ΔA||/||A||+||ΔB||/||B||);问题与实验1:试从几何的角度解释矩阵出现病态的原因，并用‘有说服力’的例子来支持你的观点；线性方程组解的敏感性的几何解释(2x2矩阵)线性方程组求解：两直线求交点下面两图分别反映了良态问题和病态问题两种情况。
《病态矩阵的例子》课件
03
CHAPTER
病态矩阵的影响
对计算结果的影响
不正确的解
病态矩阵可能导致数值计算的解完全 错误，即使算法本身是正确的。
解的不稳定性
即使使用相同的算法和相同的数据， 由于舍入误差，不同的计算可能会得 到不同的解。
对算法稳定性的影响
病态矩阵可能导致算法不稳定
在迭代算法中，病态矩阵可能导致算法发散或收敛到错误的 解。
21, 22, 23, 24, 25]
数值例子
``` 计算其条件数为inf，说明这是一个病态矩阵。
数值例子2：给定一个3x3的矩阵B，其元素如下
数值例子
``` [1000, 0.0001, 0.000001; 0.0001, 1000, 0.0001;
数值例子
• 0.000001, 0.0001, 1000]
迭代收敛性检查
在迭代过程中，检查算法 的收敛性，如果发现不收 敛的情况，则采取相应的 措施进行修正。
误差估计
在迭代过程中，估计误差 的大小，并根据误差的大 小调整算法的参数，以提 高数值精度。
病态矩阵的避免方法
避免零元素的出现
在构建矩阵时，尽量避免出现零元素，因为零元素可能会导致算 法的不稳定。
选择合适的数值方法
病态矩阵的识别
高条件数
可以通过计算矩阵的条件数来判断一个矩阵是否为病态矩 阵。如果条件数很大，则该矩阵可能是病态的。
接近奇异或退化
通过检查矩阵的特征值和奇异值，可以判断一个矩阵是否 接近奇异或退化。如果特征值或奇异值很小，则该矩阵可 能是病态的。
元素变化对结果影响
可以通过比较不同初始条件下计算结果的差异来判断一个 矩阵是否为病态矩阵。如果不同初始条件下的结果差异很 大，则该矩阵可能是病态的。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
从一个数值结果谈起
2 n 以 1, x, x , x 为基底， 构造不超过 n 次最佳平方逼近多项
式时，其法方程的系数矩阵会出现所谓 Hilbert 矩阵。 n=5 >> A=hilb(5) A= 1.0000 0.5000 0.3333 0.2500 0.2000 0.5000 0.3333 0.2500 0.2000 0.1667 0.3333 0.2500 0.2000 0.1667 0.1429 0.2500 0.2000 0.1667 0.1429 0.1250 0.2000 0.1667 0.1429 0.1250 0.1111
>> A=hilb(100); >> rank(A) ans = 18 >> A=hilb(1000); >> rank(A) ans = 24
>> rank(hilb(2000)) ans = 25 将系数矩阵作微调为 Hilbert(2000)+7*10(-7)* I 2000 >> A=hilb(2000)+0.00000007*eye(2000); >> rank(A) ans = 2000
>> Eigenvalue=eig(A) Eigenvalue = 0.0000 0.0003 0.0114 0.2085 1.5671 >> cond(A) ans = 4.7661e+005 >> b=A*ones(5,1);
>> x=inv(A)*b; >> x solution =x 1.0000 1.0000 1.0000 1.0000 1.0000 >> norm(ones(5,1)-x) ans = 3.0055e-011 >> norm(b-A*x) ans = 6.9024e-012
n=10; >> A=hilb(10); >> cond(A) ans = 1.6025e+013
> b=A*ones(10,1); >> solution=inv(A)*b; >> solution' ans = 1.0000 1.0003 1.0000 1.0007 1.0000 0.9996 1.0000 1.0001 0.9999
n=8 >> A=hilb(8); >> Eigenvalue=eig(A)' Eigenvalue = 0.0000 0.0262 >> cond(A) ans = 1.5258e+010 0.0000 0.2981 0.0000 1.6959 0.0001 0.0015
>> b=A*ones(8,1); >> x=(inv(A)*b)' x= 1.0000 1.0000 1Байду номын сангаас0000 1.0000 1.0000 1.0000 1.0000 1.0000
各类数学软件的广泛使用， 使得大量以前不能计算或 者计算非常困难的问题成 为可能或容易计算的问题。 也就是说， 现代数值方法丰 硕成果，譬如上世纪 70 年 代以来，样条函数、FFT、 有限元方法、 区域分解方法 等等都是非常有效的方法。
实际应用无止境， 需要解决 的问题也无止境。譬如 1. 我们难以验证通过这些 数以亿计的计算得到的 结果是否可信， 许多数值 分析家为此在努力， 图灵 奖 获 得 者 Wilkson 、 Kuhan 等； 2. 理论和实际的差异； 3. -----
0.9991
>> norm(ones(10,1)-solution) ans = 0.0012 >> norm(A*solution-b) ans = 7.6071e-005
n=12; >> A=hilb(12); >> cond(A) ans = 1.7945e+016 >> x=inv(A)*(A*ones(12,1)); Warning: Matrix is close to singular or badly scaled. Results 2.632091e-017. >> rank(A) ans = 11 may be inaccurate. RCOND =
从这里至少可以看出三点，一理论上的满 秩不能保证实际计算的满秩，有时甚至相 差很大；二数学软件只能解决常态的问题， 难以解决大量存在的有特殊性态的问题； 三要研究非方阵、非满秩的问题。
>> norm(x'-ones(8,1)) ans = 9.3866e-007 >> norm(A*x'-b) ans = 7.0278e-008
就是说，对于系数矩阵为 8 阶的 Hilbert 矩阵，用 Matlab 中 的”inv”语句，大致能达到 6-7 位有效数字，这时条件数为 10 的 10 次方量级。