Octopus
solvers.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19#include "global.h"
20
22module solvers_oct_m
23 use blas_oct_m
24 use debug_oct_m
25 use global_oct_m
26 use iso_c_binding
27 use, intrinsic :: iso_fortran_env
29 use loct_oct_m
32
33 implicit none
34
35 private
36 public :: &
43 didrs, &
44 zidrs
45
46 abstract interface
47 subroutine doperator_i(x, y, userdata)
48 use, intrinsic :: iso_fortran_env
49 use, intrinsic :: iso_c_binding, only: c_ptr
50 implicit none
51 real(real64), contiguous, intent(in) :: x(:)
52 real(real64), contiguous, intent(out) :: y(:)
53 type(c_ptr), intent(in) :: userdata(:)
54 end subroutine doperator_i
55 subroutine zoperator_i(x, y, userdata)
56 use, intrinsic :: iso_fortran_env
57 use, intrinsic :: iso_c_binding, only: c_ptr
58 implicit none
59 complex(real64), contiguous, intent(in) :: x(:)
60 complex(real64), contiguous, intent(out) :: y(:)
61 type(c_ptr), intent(in) :: userdata(:)
62 end subroutine zoperator_i
63 function ddotp_i(x, y)
64 import real64
65 implicit none
66 real(real64) :: ddotp_i
67 real(real64), intent(in) :: x(:)
68 real(real64), intent(in) :: y(:)
69 end function ddotp_i
70 function dnrm_i(x)
71 import real64
72 implicit none
73 real(real64), intent(in) :: x(:)
74 real(real64) :: dnrm_i
75 end function dnrm_i
76 function zdotp_i(x, y)
77 import real64
78 implicit none
79 complex(real64), intent(in) :: x(:)
80 complex(real64), intent(in) :: y(:)
81 complex(real64) :: zdotp_i
82 end function zdotp_i
83 function znrm_i(x)
84 import real64
85 implicit none
86 complex(real64), intent(in) :: x(:)
87 real(real64) :: znrm_i
88 end function znrm_i
89 end interface
90
91
99
105
106 interface dconjugate_gradients
108 end interface dconjugate_gradients
109
110 interface zconjugate_gradients
112 end interface zconjugate_gradients
113
114contains
115
116#include "undef.F90"
117#include "complex.F90"
118#include "solvers_inc.F90"
119
120#include "undef.F90"
121#include "real.F90"
122#include "solvers_inc.F90"
123
124end module solvers_oct_m
125
126!! Local Variables:
127!! mode: f90
128!! coding: utf-8
129!! End:
This module contains interfaces for BLAS routines You should not use these routines directly....
Definition: blas.F90:120
This module is intended to contain "only mathematical" functions and procedures.
Definition: solvers.F90:117
subroutine zsym_conjugate_gradients(np, x, b, op, dotp, iter, residue, threshold, userdata)
The two following subroutines, sym_conjugate_gradients, and bi_conjugate_gradients,...
Definition: solvers.F90:338
subroutine dbi_conjugate_gradients(np, x, b, op, opt, dotp, iter, residue, threshold, userdata)
Definition: solvers.F90:1690
complex(real64) function, dimension(size(b, 1), size(b, 2)), public zidrs(b, s, preconditioner, matrixvector, ddotprod, zdotprod, tolerance, maximum_iterations, variant, flag, relres, iterations, x0, U0, omega, resvec, H)
This is the "Induced Dimension Reduction", IDR(s) (for s=4). IDR(s) is a robust and efficient short r...
Definition: solvers.F90:963
real(real64) function, dimension(size(b, 1), size(b, 2)), public didrs(b, s, preconditioner, matrixvector, ddotprod, zdotprod, tolerance, maximum_iterations, variant, flag, relres, iterations, x0, U0, omega, resvec, H)
This is the "Induced Dimension Reduction", IDR(s) (for s=4). IDR(s) is a robust and efficient short r...
Definition: solvers.F90:2218
subroutine, public dqmr_sym_gen_dotu(np, x, b, op, dotu, nrm2, prec, iter, userdata, residue, threshold, showprogress, converged, use_initial_guess)
for complex symmetric matrices W Chen and B Poirier, J Comput Phys 219, 198-209 (2006)
Definition: solvers.F90:1774
subroutine, public zqmr_sym_gen_dotu(np, x, b, op, dotu, nrm2, prec, iter, userdata, residue, threshold, showprogress, converged, use_initial_guess)
for complex symmetric matrices W Chen and B Poirier, J Comput Phys 219, 198-209 (2006)
Definition: solvers.F90:519
subroutine, public zqmr_gen_dotu(np, x, b, op, opt, dotu, nrm2, prec, prect, iter, userdata, residue, threshold, showprogress, converged)
for general complex matrices taken from 'An Implementation of the QMR Method based on Coupled Two-Ter...
Definition: solvers.F90:744
subroutine, public dqmr_gen_dotu(np, x, b, op, opt, dotu, nrm2, prec, prect, iter, userdata, residue, threshold, showprogress, converged)
for general complex matrices taken from 'An Implementation of the QMR Method based on Coupled Two-Ter...
Definition: solvers.F90:1999
subroutine dsym_conjugate_gradients(np, x, b, op, dotp, iter, residue, threshold, userdata)
The two following subroutines, sym_conjugate_gradients, and bi_conjugate_gradients,...
Definition: solvers.F90:1593
subroutine zbi_conjugate_gradients(np, x, b, op, opt, dotp, iter, residue, threshold, userdata)
Definition: solvers.F90:435