Octopus
|
QMR (quasi-minimal residual) algorithm for complex symmetric matrices algorithm taken from: Parallel implementation of efficient preconditioned linear solver for grid-based applications in chemical physics. II: QMR linear solver Appendix A. Simplified QMR algorithm W Chen and B Poirier, J Comput Phys 219, 198-209 (2006)
QMR (quasi-minimal residual) algorithm for complex matrices algorithm taken from: An Implementation of the QMR Method based on Coupled Two-Term Recurrences by R. W. Freund and N. M. Nachtigal (page 25) http:
Definition at line 194 of file solvers.F90.
Private Member Functions | |
subroutine | dsym_conjugate_gradients (np, x, b, op, dotp, iter, residue, threshold) |
The two following subroutines, sym_conjugate_gradients, and bi_conjugate_gradients, must be called under a common interface: conjugate_gradients. It provides an approximate solution to the linear system problem Ax = b. Solving a symmetric linear system, which is either real or complex symmetric or Hermitian the best choice is sym_conjugate_gradients (does not need \( A^\dagger \)) Solving a real unsymmetric or a complex non-Hermitian system bi_conjugate_gradients has to be chosen (where one does need \( A^\dagger \)). More... | |
subroutine | dbi_conjugate_gradients (np, x, b, op, opt, dotp, iter, residue, threshold) |
|
private |
The two following subroutines, sym_conjugate_gradients, and bi_conjugate_gradients, must be called under a common interface: conjugate_gradients. It provides an approximate solution to the linear system problem Ax = b. Solving a symmetric linear system, which is either real or complex symmetric or Hermitian the best choice is sym_conjugate_gradients (does not need \( A^\dagger \)) Solving a real unsymmetric or a complex non-Hermitian system bi_conjugate_gradients has to be chosen (where one does need \( A^\dagger \)).
Note: the complex-valued versions (both CG and BiCG) work only with symmetric operators but they may be non-Hermitian. This is a property of the CG algorithm. A comment on this can be found in chapter 4.1 of ftp:
subroutine conjugate_gradients(np, x, b, op, [opt,] dotp, iter [, residue] [, threshold] ) integer, intent(in) :: np => The dimension of the problem. real(real64), intent(inout) :: x => On input, an estimate to the solution. => On output, the approximate solution. real(real64), intent(in) :: b => The inhomogeneous term of the equation. interface subroutine op(x, y) real(real64), intent(in) :: x(:) real(real64), intent(out) :: y(:) end subroutine op end interface => This should be a procedure that computes \( Ax = y \). interface subroutine opt(x, y) real(real64), intent(in) :: x(:) real(real64), intent(out) :: y(:) end subroutine opt end interface => If present, this should be a procedure that computes \( A^\dagger x = y \). Only useful for non-Hermitian operators. interface real(real64) function dotp(x, y) real(real64), intent(inout) :: x(:) real(real64), intent(in) :: y(:) end function dotp => Calculates the \( <x | y> \). Depending on the matrix A one should choose: complex symmetric: \( <x | y> = x^T * y \) hermitian: \( <x | y> = x^\dagger * y \) general: \( <x | y> = x^\dagger * y \) integer, intent(inout) :: iter => On input, the maximum number of iterations that the procedure is allowed to take. On output, the iterations it actually did. real(real64), intent(out) :: residue => If present, it measures the final error: residue = \( <Ax - b | Ax - b> \). real(real64), intent(in) :: threshold => If present, it sets the required accuracy threshold for the algorithm to stop. If not present, this is set to 1.0e-6. [The algorithm stops when \( <Ax - b | Ax - b> <= threshold \), or iter iterations are reached.] end subroutine conjugate_gradients
(*) NOTE: The algorithm assumes that the vectors are given in an orthonormal basis.
Definition at line 1559 of file solvers.F90.
|
private |
Definition at line 1628 of file solvers.F90.