求助householde矩阵QR分解

2020-11-25 11:24:02 字数 7158 阅读 2879

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阵。...