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