一、算法设计方案:
1. 求出 A 的压缩后的矩阵 C,并对 C 作 LU 分
解;
2. 根据反幂法求出 λ
s
的值;
若矩阵 A 存在相等的模最小的特征值,则
当这些特征值完全相等时,反幂法可以求
出此特征值和此特征值对应的特征向量。
当这些特征值互为相反数时,由反幂法的
推导过程可知:此时亦能求出这些特征值
的模,但求不出与这些特征值对应的特征
向量。综上所述,反幂法可以求出矩阵 A
模最小的特征值。同理可知:幂法可以求
出矩阵 A 模最大的特征值。
3. 根据幂法求出矩阵 A 模最大的特征值 λ
m
;
4. 计算出矩阵 B = A + |λ
m
|I 的模最大的特征值
μ
m
和模最小的特征值 μ
s
。则有:
λ
501
= μ
m
- |λ
m
|
λ
1
= μ
s
- |λ
m
|
5. 求出 μ
k
;
μ
k
=λ
m
+ k(λ
501
–λ
1
)/40
用反幂法计算出矩阵 B
k
的模最小的特征值。
B
k
= A -μ
k
I
(其中,k = 1,2,…,39)。
6. 求出 A 的条件数 cond(A)
2
以及 A 的行列式
detA,其中:
cond(A)
2
= |λ
m
|/|λ
s
|;
detA = det(LU) = detU。
二、全部源程序:
//1.主程序
#include "stdafx.h"
#include "iostream.h"
#include "iomanip.h"
#include "head.h"
#include "math.h"
void main( )
{
double *a_D, *u0, *y; //u0 为迭代初始值
double **a_L, **a_U, **c;
//下面的 n 代表 A 阵维数,s 代表上半带宽,r
代表下半带宽,可以随意设置
const int n = 501;
const int s = 2;
const int r = 2;
int i = 0, j = 0;
a_D = new double [n];
u0 = new double [n];
y = new double [n];
a_U = new double *[s];
a_L = new double *[r];
c = new double *[r+s+1];
//生成压缩矩阵
for ( i = 0; i != r+s+1; ++i )
c[i] = new double [n];
for ( i = 0; i < s; ++i )
{
a_U[i] = new double [n-1-i];
}
//求 A 的行列式
double u_max = 0, u_min = 0;//u_max 为 A
模最大的特征值,u_min 为 A 模最小的特征值
//生成矩阵 A 的非零元素
arr_ger( a_U, a_D, a_L, s, r, n, u_max);
//将带状矩阵 A 压缩为矩阵 C
arr_compress( c, a_D, a_U, a_L, s, r, n);
LU_decomp( c, s, r, n );
cout<<"A 的 行 列 式
为:"<<setiosflags( ios::scientific )<<setprecision(1
2)<<deta( c, s, r, n )<<endl;
cout<<endl;
//求出 A 的模最大的特征值
u_max = m_max( c, a_U, a_D, a_L, s, r, n );
cout<<"A 的 模 最 大 特 征 值
为:"<<u_max<<endl;
cout<<endl;
// 求出 A 的模最小的特征值
u_min = m_min( c, a_U, a_D, y, a_L, s, r, n );
cout<<"A 的 模 最 小 特 征 值
为:"<<u_min<<endl;
cout<<endl;
//求出 A 最大的特征值