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 : !> \brief Bounded geometry-only cache for native atom-grid interpolation.
8 : MODULE qs_native_grid_cache
9 : USE cell_types, ONLY: cell_type
10 : USE kinds, ONLY: dp,&
11 : int_8
12 : USE pw_grid_types, ONLY: pw_grid_type
13 :
14 : IMPLICIT NONE
15 : PRIVATE
16 :
17 : INTEGER, PARAMETER :: native_grid_interp_offset_min = -5, native_grid_interp_offset_max = 6
18 : INTEGER, PARAMETER :: native_grid_interp_npts = native_grid_interp_offset_max - native_grid_interp_offset_min + 1
19 : INTEGER(KIND=int_8), PARAMETER :: cache_budget = 128_int_8*1024*1024
20 :
21 : TYPE native_grid_interpolation_stencil_type
22 : INTEGER :: relative_index(native_grid_interp_npts, 3) = 0
23 : REAL(dp) :: weight(native_grid_interp_npts, 3) = 0.0_dp
24 : LOGICAL :: valid(native_grid_interp_npts, 3) = .FALSE.
25 : LOGICAL :: active = .FALSE.
26 : END TYPE native_grid_interpolation_stencil_type
27 :
28 : TYPE stencil_entry_type
29 : TYPE(native_grid_interpolation_stencil_type) :: stencil
30 : REAL(dp) :: point(3) = 0.0_dp
31 : LOGICAL :: ready = .FALSE.
32 : END TYPE stencil_entry_type
33 :
34 : TYPE native_grid_cache_type
35 : PRIVATE
36 : TYPE(stencil_entry_type), ALLOCATABLE :: entries(:)
37 : REAL(dp) :: dh_inv(3, 3) = 0.0_dp, hmat(3, 3) = 0.0_dp
38 : INTEGER :: npts(3) = 0, periodic(3) = 0, nrows = -1
39 : LOGICAL :: wrap = .FALSE., prepared = .FALSE.
40 : END TYPE native_grid_cache_type
41 :
42 : PUBLIC :: native_grid_cache_type, native_grid_interpolation_stencil_type
43 : PUBLIC :: native_grid_interp_offset_min, native_grid_interp_offset_max, native_grid_interp_npts
44 : PUBLIC :: prepare_native_grid_cache, release_native_grid_cache, fetch_native_grid_stencil, store_native_grid_stencil
45 : CONTAINS
46 : ! **************************************************************************************************
47 : !> \brief Prepare outside parallel regions. Cache a bounded prefix and compute the rest normally.
48 : !> \param cache per-QS-environment cache
49 : !> \param grid current auxiliary grid
50 : !> \param cell current cell and periodicity
51 : !> \param nrows local atom-grid row count
52 : !> \param wrap wrap the auxiliary cell in nonperiodic directions
53 : !> \param max_bytes optional memory budget for testing, in bytes
54 : ! **************************************************************************************************
55 264 : SUBROUTINE prepare_native_grid_cache(cache, grid, cell, nrows, wrap, max_bytes)
56 : TYPE(native_grid_cache_type), INTENT(INOUT) :: cache
57 : TYPE(pw_grid_type), INTENT(IN) :: grid
58 : TYPE(cell_type), INTENT(IN) :: cell
59 : INTEGER, INTENT(IN) :: nrows
60 : LOGICAL, INTENT(IN) :: wrap
61 : INTEGER(KIND=int_8), INTENT(IN), OPTIONAL :: max_bytes
62 :
63 : INTEGER :: count, ierr
64 : INTEGER(KIND=int_8) :: budget, bytes_per_entry
65 : LOGICAL :: unchanged
66 : TYPE(stencil_entry_type) :: entry
67 :
68 264 : budget = cache_budget
69 0 : IF (PRESENT(max_bytes)) budget = MAX(0_int_8, max_bytes)
70 264 : bytes_per_entry = INT((STORAGE_SIZE(ENTRY) + 7)/8, int_8)
71 264 : count = INT(MIN(INT(MAX(0, nrows), int_8), budget/bytes_per_entry))
72 264 : unchanged = .FALSE.
73 264 : IF (ALLOCATED(cache%entries)) THEN
74 196 : unchanged = SIZE(cache%entries) == count
75 196 : IF (.NOT. unchanged) CALL release_native_grid_cache(cache)
76 : END IF
77 264 : IF (.NOT. ALLOCATED(cache%entries)) THEN
78 140893 : ALLOCATE (cache%entries(count), STAT=ierr)
79 68 : IF (ierr /= 0) RETURN
80 : END IF
81 332 : unchanged = unchanged .AND. cache%prepared .AND. cache%nrows == nrows
82 3348 : unchanged = unchanged .AND. ALL(cache%dh_inv == grid%dh_inv) .AND. ALL(cache%npts == grid%npts)
83 3270 : unchanged = unchanged .AND. ALL(cache%hmat == cell%hmat) .AND. ALL(cache%periodic == cell%perd)
84 168 : unchanged = unchanged .AND. (cache%wrap .EQV. wrap)
85 140476 : IF (.NOT. unchanged) cache%entries%ready = .FALSE.
86 3432 : cache%dh_inv = grid%dh_inv
87 1056 : cache%npts = grid%npts
88 3432 : cache%hmat = cell%hmat
89 1056 : cache%periodic = cell%perd
90 264 : cache%nrows = nrows
91 264 : cache%wrap = wrap
92 264 : cache%prepared = .TRUE.
93 31944 : END SUBROUTINE prepare_native_grid_cache
94 :
95 : ! **************************************************************************************************
96 : !> \brief Release geometry storage, including when the QS environment is only partly released.
97 : !> \param cache ...
98 : ! **************************************************************************************************
99 9201 : SUBROUTINE release_native_grid_cache(cache)
100 : TYPE(native_grid_cache_type), INTENT(INOUT) :: cache
101 :
102 9201 : IF (ALLOCATED(cache%entries)) DEALLOCATE (cache%entries)
103 9201 : cache%prepared = .FALSE.
104 9201 : END SUBROUTINE release_native_grid_cache
105 :
106 : ! **************************************************************************************************
107 : !> \brief Read a stencil only on an exact coordinate match. Safe for concurrent readers.
108 : !> \param cache ...
109 : !> \param row local grid row
110 : !> \param point current Cartesian coordinates, including atom motion and quadrature changes
111 : !> \param stencil cached interpolation stencil on a hit
112 : !> \return whether the cached stencil matches
113 : ! **************************************************************************************************
114 108149360 : FUNCTION fetch_native_grid_stencil(cache, row, point, stencil) RESULT(hit)
115 : TYPE(native_grid_cache_type), INTENT(IN) :: cache
116 : INTEGER, INTENT(IN) :: row
117 : REAL(dp), INTENT(IN) :: point(3)
118 : TYPE(native_grid_interpolation_stencil_type), &
119 : INTENT(OUT) :: stencil
120 : LOGICAL :: hit
121 :
122 916520 : hit = .FALSE.
123 916520 : IF (.NOT. cache%prepared .OR. .NOT. ALLOCATED(cache%entries)) RETURN
124 916520 : IF (row < 1 .OR. row > SIZE(cache%entries)) RETURN
125 916520 : IF (.NOT. cache%entries(row)%ready) RETURN
126 3094560 : IF (ANY(cache%entries(row)%point /= point)) RETURN
127 766140 : stencil = cache%entries(row)%stencil
128 766140 : hit = .TRUE.
129 766140 : END FUNCTION fetch_native_grid_stencil
130 :
131 : ! **************************************************************************************************
132 : !> \brief Store from the unique forward owner of a row, never from concurrent adjoint readers.
133 : !> \param cache ...
134 : !> \param row local grid row
135 : !> \param point current Cartesian coordinates
136 : !> \param stencil complete value stencil, not an indices-only stencil
137 : ! **************************************************************************************************
138 150380 : SUBROUTINE store_native_grid_stencil(cache, row, point, stencil)
139 : TYPE(native_grid_cache_type), INTENT(INOUT) :: cache
140 : INTEGER, INTENT(IN) :: row
141 : REAL(dp), INTENT(IN) :: point(3)
142 : TYPE(native_grid_interpolation_stencil_type), &
143 : INTENT(IN) :: stencil
144 :
145 150380 : IF (.NOT. cache%prepared .OR. .NOT. ALLOCATED(cache%entries)) RETURN
146 150380 : IF (row < 1 .OR. row > SIZE(cache%entries)) RETURN
147 150380 : cache%entries(row)%stencil = stencil
148 601520 : cache%entries(row)%point = point
149 150380 : cache%entries(row)%ready = .TRUE.
150 : END SUBROUTINE store_native_grid_stencil
151 0 : END MODULE qs_native_grid_cache
|
