Skip to content

Commit

Permalink
Block-Krylov processes
Browse files Browse the repository at this point in the history
  • Loading branch information
amontoison committed Oct 20, 2023
1 parent 9b536d9 commit 03f5f95
Show file tree
Hide file tree
Showing 18 changed files with 1,401 additions and 30 deletions.
2 changes: 1 addition & 1 deletion README.md
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
# Krylov.jl: A Julia basket of hand-picked Krylov methods

| **Documentation** | **CI** | **Coverage** | **DOI** | **Downloads** |
|:-----------------:|:-------------------------------:|:------------:|:-------:|:-------------:|
|:-----------------:|:------:|:------------:|:-------:|:-------------:|
| [![docs-stable][docs-stable-img]][docs-stable-url] [![docs-dev][docs-dev-img]][docs-dev-url] | [![build-gh][build-gh-img]][build-gh-url] [![build-cirrus][build-cirrus-img]][build-cirrus-url] | [![codecov][codecov-img]][codecov-url] | [![doi][doi-img]][doi-url] | [![downloads][downloads-img]][downloads-url] |

[docs-stable-img]: https://img.shields.io/badge/docs-stable-blue.svg
Expand Down
1 change: 1 addition & 0 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@ makedocs(
pages = ["Home" => "index.md",
"API" => "api.md",
"Krylov processes" => "processes.md",
"Block Krylov processes" => "block_processes.md",
"Krylov methods" => ["Hermitian positive definite linear systems" => "solvers/spd.md",
"Hermitian indefinite linear systems" => "solvers/sid.md",
"Non-Hermitian square linear systems" => "solvers/unsymmetric.md",
Expand Down
227 changes: 227 additions & 0 deletions docs/src/block_processes.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,227 @@
# [Block Krylov processes](@id block-krylov-processes)

## [Block Hermitian Lanczos](@id block-hermitian-lanczos)

If the vector $b$ in the [Hermitian Lanczos](@ref hermitian-lanczos) process is replaced by a matrix $B$ with $p$ columns, we can derive the block Hermitian Lanczos process.

![block_hermitian_lanczos](./graphics/block_hermitian_lanczos.png)

After $k$ iterations of the block Hermitian Lanczos process, the situation may be summarized as
```math
\begin{align*}
A V_k &= V_{k+1} T_{k+1,k}, \\
V_k^H V_k &= I_{pk},
\end{align*}
```
where $V_k$ is an orthonormal basis of the block Krylov subspace $\mathcal{K}_k^{\square}(A,B)$,
```math
T_{k+1,k} =
\begin{bmatrix}
\Omega_1 & \Psi_2^H & & \\
\Psi_2 & \Omega_2 & \ddots & \\
& \ddots & \ddots & \Psi_k^H \\
& & \Psi_k & \Omega_k \\
& & & \Psi_{k+1}
\end{bmatrix}.
```

The function [`hermitian_lanczos`](@ref hermitian_lanczos(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $\Psi_1$ and $T_{k+1,k}$.

```@docs
hermitian_lanczos(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```

## [Block Non-Hermitian Lanczos](@id block-nonhermitian-lanczos)

If the vectors $b$ and $c$ in the [non-Hermitian Lanczos](@ref nonhermitian-lanczos) process are replaced by matrices $B$ and $C$ with both $p$ columns, we can derive the block non-Hermitian Lanczos process.

![block_nonhermitian_lanczos](./graphics/block_nonhermitian_lanczos.png)

After $k$ iterations of the block non-Hermitian Lanczos process, the situation may be summarized as
```math
\begin{align*}
A V_k &= V_{k+1} T_{k+1,k}, \\
A^H U_k &= U_{k+1} T_{k,k+1}^H, \\
V_k^H U_k &= U_k^H V_k = I_{pk},
\end{align*}
```
where $V_k$ and $U_k$ are bases of the block Krylov subspaces $\mathcal{K}^{\square}_k(A,B)$ and $\mathcal{K}^{\square}_k (A^H,C)$, respectively,
```math
T_{k+1,k} =
\begin{bmatrix}
\Omega_1 & \Phi_2 & & \\
\Psi_2 & \Omega_2 & \ddots & \\
& \ddots & \ddots & \Phi_k \\
& & \Psi_k & \Omega_k \\
& & & \Psi_{k+1}
\end{bmatrix}
, \qquad
T_{k,k+1}^H =
\begin{bmatrix}
\Omega_1^H & \Psi_2^H & & \\
\Phi_2^H & \Omega_2^H & \ddots & \\
& \ddots & \ddots & \Psi_k^H \\
& & \Phi_k^H & \Omega_k^H \\
& & & \Phi_{k+1}^H
\end{bmatrix}.
```

The function [`nonhermitian_lanczos`](@ref nonhermitian_lanczos(::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $\Psi_1$, $T_{k+1,k}$, $U_{k+1}$ $\Phi_1^H$ and $T_{k,k+1}^H$.

```@docs
nonhermitian_lanczos(::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```

## [Block Arnoldi](@id block-arnoldi)

If the vector $b$ in the [Arnoldi](@ref arnoldi) process is replaced by a matrix $B$ with $p$ columns, we can derive the block Arnoldi process.

![block_arnoldi](./graphics/block_arnoldi.png)

After $k$ iterations of the block Arnoldi process, the situation may be summarized as
```math
\begin{align*}
A V_k &= V_{k+1} H_{k+1,k}, \\
V_k^H V_k &= I_{pk},
\end{align*}
```
where $V_k$ is an orthonormal basis of the block Krylov subspace $\mathcal{K}_k^{\square}(A,B)$,
```math
H_{k+1,k} =
\begin{bmatrix}
\Psi_{1,1}~ & \Psi_{1,2}~ & \ldots & \Psi_{1,k} \\
\Psi_{2,1}~ & \ddots~ & \ddots & \vdots \\
& \ddots~ & \ddots & \Psi_{k-1,k} \\
& & \Psi_{k,k-1} & \Psi_{k,k} \\
& & & \Psi_{k+1,k}
\end{bmatrix}.
```

The function [`arnoldi`](@ref arnoldi(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $\Gamma$, and $H_{k+1,k}$.

```@docs
arnoldi(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```

## [Block Golub-Kahan](@id block-golub-kahan)

If the vector $b$ in the [Golub-Kahan](@ref golub-kahan) process is replaced by a matrix $B$ with $p$ columns, we can derive the block Golub-Kahan process.

![block_golub_kahan](./graphics/block_golub_kahan.png)

After $k$ iterations of the block Golub-Kahan process, the situation may be summarized as
```math
\begin{align*}
A V_k &= U_{k+1} B_k, \\
A^H U_{k+1} &= V_{k+1} L_{k+1}^H, \\
V_k^H V_k &= U_k^H U_k = I_{pk},
\end{align*}
```
where $V_k$ and $U_k$ are bases of the block Krylov subspaces $\mathcal{K}_k^{\square}(A^HA,A^HB)$ and $\mathcal{K}_k^{\square}(AA^H,B)$, respectively,
```math
B_k =
\begin{bmatrix}
\Omega_1 & & & \\
\Psi_2 & \Omega_2 & & \\
& \ddots & \ddots & \\
& & \Psi_k & \Omega_k \\
& & & \Psi_{k+1} \\
\end{bmatrix}
, \qquad
L_{k+1}^H =
\begin{bmatrix}
\Omega_1^H & \Psi_2^H & & & \\
& \Omega_2^H & \ddots & & \\
& & \ddots & \Psi_k^H & \\
& & & \Omega_k^H & \Psi_{k+1}^H \\
& & & & \Omega_{k+1}^H \\
\end{bmatrix}.
```

The function [`golub_kahan`](@ref golub_kahan(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $U_{k+1}$, $\Psi_1$ and $L_{k+1}$.

```@docs
golub_kahan(::Any, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```

## [Block Saunders-Simon-Yip](@id block-saunders-simon-yip)

If the vectors $b$ and $c$ in the [Saunders-Simon-Yip](@ref saunders-simon-yip) process are replaced by matrices $B$ and $C$ with both $p$ columns, we can derive the block Saunders-Simon-Yip process.

![block_saunders_simon_yip](./graphics/block_saunders_simon_yip.png)

After $k$ iterations of the block Saunders-Simon-Yip process, the situation may be summarized as
```math
\begin{align*}
A U_k &= V_{k+1} T_{k+1,k}, \\
A^H V_k &= U_{k+1} T_{k,k+1}^H, \\
V_k^H V_k &= U_k^H U_k = I_{pk},
\end{align*}
```
where $\begin{bmatrix} V_k & 0 \\ 0 & U_k \end{bmatrix}$ is an orthonormal basis of the block Krylov subspace $\mathcal{K}^{\square}_k \left(\begin{bmatrix} 0 & A \\ A^H & 0 \end{bmatrix}, \begin{bmatrix} B & 0 \\ 0 & C \end{bmatrix}\right)$,
```math
T_{k+1,k} =
\begin{bmatrix}
\Omega_1 & \Phi_2 & & \\
\Psi_2 & \Omega_2 & \ddots & \\
& \ddots & \ddots & \Phi_k \\
& & \Psi_k & \Omega_k \\
& & & \Psi_{k+1}
\end{bmatrix}
, \qquad
T_{k,k+1}^H =
\begin{bmatrix}
\Omega_1^H & \Psi_2^H & & \\
\Phi_2^H & \Omega_2^H & \ddots & \\
& \ddots & \ddots & \Psi_k^H \\
& & \Phi_k^H & \Omega_k^H \\
& & & \Phi_{k+1}^H
\end{bmatrix}.
```

The function [`saunders_simon_yip`](@ref saunders_simon_yip(::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $\Psi_1$, $T_{k+1,k}$, $U_{k+1}$, $\Phi_1^H$ and $T_{k,k+1}^H$.

```@docs
saunders_simon_yip(::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```

## [Block Montoison-Orban](@id block-montoison-orban)

If the vectors $b$ and $c$ in the [Montoison-Orban](@ref montoison-orban) process are replaced by matrices $D$ and $C$ with both $p$ columns, we can derive the block Montoison-Orban process.

![block_montoison_orban](./graphics/block_montoison_orban.png)

After $k$ iterations of the block Montoison-Orban process, the situation may be summarized as
```math
\begin{align*}
A U_k &= V_{k+1} H_{k+1,k}, \\
B V_k &= U_{k+1} F_{k+1,k}, \\
V_k^H V_k &= U_k^H U_k = I_{pk},
\end{align*}
```
where $\begin{bmatrix} V_k & 0 \\ 0 & U_k \end{bmatrix}$ is an orthonormal basis of the block Krylov subspace $\mathcal{K}^{\square}_k \left(\begin{bmatrix} 0 & A \\ B & 0 \end{bmatrix}, \begin{bmatrix} D & 0 \\ 0 & C \end{bmatrix}\right)$,
```math
H_{k+1,k} =
\begin{bmatrix}
\Psi_{1,1}~ & \Psi_{1,2}~ & \ldots & \Psi_{1,k} \\
\Psi_{2,1}~ & \ddots~ & \ddots & \vdots \\
& \ddots~ & \ddots & \Psi_{k-1,k} \\
& & \Psi_{k,k-1} & \Psi_{k,k} \\
& & & \Psi_{k+1,k}
\end{bmatrix}
, \qquad
F_{k+1,k} =
\begin{bmatrix}
\Phi_{1,1}~ & \Phi_{1,2}~ & \ldots & \Phi_{1,k} \\
\Phi_{2,1}~ & \ddots~ & \ddots & \vdots \\
& \ddots~ & \ddots & \Phi_{k-1,k} \\
& & \Phi_{k,k-1} & \Phi_{k,k} \\
& & & \Phi_{k+1,k}
\end{bmatrix}.
```

The function [`montoison_orban`](@ref montoison_orban(::Any, ::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)) returns $V_{k+1}$, $\Gamma$, $H_{k+1,k}$, $U_{k+1}$, $\Lambda$, and $F_{k+1,k}$.

```@docs
montoison_orban(::Any, ::Any, ::AbstractMatrix{FC}, ::AbstractMatrix{FC}, ::Int) where FC <: (Union{Complex{T}, T} where T <: AbstractFloat)
```
Binary file added docs/src/graphics/block_arnoldi.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added docs/src/graphics/block_golub_kahan.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added docs/src/graphics/block_hermitian_lanczos.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added docs/src/graphics/block_montoison_orban.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added docs/src/graphics/block_nonhermitian_lanczos.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added docs/src/graphics/block_saunders_simon_yip.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
3 changes: 3 additions & 0 deletions src/Krylov.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,10 @@ using LinearAlgebra, SparseArrays, Printf
include("krylov_utils.jl")
include("krylov_stats.jl")
include("krylov_solvers.jl")

include("processes_utils.jl")
include("krylov_processes.jl")
include("block_krylov_processes.jl")

include("cg.jl")
include("cr.jl")
Expand Down
Loading

0 comments on commit 03f5f95

Please sign in to comment.