Preconditioner: Chebyshev
This smoother provides Chekyshev polynomial smoothing of a runtime specified order.
Chebyshev Type
Different Chebyshev smoother types may be selected. See James Lottes (2022) and Malachi Phillips, Paul Fischer (2022) for discussion of Chebyshev smoother types for multigrid V-cycles.
LFAToolkit.ChebyshevType.ChebyType
— TypeChebyshev types
Types:
first
: 1st-kind Chebyshevfourth
: 4th-kind Chebyshevopt_fourth
: optimized 4th-kind Chebyshev
Example:
LFAToolkit.ChebyshevType.ChebyshevType
# output
Enum LFAToolkit.ChebyshevType.ChebyType:
first = 0
fourth = 1
opt_fourth = 2
Example
This is an example of a Chebyshev smoother.
# ------------------------------------------------------------------------------
# Chebyshev smoother example
# ------------------------------------------------------------------------------
using LFAToolkit
using LinearAlgebra
# setup
mesh = Mesh2D(1.0, 1.0)
p = 3
# diffusion operator
diffusion = GalleryOperator("diffusion", p + 1, p + 1, mesh)
# Chebyshev smoother
chebyshev = Chebyshev(diffusion)
# compute operator symbols
A = computesymbols(chebyshev, [3], [π, π])
eigenvalues = real(eigvals(A))
# ------------------------------------------------------------------------------
Plot for the symbol a cubic Chebyshev smoother for the 2D scalar diffusion problem with cubic basis.
Documentation
LFAToolkit.Chebyshev
— TypeChebyshev(operator)
Chebyshev polynomial preconditioner for finite element operators. The Chebyshev semi-iterative method is applied to the matrix $D^{-1} A$, where $D^{-1}$ is the inverse of the operator diagonal.
Arguments:
operator::Operator
: finite element operator to preconditionchebyshevtype::ChebyshevType
: Chebyshev type, first, fourth, or opt. fourth kind
Returns:
- Chebyshev preconditioner object
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
mass = GalleryOperator("mass", 4, 4, mesh);
# preconditioner
chebyshev = Chebyshev(mass);
# verify
println(chebyshev)
# output
chebyshev preconditioner:
1st-kind Chebyshev
eigenvalue estimates:
estimated minimum 0.2500
estimated maximum 1.3611
estimate scaling:
λ_min = a * estimated min + b * estimated max
λ_max = c * estimated min + d * estimated max
a = 0.0000
b = 0.1000
c = 0.0000
d = 1.0000
LFAToolkit.computesymbols
— Methodcomputesymbols(chebyshev, ω, θ)
Compute or retrieve the symbol matrix for a Chebyshev preconditioned operator
Arguments:
chebyshev::Chebyshev
: Chebyshev preconditioner to compute symbol matrix forω::Array{Real}
: smoothing parameter array [degree], [degree, $\lambda_{\text{max}}$], or [degree, $\lambda_{\text{min}}$, $\lambda_{\text{max}}$]θ::Array{Real}
: Fourier mode frequency array (one frequency per dimension)
Returns:
- symbol matrix for the Chebyshev preconditioned operator
Example:
using LinearAlgebra
for dimension = 1:3
# setup
mesh = []
if dimension == 1
mesh = Mesh1D(1.0)
elseif dimension == 2
mesh = Mesh2D(1.0, 1.0)
elseif dimension == 3
mesh = Mesh3D(1.0, 1.0, 1.0)
end
diffusion = GalleryOperator("diffusion", 3, 3, mesh)
# preconditioner
chebyshev = Chebyshev(diffusion)
# compute symbols
A = computesymbols(chebyshev, [1], π * ones(dimension))
# verify
eigenvalues = real(eigvals(A))
if dimension == 1
@assert minimum(eigenvalues) ≈ 0.15151515151515105
@assert maximum(eigenvalues) ≈ 0.27272727272727226
elseif dimension == 2
@assert minimum(eigenvalues) ≈ -0.25495098334134725
@assert maximum(eigenvalues) ≈ -0.17128758445192374
elseif dimension == 3
@assert minimum(eigenvalues) ≈ -0.8181818181818181
@assert maximum(eigenvalues) ≈ -0.357575757575757
end
end
# output