Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2022 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 : MODULE fft_lib
8 :
9 : USE fft_kinds, ONLY: dp
10 : USE fft_plan, ONLY: fft_plan_type
11 : USE fftsg_lib, ONLY: fftsg1dm,&
12 : fftsg3d,&
13 : fftsg_do_cleanup,&
14 : fftsg_do_init,&
15 : fftsg_get_lengths
16 : USE fftw3_lib, ONLY: fftw31dm,&
17 : fftw33d,&
18 : fftw3_create_plan_1dm,&
19 : fftw3_create_plan_3d,&
20 : fftw3_destroy_plan,&
21 : fftw3_do_cleanup,&
22 : fftw3_do_init,&
23 : fftw3_get_lengths
24 : #include "../../base/base_uses.f90"
25 :
26 : IMPLICIT NONE
27 : PRIVATE
28 :
29 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'fft_lib'
30 :
31 : PUBLIC :: fft_do_cleanup, fft_do_init, fft_get_lengths, fft_create_plan_3d
32 : PUBLIC :: fft_create_plan_1dm, fft_1dm, fft_library, fft_3d, fft_destroy_plan
33 :
34 : CONTAINS
35 : ! **************************************************************************************************
36 : !> \brief Interface to FFT libraries
37 : !> \param fftlib ...
38 : !> \return ...
39 : !> \par History
40 : !> IAB 09-Jan-2009 : Modified to use fft_plan_type
41 : !> (c) The Numerical Algorithms Group (NAG) Ltd, 2009 on behalf of the HECToR project
42 : !> \author JGH
43 : ! **************************************************************************************************
44 8213 : FUNCTION fft_library(fftlib) RESULT(flib)
45 :
46 : CHARACTER(len=*), INTENT(IN) :: fftlib
47 : INTEGER :: flib
48 :
49 : SELECT CASE (fftlib)
50 : CASE DEFAULT
51 10 : flib = -1
52 : CASE ("FFTSG")
53 10 : flib = 1
54 : CASE ("FFTW3")
55 8213 : flib = 3
56 : END SELECT
57 :
58 8213 : END FUNCTION fft_library
59 :
60 : ! **************************************************************************************************
61 : !> \brief ...
62 : !> \param fft_type ...
63 : !> \param wisdom_file ...
64 : ! **************************************************************************************************
65 8213 : SUBROUTINE fft_do_init(fft_type, wisdom_file)
66 : INTEGER, INTENT(IN) :: fft_type
67 : CHARACTER(LEN=*), INTENT(IN) :: wisdom_file
68 :
69 8213 : SELECT CASE (fft_type)
70 : CASE DEFAULT
71 0 : CPABORT("fft_do_init")
72 : CASE (1)
73 10 : CALL fftsg_do_init()
74 : CASE (3)
75 8213 : CALL fftw3_do_init(wisdom_file)
76 : END SELECT
77 :
78 8213 : END SUBROUTINE
79 :
80 : ! **************************************************************************************************
81 : !> \brief ...
82 : !> \param fft_type ...
83 : !> \param wisdom_file ...
84 : !> \param ionode ...
85 : ! **************************************************************************************************
86 8003 : SUBROUTINE fft_do_cleanup(fft_type, wisdom_file, ionode)
87 : INTEGER, INTENT(IN) :: fft_type
88 : CHARACTER(LEN=*), INTENT(IN) :: wisdom_file
89 : LOGICAL, INTENT(IN) :: ionode
90 :
91 8003 : SELECT CASE (fft_type)
92 : CASE DEFAULT
93 0 : CPABORT("fft_do_cleanup")
94 : CASE (1)
95 10 : CALL fftsg_do_cleanup()
96 : CASE (3)
97 8003 : CALL fftw3_do_cleanup(wisdom_file, ionode)
98 : END SELECT
99 :
100 8003 : END SUBROUTINE
101 :
102 : ! **************************************************************************************************
103 : !> \brief ...
104 : !> \param fft_type ...
105 : !> \param DATA ...
106 : !> \param max_length ...
107 : ! **************************************************************************************************
108 93374 : SUBROUTINE fft_get_lengths(fft_type, DATA, max_length)
109 :
110 : INTEGER, INTENT(IN) :: fft_type
111 : INTEGER, DIMENSION(*) :: DATA
112 : INTEGER, INTENT(INOUT) :: max_length
113 :
114 93374 : SELECT CASE (fft_type)
115 : CASE DEFAULT
116 0 : CPABORT("fft_get_lengths")
117 : CASE (1)
118 92870 : CALL fftsg_get_lengths(DATA, max_length)
119 : CASE (3)
120 93374 : CALL fftw3_get_lengths(DATA, max_length)
121 : END SELECT
122 :
123 93374 : END SUBROUTINE fft_get_lengths
124 :
125 : ! **************************************************************************************************
126 :
127 : ! **************************************************************************************************
128 : !> \brief ...
129 : !> \param plan ...
130 : !> \param fft_type ...
131 : !> \param fft_in_place ...
132 : !> \param fsign ...
133 : !> \param n ...
134 : !> \param zin ...
135 : !> \param zout ...
136 : !> \param plan_style ...
137 : ! **************************************************************************************************
138 37564 : SUBROUTINE fft_create_plan_3d(plan, fft_type, fft_in_place, fsign, n, zin, zout, plan_style)
139 :
140 : TYPE(fft_plan_type), INTENT(INOUT) :: plan
141 : INTEGER, INTENT(IN) :: fft_type
142 : LOGICAL, INTENT(IN) :: fft_in_place
143 : INTEGER, INTENT(IN) :: fsign
144 : INTEGER, DIMENSION(3), INTENT(IN) :: n
145 : COMPLEX(KIND=dp), DIMENSION(*), INTENT(INOUT) :: zin, zout
146 : INTEGER, INTENT(IN) :: plan_style
147 :
148 37564 : plan%fft_type = fft_type
149 37564 : plan%fsign = fsign
150 37564 : plan%fft_in_place = fft_in_place
151 150256 : plan%n_3d = n
152 37564 : !$ plan%need_alt_plan = .FALSE.
153 :
154 : ! Planning only needed for FFTW3
155 37564 : IF (fft_type .EQ. 3) THEN
156 37492 : CALL fftw3_create_plan_3d(plan, zin, zout, plan_style)
157 37492 : plan%valid = .TRUE.
158 : END IF
159 :
160 37564 : END SUBROUTINE fft_create_plan_3d
161 :
162 : !
163 : ! really ugly, plan is intent out, because plan%fsign is also a status flag
164 : ! if something goes wrong, plan%fsign is set to zero, and the plan becomes invalid
165 : !
166 : ! **************************************************************************************************
167 : !> \brief ...
168 : !> \param plan ...
169 : !> \param scale ...
170 : !> \param zin ...
171 : !> \param zout ...
172 : !> \param stat ...
173 : ! **************************************************************************************************
174 527662 : SUBROUTINE fft_3d(plan, scale, zin, zout, stat)
175 : TYPE(fft_plan_type), INTENT(IN) :: plan
176 : REAL(KIND=dp), INTENT(IN) :: scale
177 : COMPLEX(KIND=dp), DIMENSION(*), INTENT(INOUT) :: zin, zout
178 : INTEGER, INTENT(OUT) :: stat
179 :
180 527662 : stat = plan%fsign
181 527662 : IF (plan%n_3d(1)*plan%n_3d(2)*plan%n_3d(3) > 0) THEN
182 527662 : SELECT CASE (plan%fft_type)
183 : CASE DEFAULT
184 0 : CPABORT("fft_3d")
185 : CASE (1)
186 886 : CALL fftsg3d(plan%fft_in_place, stat, scale, plan%n_3d, zin, zout)
187 : CASE (3)
188 527662 : CALL fftw33d(plan, scale, zin, zout, stat)
189 : END SELECT
190 : END IF
191 : ! stat is set to zero on error, -1,+1 are OK
192 527662 : IF (stat .EQ. 0) THEN
193 0 : stat = 1
194 : ELSE
195 527662 : stat = 0
196 : END IF
197 :
198 527662 : END SUBROUTINE fft_3d
199 :
200 : ! **************************************************************************************************
201 :
202 : ! **************************************************************************************************
203 : !> \brief ...
204 : !> \param plan ...
205 : !> \param fft_type ...
206 : !> \param fsign ...
207 : !> \param trans ...
208 : !> \param n ...
209 : !> \param m ...
210 : !> \param zin ...
211 : !> \param zout ...
212 : !> \param plan_style ...
213 : ! **************************************************************************************************
214 125152 : SUBROUTINE fft_create_plan_1dm(plan, fft_type, fsign, trans, n, m, zin, zout, plan_style)
215 : TYPE(fft_plan_type), INTENT(INOUT) :: plan
216 : INTEGER, INTENT(IN) :: fft_type, fsign
217 : LOGICAL, INTENT(IN) :: trans
218 : INTEGER, INTENT(IN) :: n, m
219 : COMPLEX(KIND=dp), DIMENSION(*), INTENT(IN) :: zin, zout
220 : INTEGER, INTENT(IN) :: plan_style
221 :
222 125152 : plan%fft_type = fft_type
223 125152 : plan%fsign = fsign
224 125152 : plan%trans = trans
225 125152 : plan%n = n
226 125152 : plan%m = m
227 125152 : !$ plan%need_alt_plan = .FALSE.
228 :
229 : ! Planning only needed for FFTW3
230 125152 : IF ((fft_type .EQ. 3) .AND. (n*m .NE. 0)) THEN
231 124984 : CALL fftw3_create_plan_1dm(plan, zin, zout, plan_style)
232 124984 : plan%valid = .TRUE.
233 : ELSE
234 168 : plan%valid = .FALSE.
235 : END IF
236 :
237 125152 : END SUBROUTINE fft_create_plan_1dm
238 :
239 : ! **************************************************************************************************
240 : !> \brief ...
241 : !> \param plan ...
242 : ! **************************************************************************************************
243 230776 : SUBROUTINE fft_destroy_plan(plan)
244 : TYPE(fft_plan_type), INTENT(INOUT) :: plan
245 :
246 : ! Planning only needed for FFTW3
247 :
248 230776 : IF (plan%valid .AND. plan%fft_type .EQ. 3) THEN
249 162476 : CALL fftw3_destroy_plan(plan)
250 162476 : plan%valid = .FALSE.
251 : END IF
252 :
253 230776 : END SUBROUTINE
254 :
255 : ! **************************************************************************************************
256 : !> \brief ...
257 : !> \param plan ...
258 : !> \param zin ...
259 : !> \param zout ...
260 : !> \param scale ...
261 : !> \param stat ...
262 : ! **************************************************************************************************
263 5922796 : SUBROUTINE fft_1dm(plan, zin, zout, scale, stat)
264 : TYPE(fft_plan_type), INTENT(IN) :: plan
265 : COMPLEX(KIND=dp), DIMENSION(*), INTENT(INOUT) :: zin, zout
266 : REAL(KIND=dp), INTENT(IN) :: scale
267 : INTEGER, INTENT(OUT) :: stat
268 :
269 5922796 : stat = plan%fsign
270 5922796 : IF (plan%n*plan%m > 0) THEN
271 5922796 : SELECT CASE (plan%fft_type)
272 : CASE DEFAULT
273 0 : CPABORT("fft_1dm")
274 : CASE (1)
275 15402 : CALL fftsg1dm(stat, plan%trans, plan%n, plan%m, zin, zout, scale)
276 : CASE (3)
277 5922796 : CALL fftw31dm(plan, zin, zout, scale, stat)
278 : END SELECT
279 : END IF
280 : ! stat is set to zero on error, -1,+1 are OK
281 5922796 : IF (stat .EQ. 0) THEN
282 0 : stat = 1
283 : ELSE
284 5922796 : stat = 0
285 : END IF
286 :
287 5922796 : END SUBROUTINE fft_1dm
288 :
289 : END MODULE
290 :
|