Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief pw_types
10 : !> \author CJM
11 : ! **************************************************************************************************
12 : MODULE ewald_pw_types
13 : USE ao_util, ONLY: exp_radius
14 : USE cell_types, ONLY: cell_type
15 : USE cp_log_handling, ONLY: cp_get_default_logger,&
16 : cp_logger_type
17 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
18 : cp_print_key_unit_nr
19 : USE cp_realspace_grid_init, ONLY: init_input_type
20 : USE dg_types, ONLY: dg_create,&
21 : dg_release,&
22 : dg_type
23 : USE dgs, ONLY: dg_pme_grid_setup
24 : USE ewald_environment_types, ONLY: ewald_env_get,&
25 : ewald_environment_type
26 : USE input_section_types, ONLY: section_vals_get_subs_vals,&
27 : section_vals_type
28 : USE kinds, ONLY: dp
29 : USE mathconstants, ONLY: pi
30 : USE message_passing, ONLY: mp_comm_self,&
31 : mp_para_env_type
32 : USE pw_grid_types, ONLY: HALFSPACE,&
33 : pw_grid_type
34 : USE pw_grids, ONLY: pw_grid_create,&
35 : pw_grid_release,&
36 : pw_grid_setup
37 : USE pw_poisson_methods, ONLY: pw_poisson_set
38 : USE pw_poisson_read_input, ONLY: pw_poisson_read_parameters
39 : USE pw_poisson_types, ONLY: do_ewald_ewald,&
40 : do_ewald_none,&
41 : do_ewald_pme,&
42 : do_ewald_spme,&
43 : pw_poisson_parameter_type,&
44 : pw_poisson_type
45 : USE pw_pool_types, ONLY: pw_pool_create,&
46 : pw_pool_p_type,&
47 : pw_pool_release,&
48 : pw_pool_type
49 : USE realspace_grid_types, ONLY: &
50 : realspace_grid_desc_type, realspace_grid_input_type, realspace_grid_type, rs_grid_create, &
51 : rs_grid_create_descriptor, rs_grid_print, rs_grid_release, rs_grid_release_descriptor, &
52 : rs_grid_retain_descriptor
53 : #include "./base/base_uses.f90"
54 :
55 : IMPLICIT NONE
56 :
57 : PRIVATE
58 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'ewald_pw_types'
59 : PUBLIC :: ewald_pw_type, ewald_pw_release, &
60 : ewald_pw_create, &
61 : ewald_pw_get, ewald_pw_set
62 :
63 : ! **************************************************************************************************
64 : TYPE ewald_pw_type
65 : PRIVATE
66 : TYPE(pw_pool_type), POINTER :: pw_small_pool => NULL()
67 : TYPE(pw_pool_type), POINTER :: pw_big_pool => NULL()
68 : TYPE(realspace_grid_desc_type), POINTER :: rs_desc => NULL()
69 : TYPE(pw_poisson_type), POINTER :: poisson_env => NULL()
70 : TYPE(dg_type), POINTER :: dg => NULL()
71 : END TYPE ewald_pw_type
72 :
73 : CONTAINS
74 :
75 : ! **************************************************************************************************
76 : !> \brief creates the structure ewald_pw_type
77 : !> \param ewald_pw ...
78 : !> \param ewald_env ...
79 : !> \param cell ...
80 : !> \param cell_ref ...
81 : !> \param print_section ...
82 : ! **************************************************************************************************
83 4169 : SUBROUTINE ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section)
84 : TYPE(ewald_pw_type), INTENT(OUT) :: ewald_pw
85 : TYPE(ewald_environment_type), POINTER :: ewald_env
86 : TYPE(cell_type), POINTER :: cell, cell_ref
87 : TYPE(section_vals_type), POINTER :: print_section
88 :
89 : NULLIFY (ewald_pw%pw_big_pool)
90 : NULLIFY (ewald_pw%pw_small_pool)
91 : NULLIFY (ewald_pw%rs_desc)
92 : NULLIFY (ewald_pw%poisson_env)
93 4169 : ALLOCATE (ewald_pw%dg)
94 4169 : CALL dg_create(ewald_pw%dg)
95 4169 : CALL ewald_pw_init(ewald_pw, ewald_env, cell, cell_ref, print_section)
96 4169 : END SUBROUTINE ewald_pw_create
97 :
98 : ! **************************************************************************************************
99 : !> \brief releases the memory used by the ewald_pw
100 : !> \param ewald_pw ...
101 : ! **************************************************************************************************
102 4169 : SUBROUTINE ewald_pw_release(ewald_pw)
103 : TYPE(ewald_pw_type), INTENT(INOUT) :: ewald_pw
104 :
105 4169 : CALL pw_pool_release(ewald_pw%pw_small_pool)
106 4169 : CALL pw_pool_release(ewald_pw%pw_big_pool)
107 4169 : CALL rs_grid_release_descriptor(ewald_pw%rs_desc)
108 4169 : IF (ASSOCIATED(ewald_pw%poisson_env)) THEN
109 2302 : CALL ewald_pw%poisson_env%release()
110 2302 : DEALLOCATE (ewald_pw%poisson_env)
111 : END IF
112 4169 : CALL dg_release(ewald_pw%dg)
113 4169 : DEALLOCATE (ewald_pw%dg)
114 :
115 4169 : END SUBROUTINE ewald_pw_release
116 :
117 : ! **************************************************************************************************
118 : !> \brief ...
119 : !> \param ewald_pw ...
120 : !> \param ewald_env ...
121 : !> \param cell ...
122 : !> \param cell_ref ...
123 : !> \param print_section ...
124 : !> \par History
125 : !> JGH (12-Jan-2001): Added SPME part
126 : !> JGH (15-Mar-2001): Work newly distributed between initialize, setup,
127 : !> and force routine
128 : !> \author CJM
129 : ! **************************************************************************************************
130 4169 : SUBROUTINE ewald_pw_init(ewald_pw, ewald_env, cell, cell_ref, print_section)
131 : TYPE(ewald_pw_type), INTENT(INOUT) :: ewald_pw
132 : TYPE(ewald_environment_type), POINTER :: ewald_env
133 : TYPE(cell_type), POINTER :: cell, cell_ref
134 : TYPE(section_vals_type), POINTER :: print_section
135 :
136 : CHARACTER(len=*), PARAMETER :: routineN = 'ewald_pw_init'
137 :
138 : INTEGER :: bo(2, 3), ewald_type, gmax(3), handle, &
139 : npts_s(3), ns_max, o_spline, &
140 : output_unit
141 : REAL(KIND=dp) :: alpha, alphasq, cutoff_radius, epsilon, &
142 : norm
143 : TYPE(cp_logger_type), POINTER :: logger
144 : TYPE(mp_para_env_type), POINTER :: para_env
145 : TYPE(pw_grid_type), POINTER :: pw_big_grid, pw_small_grid
146 : TYPE(pw_poisson_parameter_type) :: poisson_params
147 4169 : TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pw_pools
148 : TYPE(pw_pool_type), POINTER :: pw_pool
149 : TYPE(realspace_grid_desc_type), POINTER :: rs_desc
150 : TYPE(realspace_grid_input_type) :: input_settings
151 : TYPE(section_vals_type), POINTER :: poisson_section, rs_grid_section
152 :
153 4169 : CALL timeset(routineN, handle)
154 :
155 4169 : NULLIFY (pw_big_grid)
156 4169 : NULLIFY (pw_small_grid, poisson_section)
157 :
158 4169 : CPASSERT(ASSOCIATED(ewald_env))
159 4169 : CPASSERT(ASSOCIATED(cell))
160 : CALL ewald_env_get(ewald_env=ewald_env, &
161 : para_env=para_env, &
162 : gmax=gmax, alpha=alpha, &
163 : ns_max=ns_max, &
164 : ewald_type=ewald_type, &
165 : o_spline=o_spline, &
166 : poisson_section=poisson_section, &
167 4169 : epsilon=epsilon)
168 :
169 4169 : rs_grid_section => section_vals_get_subs_vals(poisson_section, "EWALD%RS_GRID")
170 :
171 821 : SELECT CASE (ewald_type)
172 : CASE (do_ewald_ewald)
173 : ! set up Classic EWALD sum
174 821 : logger => cp_get_default_logger()
175 821 : output_unit = cp_print_key_unit_nr(logger, print_section, "", extension=".Log")
176 821 : CALL pw_grid_create(pw_big_grid, mp_comm_self)
177 :
178 3284 : IF (ANY(gmax == 2*(gmax/2))) THEN
179 0 : CPABORT("gmax has to be odd.")
180 : END IF
181 3284 : bo(1, :) = -gmax/2
182 3284 : bo(2, :) = +gmax/2
183 : CALL pw_grid_setup(cell_ref%hmat, pw_big_grid, grid_span=HALFSPACE, bounds=bo, spherical=.TRUE., &
184 821 : fft_usage=.FALSE., iounit=output_unit)
185 821 : NULLIFY (pw_pool)
186 821 : CALL pw_pool_create(pw_pool, pw_grid=pw_big_grid)
187 821 : ewald_pw%pw_big_pool => pw_pool
188 821 : CALL pw_grid_release(pw_big_grid)
189 821 : CALL cp_print_key_finished_output(output_unit, logger, print_section, "")
190 :
191 : CASE (do_ewald_pme)
192 : ! set up Particle-Mesh EWALD sum
193 86 : logger => cp_get_default_logger()
194 86 : output_unit = cp_print_key_unit_nr(logger, print_section, "", extension=".Log")
195 86 : IF (.NOT. ASSOCIATED(ewald_pw%poisson_env)) THEN
196 1376 : ALLOCATE (ewald_pw%poisson_env)
197 86 : CALL ewald_pw%poisson_env%create()
198 : END IF
199 86 : CALL pw_grid_create(pw_small_grid, mp_comm_self)
200 86 : CALL pw_grid_create(pw_big_grid, para_env)
201 86 : IF (ns_max == 2*(ns_max/2)) THEN
202 0 : CPABORT("ns_max has to be odd.")
203 : END IF
204 344 : npts_s(:) = ns_max
205 : ! compute cut-off radius
206 86 : alphasq = alpha**2
207 86 : norm = (2.0_dp*alphasq/pi)**(1.5_dp)
208 86 : cutoff_radius = exp_radius(0, 2.0_dp*alphasq, epsilon, norm)
209 :
210 : CALL dg_pme_grid_setup(cell_ref%hmat, npts_s, cutoff_radius, &
211 : pw_small_grid, pw_big_grid, rs_dims=(/para_env%num_pe, 1/), &
212 258 : iounit=output_unit, fft_usage=.TRUE.)
213 : ! Write some useful info
214 86 : IF (output_unit > 0) THEN
215 : WRITE (output_unit, '( A,T71,E10.4 )') &
216 43 : ' EWALD| Gaussian tolerance (effective) ', epsilon
217 : WRITE (output_unit, '( A,T63,3I6 )') &
218 43 : ' EWALD| Small box grid ', pw_small_grid%npts
219 : WRITE (output_unit, '( A,T63,3I6 )') &
220 43 : ' EWALD| Full box grid ', pw_big_grid%npts
221 : END IF
222 :
223 : ! pw pools initialized
224 86 : NULLIFY (pw_pool)
225 86 : CALL pw_pool_create(pw_pool, pw_grid=pw_big_grid)
226 86 : ewald_pw%pw_big_pool => pw_pool
227 :
228 86 : NULLIFY (pw_pool)
229 86 : CALL pw_pool_create(pw_pool, pw_grid=pw_small_grid)
230 86 : ewald_pw%pw_small_pool => pw_pool
231 :
232 86 : NULLIFY (rs_desc)
233 : CALL init_input_type(input_settings, nsmax=MAXVAL(pw_small_grid%npts(1:3)), &
234 : rs_grid_section=rs_grid_section, ilevel=1, &
235 344 : higher_grid_layout=(/-1, -1, -1/))
236 86 : CALL rs_grid_create_descriptor(rs_desc, pw_big_grid, input_settings)
237 :
238 258 : BLOCK
239 1376 : TYPE(realspace_grid_type) :: rs
240 86 : CALL rs_grid_create(rs, rs_desc)
241 86 : CALL rs_grid_print(rs, output_unit)
242 258 : CALL rs_grid_release(rs)
243 : END BLOCK
244 :
245 86 : CALL cp_print_key_finished_output(output_unit, logger, print_section, "")
246 :
247 86 : ewald_pw%rs_desc => rs_desc
248 :
249 86 : CALL rs_grid_retain_descriptor(ewald_pw%rs_desc)
250 86 : CALL rs_grid_release_descriptor(rs_desc)
251 :
252 86 : CALL pw_grid_release(pw_small_grid)
253 86 : CALL pw_grid_release(pw_big_grid)
254 :
255 : CASE (do_ewald_spme)
256 : ! set up the Smooth-Particle-Mesh EWALD sum
257 2216 : logger => cp_get_default_logger()
258 2216 : output_unit = cp_print_key_unit_nr(logger, print_section, "", extension=".Log")
259 2216 : IF (.NOT. ASSOCIATED(ewald_pw%poisson_env)) THEN
260 35456 : ALLOCATE (ewald_pw%poisson_env)
261 2216 : CALL ewald_pw%poisson_env%create()
262 : END IF
263 2216 : CALL pw_grid_create(pw_big_grid, para_env)
264 2216 : npts_s = gmax
265 : CALL pw_grid_setup(cell_ref%hmat, pw_big_grid, grid_span=HALFSPACE, npts=npts_s, spherical=.TRUE., &
266 6648 : rs_dims=(/para_env%num_pe, 1/), iounit=output_unit, fft_usage=.TRUE.)
267 :
268 : ! pw pools initialized
269 2216 : NULLIFY (pw_pool)
270 2216 : CALL pw_pool_create(pw_pool, pw_grid=pw_big_grid)
271 2216 : ewald_pw%pw_big_pool => pw_pool
272 :
273 2216 : NULLIFY (rs_desc)
274 : CALL init_input_type(input_settings, nsmax=o_spline, &
275 : rs_grid_section=rs_grid_section, ilevel=1, &
276 2216 : higher_grid_layout=(/-1, -1, -1/))
277 2216 : CALL rs_grid_create_descriptor(rs_desc, pw_big_grid, input_settings)
278 :
279 6648 : BLOCK
280 35456 : TYPE(realspace_grid_type) :: rs
281 :
282 2216 : CALL rs_grid_create(rs, rs_desc)
283 2216 : CALL rs_grid_print(rs, output_unit)
284 6648 : CALL rs_grid_release(rs)
285 : END BLOCK
286 2216 : CALL cp_print_key_finished_output(output_unit, logger, print_section, "")
287 :
288 2216 : ewald_pw%rs_desc => rs_desc
289 :
290 2216 : CALL rs_grid_retain_descriptor(ewald_pw%rs_desc)
291 2216 : CALL rs_grid_release_descriptor(rs_desc)
292 :
293 2216 : CALL pw_grid_release(pw_big_grid)
294 : CASE (do_ewald_none)
295 : ! No EWALD sums..
296 : CASE default
297 4169 : CPABORT("")
298 : END SELECT
299 : ! Poisson Environment
300 4169 : IF (ASSOCIATED(ewald_pw%poisson_env)) THEN
301 4604 : ALLOCATE (pw_pools(1))
302 2302 : pw_pools(1)%pool => ewald_pw%pw_big_pool
303 2302 : CALL pw_poisson_read_parameters(poisson_section, poisson_params)
304 2302 : poisson_params%ewald_type = ewald_type
305 2302 : poisson_params%ewald_o_spline = o_spline
306 2302 : poisson_params%ewald_alpha = alpha
307 : CALL pw_poisson_set(ewald_pw%poisson_env, cell_hmat=cell%hmat, parameters=poisson_params, &
308 2302 : use_level=1, pw_pools=pw_pools)
309 2302 : DEALLOCATE (pw_pools)
310 : END IF
311 4169 : CALL timestop(handle)
312 29183 : END SUBROUTINE ewald_pw_init
313 :
314 : ! **************************************************************************************************
315 : !> \brief get the ewald_pw environment to the correct program.
316 : !> \param ewald_pw ...
317 : !> \param pw_big_pool ...
318 : !> \param pw_small_pool ...
319 : !> \param rs_desc ...
320 : !> \param poisson_env ...
321 : !> \param dg ...
322 : !> \author CJM
323 : ! **************************************************************************************************
324 94318 : SUBROUTINE ewald_pw_get(ewald_pw, pw_big_pool, pw_small_pool, rs_desc, poisson_env, dg)
325 :
326 : TYPE(ewald_pw_type), INTENT(IN) :: ewald_pw
327 : TYPE(pw_pool_type), OPTIONAL, POINTER :: pw_big_pool, pw_small_pool
328 : TYPE(realspace_grid_desc_type), OPTIONAL, POINTER :: rs_desc
329 : TYPE(pw_poisson_type), OPTIONAL, POINTER :: poisson_env
330 : TYPE(dg_type), OPTIONAL, POINTER :: dg
331 :
332 94318 : IF (PRESENT(poisson_env)) poisson_env => ewald_pw%poisson_env
333 94318 : IF (PRESENT(pw_big_pool)) pw_big_pool => ewald_pw%pw_big_pool
334 94318 : IF (PRESENT(pw_small_pool)) pw_small_pool => ewald_pw%pw_small_pool
335 94318 : IF (PRESENT(rs_desc)) rs_desc => ewald_pw%rs_desc
336 94318 : IF (PRESENT(dg)) dg => ewald_pw%dg
337 :
338 94318 : END SUBROUTINE ewald_pw_get
339 :
340 : ! **************************************************************************************************
341 : !> \brief set the ewald_pw environment to the correct program.
342 : !> \param ewald_pw ...
343 : !> \param pw_big_pool ...
344 : !> \param pw_small_pool ...
345 : !> \param rs_desc ...
346 : !> \param dg ...
347 : !> \param poisson_env ...
348 : !> \author CJM
349 : ! **************************************************************************************************
350 10726 : SUBROUTINE ewald_pw_set(ewald_pw, pw_big_pool, pw_small_pool, rs_desc, dg, &
351 : poisson_env)
352 :
353 : TYPE(ewald_pw_type), INTENT(INOUT) :: ewald_pw
354 : TYPE(pw_pool_type), OPTIONAL, POINTER :: pw_big_pool, pw_small_pool
355 : TYPE(realspace_grid_desc_type), OPTIONAL, POINTER :: rs_desc
356 : TYPE(dg_type), OPTIONAL, POINTER :: dg
357 : TYPE(pw_poisson_type), OPTIONAL, POINTER :: poisson_env
358 :
359 10726 : IF (PRESENT(pw_big_pool)) THEN
360 10726 : CALL pw_big_pool%retain()
361 10726 : CALL pw_pool_release(ewald_pw%pw_big_pool)
362 10726 : ewald_pw%pw_big_pool => pw_big_pool
363 : END IF
364 10726 : IF (PRESENT(pw_small_pool)) THEN
365 86 : CALL pw_small_pool%retain()
366 86 : CALL pw_pool_release(ewald_pw%pw_small_pool)
367 86 : ewald_pw%pw_small_pool => pw_small_pool
368 : END IF
369 10726 : IF (PRESENT(rs_desc)) THEN
370 0 : CALL rs_grid_retain_descriptor(rs_desc)
371 0 : CALL rs_grid_release_descriptor(ewald_pw%rs_desc)
372 0 : ewald_pw%rs_desc => rs_desc
373 : END IF
374 10726 : IF (PRESENT(dg)) THEN
375 0 : CALL dg_release(ewald_pw%dg)
376 0 : ewald_pw%dg => dg
377 : END IF
378 10726 : IF (PRESENT(poisson_env)) THEN
379 10726 : IF (ASSOCIATED(ewald_pw%poisson_env)) THEN
380 8876 : IF (.NOT. ASSOCIATED(ewald_pw%poisson_env, poisson_env)) THEN
381 0 : CALL ewald_pw%poisson_env%release()
382 0 : DEALLOCATE (ewald_pw%poisson_env)
383 : END IF
384 : END IF
385 10726 : ewald_pw%poisson_env => poisson_env
386 : END IF
387 :
388 10726 : END SUBROUTINE ewald_pw_set
389 :
390 0 : END MODULE ewald_pw_types
|