Preconditioner: Chebyshev
LFAToolkit.getoperatordiagonalinverse
— Methodgetoperatordiagonalinverse(chebyshev)
Compute or retrieve the inverse of the symbol matrix diagonal for a Chebyshev preconditioner
Returns:
- symbol matrix diagonal inverse for the operator
Example:
# setup
mesh = Mesh1D(1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
# preconditioner
chebyshev = Chebyshev(diffusion);
# verify operator diagonal inverse
@assert chebyshev.operatordiagonalinverse ≈ [6/7 0; 0 3/4]
# output
LFAToolkit.seteigenvalueestimatescaling
— Functionseteigenvalueestimatescaling(chebyshev, eigenvaluescaling)
Set the scaling of the eigenvalue estimates for a Chebyshev preconditioner
Arguments:
chebyshev::Chebyshev
: preconditioner to set eigenvalue estimate scalingeigenvaluescaling::Array{Float64,1}
: array of 4 scaling factors to use when setting $\lambda_{\text{min}}$ and $\lambda_{\text{max}}$ based on eigenvalue estimates
$\lambda_{\text{min}}$ = a * estimated min + b * estimated max
$\lambda_{\text{max}}$ = c * estimated min + d * estimated max
Example:
# setup
mesh = Mesh1D(1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
# preconditioner
chebyshev = Chebyshev(diffusion)
# set eigenvalue estimate scaling
# PETSc default is to use 1.1, 0.1 of max eigenvalue estimate
# https://www.mcs.anl.gov/petsc/petsc-3.3/docs/manualpages/KSP/KSPChebyshevSetEstimateEigenvalues.html
chebyshev.eigenvaluescaling = [0.0, 0.1, 0.0, 1.1];
println(chebyshev)
# output
chebyshev preconditioner:
1st-kind Chebyshev
eigenvalue estimates:
estimated minimum 0.0000
estimated maximum 2.1429
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.1000
LFAToolkit.geteigenvalueestimates
— Methodgeteigenvalueestimates(chebyshev)
Compute or retrieve the eigenvalue estimates for a Chebyshev preconditioner
Returns:
- eigenvalue estimates for the operator
Example:
# setup
mesh = Mesh1D(1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
# preconditioner
chebyshev = Chebyshev(diffusion);
# verify eigenvalue estimates
@assert chebyshev.eigenvalueestimates ≈ [0, 15 / 7]
# output
LFAToolkit.getbetas
— Methodgetbetas(k)
Compute beta coefficients associated with fourth/optimal fourth-kind Chebyshev polynomials
Returns:
- beta coefficients associated with fourth/optimal fourth-kind Chebyshev polynomials