Jacobi迭代求矩阵特征值和特征向量+C代码

Jacobi计算过程如下：

1. 选择矩阵A非对角元中最大值A[i][j]，运用公式 tan 2O = 2*A[i][j] / (A[i][i] -A[j][j]) 获得选择平面矩阵J，使J * A *J'=A1

2. 计算得到A1后，重复第1歩，计算A2 =J2 *A1 *J2'，连续这个过程，直到非对角线元素全化为充分小，得到特征值矩阵Am。

3. 特征向量的计算为 D = J1' *J2'* ... *Jm' ，D 的列向量就是对应特征值矩阵 Am 对角线的特征向量。

下面是一个课本上的例子：

计算结果：（A的对角线元素为矩阵的特征值，Q为对应的特征向量(列对应)）

实现代码：（数据为书上的另一个例子）

#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "vector"
using namespace std;

#define N 5
#define E 0.000001 //非对角线数充分小
#define PI 3.141592657
#define MAXITER 10000 //最大迭代次数

typedef vector<double> dim1Vector;
typedef vector<dim1Vector> dim2Vector;
typedef vector<dim2Vector> dim3Vector;

bool QueryArray(dim2Vector Array); //检查是否满足
dim2Vector matTran(dim2Vector Array); //矩阵转置
dim2Vector matMul(dim2Vector mat1, dim2Vector mat2); //矩阵相乘

void main()
{
int i, j;
int count;
bool flag = false;
double A[N][N] = {1267.9, -307.23, 0, 0, 0,
-307.23, 710.24, -403.01, 0, 0,
0, -403.01, 927.21, -524.2, 0,
0, 0, -524.2, 770.07, -245.88,
0, 0, 0, -245.88, 245.88};
dim1Vector tempArray(N, 0);
dim2Vector Array;
dim2Vector charatMat(N, tempArray);
dim2Vector dim2Jac;
dim2Vector dim2JacT;
dim3Vector dim3Jac;
double maxArrayNum;
int laber_j, laber_i;
double theta;

for(i=0; i<N; i++)
{
tempArray.clear();
for(j=0; j<N; j++)
tempArray.push_back(A[i][j]);
Array.push_back(tempArray);
}

//开始迭代
count = 0;
tempArray.clear();
tempArray.resize(N, 0);
while(count<MAXITER && !flag)
{
count++;
dim2Jac.clear();
dim2Jac.resize(N, tempArray);
maxArrayNum = 0;
laber_i = laber_j = 0;

//寻找非对角元中绝对值最大的A[i][j]
for(i=0; i<N; i++)
for(j=0; j<N; j++)
{
if(i==j)
continue;

if(maxArrayNum<fabs(Array[i][j]))
{
maxArrayNum = fabs(Array[i][j]);
laber_i = i;
laber_j = j;
}
}

theta = atanf(Array[laber_i][laber_j]*2/(Array[laber_i][laber_i]-Array[laber_j][laber_j]+E));

//构造雅克比矩阵
for(i=0; i<N; i++)
dim2Jac[i][i] = 1;

dim2Jac[laber_i][laber_i] = dim2Jac[laber_j][laber_j] = cosf(theta/2);
dim2Jac[laber_i][laber_j] = sinf(theta/2);
dim2Jac[laber_j][laber_i] = -sinf(theta/2);

dim2JacT = matTran(dim2Jac); //矩阵转置
dim3Jac.push_back(dim2JacT); //保存矩阵

Array = matMul(matMul(dim2Jac, Array), dim2JacT);

if(QueryArray(Array))
flag = true;

}

//初始化特征矩阵
for(i=0; i<N; i++)
charatMat[i][i] = 1;

//计算特征矩阵
for(i=0; i<dim3Jac.size(); i++)
charatMat = matMul(charatMat, dim3Jac[i]);

printf("迭代次数：%d\n", count);

printf("\n特征值：\n");
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
printf("%lf ", Array[i][j]);
printf("\n");
}

printf("\n特征矩阵：\n");
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
printf("%lf ", charatMat[i][j]);
printf("\n");
}

}

//检查是否满足
bool QueryArray(dim2Vector Array)
{
int i, j;

for(i=0; i<N; i++)
for(j=0; j<N; j++)
{
if(i==j)
continue;

if(fabs(Array[i][j])>E)
return false;
}

return true;
}

//矩阵转置
dim2Vector matTran(dim2Vector Array)
{
int i, j;
dim1Vector temp(N, 0);
dim2Vector dst(N, temp);

for(i=0; i<N; i++)
for(j=0; j<N; j++)
dst[j][i] = Array[i][j];

return dst;
}

//矩阵相乘
dim2Vector matMul(dim2Vector mat1, dim2Vector mat2)
{
int i, j, k;
dim1Vector temp(N, 0);
dim2Vector dst(N, temp);

for(i=0; i<mat1.size(); i++)
for(j=0; j<mat2[0].size(); j++)
for(k=0; k<mat2.size(); k++)
dst[i][j] += mat1[i][k]*mat2[k][j];

return dst;
}

运行结果：

Jacobi迭代求矩阵特征值和特征向量+C代码相关推荐

求矩阵特征值和特征向量
求矩阵特征值和特征向量的一个小程序代码较长,如果不能执行,就是要建立结构体,大家试试吧,希望能用. // // 实对称三对角阵的全部特征值与特征向量的计算 // // 参数: // 1. doubl ...
QR分解求矩阵特征值、特征向量 C语言
最近在看一个高光谱图像压缩算法,其中涉及到正交变换,计算正交变换时,需要对普通矩阵求其特征向量.想要在网上找一个现成的程序,可能是我百度的能力不强吧,居然真的没找见.好了废话不多说,下面进入正题. 计 ...
matlab求矩阵特征值和特征向量、行列式
如果A为方阵,满足AX=λX的λ称为A的特征值,X称为A的特征向量. 计算A的特征值用eig(A). 例:A=[1 2 3;4 5 6;7 8 9]; Z=eig(A) Z = 16.1168 ...
python求矩阵特征值和特征向量_特征值，特征向量，标准正交向量组与numpy
使用python的数值计算库numpy来计算矩阵的特征值,特征向量与标准正交向量组 import numpy as np 1.求矩阵的特征值和各特征值所对应的特征向量 x = np.array([[ ...
MATLAB中求矩阵特征值和特征向量
格式:[X,B]=eig(A) %求矩阵A的特征值和特征向量,其中B的对角线元素是特征值 eg: 特征值0.8105对应的特征向量是[0.3365 -0.9417]: 特征值21.2895对应的特征向 ...
移位取逆迭代(shifted inverse iteration)求最近特征值和特征向量(python,数值积分)
第二十七篇移位取逆迭代求最近特征值和特征向量移位逆迭代一种比"最大"特征值法更直接实现向量迭代收敛的的特征值方法是将移位向量迭代法改写为下面形式其中p是一个标量" ...
雅可比算法求矩阵特征值C语言源代码,雅可比(Jacobi)计算特征值和特征向量
雅可比迭代法法在图形图像中不少地方用到求矩阵的特征值和特征向量,好比主成分分析.OBB包围盒等.编程时通常都是用数值分析的方法来计算,这里介绍一下雅可比迭代法求解特征值和特征向量.雅可比迭代法的原理 ...
pythonnumpy库求特征向量_Python使用numpy计算矩阵特征值、特征向量与逆矩阵
Python扩展库numpy.linalg的eig()函数可以用来计算矩阵的特征值与特征向量,而numpy.linalg.inv()函数用来计算可逆矩阵的逆矩阵. >>> impor ...
使用MTL库求解矩阵特征值和特征向量
关于矩阵的特征值和特征向量求解,大部分的数学运算库都进行了提供,下面是使用MTL库的接口进行封装. #include <mtl/matrix.h> #include <mtl/mtl ...

Jacobi迭代求矩阵特征值和特征向量+C代码

Jacobi迭代求矩阵特征值和特征向量+C代码相关推荐

最新文章

热门文章