求矩阵特征值算法及程序
求矩阵特征值算法及程序
求矩阵特征值算法及程序简介1.幂法1、幂法规范化算法（1）输入矩阵A 、初始向量)0(μ，误差eps ；（2）1⇐k ； （3）计算)1()(-⇐k k A Vμ；（4）)max (,)max ()1(1)(--⇐⇐k k k k V m V m ；（5）k k k m V /)()(⇐μ；（6）如果eps m m k k <--1,则显示特征值1λ和对应的特征向量)1(x ),终止；（7）1+⇐k k ,转（3）注：如上算法中的符号)max(V 表示取向量V 中绝对值最大的分量。
本算法使用了数据规范化处理技术以防止计算过程中出现益出错误。
2、规范化幂法程序 Clear[a,u,x];a=Input["系数矩阵A="];u=Input["初始迭代向量u(0)="]; n=Length[u];eps=Input["误差精度eps ="];nmax=Input["迭代允许最大次数nmax="]; fmax[x_]:=Module[{m=0,m1,m2}, Do[m1=Abs[x[[k]]];If[m1>m,m2=x[[k]];m=m1], {k,1,Length[x]}]; m2] v=a.u;m0=fmax[u]; m1=fmax[v];t=Abs[m1-m0]//N; k=0;While[t>eps&&k<nmax, u=v/m1; v=a.u; k=k+1;m0=m1;m1=fmax[v];t=Abs[m1-m0]//N;Print["k=",k," 特征值=",N[m1,10]," 误差=",N[t,10]]; Print[" 特征向量=",N[u,10]]]; If[k ≥nmax,Print["迭代超限"]]说明：本程序用于求矩阵A 按模最大的特征值及其相应特征向量。
计算方法(5)第四章 矩阵特征值和特征向量的计算
n
使得u 0
i xi
i 1
n
n
uk Auk1 Aku0 Ak (i xi ) iik xi
i 1
i 1
1k [1x1
n i2
( i 1
)k i xi ]
由1 0, 1 i (i 2, 3,L , n) 得
lim(
对矩阵A1用乘幂法得 uk
A-1u
k
，
1
因为A1 的计算
比较麻烦，而且往往不能保持矩阵A 的一些好性质
（如稀疏性），因此，反幂法在实际计算时以求解
方程组
Auk
u
k
，代替迭代
1
uk
A-1uk1求得uk，每
迭代一次要解一线性方程组。 由于矩阵在迭代过
程中不变，故可对A 先进行三角分解，每次迭代只 要解两个三角形方程组。
且
2 p 2 n
2 n
2 n 2
1 p 21 2 n 1 n 1 2 1 n 1
因此，用原点平移法求1可使收敛速度加快。
三、反幂法
反幂法是计算矩阵按模最小的特征值及特征向 量的方法，也是修正特征值、求相应特征向量的最 有效的方法。
0
0.226
0.975
做正交相似变换后得到
3.366
A3 =R2 AR2T
0.0735
0.317
0.0735 1.780
0
0.317
0
1.145
雅可比方法是一个迭代过程，它生成的是一个矩阵的
序列 Ak，当k越大时Ak就越接近于对角矩阵，从而
矩阵特征值快速求法
矩阵特征值快速求法矩阵特征值是矩阵分析中十分重要的概念。
它在物理、工程、数学等许多领域都有着广泛的应用。
矩阵特征值是指矩阵运动时特殊的运动状态，是一种宏观量度矩阵运动的指标。
求解矩阵特征值是一项复杂的任务，通常需要使用高级算法来完成。
本文将介绍几种常用的求解矩阵特征值的算法，其中包括幂法、反幂法、QR算法、分裂Broyden算法等。
一、幂法幂法是求解矩阵特征值的一种基础算法，其基本思想是通过迭代来逐步逼近矩阵的最大特征值。
幂法的核心公式如下：x_(k+1)=A*x_k/||A*x_k||其中，x_k表示第k次迭代中得到的特征向量，A表示原始矩阵。
幂法通过不断的迭代来逼近A的最大特征值，当迭代次数趋近于无限大时，得到的特征向量就是A的最大特征值所对应的特征向量。
幂法的运算量较小，适用于比较简单的矩阵。
反幂法与幂法类似，不同之处在于每次迭代时采用的是A的逆矩阵来进行计算。
其核心公式如下：x_(k+1)=(A-λI)^(-1)*x_k其中，λ表示要求解的特征值。
反幂法能够求解非常接近于特征值λ的特征向量，并且对于奇异矩阵同样适用。
需要注意的是，在实际计算中，如果A-λI的秩不满，那么反幂法就无法使用。
三、QR算法1. 将原矩阵A进行QR分解，得到A=Q*R。
2. 计算A的近似特征矩阵A1=R*Q。
5. 重复步骤3-4，直到A的对角线元素全部趋近于所求特征值为止。
QR算法的计算量较大，但其具有收敛速度快、精度高等优点，广泛应用于科学计算中。
四、分裂Broyden算法分裂Broyden算法是QR算法的一种改进算法，其基本思想是将矩阵分解成上下三角形式，然后再对其进行QR分解，以减少QR算法中的乘法运算量。
具体实现过程如下：2. 构造一个倒数矩阵B=U^(-1)*L^(-1)。
4. 计算A的近似特征矩阵A1=Q^(-1)*L^(-1)*A*R^(-1)*U^(-1)*Q。
分裂Broyden算法的计算量较小，能够有效地解决QR算法中的乘法运算量过大的问题。
和法求矩阵最大特征值
和法求矩阵最大特征值是一种常用的算法，用于求解大规模矩阵的最大特征值及其对应的特征向量。
这种方法的基本思想是将矩阵分解为若干个子矩阵，并将这些子矩阵的和作为输入，通过求解子矩阵的最大特征值来逐步逼近原矩阵的最大特征值。
具体来说，和法求矩阵最大特征值的基本步骤如下：
1. 将原矩阵分解为若干个子矩阵，通常选择主对角线上的子矩阵作为基本子矩阵。
2. 对于每个基本子矩阵，使用适当的算法（如Jacobi方法、SOR方法等）求解其最大特征值。
3. 将所有基本子矩阵的最大特征值相加得到一个近似值，该值即为原矩阵的最大特征值的近似值。
4. 重复步骤2和3，直到达到预定的精度要求或达到最大迭代次数。
在求解过程中，需要注意以下几点：
* 算法的收敛性：和法求矩阵最大特征值算法需要保证收敛到真实解，因此需要选择合适的算法和参数设置。
* 算法的稳定性：算法需要保证稳定运行，避免出现数值不稳定的情况。
* 矩阵分解的精度：分解的子矩阵大小和数量会影响到求解的精度和速度，需要根据实际情况进行选择。
总体来说，和法求矩阵最大特征值是一种常用的算法，适用于大规模矩阵的特征值求解问题。
通过将矩阵分解为若干个子矩阵并逐步逼近原矩阵的最大特征值，可以获得相对准确的解。
但是，该方法需要较长时间和计算资源，因此在处理大规模问题时需要权衡精度和效率。
矩阵特征值的计算
物理、力学和工程技术中的许多问题在数学上都归结为求矩 阵的特征值和特征向量问题。
� 计算方阵 A 的特征值，就是求特征多项式方程：
| A − λI |= 0 即 λn + p1λn−1 + ⋅ ⋅ ⋅ + pn−1λ + pn = 0
的根。求出特征值 λ 后，再求相应的齐次线性方程组：
(13)
为了防止溢出，计算公式为
⎧ Ay k = xk −1
⎪ ⎨
m
k
=
max(
yk )
( k = 1, 2, ⋅ ⋅⋅)
⎪ ⎩
x
k
=
yk
/ mk
(14)
相应地取
⎧ ⎪
λ
n
⎨
≈
1 mk
⎪⎩ v n ≈ y k ( 或 x k )
(15)
9
(13)式中方程组有相同的系数矩阵 A ，为了节省工作量，可先对
11
11
≤ ≤ ⋅⋅⋅ ≤
<
λ1 λ2
λn −1
λn
对应的特征向量仍然为 v1, v2 ,⋅⋅⋅, vn 。因此，计算矩阵 A 的按模
最小的特征值，就是计算 A−1 的按模最大的特征值。
� 反幂法的基本思想：把幂法用到 A−1 上。
任取一个非零的初始向量 x0 ，由矩阵 A−1 构造向量序列：
xk = A−1xk−1 , k = 1, 2, ⋅⋅⋅
如果 p 是矩阵 A 的特征值 λi 的一个近似值，且
| λi − p |<| λ j − p | , i ≠ j
1 则 λ i − p 是矩阵 ( A − pI )−1 的按模最大的特征值。因此，当给
求矩阵特征值算法及程序
求矩阵特征值算法及程序简介1.幂法1、幂法规范化算法（1）输入矩阵A 、初始向量)0(μ，误差eps ； （2）1⇐k ；（3）计算)1()(-⇐k k A V μ；（4）)max (,)max ()1(1)(--⇐⇐k k k k V m V m ； （5）k k k m V /)()(⇐μ；（6）如果eps m m k k <--1,则显示特征值1λ和对应的特征向量)1(x),终止；（7）1+⇐k k ,转（3） 注：如上算法中的符号)max(V 表示取向量V 中绝对值最大的分量。
本算法使用了数据规范化处理技术以防止计算过程中出现益出错误。
2、规范化幂法程序Clear[a,u,x];a=Input["系数矩阵A="];u=Input["初始迭代向量u(0)="];n=Length[u];eps=Input["误差精度eps ="];nmax=Input["迭代允许最大次数nmax="];fmax[x_]:=Module[{m=0,m1,m2},Do[m1=Abs[x[[k]]];If[m1>m,m2=x[[k]];m=m1],{k,1,Length[x]}];m2]v=a.u;m0=fmax[u];m1=fmax[v];t=Abs[m1-m0]//N;k=0;While[t>eps&&k<nmax,u=v/m1;v=a.u;k=k+1;m0=m1;m1=fmax[v];t=Abs[m1-m0]//N;Print["k=",k," 特征值=",N[m1,10]," 误差=",N[t,10]];Print[" 特征向量=",N[u,10]]];If[k ≥nmax,Print["迭代超限"]]说明：本程序用于求矩阵A 按模最大的特征值及其相应特征向量。
特征值分解 算法流程
特征值分解算法流程特征值分解算法流程简介特征值分解是一种重要的矩阵分解方法，广泛应用于数学、物理、工程和计算机科学等领域。
它能将一个方阵分解为特征向量和特征值的乘积形式，揭示出矩阵的结构和性质。
算法流程特征值分解的算法流程可以概括为以下几个步骤：1.输入方阵：将待分解的方阵记为A。
2.求解特征方程：通过求解特征方程det(A - λI) =0，得到特征值λ。
其中，A是输入方阵，λ是待求的特征值，I是单位矩阵。
3.求解特征向量：对于每一个特征值，利用A - λI，解出特征向量。
4.归一化特征向量：对每一个特征向量进行归一化处理，使其模长为1。
5.组成特征值矩阵和特征向量矩阵：将求得的特征值按对角线排列形成特征值矩阵Λ，将求得的特征向量按列形成特征向量矩阵P。
6.完成分解：得到特征值矩阵Λ和特征向量矩阵P后，即完成了对输入方阵A的特征值分解。
示例为了更好地理解上述算法流程，以下给出一个简单的示例：输入方阵A为：[ 4 2 ][ 1 3 ]1.求解特征方程：det(A - λI) = 0=> det([4-λ 2 | 1 0 ])[ 1 3-λ | 0 1 ]= (4-λ)(3-λ) - 2·1 = λ^2 - 7λ + 10 = (λ - 2)(λ - 5) = 0得到特征值λ1 = 2 和λ2 = 5。
2.求解特征向量：对于λ1 = 2：[2-2 2 | 1 0 ][ 1 3-2 | 0 1 ]解得特征向量v1 = [1, -1]。
对于λ2 = 5：[4-5 2 | 1 0 ][ 1 3-5 | 0 1 ]解得特征向量v2 = [1, 1]。
3.归一化特征向量：归一化v1 = [1/√2, -1/√2]，归一化v2 = [1/√2, 1/√2]。
4.组成特征值矩阵和特征向量矩阵：特征值矩阵Λ为：[ 2 0 ][ 0 5 ]特征向量矩阵P为：[ 1/√2 1/√2 ][ -1/√2 1/√2 ]5.完成分解：得到特征值矩阵Λ和特征向量矩阵P后，即完成了对输入方阵A的特征值分解。
QR迭代法求矩阵特征值
本文用QR迭代法求解矩阵A的特征值:第一步先用豪斯荷尔德变换将矩阵A化为上海 森伯格矩阵AH ，第二步再对AH 进行QR迭代(使用吉文斯变换)，当迭代满足精度要求时 输出n个特征值。
1
问题描述
对于给定的矩阵A(存在n个不相等的实特征值),采用QR方法求其n个特征值。本文中的矩 阵A由文件”gr9 009 00c rg.mm”给出。并以.mm格式描述。
2 AHi,i 2 + AHi +1,i
si = 于是，第i行和第i+1行的元素变换为：
AHi+1,i
2 + AH 2 AHi,i i+1,i
AHi,j = ci ∗ AHi,j + si ∗ AHi+1,j AHi+1,j = ci ∗ AHi+1,j − si ∗ AHi,j j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效实现QR分解，并得到上三角矩阵R。 在第二个过程中，R矩阵右乘上Q矩阵，即： AHnew = RQ = RG(1, 2; θ1 )T · · · · · · G(n − 2, n − 1; θn−1 )T G(n − 1, n; θn−1 )T 等效于依次对R矩阵的第i列和第i+1列进行列变换。变换后的元素为： Rj,i = ci ∗ Rj,i + si ∗ Rj,i+1 Rj,i+1 = ci ∗ Rj,i+1 − si ∗ Rj,i j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效地实现AHnew = RQ过程。 至此，完成了一次QR迭代。 本函数(AHnew=QR-iteration-once(AH,n))的流程图如图三所示。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
求矩阵特征值算法及程序简介1.幂法1、幂法规范化算法(1)输入矩阵A、初始向量( 0)，误差eps；(2) k 1；(3)计算V(k)A(k 1)；(4)m k max(V(k)) ,m k1max( V ( k 1))；(5) (k)V(k)/m k；(6)如果m k m k 1eps,则显示特征值1和对应的特征向量x(1) ),终止；(7)k k 1, 转(3)注：如上算法中的符号max(V )表示取向量V 中绝对值最大的分量。
本算法使用了数据规范化处理技术以防止计算过程中出现益出错误。
2、规范化幂法程序Clear[a,u,x];a=Input[" 系数矩阵A="];u=Input[" 初始迭代向量u(0)="]; n=Length[u];eps=Input[" 误差精度eps ="]; nmax=Input[" 迭代允许最大次数nmax="];fmax[x_]:=Module[{m=0,m1,m2},Do[m1=Abs[x[[k]]];If[m1>m,m2=x[[k]];m=m1], {k,1,Length[x]}];m2]v=a.u;m0=fmax[u];m1=fmax[v]; t=Abs[m1-m0]//N;k=0;While[t>eps&&k<nmax,u=v/m1;v=a.u;k=k+1;m0=m1;m1=fmax[v];t=Abs[m1-m0]//N;Print["k=",k," 特征值=",N[m1,10]," 误差=",N[t,10]];Print[" 特征向量=",N[u,10]]];If[k nmax,Print[" 迭代超限"]]说明：本程序用于求矩阵A按模最大的特征值及其相应特征向量。
程序执行后，先通过键盘输入矩阵A 、迭代初值向量(0)、精度控制eps和迭代允许最大次数n max ，程序即可给出每次迭代的次数和对应的迭代特征值、特征向量及误差序列，它们都按10 位有效数输出。
其中最后输出的结果即为所求的特征值和特征向量序列。
如果迭代超出nmax 次还没有求出满足精度的根则输出迭代超限提示，此时可以根据输出序列判别收敛情况。
程序中变量说明a: 存放矩阵A ；u: 初始向量(0)和迭代过程中的向量(k)及所求特征向量；v: 存放迭代过程中的向量V (k )；m1:存放所求特征值和迭代过程中的近似特征值；nmax:存放迭代允许的最大次数；eps: 存放误差精度；fmax[x]: 给出向量x 中绝对值最大的分量；k: 记录迭代次数；t1: 临时变量；注：迭代最大次数可以修改为其他数字。
3、例题与实验133 6 135例1. 用幂法求矩阵A 44 5 4688 6 90的按模最大的特征值及其相应特征向量，要求误差eps 10解：执行幂法程序后在输入的四个窗口中按提示分别输入:{{133,6,135},{44,5,46},{-88,-6,-90}} 、{1,1,1} 、0.0001 、20每次输入后用鼠标点击窗口的“ OK”按扭，得如下输出结果：此结果说明迭代6 次，求得误差为0.0000101442的按模最大的特征值为44.99999952,及其对应的一个特征向量: {1.000000000,0.3333333371,-0.6666666704}2.反幂法1、反幂法规范化算法(1)输入矩阵A 、初始向量(0)，误差eps；(2) k 1；(3)计算AV(k) (k 1)求出解V(k)；(4)m k max(V (k)) ,m k 1max(V (k 1))；(5) (k)V(k)/m k；(6)如果m k m k 1eps,则显示特征值1和对应的特征向量x(1) ),终止；7)k k 1, 转(3) 注：如上算法中解方程AV(k) (k 1)可以使用Dololittle 分解法。
本算法使用了数据规范化处理技术以防止计算过程中出现益出错误。
2、规范化反幂法程序Clear[a,u,x];a=Input[" 系数矩阵A="];u=Input[" 初始迭代向量u(0)="];n=Length[u];eps=Input[" 误差精度eps ="];nmax=Input[" 迭代允许最大次数nmax="]; fmax[x_]:=Module[{m=0,m1,m2},Do[m1=Abs[x[[k]]];If[m1>m,m2=x[[k]];m=m1],{k,1,Length[x]}];m2];v=a.u;a1=Inverse[a];m0=fmax[u];m1=fmax[v];t=Abs[m1-m0]//N;k=0;While[t>eps&&k<nmax,u=v/m1;v=a1.u;k=k+1;m0=m1;m1=fmax[v];t=Abs[m1-m0]//N;t1=Abs[1/m1-1/m0]//N;Print["k=",k," 特征值=",N[1/m1,10]," 误差=",N[t1,10]];Print[" 特征向量=",N[u,10]]];If[k nmax,Print[" 迭代超限"]]说明：本程序用于求矩阵A按模最小的特征值及其相应特征向量。
程序执行后，先通过键盘输入矩阵A 、迭代初值向量(0)、精度控制eps和迭代允许最大次数n max ，程序即可给出每次迭代的次数和对应的迭代特征值、特征向量及误差序列，它们都按10 位有效数输出。
其中最后输出的结果即为所求的特征值和特征向量序。
如果迭代超出nmax 次还没有求出满足精度的根则输出迭代超限提示，此时可以根据输出序列判别收敛情况。
程序中变量说明a：存放矩阵Au：初始向量(0)和迭代过程中的向量(k)及所求特征向量v: 存放迭代过程中的向量V(k )a1：存放逆矩阵A 1m1: 存放所求特征值和迭代过程中的近似特征值nmax:存放迭代允许的最大次数eps: 存放误差精度fmax[x]给出向量x 中绝对值最大k: 记录迭代次数t1: 临时变量注：迭代最大次数可以修改为其他数字3、例题与实验2 2 3例3. 用反幂法求矩阵A 1 1 11 3 1的按模最小的特征值及其相应特征向量，要求误差eps 10 解：执行幂法程序后在输入的四个窗口中按提示分别输入:{{2.,-2.,3.},{1,1.,1},{1.,3,-1}} 、{1,0,1} 、0.00001 、100每次输入后用鼠标点击窗口的“ OK”按扭，得如下输出结果：注意到本题按模最小的特征值为1，因此求解效果较满意。
3.J acobi 方法1、Jacobi 旋转法算法1、输入矩阵A ，误差；2、For k 1,2, ；(2.1 )选取a(p k q1)max ai(j k 1)记录p, q2a(p k q1)(2.2 )由tan 2 (k21)a pq(k1), 确定旋转角，获得旋转矩阵J p,q,a pp a qq 42.3 ) a (p k j) pj a (p k j1) pj cos a q(k j1)sin , a(jp k) aqj jp (k) pjj p,q；2.4 )a(k) a qj a(k a pj 1)sin(k 1) a qj c os ,a(jp k)a(k) a qj jp,q；2.5 )a i(j k) ij a i(j k 1)iji, j p, q2.6 )a(k) a pp a(k 1)a ppcos2a q(k q1)sin22a(p kq 1)sin cos2.7 ) (k) a qq (k 1)a pp2 sin a(qq k 1)cos2(k2a pq 1)sin cos2.8 ) 计算E(A k ) ( k) 2 a ij2.9 )如果E(A k ) ，输出对角矩阵D diag( 1, 2注：如上算法中A k (a i(j k))，A0 (a i(j0)) (a ij )n )和特征向量矩阵J ，停止ij2、Jacobi 算法程序Clear[a,bb]; a=Input[" 矩阵A="]; n=Input[" 矩阵阶数n="];eps=Input[" 误差精度eps ="];nmax=Input[" 迭代允许最大次数nmax="]; k=0;bb=IdentityMatrix[n]; ea=Sum[a[[i,j]]^2,{i,1,n},{j,1,n}]-Sum[a[[i,i]]^2,{i,1,n}]//N;While[ea>eps&&k<nmax,m=0;Print[" 迭代次数k=",k];Do[If[Abs[a[[i,j]]]>m,m=Abs[a[[i,j]]];p=i;q=j],{i,1,n},{j,i+1,n}]; mu=a[[p,p]]-a[[q,q]];If[mu 0,thi=Pi/4,thi=ArcTan[2*a[[p,q]]/mu]/2]; s=Sin[thi]//N;c=Sqrt[1-s^2]; a1=bb[[p]];bb[[p]]=c*bb[[p]]+s*bb[[q]]; bb[[q]]=-s*a1+c*bb[[q]]; pp=a[[p,p]]*c*c+a[[q,q]]*s*s+2a[[p,q]]*s*c;qq=a[[p,p]]*s*s+a[[q,q]]*c*c-2a[[p,q]]*s*c;Do[a1=a[[p,j]]; a[[p,j]]=c*a[[p,j]]+s*a[[q,j]];a[[j,p]]=a[[p,j]]; a[[q,j]]=c*a[[q,j]]-s*a1;a[[j,q]]=a[[q,j]],{j,1,n}];a[[p,p]]=pp; a[[q,q]]=qq;a[[p,q]]=0; a[[q,p]]=0; ea=Sum[a[[i,j]]^2,{i,1,n},{j,1,n}]-Sum[a[[i,i]]^2,{i,1,n}]//N; k=k+1;Print[" 误差=",ea];Print[" 相似矩阵A="];Print[MatrixForm[a]];Print[" 特征向量J"]; Print[MatrixForm[Transpose[bb]]]];说明本程序用于求对称矩阵A 的所有特征值及其相应特征向量。