1!--------------------------------------------------------------------------------------------------!
2!  DFTB+: general package for performing fast atomistic simulations                                !
3!  Copyright (C) 2006 - 2019  DFTB+ developers group                                               !
4!                                                                                                  !
5!  See the LICENSE file for terms of usage and distribution.                                       !
6!--------------------------------------------------------------------------------------------------!
7
8#:include 'common.fypp'
9
10!> Interfaces for the ARPACK routines needed in DFTB+ (currently for the linear response excited
11!> state calculations).
12module dftbp_arpack
13  use dftbp_accuracy, only : rsp, rdp
14  implicit none
15  private
16
17  public :: withArpack
18
19#:if WITH_ARPACK
20
21  public :: saupd, seupd
22
23
24  !> Whether code was built with Arpack support
25  logical, parameter :: withArpack = .true.
26
27
28  !> Wrapper around ARPACK routines ssaupd/dsaupd.
29  interface saupd
30
31    !> single precision Arnoldi solver call
32    subroutine ssaupd(ido, bmat, n, which, nev, tol, resid, ncv, v, ldv,&
33        & iparam, ipntr, workd, workl, lworkl, info)
34      import :: rsp
35      implicit none
36      integer, intent(inout) :: ido
37      character, intent(in) :: bmat
38      integer, intent(in) :: n
39      character(2), intent(in) :: which
40      integer, intent(in) :: nev
41      real(rsp), intent(in) :: tol
42      real(rsp), intent(inout) :: resid(n)
43      integer, intent(in) :: ncv
44      integer, intent(in) :: ldv
45      real(rsp), intent(out) :: v(ldv, ncv)
46      integer, intent(inout) :: iparam(11)
47      integer, intent(out) :: ipntr(11)
48      real(rsp), intent(inout) :: workd(3 * n)
49      integer, intent(in) :: lworkl
50      real(rsp), intent(inout) :: workl(lworkl)
51      integer, intent(inout) :: info
52    end subroutine ssaupd
53
54
55    !> double precision Arnoldi solver call
56    subroutine dsaupd(ido, bmat, n, which, nev, tol, resid, ncv, v, ldv,&
57        & iparam, ipntr, workd, workl, lworkl, info)
58      import :: rdp
59      implicit none
60      integer, intent(inout) :: ido
61      character, intent(in) :: bmat
62      integer, intent(in) :: n
63      character(2), intent(in) :: which
64      integer, intent(in) :: nev
65      real(rdp), intent(in) :: tol
66      real(rdp), intent(inout) :: resid(n)
67      integer, intent(in) :: ncv
68      integer, intent(in) :: ldv
69      real(rdp), intent(out) :: v(ldv, ncv)
70      integer, intent(inout) :: iparam(11)
71      integer, intent(out) :: ipntr(11)
72      real(rdp), intent(inout) :: workd(3 * n)
73      integer, intent(in) :: lworkl
74      real(rdp), intent(inout) :: workl(lworkl)
75      integer, intent(inout) :: info
76    end subroutine dsaupd
77  end interface saupd
78
79
80  !> Wrapper around ARPACK routines sseupd/dseupd.
81  interface seupd
82
83
84    !> single precision return from the results of the solver
85    subroutine sseupd(rvec, howmny, sel, d, z, ldz, sigma, bmat, n, which, nev,&
86        & tol, resid, ncv, v, ldv, iparam, ipntr, workd, workl, lworkl, info)
87      import :: rsp
88      logical, intent(in) :: rvec
89      character, intent(in) :: howmny
90      integer, intent(in) :: ncv
91      logical, intent(in) :: sel(ncv)
92      integer, intent(in) :: nev
93      real(rsp), intent(out) :: d(nev)
94      integer, intent(in) :: ldz
95      real(rsp), intent(out) :: z(ldz, nev)
96      real(rsp), intent(in) :: sigma
97      character, intent(in) :: bmat
98      integer, intent(in) :: n
99      character(2), intent(in) :: which
100      real(rsp), intent(in) :: tol
101      real(rsp), intent(in) :: resid(n)
102      integer, intent(in) :: ldv
103      real(rsp), intent(inout) :: v(ldv, ncv)
104      integer, intent(in) :: iparam(7)
105      integer, intent(inout) :: ipntr(11)
106      real(rsp), intent(inout) :: workd(2 * n)
107      integer, intent(in) :: lworkl
108      real(rsp), intent(inout) :: workl(lworkl)
109      integer, intent(inout) :: info
110    end subroutine sseupd
111
112
113    !> double precision return from the results of the solver
114    subroutine dseupd(rvec, howmny, sel, d, z, ldz, sigma, bmat, n, which, nev,&
115        & tol, resid, ncv, v, ldv, iparam, ipntr, workd, workl, lworkl, info)
116      import :: rdp
117      logical, intent(in) :: rvec
118      character, intent(in) :: howmny
119      integer, intent(in) :: ncv
120      logical, intent(in) :: sel(ncv)
121      integer, intent(in) :: nev
122      real(rdp), intent(out) :: d(nev)
123      integer, intent(in) :: ldz
124      real(rdp), intent(out) :: z(ldz, nev)
125      real(rdp), intent(in) :: sigma
126      character, intent(in) :: bmat
127      integer, intent(in) :: n
128      character(2), intent(in) :: which
129      real(rdp), intent(in) :: tol
130      real(rdp), intent(in) :: resid(n)
131      integer, intent(in) :: ldv
132      real(rdp), intent(inout) :: v(ldv, ncv)
133      integer, intent(in) :: iparam(7)
134      integer, intent(inout) :: ipntr(11)
135      real(rdp), intent(inout) :: workd(2 * n)
136      integer, intent(in) :: lworkl
137      real(rdp), intent(inout) :: workl(lworkl)
138      integer, intent(inout) :: info
139    end subroutine dseupd
140  end interface seupd
141
142#:else
143
144
145  !> Whether code was built with ARPACK support
146  logical, parameter :: withArpack = .false.
147
148#:endif
149
150end module dftbp_arpack
151