1楼:匿名用户
|[householder阵]
(1) 设a rn, = ||a||2,通常取 与a1同号,记h=i-2vvt,(v= ),
则ha= - e1. h=i -2vvt称为householder阵。
(2) 更一般地,对a=(a1,a2,…am,am+1,…,an)t,记 = ,可求出h,使
ha=(a1,a2,…am, ,0,…,0)t。
为此,先在rn-m中求 使 满足
=(am+1,…,an)t=(- ,0,…,0,0)t,
再作h= ,则ha= (a1,a2,…am,am+1,…,an)t =( a1,a2,…am,- ,0,…,0,0)t
[用householder方法求矩阵的qr分解]
记a=(aij)n*n,由1可知,存在h1=i -2v1v1t,使
h1(a11,a21,…,an1)t=(a11(1),0,…,0)t,
于是 h1a=
又由1知,存在h2= ,使 ,于是
h1a= =
类似地依次进行n-1次,得出
hn-1hn-2…h1a= 。
记r=hn-1hn-2…h1a,q=hnhn-1…h1,得a=q*r
如何用householder变换求矩阵的qr分解 例子
2楼:匿名用户
||[householder阵]
(1) 设a rn, = ||a||2,通常取 与a1同号,记h=i-2vvt,(v= ),
则ha= - e1. h=i -2vvt称为householder阵。
(2) 更一般地,对a=(a1,a2,…am,am+1,…,an)t,记 = ,可求出h,使
ha=(a1,a2,…am, ,0,…,0)t。
为此,先在rn-m中求 使 满足
=(am+1,…,an)t=(- ,0,…,0,0)t,
再作h= ,则ha= (a1,a2,…am,am+1,…,an)t =( a1,a2,…am,- ,0,…,0,0)t
[用householder方法求矩阵的qr分解]
记a=(aij)n*n,由1可知,存在h1=i -2v1v1t,使
h1(a11,a21,…,an1)t=(a11(1),0,…,0)t,
于是 h1a=
又由1知,存在h2= ,使 ,于是
h1a= =
类似地依次进行n-1次,得出
hn-1hn-2…h1a= 。
记r=hn-1hn-2…h1a,q=hnhn-1…h1,得a=q*r
qr分解采用的household什么算法
3楼:电灯剑客
也没有什么特别的名字,就是一个用镜像变换做消去(或者说上三角化)的算法
求助!怎样求矩阵特征值,要用household方法(有c语言程序最好)
4楼:饮水思春
//数值计算程序-特征值和特征向量
//约化对称矩阵为三对角对称矩阵
//利用householder变换将n阶实对称矩阵约化为对称三对角矩阵
//a-长度为n*n的数组,存放n阶实对称矩阵
//n-矩阵的阶数
//q-长度为n*n的数组,返回时存放householder变换矩阵
//b-长度为n的数组,返回时存放三对角阵的主对角线元素
//c-长度为n的数组,返回时前n-1个元素存放次对角线元素
void eastrq(double a,int n,double q,double b,double c);
//求实对称三对角对称矩阵的全部特征值及特征向量
//利用变型qr方法计算实对称三对角矩阵全部特征值及特征向量
//n-矩阵的阶数
//b-长度为n的数组,返回时存放三对角阵的主对角线元素
//c-长度为n的数组,返回时前n-1个元素存放次对角线元素
//q-长度为n*n的数组,若存放单位矩阵,则返回实对称三对角矩阵的特征向量组
// 若存放householder变换矩阵,则返回实对称矩阵a的特征向量组
//a-长度为n*n的数组,存放n阶实对称矩阵
int ebstq(int n,double b,double c,double q,double eps,int l);
//约化实矩阵为赫申伯格(hessen berg)矩阵
//利用初等相似变换将n阶实矩阵约化为上h矩阵
//a-长度为n*n的数组,存放n阶实矩阵,返回时存放上h矩阵
//n-矩阵的阶数
void echbg(double a,int n);
//求赫申伯格(hessen berg)矩阵的全部特征值
//返回值小于0表示超过迭代jt次仍未达到精度要求
//返回值大于0表示正常返回
//利用带原点位移的双重步qr方法求上h矩阵的全部特征值
//a-长度为n*n的数组,存放上h矩阵
//n-矩阵的阶数
//u-长度为n的数组,返回n个特征值的实部
//v-长度为n的数组,返回n个特征值的虚部
//eps-控制精度要求
//jt-整型变量,控制最大迭代次数
int edqr(double a,int n,double u,double v,double eps,int jt);
//求实对称矩阵的特征值及特征向量的雅格比法
//利用雅格比(jacobi)方法求实对称矩阵的全部特征值及特征向量
//返回值小于0表示超过迭代jt次仍未达到精度要求
//返回值大于0表示正常返回
//a-长度为n*n的数组,存放实对称矩阵,返回时对角线存放n个特征值
//n-矩阵的阶数
//u-长度为n*n的数组,返回特征向量(按列存储)
//eps-控制精度要求
//jt-整型变量,控制最大迭代次数
int eejcb(double a,int n,double v,double eps,int jt);
选自《徐世良数值计算程序集(c)>>
每个程序都加上了适当地注释,陆陆续续干了几个月才整理出来的啊。
今天都给贴出来了
#include "stdio.h"
#include "math.h"
//约化对称矩阵为三对角对称矩阵
//利用householder变换将n阶实对称矩阵约化为对称三对角矩阵
//a-长度为n*n的数组,存放n阶实对称矩阵
//n-矩阵的阶数
//q-长度为n*n的数组,返回时存放householder变换矩阵
//b-长度为n的数组,返回时存放三对角阵的主对角线元素
//c-长度为n的数组,返回时前n-1个元素存放次对角线元素
void eastrq(double a,int n,double q,double b,double c)
}for (i=n-1; i>=1; i--)
}if (h+1.0==1.0)
b[i]=0.0;
}else
h=h-q[u]*c[i-1];
q[u]=q[u]-c[i-1];
f=0.0;
for (j=0; j<=i-1; j++)
if (j+1<=i-1)
}c[j-1]=g/h;
f=f+g*q[j*n+i];
}h2=f/(h+h);
for (j=0; j<=i-1; j++)
}b[i]=h;}}
b[0]=0.0;
for (i=0; i<=n-1; i++)
for (k=0; k<=i-1; k++)}}
u=i*n+i;
b[i]=q[u];
q[u]=1.0;
if (i-1>=0)}}
return;
//求实对称三对角对称矩阵的全部特征值及特征向量
//利用变型qr方法计算实对称三对角矩阵全部特征值及特征向量
//n-矩阵的阶数
//b-长度为n的数组,返回时存放三对角阵的主对角线元素
//c-长度为n的数组,返回时前n-1个元素存放次对角线元素
//q-长度为n*n的数组,若存放单位矩阵,则返回实对称三对角矩阵的特征向量组
// 若存放householder变换矩阵,则返回实对称矩阵a的特征向量组
//a-长度为n*n的数组,存放n阶实对称矩阵
int ebstq(int n,double b,double c,double q,double eps,int l)
m=j;
while ((m<=n-1)&&(fabs(c[m])>d))
if (m!=j)
it=it+1;
g=b[j];
p=(b[j+1]-g)/(2.0*c[j]);
r=sqrt(p*p+1.0);
if (p>=0.0)
else
h=g-b[j];
for (i=j+1; i<=n-1; i++)
f=f+h;
p=b[m];
e=1.0;
s=0.0;
for (i=m-1; i>=j; i--)
else
p=e*b[i]-s*g;
b[i+1]=h+s*(e*g+s*b[i]);
for (k=0; k<=n-1; k++)
}c[j]=s*p;
b[j]=e*p;
}while (fabs(c[j])>d);
}b[j]=b[j]+f;
}for (i=0; i<=n-1; i++)
}if (k!=i)}}
return(1);
}//约化实矩阵为赫申伯格(hessen berg)矩阵
//利用初等相似变换将n阶实矩阵约化为上h矩阵
//a-长度为n*n的数组,存放n阶实矩阵,返回时存放上h矩阵
//n-矩阵的阶数
void echbg(double a,int n)
}if (fabs(d)+1.0!=1.0)
for (j=0; j<=n-1; j++)
}for (i=k+1; i<=n-1; i++)
for (j=0; j<=n-1; j++)}}
}return;
}//求赫申伯格(hessen berg)矩阵的全部特征值
//利用带原点位移的双重步qr方法求上h矩阵的全部特征值
//返回值小于0表示超过迭代jt次仍未达到精度要求
//返回值大于0表示正常返回
//a-长度为n*n的数组,存放上h矩阵
//n-矩阵的阶数
//u-长度为n的数组,返回n个特征值的实部
//v-长度为n的数组,返回n个特征值的虚部
//eps-控制精度要求
//jt-整型变量,控制最大迭代次数
int edqr(double a,int n,double u,double v,double eps,int jt)
ii=(m-1)*n+m-1;
jj=(m-1)*n+m-2;
kk=(m-2)*n+m-1;
ll=(m-2)*n+m-2;
if (l==m-1)
else if (l==m-2)
u[m-1]=(-b-xy*y)/2.0;
u[m-2]=c/u[m-1];
v[m-1]=0.0; v[m-2]=0.0;
}else
m=m-2;
it=0;
}else
it=it+1;
for (j=l+2; j<=m-1; j++)
for (j=l+3; j<=m-1; j++)
for (k=l; k<=m-2; k++)
}else
if ((fabs(p)+fabs(q)+fabs(r))!=0.0)
s=xy*sqrt(p*p+q*q+r*r);
if (k!=l)
e=-q/s;
f=-r/s;
x=-p/s;
y=-x-f*r/(p+s);
g=e*r/(p+s);
z=-x-e*q/(p+s);
for (j=k; j<=m-1; j++)
a[jj]=q;
a[ii]=p;
}j=k+3;
if (j>=m-1)
for (i=l; i<=j; i++)
a[jj]=q;
a[ii]=p;}}
}}}return(1);
}//求实对称矩阵的特征值及特征向量的雅格比法
//利用雅格比(jacobi)方法求实对称矩阵的全部特征值及特征向量
//返回值小于0表示超过迭代jt次仍未达到精度要求
//返回值大于0表示正常返回
//a-长度为n*n的数组,存放实对称矩阵,返回时对角线存放n个特征值
//n-矩阵的阶数
//u-长度为n*n的数组,返回特征向量(按列存储)
//eps-控制精度要求
//jt-整型变量,控制最大迭代次数
int eejcb(double a,int n,double v,double eps,int jt)}}
while (1==1)}}
if (fmjt)
l=l+1;
u=p*n+q;
w=p*n+p;
t=q*n+p;
s=q*n+q;
x=-a[u];
y=(a[s]-a[w])/2.0;
omega=x/sqrt(x*x+y*y);
if (y<0.0)
sn=1.0+sqrt(1.0-omega*omega);
sn=omega/sqrt(2.0*sn);
**=sqrt(1.0-sn*sn);
fm=a[w];
a[w]=fm******+a[s]*sn*sn+a[u]*omega;
a[s]=fm*sn*sn+a[s]******-a[u]*omega;
a[u]=0.0;
a[t]=0.0;
for (j=0; j<=n-1; j++)
}for (i=0; i<=n-1; i++)
}for (i=0; i<=n-1; i++)
}return(1);}
对矩阵A进行QR分解,只有唯一的一种情况吗
1楼 匿名用户 是唯一的! 矩阵a进行qr分解后,r中对角元素一定是实数,由于这个条件,使得它唯一! 如何证明矩阵a的qr分解的唯一性 2楼 电灯剑客 a qr a a r r 用cholesky分解的唯一性得到r的唯一性,从而q ar 也唯一 3楼 挚爱红军 假设a是实的非奇异阵,则在要求r的对角...
关于矩阵的QR分解,我不明白下图中的R是怎么来的
1楼 匿名用户 a qr r q a q t a 因为q是正交矩阵 r11 q11 a11 q21 a12 a1 t q1 以此类推 矩阵qr分解的证明题 2楼 电灯剑客 r中所有对角元素非零 rank r n rank r hr n rank a ha n rank a n 至于第二个问题,这个没...
qr分解采用的household什么算法
1楼 电灯剑客 也没有什么特别的名字,就是一个用镜像变换做消去 或者说上三角化 的算法 求助householde矩阵qr分解 2楼 匿名用户 householder阵 1 设a rn, a 2,通常取 与a1同号,记h i 2vvt, v 则ha e1 h i 2vvt称为householder阵。...