The Open FUSION Toolkit 1.0.0-8905cc5
Modeling tools for plasma and fusion research and engineering
|
Implicity Restarted Lanczos Method.
Public Member Functions | |
procedure | apply (self, u, alam) |
Solve system. | |
procedure(eigsolver_apply), deferred | apply (self, u, alam) |
Solve Eigen-system. | |
procedure | delete (self) |
Clean-up internal storage. | |
procedure(eigsolver_delete), deferred | delete (self) |
Clean-up internal storage. | |
procedure | max (self, u, alam) |
Solve system. | |
Public Attributes | |
class(oft_matrix), pointer | a => NULL() |
LHS matrix. | |
real(r8) | atol = 1.d-14 |
Absolute convergence tolerance \( |res| < atol \). | |
class(oft_solver_bc), pointer | bc => NULL() |
Boundary condition. | |
integer(i4) | cits = 0 |
Number of iteractions to convergence. | |
real(r8), dimension(:,:), pointer | eig_val => NULL() |
Eigenvalues. | |
real(r8), dimension(:,:), pointer | eig_vec => NULL() |
Eigenvectors. | |
integer(i4) | info = 0 |
Solver status/return code. | |
logical | initialized = .FALSE. |
Solver has been constructed. | |
integer(i4) | itplot =10 |
Output frequency for iterative solvers when pm=.TRUE. | |
integer(i4) | its = -1 |
Maximum iteration count. | |
class(oft_matrix), pointer | m => NULL() |
RHS matrix. | |
class(oft_solver), pointer | minv => NULL() |
RHS inversion operator. | |
integer(i4) | mode = 1 |
Operational mode. | |
integer(i4) | ncv = 0 |
Size of Ritz space (determined in call to apply ) | |
integer(i4) | nev = 1 |
Number of eigenvalues to compute. | |
class(oft_solver_bc), pointer | orthog => NULL() |
Orthogonalization. | |
logical | pm = .FALSE. |
Performance monitor override. | |
class(oft_solver), pointer | pre => NULL() |
Preconditioner. | |
real(r8) | rtol = 1.d-14 |
Relative convergence tolerance \( |res|/|res_0| < rtol \). | |
real(r8) | tol = 1.E-10_r8 |
Solver tolerance. | |
character(len=2) | which = 'LM' |
Spectrum search flag. | |
procedure apply | ( | class(oft_irlm_eigsolver), intent(inout) | self, |
class(oft_vector), intent(inout) | u, | ||
real(r8), intent(inout) | alam | ||
) |
Solve system.
Solver employs the ARPACK symmetric driver routine dsaupd
. The location of the eigenvalue is set in the calling class.
[in,out] | u | Guess field/Eigenvector |
[in,out] | alam | Eigenvalue |
|
pure virtualinherited |
Solve Eigen-system.
[in,out] | self | Solver object |
[in,out] | u | Guess/Solution field |
[in,out] | alam | Eigenvalue |
procedure delete | ( | class(oft_irlm_eigsolver), intent(inout) | self | ) |
Clean-up internal storage.
|
pure virtualinherited |
Clean-up internal storage.
[in,out] | self | Solver object |
procedure max | ( | class(oft_irlm_eigsolver), intent(inout) | self, |
class(oft_vector), intent(inout) | u, | ||
real(r8), intent(inout) | alam | ||
) |
Solve system.
Solver employs the ARPACK symmetric driver routine. Currently the guess field is only used to create work vectors and a random initialization is used for the solver.
[in,out] | u | Guess field/Eigenvector |
[in,out] | alam | Eigenvalue |
|
inherited |
LHS matrix.
|
inherited |
Absolute convergence tolerance \( |res| < atol \).
class(oft_solver_bc), pointer bc => NULL() |
Boundary condition.
|
inherited |
Number of iteractions to convergence.
real(r8), dimension(:,:), pointer eig_val => NULL() |
Eigenvalues.
real(r8), dimension(:,:), pointer eig_vec => NULL() |
Eigenvectors.
integer(i4) info = 0 |
Solver status/return code.
|
inherited |
Solver has been constructed.
|
inherited |
Output frequency for iterative solvers when pm=.TRUE.
|
inherited |
Maximum iteration count.
|
inherited |
RHS matrix.
class(oft_solver), pointer minv => NULL() |
RHS inversion operator.
integer(i4) mode = 1 |
Operational mode.
integer(i4) ncv = 0 |
Size of Ritz space (determined in call to apply
)
integer(i4) nev = 1 |
Number of eigenvalues to compute.
class(oft_solver_bc), pointer orthog => NULL() |
Orthogonalization.
|
inherited |
Performance monitor override.
|
inherited |
Preconditioner.
|
inherited |
Relative convergence tolerance \( |res|/|res_0| < rtol \).
real(r8) tol = 1.E-10_r8 |
Solver tolerance.
character(len=2) which = 'LM' |
Spectrum search flag.