Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Per-nnp persistent neighbour-interface state for the NNP hot path.
10 : !> Separates neighbour bookkeeping from the ACSF and network loops:
11 : !> species-pair routing is precomputed once and per-element work
12 : !> buffers are reused. State lives on the parent nnp_type, so each
13 : !> &NNP force_eval keeps its own copy and the lifetime tracks
14 : !> nnp_env_release.
15 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
16 : !> \author Christoph Schran (christoph.schran@rub.de)
17 : !> \date 2026-05-21
18 : ! **************************************************************************************************
19 : MODULE nnp_neighbor_interface
20 :
21 : USE kinds, ONLY: dp
22 : USE nnp_environment_types, ONLY: nnp_dGdr_grp_type,&
23 : nnp_neigh_grp_type,&
24 : nnp_neighbor_interface_state_release,&
25 : nnp_neighbor_workspace_type,&
26 : nnp_type
27 :
28 : !$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num
29 : #include "./base/base_uses.f90"
30 :
31 : IMPLICIT NONE
32 :
33 : PRIVATE
34 :
35 : PUBLIC :: nnp_neighbor_interface_prepare, &
36 : nnp_neighbor_interface_reset_neighbor, &
37 : nnp_grp_grow_dGdr, &
38 : nnp_neigh_grp_grow, &
39 : nnp_workspace_grow_caches
40 :
41 : CONTAINS
42 :
43 : ! **************************************************************************************************
44 : !> \brief Ensure pair-routing metadata and reusable workspaces are ready for the current NNP model.
45 : !> \param nnp NNP environment whose neighbor_interface_state will be (re)built if needed.
46 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
47 : ! **************************************************************************************************
48 55492 : SUBROUTINE nnp_neighbor_interface_prepare(nnp)
49 :
50 : TYPE(nnp_type), INTENT(INOUT) :: nnp
51 :
52 : CHARACTER(len=*), PARAMETER :: routineN = 'nnp_neighbor_interface_prepare'
53 :
54 : INTEGER :: handle
55 : LOGICAL :: rebuild
56 :
57 55492 : CALL timeset(routineN, handle)
58 :
59 55492 : rebuild = .NOT. nnp%neighbor_interface_state%initialized
60 55492 : IF (.NOT. rebuild) CALL nnp_neighbor_interface_needs_rebuild(nnp, rebuild)
61 :
62 55492 : IF (rebuild) THEN
63 17 : CALL nnp_neighbor_interface_state_release(nnp%neighbor_interface_state)
64 17 : CALL nnp_neighbor_interface_build_pair_maps(nnp)
65 : END IF
66 :
67 55492 : CALL timestop(handle)
68 :
69 55492 : END SUBROUTINE nnp_neighbor_interface_prepare
70 :
71 : ! **************************************************************************************************
72 : !> \brief Reset per-group neighbour counters for one central element before refilling
73 : !> the reusable buffers. Zeroes n_rad/n_ang1/n_ang2; the (ind, dist) slabs
74 : !> are kept allocated and overwritten in place by the next push.
75 : !> \param nnp NNP environment whose neighbor_interface_state will have its counters cleared.
76 : !> \param ind central-element index (1..nnp%n_ele) whose workspace is reset.
77 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
78 : ! **************************************************************************************************
79 252785 : SUBROUTINE nnp_neighbor_interface_reset_neighbor(nnp, ind)
80 :
81 : TYPE(nnp_type), INTENT(INOUT) :: nnp
82 : INTEGER, INTENT(IN) :: ind
83 :
84 : INTEGER :: tid
85 :
86 252785 : tid = 1
87 252785 : !$ tid = omp_get_thread_num() + 1
88 868699 : nnp%neighbor_interface_state%workspace(ind, tid)%neighbor%n_rad(:) = 0
89 769054 : nnp%neighbor_interface_state%workspace(ind, tid)%neighbor%n_ang1(:) = 0
90 769054 : nnp%neighbor_interface_state%workspace(ind, tid)%neighbor%n_ang2(:) = 0
91 :
92 252785 : END SUBROUTINE nnp_neighbor_interface_reset_neighbor
93 :
94 : ! **************************************************************************************************
95 : !> \brief Check whether the persistent state matches the current NNP model. Assumes
96 : !> the per-group routing (ele_ind, cutoff) is fixed after nnp_init_acsf_groups,
97 : !> so only element and per-element group counts are compared.
98 : !> \param nnp NNP environment whose persistent neighbour-interface state will be inspected
99 : !> \param rebuild (out) .TRUE. if element count or per-element group sizes have changed
100 : ! **************************************************************************************************
101 55475 : SUBROUTINE nnp_neighbor_interface_needs_rebuild(nnp, rebuild)
102 :
103 : TYPE(nnp_type), INTENT(IN) :: nnp
104 : LOGICAL, INTENT(OUT) :: rebuild
105 :
106 : INTEGER :: i
107 :
108 55475 : rebuild = .FALSE.
109 55475 : IF (nnp%neighbor_interface_state%n_ele /= nnp%n_ele) THEN
110 0 : rebuild = .TRUE.
111 0 : RETURN
112 : END IF
113 55475 : IF (.NOT. ALLOCATED(nnp%neighbor_interface_state%n_rad)) THEN
114 0 : rebuild = .TRUE.
115 0 : RETURN
116 : END IF
117 :
118 221596 : DO i = 1, nnp%n_ele
119 : IF (nnp%neighbor_interface_state%n_rad(i) /= nnp%n_rad(i) .OR. &
120 : nnp%neighbor_interface_state%n_ang(i) /= nnp%n_ang(i) .OR. &
121 166121 : nnp%neighbor_interface_state%n_radgrp(i) /= nnp%rad(i)%n_symfgrp .OR. &
122 55475 : nnp%neighbor_interface_state%n_anggrp(i) /= nnp%ang(i)%n_symfgrp) THEN
123 0 : rebuild = .TRUE.
124 0 : RETURN
125 : END IF
126 : END DO
127 :
128 : END SUBROUTINE nnp_neighbor_interface_needs_rebuild
129 :
130 : ! **************************************************************************************************
131 : !> \brief Build species-pair routing tables and initialize reusable workspaces.
132 : !> \param nnp NNP environment; pair_map and workspace arrays are (re-)allocated on its neighbor_interface_state
133 : ! **************************************************************************************************
134 17 : SUBROUTINE nnp_neighbor_interface_build_pair_maps(nnp)
135 :
136 : TYPE(nnp_type), INTENT(INOUT) :: nnp
137 :
138 : INTEGER :: i, nthreads, t
139 :
140 : ASSOCIATE (state => nnp%neighbor_interface_state)
141 17 : state%n_ele = nnp%n_ele
142 51 : ALLOCATE (state%n_rad(nnp%n_ele))
143 34 : ALLOCATE (state%n_ang(nnp%n_ele))
144 34 : ALLOCATE (state%n_radgrp(nnp%n_ele))
145 34 : ALLOCATE (state%n_anggrp(nnp%n_ele))
146 176 : ALLOCATE (state%pair_map(nnp%n_ele, nnp%n_ele))
147 : ! one workspace column per thread of the atom loop in nnp_calc_energy_force
148 17 : nthreads = 1
149 17 : !$ nthreads = omp_get_max_threads()
150 171 : ALLOCATE (state%workspace(nnp%n_ele, nthreads))
151 :
152 52 : DO i = 1, nnp%n_ele
153 35 : state%n_rad(i) = nnp%n_rad(i)
154 35 : state%n_ang(i) = nnp%n_ang(i)
155 35 : state%n_radgrp(i) = nnp%rad(i)%n_symfgrp
156 52 : state%n_anggrp(i) = nnp%ang(i)%n_symfgrp
157 : END DO
158 :
159 52 : DO i = 1, nnp%n_ele
160 35 : CALL nnp_neighbor_interface_build_pair_map_for_element(nnp, i)
161 87 : DO t = 1, nthreads
162 70 : CALL nnp_neighbor_interface_init_workspace_metadata(nnp, i, t)
163 : END DO
164 : END DO
165 :
166 17 : state%initialized = .TRUE.
167 : END ASSOCIATE
168 :
169 17 : END SUBROUTINE nnp_neighbor_interface_build_pair_maps
170 :
171 : ! **************************************************************************************************
172 : !> \brief Build all species-pair routes for one central element.
173 : !> \param nnp NNP environment providing the radial/angular SF groups
174 : !> \param ind central-element index whose pair-map row is built
175 : ! **************************************************************************************************
176 35 : SUBROUTINE nnp_neighbor_interface_build_pair_map_for_element(nnp, ind)
177 :
178 : TYPE(nnp_type), INTENT(INOUT) :: nnp
179 : INTEGER, INTENT(IN) :: ind
180 :
181 : INTEGER :: idx, neighbor_ind, s
182 :
183 108 : DO neighbor_ind = 1, nnp%n_ele
184 35 : ASSOCIATE (pair_map => nnp%neighbor_interface_state%pair_map(ind, neighbor_ind))
185 73 : pair_map%n_rad = 0
186 73 : pair_map%n_ang1 = 0
187 73 : pair_map%n_ang2 = 0
188 73 : pair_map%max_relevant_cutoff = 0.0_dp
189 :
190 222 : DO s = 1, nnp%rad(ind)%n_symfgrp
191 222 : IF (nnp%rad(ind)%symfgrp(s)%ele_ind(1) == neighbor_ind) pair_map%n_rad = pair_map%n_rad + 1
192 : END DO
193 251 : DO s = 1, nnp%ang(ind)%n_symfgrp
194 178 : IF (nnp%ang(ind)%symfgrp(s)%ele_ind(1) == neighbor_ind) pair_map%n_ang1 = pair_map%n_ang1 + 1
195 251 : IF (nnp%ang(ind)%symfgrp(s)%ele_ind(2) == neighbor_ind) pair_map%n_ang2 = pair_map%n_ang2 + 1
196 : END DO
197 :
198 219 : ALLOCATE (pair_map%rad_groups(MAX(1, pair_map%n_rad)))
199 219 : ALLOCATE (pair_map%ang1_groups(MAX(1, pair_map%n_ang1)))
200 219 : ALLOCATE (pair_map%ang2_groups(MAX(1, pair_map%n_ang2)))
201 :
202 73 : idx = 0
203 222 : DO s = 1, nnp%rad(ind)%n_symfgrp
204 222 : IF (nnp%rad(ind)%symfgrp(s)%ele_ind(1) == neighbor_ind) THEN
205 71 : idx = idx + 1
206 71 : pair_map%rad_groups(idx) = s
207 71 : pair_map%max_relevant_cutoff = MAX(pair_map%max_relevant_cutoff, nnp%rad(ind)%symfgrp(s)%cutoff)
208 : END IF
209 : END DO
210 :
211 73 : idx = 0
212 251 : DO s = 1, nnp%ang(ind)%n_symfgrp
213 251 : IF (nnp%ang(ind)%symfgrp(s)%ele_ind(1) == neighbor_ind) THEN
214 86 : idx = idx + 1
215 86 : pair_map%ang1_groups(idx) = s
216 86 : pair_map%max_relevant_cutoff = MAX(pair_map%max_relevant_cutoff, nnp%ang(ind)%symfgrp(s)%cutoff)
217 : END IF
218 : END DO
219 :
220 73 : idx = 0
221 324 : DO s = 1, nnp%ang(ind)%n_symfgrp
222 251 : IF (nnp%ang(ind)%symfgrp(s)%ele_ind(2) == neighbor_ind) THEN
223 86 : idx = idx + 1
224 86 : pair_map%ang2_groups(idx) = s
225 86 : pair_map%max_relevant_cutoff = MAX(pair_map%max_relevant_cutoff, nnp%ang(ind)%symfgrp(s)%cutoff)
226 : END IF
227 : END DO
228 : END ASSOCIATE
229 : END DO
230 :
231 35 : END SUBROUTINE nnp_neighbor_interface_build_pair_map_for_element
232 :
233 : ! **************************************************************************************************
234 : !> \brief Cache max scratch sizes for one central element.
235 : !> \param nnp NNP environment providing per-element SF group definitions
236 : !> \param ind central-element index whose scratch workspace metadata is cached
237 : !> \param tid thread index whose workspace column is initialised
238 : ! **************************************************************************************************
239 35 : SUBROUTINE nnp_neighbor_interface_init_workspace_metadata(nnp, ind, tid)
240 :
241 : TYPE(nnp_type), INTENT(INOUT) :: nnp
242 : INTEGER, INTENT(IN) :: ind, tid
243 :
244 : INTEGER :: s
245 :
246 : ASSOCIATE (workspace => nnp%neighbor_interface_state%workspace(ind, tid))
247 35 : workspace%max_rad_symf = 0
248 35 : workspace%max_ang_symf = 0
249 35 : workspace%n_input_nodes = nnp%n_rad(ind) + nnp%n_ang(ind)
250 :
251 106 : DO s = 1, nnp%rad(ind)%n_symfgrp
252 106 : workspace%max_rad_symf = MAX(workspace%max_rad_symf, nnp%rad(ind)%symfgrp(s)%n_symf)
253 : END DO
254 121 : DO s = 1, nnp%ang(ind)%n_symfgrp
255 121 : workspace%max_ang_symf = MAX(workspace%max_ang_symf, nnp%ang(ind)%symfgrp(s)%n_symf)
256 : END DO
257 :
258 : ! Per-SF scratch reused inside nnp_calc_rad / nnp_calc_ang, sized by max_*_symf.
259 35 : IF (ALLOCATED(workspace%radial_sym)) DEALLOCATE (workspace%radial_sym)
260 35 : IF (ALLOCATED(workspace%radial_force)) DEALLOCATE (workspace%radial_force)
261 35 : IF (ALLOCATED(workspace%angular_sym)) DEALLOCATE (workspace%angular_sym)
262 35 : IF (ALLOCATED(workspace%angular_force)) DEALLOCATE (workspace%angular_force)
263 105 : ALLOCATE (workspace%radial_sym(MAX(1, workspace%max_rad_symf)))
264 105 : ALLOCATE (workspace%radial_force(3, MAX(1, workspace%max_rad_symf)))
265 105 : ALLOCATE (workspace%angular_sym(MAX(1, workspace%max_ang_symf)))
266 105 : ALLOCATE (workspace%angular_force(3, 3, MAX(1, workspace%max_ang_symf)))
267 :
268 : ! 1D angular cutoff caches; start small and grow lazily to the peak
269 : ! per-element angular neighbour count (nnp_workspace_grow_caches).
270 35 : IF (ALLOCATED(workspace%fc_cache1)) DEALLOCATE (workspace%fc_cache1)
271 35 : IF (ALLOCATED(workspace%dfc_cache1)) DEALLOCATE (workspace%dfc_cache1)
272 35 : IF (ALLOCATED(workspace%fc_cache2)) DEALLOCATE (workspace%fc_cache2)
273 35 : IF (ALLOCATED(workspace%dfc_cache2)) DEALLOCATE (workspace%dfc_cache2)
274 35 : workspace%cache_cap = 8
275 35 : ALLOCATE (workspace%fc_cache1(workspace%cache_cap))
276 35 : ALLOCATE (workspace%dfc_cache1(workspace%cache_cap))
277 35 : ALLOCATE (workspace%fc_cache2(workspace%cache_cap))
278 35 : ALLOCATE (workspace%dfc_cache2(workspace%cache_cap))
279 :
280 : ! self_dGdr sized once here; independent of the candidate pool.
281 35 : IF (ALLOCATED(workspace%self_dGdr)) DEALLOCATE (workspace%self_dGdr)
282 105 : ALLOCATE (workspace%self_dGdr(3, MAX(1, workspace%n_input_nodes)))
283 :
284 : ! Per-group neighbour counters and dense (ind, dist) containers, sized to
285 : ! n_symfgrp here; the per-group slabs grow lazily via nnp_neigh_grp_grow.
286 35 : CALL nnp_release_neighbor_local(workspace)
287 105 : ALLOCATE (workspace%neighbor%n_rad(MAX(1, nnp%rad(ind)%n_symfgrp)))
288 105 : ALLOCATE (workspace%neighbor%n_ang1(MAX(1, nnp%ang(ind)%n_symfgrp)))
289 105 : ALLOCATE (workspace%neighbor%n_ang2(MAX(1, nnp%ang(ind)%n_symfgrp)))
290 176 : ALLOCATE (workspace%neighbor%rad(MAX(1, nnp%rad(ind)%n_symfgrp)))
291 191 : ALLOCATE (workspace%neighbor%ang1(MAX(1, nnp%ang(ind)%n_symfgrp)))
292 191 : ALLOCATE (workspace%neighbor%ang2(MAX(1, nnp%ang(ind)%n_symfgrp)))
293 106 : workspace%neighbor%n_rad(:) = 0
294 121 : workspace%neighbor%n_ang1(:) = 0
295 121 : workspace%neighbor%n_ang2(:) = 0
296 140 : workspace%neighbor%pbc_copies = 0
297 :
298 : ! Per-group dG/dr buffers: container arrays sized to n_symfgrp and each
299 : ! entry's n_symf recorded. The %data slab stays unallocated until
300 : ! nnp_grp_grow_dGdr sees a real neighbour count.
301 35 : CALL nnp_release_dGdr_grp_array(workspace%dGdr_rad)
302 35 : CALL nnp_release_dGdr_grp_array(workspace%dGdr_ang_jj)
303 35 : CALL nnp_release_dGdr_grp_array(workspace%dGdr_ang_kk)
304 176 : ALLOCATE (workspace%dGdr_rad(MAX(1, nnp%rad(ind)%n_symfgrp)))
305 191 : ALLOCATE (workspace%dGdr_ang_jj(MAX(1, nnp%ang(ind)%n_symfgrp)))
306 191 : ALLOCATE (workspace%dGdr_ang_kk(MAX(1, nnp%ang(ind)%n_symfgrp)))
307 106 : DO s = 1, nnp%rad(ind)%n_symfgrp
308 71 : workspace%dGdr_rad(s)%n_symf = nnp%rad(ind)%symfgrp(s)%n_symf
309 106 : workspace%dGdr_rad(s)%cap = 0
310 : END DO
311 156 : DO s = 1, nnp%ang(ind)%n_symfgrp
312 86 : workspace%dGdr_ang_jj(s)%n_symf = nnp%ang(ind)%symfgrp(s)%n_symf
313 86 : workspace%dGdr_ang_jj(s)%cap = 0
314 86 : workspace%dGdr_ang_kk(s)%n_symf = nnp%ang(ind)%symfgrp(s)%n_symf
315 121 : workspace%dGdr_ang_kk(s)%cap = 0
316 : END DO
317 : END ASSOCIATE
318 :
319 35 : END SUBROUTINE nnp_neighbor_interface_init_workspace_metadata
320 :
321 : ! **************************************************************************************************
322 : !> \brief Release the workspace's neighbor%(rad,ang1,ang2) per-group slabs before
323 : !> re-allocation, keeping the persistent caches that survive a rebuild.
324 : !> \param workspace per-element scratch workspace whose neighbor%(rad,ang1,ang2) slabs will be released
325 : ! **************************************************************************************************
326 35 : SUBROUTINE nnp_release_neighbor_local(workspace)
327 :
328 : TYPE(nnp_neighbor_workspace_type), INTENT(INOUT) :: workspace
329 :
330 : INTEGER :: s
331 :
332 35 : IF (ALLOCATED(workspace%neighbor%rad)) THEN
333 0 : DO s = 1, SIZE(workspace%neighbor%rad)
334 0 : IF (ALLOCATED(workspace%neighbor%rad(s)%ind)) DEALLOCATE (workspace%neighbor%rad(s)%ind)
335 0 : IF (ALLOCATED(workspace%neighbor%rad(s)%dist)) DEALLOCATE (workspace%neighbor%rad(s)%dist)
336 : END DO
337 0 : DEALLOCATE (workspace%neighbor%rad)
338 : END IF
339 35 : IF (ALLOCATED(workspace%neighbor%ang1)) THEN
340 0 : DO s = 1, SIZE(workspace%neighbor%ang1)
341 0 : IF (ALLOCATED(workspace%neighbor%ang1(s)%ind)) DEALLOCATE (workspace%neighbor%ang1(s)%ind)
342 0 : IF (ALLOCATED(workspace%neighbor%ang1(s)%dist)) DEALLOCATE (workspace%neighbor%ang1(s)%dist)
343 : END DO
344 0 : DEALLOCATE (workspace%neighbor%ang1)
345 : END IF
346 35 : IF (ALLOCATED(workspace%neighbor%ang2)) THEN
347 0 : DO s = 1, SIZE(workspace%neighbor%ang2)
348 0 : IF (ALLOCATED(workspace%neighbor%ang2(s)%ind)) DEALLOCATE (workspace%neighbor%ang2(s)%ind)
349 0 : IF (ALLOCATED(workspace%neighbor%ang2(s)%dist)) DEALLOCATE (workspace%neighbor%ang2(s)%dist)
350 : END DO
351 0 : DEALLOCATE (workspace%neighbor%ang2)
352 : END IF
353 35 : IF (ALLOCATED(workspace%neighbor%n_rad)) DEALLOCATE (workspace%neighbor%n_rad)
354 35 : IF (ALLOCATED(workspace%neighbor%n_ang1)) DEALLOCATE (workspace%neighbor%n_ang1)
355 35 : IF (ALLOCATED(workspace%neighbor%n_ang2)) DEALLOCATE (workspace%neighbor%n_ang2)
356 140 : workspace%neighbor%pbc_copies = -1
357 :
358 35 : END SUBROUTINE nnp_release_neighbor_local
359 :
360 : ! **************************************************************************************************
361 : !> \brief Release every per-group dG/dr buffer in a container array, then deallocate the container.
362 : !> \param grps per-group dG/dr container array to release (each entry's %data is freed first)
363 : ! **************************************************************************************************
364 105 : SUBROUTINE nnp_release_dGdr_grp_array(grps)
365 :
366 : TYPE(nnp_dGdr_grp_type), ALLOCATABLE, &
367 : INTENT(INOUT) :: grps(:)
368 :
369 : INTEGER :: s
370 :
371 105 : IF (.NOT. ALLOCATED(grps)) RETURN
372 0 : DO s = 1, SIZE(grps)
373 0 : IF (ALLOCATED(grps(s)%data)) DEALLOCATE (grps(s)%data)
374 0 : grps(s)%cap = 0
375 0 : grps(s)%n_symf = 0
376 : END DO
377 0 : DEALLOCATE (grps)
378 :
379 : END SUBROUTINE nnp_release_dGdr_grp_array
380 :
381 : ! **************************************************************************************************
382 : !> \brief Ensure a per-group dG/dr buffer holds n_needed neighbours, growing by
383 : !> 1.5x (with an additive floor) so reallocation amortizes to O(1) over a
384 : !> trajectory. grp%data is allocated on return (cap >= 1), so callers can
385 : !> ASSOCIATE-bind it even for a zero-trip loop.
386 : !> \param grp per-element dGdr group whose %data slab is reallocated if undersized
387 : !> \param n_needed minimum required capacity (third dim of grp%data)
388 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
389 : ! **************************************************************************************************
390 157020 : SUBROUTINE nnp_grp_grow_dGdr(grp, n_needed)
391 :
392 : TYPE(nnp_dGdr_grp_type), INTENT(INOUT) :: grp
393 : INTEGER, INTENT(IN) :: n_needed
394 :
395 : INTEGER :: new_cap
396 :
397 157020 : IF (ALLOCATED(grp%data) .AND. grp%cap >= n_needed) RETURN
398 313 : IF (ALLOCATED(grp%data)) DEALLOCATE (grp%data)
399 313 : new_cap = MAX(MAX(1, n_needed), INT(grp%cap*1.5_dp) + 8)
400 1252 : ALLOCATE (grp%data(3, MAX(1, grp%n_symf), new_cap))
401 313 : grp%cap = new_cap
402 :
403 : END SUBROUTINE nnp_grp_grow_dGdr
404 :
405 : ! **************************************************************************************************
406 : !> \brief Ensure a per-group (ind, dist) neighbour buffer holds n_needed entries.
407 : !> Called from inside the linked-cell push loop, so existing entries are
408 : !> preserved on grow via MOVE_ALLOC. 1.5x growth amortizes reallocation to
409 : !> O(1) over a trajectory.
410 : !> \param grp per-species-pair neighbour group whose (ind, dist) buffers grow
411 : !> \param n_needed minimum required entry count
412 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
413 : ! **************************************************************************************************
414 645 : SUBROUTINE nnp_neigh_grp_grow(grp, n_needed)
415 :
416 : TYPE(nnp_neigh_grp_type), INTENT(INOUT) :: grp
417 : INTEGER, INTENT(IN) :: n_needed
418 :
419 : INTEGER :: n_old, new_cap
420 645 : INTEGER, ALLOCATABLE :: new_ind(:)
421 645 : REAL(KIND=dp), ALLOCATABLE :: new_dist(:, :)
422 :
423 645 : IF (ALLOCATED(grp%dist) .AND. grp%cap >= n_needed) RETURN
424 :
425 645 : new_cap = MAX(MAX(8, n_needed), INT(grp%cap*1.5_dp) + 8)
426 :
427 645 : IF (ALLOCATED(grp%dist)) THEN
428 476 : n_old = grp%cap
429 1428 : ALLOCATE (new_dist(4, new_cap))
430 1428 : ALLOCATE (new_ind(new_cap))
431 476 : IF (n_old > 0) THEN
432 64876 : new_dist(:, 1:n_old) = grp%dist(:, 1:n_old)
433 13356 : new_ind(1:n_old) = grp%ind(1:n_old)
434 : END IF
435 476 : CALL MOVE_ALLOC(new_dist, grp%dist)
436 476 : CALL MOVE_ALLOC(new_ind, grp%ind)
437 : ELSE
438 507 : ALLOCATE (grp%dist(4, new_cap))
439 507 : ALLOCATE (grp%ind(new_cap))
440 : END IF
441 645 : grp%cap = new_cap
442 :
443 : END SUBROUTINE nnp_neigh_grp_grow
444 :
445 : ! **************************************************************************************************
446 : !> \brief Ensure the four per-element angular cutoff caches hold at least n_needed
447 : !> entries. These 1D scratch arrays are reused across angular groups within
448 : !> one nnp_calc_acsf call, so the peak is MAX_s(n_ang1) and MAX_s(n_ang2).
449 : !> 1.5x lazy growth.
450 : !> \param workspace per-element workspace whose fc_cache1/dfc_cache1/fc_cache2/dfc_cache2 grow
451 : !> \param n_needed minimum required cache length
452 : !> \author Dhruv Sharma (ds2173@cam.ac.uk)
453 : ! **************************************************************************************************
454 204699 : SUBROUTINE nnp_workspace_grow_caches(workspace, n_needed)
455 :
456 : TYPE(nnp_neighbor_workspace_type), INTENT(INOUT) :: workspace
457 : INTEGER, INTENT(IN) :: n_needed
458 :
459 : INTEGER :: new_cap
460 :
461 204699 : IF (workspace%cache_cap >= n_needed) RETURN
462 :
463 42 : new_cap = MAX(MAX(8, n_needed), INT(workspace%cache_cap*1.5_dp) + 8)
464 : ! Contents not preserved: nnp_fill_fc_dfc_cache overwrites fully before any read.
465 42 : IF (ALLOCATED(workspace%fc_cache1)) DEALLOCATE (workspace%fc_cache1)
466 42 : IF (ALLOCATED(workspace%dfc_cache1)) DEALLOCATE (workspace%dfc_cache1)
467 42 : IF (ALLOCATED(workspace%fc_cache2)) DEALLOCATE (workspace%fc_cache2)
468 42 : IF (ALLOCATED(workspace%dfc_cache2)) DEALLOCATE (workspace%dfc_cache2)
469 126 : ALLOCATE (workspace%fc_cache1(new_cap))
470 84 : ALLOCATE (workspace%dfc_cache1(new_cap))
471 84 : ALLOCATE (workspace%fc_cache2(new_cap))
472 84 : ALLOCATE (workspace%dfc_cache2(new_cap))
473 42 : workspace%cache_cap = new_cap
474 :
475 : END SUBROUTINE nnp_workspace_grow_caches
476 :
477 : END MODULE nnp_neighbor_interface
|
