大规模矩阵问题的Krylov子空间方法综述
Krylov子空间方法II
max {q (λi )2 }
2 yi λi
min
max {q (λi ) } y Λy
5/76
= = 即
q ∈Pk , q (0)=1 1≤i≤n q ∈Pk , q (0)=1 1≤i≤n
min
max {q (λi )2 } ϵ0 Aϵ0 max {q (λi )2 } ∥ϵ0 ∥2 A,
⊺
min
11/76
(4.6)
又 Λ 是对角矩阵, 所以 ∥q (Λ)∥2 = max |q (λi )|.
1≤i≤n
设 x(k) 是由 GMRES 方法得到的近似解. 由 GMRES 方法的最优性可 知, x(k) 极小化残量的 2 范数. 因此, ∥b − Ax(k) ∥2 = = ≤
x∈x(0) +Kk (A,r0 ) q ∈Pk , q (0)=1 q ∈Pk , q (0)=1
1+ε δ
10/76
5.2 GMRES 方法的收敛性
正规矩阵情形 设 A 是正规矩阵, 即 A = U ΛU ∗ , 其中 Λ = diag(λ1 , λ2 , . . . , λn ) 的对角线元素 λi ∈ C 为 A 的特征值. 设 x ∈ x(0) + Kk (A, r0 ), 则存在多项式 p(t) ∈ Pk−1 使得 x = x(0) + p(A)r0 . 于是 b − Ax = b − Ax(0) − Ap(A)r0 = (I − Ap(A))r0 ≜ q (A)r0 , 其中 q (t) = 1 − t p(t) ∈ Pk 满足 q (0) = 1. 直接计算可知 ∥b − Ax∥2 = ∥q (A)r0 ∥2 = ∥U q (Λ)U ∗ r0 ∥2 ≤ ∥U ∥2 ∥U ∗ ∥2 ∥q (Λ)∥2 ∥r0 ∥2 = ∥q (Λ)∥2 ∥r0 ∥2 .
gmres算法范文
gmres算法范文GMRES（Generalized Minimal RESidual）算法是一种用于求解稀疏线性方程组的迭代算法。
它可以用于求解大规模稀疏矩阵的线性方程组问题，特别适用于非对称且非正定矩阵的情况。
GMRES算法能够通过矩阵向量乘法来逐步逼近线性方程组的解，从而在求解过程中保持向量的稀疏性，节约了计算和存储资源。
GMRES算法的核心思想是基于Krylov子空间的最小化残差，通过在Krylov子空间中找到一个最优的近似解向量来逼近线性方程组的解。
Krylov子空间是由矩阵A和初始向量b生成的线性空间，通过不断迭代计算可以得到Krylov子空间的一组正交基。
GMRES算法的具体步骤如下：1. 初始化：给定一个初始解向量x0和初始残差r0=b-Ax0，将正交化后的r0作为初始Krylov子空间的基向量v12. Arnoldi迭代：对于k=1到m，进行以下步骤：a. 计算w=Avk；b. 通过Gram-Schmidt过程对w与之前的基向量进行正交化，得到新的正交基向量v_k+1；c. 构造大小为(k+1)×k的Hessenberg矩阵H，其中H=Q_k^T*A*Q_k，Q_k是由正交基向量v1,v2,...,vk构成的正交矩阵；d.使用QR分解求解H的最小二乘问题，得到近似解向量y。
3.更新解向量：更新解向量为x_k=x_0+Q_k*y。
4.检测终止条件：如果达到了预定的收敛条件或者迭代次数达到了最大限制，则结束迭代；否则返回步骤2GMRES算法的核心在于利用Krylov子空间的正交基向量来构造Hessenberg矩阵，并通过最小二乘法求解近似解向量。
通过在每一步迭代中更新解向量，可以逐步逼近线性方程组的解。
当算法能够达到预定的收敛条件时，解向量可以近似地满足线性方程组。
GMRES算法的优点是可以求解大规模稀疏矩阵的线性方程组，并且能够保持向量的稀疏性。
它还可以通过调整收敛条件来控制算法的精度和计算资源的消耗。
Krylov子空间方法
由于 x ˜ ∈ x(0) + K, 因此存在向量 y ∈ Rm 使得 x ˜ = x(0) + V y 由正交性条件 (4.4) 可知 r0 − AV y ⊥ wi , i = 1, 2, . . . , m , 即 W ⊺ AV y = W ⊺ r0 .
x ˆ≜x ˜ − x(0) = V y
9/115
Arnoldi 过程: 计算 Km 的一组正交基
算法 2.1 基于 Gram-Schmidt 正交化的 Arnoldi 过程
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13:
给定非零向量 r, 计算 v1 = r/∥r∥2 for j = 1, 2, . . . , m − 1 do wj = Avj for i = 1, 2, . . . , j do hij = (wj , vi ) end for j ∑ wj = wj − hij vi hj +1,j = ∥wj ∥2 if hj +1,j = 0 then break end if vj +1 = wj /hj +1,j end for
若给定初值 x(0) ∈ Rn , 则改用仿射空间 x(0) + K, 即 find x ˜ ∈ x(0) + K such that b − Ax ˜ ⊥ L. (4.3)
好的初值一般都包含有价值 的信息
事实上, 如果将 x ˜ 写成: x ˜ = x(0) + x ˆ, 其中 x ˆ ∈ K, 则 (4.3) 就等价于 find x ˆ∈K such that r0 − Ax ˆ ⊥ L, (4.4)
定解条件
r = b − Ax ˜⊥L 其中 x ˜ 是近似解, L 是另一个 m 维子空间. 不同的 L 对应不同的投影方法 当 L = K 时, 我们称为 正交投影法 , 否则称为 斜投影法
Krylov子空间迭代法
采用IOM后，仍然需要存储v(1), v(2), …v(m),因为在第(vi)步 中仍然需要这些向量. 解决这个问题可以考虑采用H的LU分解，通过自身分解的迭代更新以减少每 一步的存储量 使xm的更新依赖于xm-1,
14
Arnoldi方法-DIOM
lower bidiagonal
banded upper triangular
15
Arnoldi方法-DIOM
16
Arnoldi方法-DIOM
17
Thanks for your time !
18
得到基于Galerkin原 理构成的算法
5
Arnoldi方法-基本算法
6
Arnoldi方法-基本算法
7
Arnoldi方法-MGS
8
Arnoldi方法-HO
9
Arnoldi方法-FOM
10
Arnoldi方法-FOM
11
Arnoldi方法-FOM（m）
12
Arnoldi方法-IOM
13
Arnoldi方法-DIOM
Krylov子空间方法
March 23, 2016
内
• Arnoldi算法
– Arnoldi过程 – Gram-Schmidt Arnoldi – HouseHolder Arnoldi
容
• 子空间和Krylov子空间
• FOM
– IOM – DIOM
2
子空间
• 空间
– 集合，元素都是向量 – 线性空间（向量空间）
• 线性空间（交换律，结合律，幺元性，零元性，可 逆性，数乘分配律等）
• 子空间
– 线性空间的非空子集
Krylov子空间、优化问题与共轭梯度法
Krylov 子空间、优化问题与共轭梯度法自动化 富晓鹏工程实践中经常需要求解大型线性系统KU=F 。
在很多情况下矩阵K 是非常稀疏的，比如来自偏微分方程的离散化等，此时矩阵中每行仅有较少的非零元素。
面临这样的问题，我们首先面对的问题是，应该采用直接消元法还是迭代方法。
对前者来说，为充分利用系数特性，节点重编号是重要的；而对后者来说，适当的预处理是关键。
本文将重点放在后一类方法中的一种进行介绍与分析，即共轭梯度法。
共轭梯度法适用于矩阵K 为对称阵的情况，算法本身简洁高效，且与一些其他的数学理论、概念相紧密联系，本文分析了共轭梯度法与Krylov 子空间，以及优化问题之间隐含的联系，并简要给出算法框架。
1. 线性方程组迭代解法与Krylov 子空间我们考虑迭代法求解线性方程组Ax=b 。
假定未采用预处理矩阵P ，或P 矩阵已经隐含在A 与b 中。
迭代法求解格式如下：1()k k P x P A x b +⋅=-⋅+ (1)为说明问题，我们考虑简单的迭代格式P=I ，并且x 1=b 。
则迭代的最初几步为：2()2x I A b b b Ab =-+=- (2)232()33x I A x b b Ab A b =-+=-+ (3) …由上面几个式子可得，以上迭代格式第j 步的解x j 是b ，Ab ，…，A j -1b 的线性组合。
当A 矩阵稀疏时，这些向量可以采用矩阵向量乘法的稀疏技巧很快得到。
以上发现自然与Krylov 子空间的概念相联系起来。
Krylov 矩阵： K j = [b Ab A 2b … A j -1b]Krylov 子空间：K j = b ，Ab ，…，A j -1b 的所有线性组合Krylov 命名了向量b ，Ab ，…，A j -1b 的全部线性组合构成的子空间，并认为在这一子空间中，有比上例中特定元素更与线性方程组的解相接近的元素。
共轭梯度法就是在这一子空间中，每一步迭代都依照某种标准寻求最优元素的线性方程组解法。
(完整版)Krylov子空间迭代法
February 10, 2020
内容
• 子空间和Krylov子空间
• Arnoldi算法
– Arnoldi过程 – Gram-Schmidt Arnoldi – HouseHolder Arnoldi
• FOM
– IOM – DIOM
2
子空间
• 空间
– 集合，元素都是向量 – 线性空间（向量空间）
根据Cayley-Hamilton定理有
������������ + ������������−1������������−1+. . . +������1������1 + ������0������0 = 0
即
VP= -������������ 其中������ = [������0, ������1, . . . , ������������−1ሿ，������ = ������0, ������1, . . . , ������������−1 ������ Krylov子空间： ������������(������, ������) = ������������������������{������, ������������, . . . , ������������−1������ሽ Krylov矩阵： ������������(������, ������) = [������, ������������, . . . , ������������−1������ሿ
• 线性空间（交换律，结合律，幺元性，零元性，可 逆性，数乘分配律等）
• 子空间
– 线性空间的非空子集
• 包含零元素，并且满足加法和乘法的封闭性
– 扩张（符合记作span）
大规模矩阵问题的Krylov
定理
回顾
如果记以v(1),v(2),……v(m)为列构成的矩阵为 Vm,由hij定义的(m+1)×m阶上Hessenberg矩 阵为 H m,删除最后一行得到的矩阵为Hm,则
AVm Vm Hm w (e ) Vm1 H m
( m) ( m) T
在Arnoldi算法中，可能有较大舍入误差，改写
K(L)
K X(0) X(1)
投影方法
L X(0)
X(1)
• 线性方程组的投影方法 方程组Ax=b,A是n×n的矩阵. 给定初始x(0),在m维空间K(右子空间）中寻 找x的近似解x(1) 满足残向量r=b-Ax(1) 与m维空 间L(左子空间）正交，即 b-Ax(1) ⊥L 此条件称为Petrov-Galerkin条件 . 当空间K=L时，称相应的投影法为正交投 影法，否则称为斜交投影法.
Krylov子空间方法解线性方程组
• 误差投影型方法 取L=K时的正交投影法 １）非对称矩阵的FOM方法 2) 非对称矩阵的IOM方法和DIOM方法
所谓不完全正交化方法(IOM),是指在正交化过 程中,v(j+1)仅与最近k个v(j-k+1),…v(j-1),v(j)正交,这 样做虽然破坏的正交性,但是降低了计算量.当 然k选得越小，对每个j对应的计算量也越小,但 可能要选更大的m才能取得满足精度要求的近 似解. IOM算法仅仅是把FOM算法的第三步改为 (iii’)对i=max(1,j-k+1),…,j, 计算hij与w(j).
• Arnoldi方法 基于Gram-Schmit正交化方法 首先，选取一个Euclid范数为1的向量v(1),对 v K m ( A, v) ,通常可取 v 在已知 ( v, v ) v(1),v(2),……v(j)的情况下,不妨设 v(1),v(2),……v(j),Av(j)线性无关(否则构造完毕), 则可求出与每个都正交的向量
krylov子空间算法
Krylov 子空间的定义:定义:令N R υ∈,由1m A υυυ-L ，，，A 所生成的子空间称之为由υ与A 所生成的m 维Krylov 子空间,并记(),m K A v 。
主要思想就是为各迭代步递归地造残差向量,即第n 步的残差向量()n r 通过系数矩阵A 的某个多项式与第一个残差向量()0r 相乘得到。
即()()()0n r p A r =。
但要注意,迭代多项式的选取应该使所构造的残差向量在某种内积意义下相互正交,从而保证某种极小性(极小残差性),达到快速收敛的目的。
Krylov 子空间方法具有两个特征:1、极小残差性,以保证收敛速度快。
2、每一迭代的计算量与存储量较少,以保证计算的高效性。
投影方法线性方程组的投影方法方程组Ax b =,A 就是n n ⨯的矩阵。
给定初始()0x ,在m 维空间K(右子空间)中寻找x 的近似解()1x 满足残向量()1r b Ax =-与m 维空间L(左子空间)正交,即()1b Ax L -⊥,此条件称为Petrov-Galerkin 条件。
当空间K=L 时,称相应的投影法为正交投影法,否则称为斜交投影法、投影方法的最优性:１、 (误差投影)设A 为对称正定矩阵,()0x 为初始近似解,且K=L,则()1x 为采用投影方法得到的新近似解的充要条件就是()()()()01min z x Kx z ϕϕ∈+=其中,()()()12,z A x z x z ϕ=--2.(残量投影)设A 为任意方阵,()0x 为初始近似解,且L AK =,则()1x 为采用投影方法得到的新近似解的充要条件就是()()()()01min z x Kx z ψψ∈+=其中()()122,z b Az b Az b Az ψ=-=--矩阵特征值的投影方法对于特征值问题Ax x λ=,其中A 就是n ×n 的矩阵,斜交投影法就是在m 维右子空间K 中寻找i x 与复数i λ满足i i i Ax x L λ-⊥,其中L 为m 维左子空间、当L=K 时,称此投影方法为正交投影法、 误差投影型方法: 取L=K 的正交投影法非对称矩阵的FOM 方法(完全正交法) 对称矩阵的IOM 方法与DIOM 方法 对称矩阵的Lanczos 方法 对称正定矩阵的CG 方法 残量投影型方法: 取L=AK 时的斜交投影法 GMERS 方法(广义最小残量法)重启型GMERS 方法、QGMERS 、DGMERSArnoldi 方法标准正交基方法:Arnoldi 方法就是求解非对称矩阵的一种正交投影方法。
《2024年解线性方程组的VRP-GMRES(m)迭代法》范文
《解线性方程组的VRP-GMRES(m)迭代法》篇一一、引言随着科技的发展，线性方程组的求解在众多领域如物理学、工程学和计算机科学等均具有重要意义。
针对大型或复杂线性方程组的求解，传统方法常常无法满足需求。
为此，一种名为VRP-GMRES(m)的迭代法逐渐被广大研究者所接受，它以高效、稳定的特性在众多领域中取得了广泛应用。
本文将详细介绍VRP-GMRES(m)迭代法的基本原理、应用场景及其实施步骤。
二、VRP-GMRES(m)迭代法的基本原理VRP-GMRES(m)迭代法是一种用于解线性方程组的迭代算法。
该算法结合了GMRES（广义最小残差法）和变体投影方法（Variational Refinement Procedure），既保证了收敛性，又提高了计算效率。
该算法主要应用于求解大规模、非对称、复杂矩阵的线性方程组。
三、VRP-GMRES(m)迭代法的应用场景1. 物理学：在求解复杂的物理系统如流体动力学、电磁场等领域，需要求解大规模的线性方程组，此时，VRP-GMRES(m)迭代法可以发挥其优势。
2. 工程学：在结构分析、优化设计等领域，需要处理大量的线性方程组，VRP-GMRES(m)迭代法能够快速、准确地求解这些问题。
3. 计算机科学：在图像处理、机器学习等领域，经常需要求解大规模的线性方程组。
VRP-GMRES(m)迭代法可以有效地处理这些问题。
四、VRP-GMRES(m)迭代法的实施步骤1. 初始化：设定初始向量和初始残差向量，并计算初始解向量。
2. 构建Krylov子空间：通过与系数矩阵相乘的方式，生成新的向量并扩充Krylov子空间。
3. 正交化过程：利用GMRES算法的正交化过程，将Krylov 子空间中的向量进行正交化处理。
4. 最小二乘问题求解：在正交化后的Krylov子空间中，求解最小二乘问题以得到下一步的搜索方向。
5. 投影与变体：利用变体投影方法对解进行投影和变体处理，以改善解的精度和收敛性。
第四讲 Krylov 子空间方法
∑j Avj = hj+1,j vj+1 + hij vi
∑j wj = wj − hij vi
i=1
vj+1 = wj /hj+1,j
i=1
=
[v1,
v2,
.
.
.
,
vj+1]
h1j
h2j ...
hj+1,j
=
[v1,
.
.
.
,
vj+1,
vj+2,
.
.
.
,
vm+1]
hjh+0......11j,j
6: end for
∑j 7: wj = wj − hij vi
i=1
8:
hj+1,j = ∥wj ∥2
9: if hj+1,j = 0 then
10:
break
11: end if
12:
vj+1 = wj /hj+1,j
13: end for
做计算时一定要注意任何可 能存在的 “中断”.
10/106
0 0 0 ···
h1,m−1 h2,m−1 h3,m−1 h4,m−1
... hm,m−1
0
h1,m
h2,m h2,m h4,m
... hm,m
∈ R(m+1)×m
hm+1,m
这里 hij 是由 Arnoldi 过程所定义的.
wj = Avj
for i = 1, 2, . . . , j do
hij = (wj , vi)
= Vm
+
Hm+1,m
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
Kyrlov子空间
K m ( A, r ) span{r , Ar , A r ，…，A
(0) (0) (0) 2 (0) m - 1 (0)
r }
定义为m维Krylov子空间
dim K m ( A, v) min{m, }
μ为v的次数，即使得q(A)v＝０的非零首一多 项式的最低次数
Krylov子空间的标准正交基
y (W AV ) W r 从而得到迭代公式
T
1
T (0)
x(1) x(0) V (W T AV )1W T r (0)
不难看出，当采用上述FOM算法时，需要 存储所有的vi，(i=1,2,…m)，当m增大时，存 储量以O(mn)量级增大.而FOM计算量是O(m2 n).可见其代价十分高昂.因此我们考虑重启 的 FOM算法
(1) zx K
其中
( z) b Az 2 (b Az, b Az)1/ 2
• 矩阵特征值的投影方法
对于特征值问题Ax=λx,其中A是n×n的矩阵， 斜交投影法是在m维右子空间K中寻找xi和复数λi 满足Axi- λixi ⊥L,其中L为m维左子空间.当L=K时， 称此投影方法为正交投影法.
( m) ( m) T
在Arnoldi算法中，可能有较大舍入误差，改写
hij ( Av h1 j v h2 j v hi1, j v
( j) (1) (2) (i 1)
V AVm Hm
T m
, v ), i 1,2, , j.
(i )
Krylov子空间方法解线性方程组
h j 1, j
得到一组标准正交基.若期间某个j使得hj+1,j=0, 则说明v的次数是j,且Kj是A的不变子空间.
•
定理
如果记以v(1),v(2),……v(m)为列构成的矩阵为 Vm,由hij定义的(m+1)×m阶上Hessenberg矩 阵为 H m ,删除最后一行得到的矩阵为Hm,则
AVm Vm Hm w (e ) Vm1 H m
• 解方程组的投影法的矩阵表示 设n×m阶矩阵V=[v(1), v(2), …v(m)]与W=[w(1), w(2), …w(m)] 的列分别构成 K 与 L 的一组基 . 记 z=x(1)-x(0),z=Vy,有 W T (r (0) AVy) 0 当 W T AV 非奇异时（通常情况成立，见定 理1.1），有
• 误差投影型方法 取L=K时的正交投影法 １）非对称矩阵的FOM方法
回顾
• 解方程组的投影法的矩阵表示 设n×m阶矩阵V=[v(1), v(2), …v(m)]与W=[w(1), w(2), …w(m)] 的列分别构成 K 与 L 的一组基 . 记 z=x(1)-x(0),z=Vy,有 W T (r (0) AVy) 0 当 W T AV 非奇异时（通常情况成立，见定 理1.1），有
K(L)
K X(0) X(1)
投影方法
L X(0)
X(1)
• 线性方程组的投影方法 方程组Ax=b,A是n×n的矩阵. 给定初始x(0),在m维空间K(右子空间）中寻 找x的近似解x(1) 满足残向量r=b-Ax(1) 与m维空 间L(左子空间）正交，即 b-Ax(1) ⊥L 此条件称为Petrov-Galerkin条件 . 当空间K=L时，称相应的投影法为正交投 影法，否则称为斜交投影法.
y (W AV ) W r 从而得到迭代公式
T
1
T (0)
x(1) x(0) V (W T AV )1W T r (0)
• 投影方法的最优性 １. (误差投影)设A为对称正定矩阵, x(0)为初始近 似解,且K=L,则x(1)为采用投影方法得到的新近 似解的充要条件是
( x(1) )
大规模矩阵问题的Krylov 子空间方法综述
Krylov子空间方法综述
• • • • • • 背景介绍 投影方法 Krylov子空间及其标准正交基 Krylov子空间方法求解线性方程组 Krylov子空间方法求解矩阵的特征值 研究热点和尚未解决的问题
背景
• 大规模线性方程组的求解 很多科学工程计算问题都转化为求解方程 组Ax=b.如偏微分方程组的差分格式，有限元 方法离散得到刚度矩阵. • 大规模矩阵特征值和特征向量的计算 工程计算领域十分常见。如量子物理中的 Kohn-Sham方程求解化为哈密顿矩阵某些关 键特征值对的计算.
(1) 1/ 2
w( j ) Av( j ) h1 j v(1) h2 j v(2)
hjj v( j )
(2.1)
而不难看出 hij ( Av( j ) , v(i ) ), i 1,2, , j ，再记 ( j) ( j ) 1/ 2 hj 1, j (w , w ) ,得到与v(1),v(2),……v(j)都 ( j) w 正交的向量 v( j 1) 重复此过程，即可
• Arnoldi方法 基于Gram-Schmit正交化方法 首先，选取一个Euclid范数为1的向量v(1),对 K m ( A, v) ,通常可取 v v 在已知 ( v, v ) v(1),v(2),……v(j)的情况下,不妨设 v(1),v(2),……v(j),Av(j)线性无关(否则构造完毕), 则可求出与每个都正交的向量
Krylov子空间方法解线性方程组
• 误差投影型方法 取L=K时的正交投影法 １）非对称矩阵的FOM方法 2) 非对称矩阵的IOM方法和DIOM方法
所谓不完全正交化方法(IOM),是指在正交化过 程中,v(j+1)仅与最近k个v(j-k+1),…v(j-1),v(j)正交,这 样做虽然破坏的正交性,但是降低了计算量.当 然k选得越小，对每个j对应的计算量也越小,但 可能要选更大的m才能取得满足精度要求的近 似解. IOM算法仅仅是把FOM算法的第三步改为 (iii’)对i=max(1,j-k+1),…,j, 计算hij与w(j).
其中, 确解.
( z) ( A( x z), x z)1/ 2
min ( z) (0) zx K
且x为(1.1)的精
２. (残量投影)设A为任意方阵, x(0)为初始近似解, 且 L=AK, 则x(1)为采用投影方法得到的新 近似解的充要条件是
( x ) min (