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 The types needed for the calculation of active space Hamiltonians
10 : !> \par History
11 : !> 04.2016 created [JGH]
12 : !> \author JGH
13 : ! **************************************************************************************************
14 : MODULE qs_active_space_types
15 : USE cp_dbcsr_api, ONLY: dbcsr_csr_destroy,&
16 : dbcsr_csr_p_type,&
17 : dbcsr_p_type
18 : USE cp_dbcsr_operations, ONLY: dbcsr_deallocate_matrix_set
19 : USE cp_fm_types, ONLY: cp_fm_release,&
20 : cp_fm_type
21 : USE input_section_types, ONLY: section_vals_type
22 : USE kinds, ONLY: default_path_length,&
23 : dp
24 : USE message_passing, ONLY: mp_comm_null,&
25 : mp_comm_type,&
26 : mp_para_env_release,&
27 : mp_para_env_type
28 : USE qs_mo_types, ONLY: deallocate_mo_set,&
29 : mo_set_type
30 : #include "./base/base_uses.f90"
31 :
32 : IMPLICIT NONE
33 : PRIVATE
34 :
35 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_active_space_types'
36 :
37 : PUBLIC :: active_space_type, eri_type, eri_type_eri_element_func
38 : PUBLIC :: create_active_space_type, release_active_space_type
39 : PUBLIC :: csr_idx_to_combined, csr_idx_from_combined
40 :
41 : ! **************************************************************************************************
42 : !> \brief Quantities needed for AS determination
43 : !> \author JGH
44 : ! **************************************************************************************************
45 : TYPE eri_gpw_type
46 : LOGICAL :: redo_poisson = .FALSE.
47 : LOGICAL :: store_wfn = .FALSE.
48 : REAL(KIND=dp) :: cutoff = 0.0_dp
49 : REAL(KIND=dp) :: rel_cutoff = 0.0_dp
50 : REAL(KIND=dp) :: eps_grid = 0.0_dp
51 : REAL(KIND=dp) :: eps_filter = 0.0_dp
52 : INTEGER :: print_level = 0
53 : INTEGER :: group_size = 0
54 : END TYPE eri_gpw_type
55 :
56 : TYPE eri_type
57 : INTEGER :: method = 0
58 : INTEGER :: operator = 0
59 : LOGICAL :: enlarge_cell = .FALSE.
60 : REAL(KIND=dp) :: omega = 0.0_dp
61 : INTEGER, DIMENSION(3) :: periodicity = 0
62 : REAL(KIND=dp), DIMENSION(3) :: eri_cell = 0
63 : REAL(KIND=dp), DIMENSION(3) :: eri_cell_angles = 0
64 : REAL(KIND=dp) :: cutoff_radius = 0.0_dp
65 : REAL(KIND=dp) :: eps_integral = 0.0_dp
66 : TYPE(eri_gpw_type) :: eri_gpw = eri_gpw_type()
67 : TYPE(dbcsr_csr_p_type), &
68 : DIMENSION(:), POINTER :: eri => NULL()
69 : INTEGER :: norb = 0
70 : TYPE(mp_para_env_type), POINTER :: para_env_sub => NULL()
71 : TYPE(mp_comm_type) :: comm_exchange = mp_comm_null
72 : CONTAINS
73 : PROCEDURE :: eri_foreach => eri_type_eri_foreach
74 : END TYPE eri_type
75 :
76 : ! **************************************************************************************************
77 : !> \brief Abstract function object for the `eri_type_eri_foreach` method
78 : ! **************************************************************************************************
79 : TYPE, ABSTRACT :: eri_type_eri_element_func
80 : CONTAINS
81 : PROCEDURE(eri_type_eri_element_func_interface), DEFERRED :: func
82 : END TYPE eri_type_eri_element_func
83 :
84 : TYPE active_space_type
85 : INTEGER :: nelec_active = 0
86 : INTEGER :: nelec_inactive = 0
87 : INTEGER :: nelec_total = 0
88 : INTEGER, POINTER, DIMENSION(:, :) :: active_orbitals => NULL()
89 : INTEGER, POINTER, DIMENSION(:, :) :: inactive_orbitals => NULL()
90 : INTEGER :: nmo_active = 0
91 : INTEGER :: nmo_inactive = 0
92 : INTEGER :: multiplicity = 0
93 : INTEGER :: nspins = 0
94 : LOGICAL :: molecule = .FALSE.
95 : INTEGER :: model = 0
96 : REAL(KIND=dp) :: energy_total = 0.0_dp
97 : REAL(KIND=dp) :: energy_ref = 0.0_dp
98 : REAL(KIND=dp) :: energy_inactive = 0.0_dp
99 : REAL(KIND=dp) :: energy_active = 0.0_dp
100 : REAL(KIND=dp) :: alpha = 0.0_dp
101 : LOGICAL :: do_scf_embedding = .FALSE.
102 : LOGICAL :: qcschema = .FALSE.
103 : LOGICAL :: fcidump = .FALSE.
104 : CHARACTER(LEN=default_path_length) :: qcschema_filename = ''
105 : TYPE(eri_type) :: eri = eri_type()
106 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_active => NULL()
107 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_inactive => NULL()
108 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: p_active => NULL()
109 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: ks_sub => NULL()
110 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: vxc_sub => NULL()
111 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: h_sub => NULL()
112 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: fock_sub => NULL()
113 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: sab_sub => NULL()
114 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: pmat_inactive => NULL()
115 : TYPE(section_vals_type), POINTER :: xc_section => NULL()
116 : END TYPE active_space_type
117 :
118 : ABSTRACT INTERFACE
119 : ! **************************************************************************************************
120 : !> \brief The function signature to be implemented by a child of `eri_type_eri_element_func`
121 : !> \param this object reference
122 : !> \param i i-index
123 : !> \param j j-index
124 : !> \param k k-index
125 : !> \param l l-index
126 : !> \param val value of the integral at (i,j,k.l)
127 : !> \return True if the ERI foreach loop should continue, false, if not
128 : ! **************************************************************************************************
129 : LOGICAL FUNCTION eri_type_eri_element_func_interface(this, i, j, k, l, val)
130 : IMPORT :: eri_type_eri_element_func, dp
131 : CLASS(eri_type_eri_element_func), INTENT(inout) :: this
132 : INTEGER, INTENT(in) :: i, j, k, l
133 : REAL(KIND=dp), INTENT(in) :: val
134 : END FUNCTION eri_type_eri_element_func_interface
135 : END INTERFACE
136 :
137 : ! **************************************************************************************************
138 :
139 : CONTAINS
140 :
141 : ! **************************************************************************************************
142 : !> \brief Creates an active space environment type, nullifying all quantities.
143 : !> \param active_space_env the active space environment to be initialized
144 : ! **************************************************************************************************
145 72 : SUBROUTINE create_active_space_type(active_space_env)
146 : TYPE(active_space_type), POINTER :: active_space_env
147 :
148 72 : IF (ASSOCIATED(active_space_env)) THEN
149 0 : CALL release_active_space_type(active_space_env)
150 : END IF
151 :
152 864 : ALLOCATE (active_space_env)
153 : NULLIFY (active_space_env%active_orbitals, active_space_env%inactive_orbitals)
154 : NULLIFY (active_space_env%mos_active, active_space_env%mos_inactive)
155 : NULLIFY (active_space_env%ks_sub, active_space_env%p_active)
156 : NULLIFY (active_space_env%vxc_sub, active_space_env%h_sub)
157 : NULLIFY (active_space_env%fock_sub, active_space_env%pmat_inactive)
158 :
159 72 : END SUBROUTINE create_active_space_type
160 :
161 : ! **************************************************************************************************
162 : !> \brief Releases all quantities in the active space environment.
163 : !> \param active_space_env the active space environment to be released
164 : ! **************************************************************************************************
165 72 : SUBROUTINE release_active_space_type(active_space_env)
166 : TYPE(active_space_type), POINTER :: active_space_env
167 :
168 : INTEGER :: imo
169 :
170 72 : IF (ASSOCIATED(active_space_env)) THEN
171 :
172 72 : IF (ASSOCIATED(active_space_env%active_orbitals)) THEN
173 72 : DEALLOCATE (active_space_env%active_orbitals)
174 : END IF
175 :
176 72 : IF (ASSOCIATED(active_space_env%inactive_orbitals)) THEN
177 72 : DEALLOCATE (active_space_env%inactive_orbitals)
178 : END IF
179 :
180 72 : IF (ASSOCIATED(active_space_env%mos_active)) THEN
181 158 : DO imo = 1, SIZE(active_space_env%mos_active)
182 158 : CALL deallocate_mo_set(active_space_env%mos_active(imo))
183 : END DO
184 72 : DEALLOCATE (active_space_env%mos_active)
185 : END IF
186 :
187 72 : IF (ASSOCIATED(active_space_env%mos_inactive)) THEN
188 158 : DO imo = 1, SIZE(active_space_env%mos_inactive)
189 158 : CALL deallocate_mo_set(active_space_env%mos_inactive(imo))
190 : END DO
191 72 : DEALLOCATE (active_space_env%mos_inactive)
192 : END IF
193 :
194 72 : CALL release_eri_type(active_space_env%eri)
195 :
196 72 : CALL cp_fm_release(active_space_env%p_active)
197 72 : CALL cp_fm_release(active_space_env%ks_sub)
198 72 : CALL cp_fm_release(active_space_env%vxc_sub)
199 72 : CALL cp_fm_release(active_space_env%h_sub)
200 72 : CALL cp_fm_release(active_space_env%fock_sub)
201 72 : CALL cp_fm_release(active_space_env%sab_sub)
202 :
203 72 : IF (ASSOCIATED(active_space_env%pmat_inactive)) &
204 72 : CALL dbcsr_deallocate_matrix_set(active_space_env%pmat_inactive)
205 :
206 72 : DEALLOCATE (active_space_env)
207 : END IF
208 :
209 72 : END SUBROUTINE release_active_space_type
210 :
211 : ! **************************************************************************************************
212 : !> \brief Releases the ERI environment type.
213 : !> \param eri_env the ERI environment to be released
214 : ! **************************************************************************************************
215 72 : SUBROUTINE release_eri_type(eri_env)
216 : TYPE(eri_type) :: eri_env
217 :
218 : INTEGER :: i
219 :
220 72 : IF (ASSOCIATED(eri_env%eri)) THEN
221 :
222 168 : DO i = 1, SIZE(eri_env%eri)
223 96 : CALL dbcsr_csr_destroy(eri_env%eri(i)%csr_mat)
224 168 : DEALLOCATE (eri_env%eri(i)%csr_mat)
225 : END DO
226 72 : CALL mp_para_env_release(eri_env%para_env_sub)
227 72 : CALL eri_env%comm_exchange%free()
228 72 : DEALLOCATE (eri_env%eri)
229 :
230 : END IF
231 :
232 72 : END SUBROUTINE release_eri_type
233 :
234 : ! **************************************************************************************************
235 : !> \brief calculates combined index (ij)
236 : !> \param i Index j
237 : !> \param j Index i
238 : !> \param n Dimension in i or j direction
239 : !> \returns The combined index
240 : !> \par History
241 : !> 04.2016 created [JGH]
242 : ! **************************************************************************************************
243 4574 : INTEGER FUNCTION csr_idx_to_combined(i, j, n) RESULT(ij)
244 : INTEGER, INTENT(IN) :: i, j, n
245 :
246 4574 : CPASSERT(i <= j)
247 4574 : CPASSERT(i <= n)
248 4574 : CPASSERT(j <= n)
249 :
250 4574 : ij = (i - 1)*n - ((i - 1)*(i - 2))/2 + (j - i + 1)
251 :
252 4574 : CPASSERT(ij <= (n*(n + 1))/2 .AND. 0 <= ij)
253 :
254 4574 : END FUNCTION csr_idx_to_combined
255 :
256 : ! **************************************************************************************************
257 : !> \brief extracts indices i and j from combined index ij
258 : !> \param ij The combined index
259 : !> \param n Dimension in i or j direction
260 : !> \param i Resulting i index
261 : !> \param j Resulting j index
262 : !> \par History
263 : !> 04.2016 created [JGH]
264 : ! **************************************************************************************************
265 5642 : SUBROUTINE csr_idx_from_combined(ij, n, i, j)
266 : INTEGER, INTENT(IN) :: ij, n
267 : INTEGER, INTENT(OUT) :: i, j
268 :
269 : INTEGER :: m, m0
270 :
271 5642 : m = MAX(ij/n, 1)
272 14872 : DO i = m, n
273 14872 : m0 = (i - 1)*n - ((i - 1)*(i - 2))/2
274 14872 : j = ij - m0 + i - 1
275 14872 : IF (j <= n) EXIT
276 : END DO
277 :
278 5642 : CPASSERT(i > 0 .AND. i <= n)
279 5642 : CPASSERT(j > 0 .AND. j <= n)
280 5642 : CPASSERT(i <= j)
281 :
282 5642 : END SUBROUTINE csr_idx_from_combined
283 :
284 : ! **************************************************************************************************
285 : !> \brief Calls the provided function for each element in the ERI
286 : !> \param this object reference
287 : !> \param nspin The spin number
288 : !> \param active_orbitals the active orbital indices
289 : !> \param fobj The function object from which to call `func(i, j, k, l, val)`
290 : !> \param spin1 the first spin value
291 : !> \param spin2 the second spin value
292 : !> \par History
293 : !> 04.2016 created [JHU]
294 : !> 06.2016 factored out from qs_a_s_methods:fcidump [TMU]
295 : !> \note Calls MPI, must be executed on all ranks.
296 : ! **************************************************************************************************
297 192 : SUBROUTINE eri_type_eri_foreach(this, nspin, active_orbitals, fobj, spin1, spin2)
298 : CLASS(eri_type), INTENT(in) :: this
299 : CLASS(eri_type_eri_element_func) :: fobj
300 : INTEGER, DIMENSION(:, :), INTENT(IN) :: active_orbitals
301 : INTEGER, OPTIONAL :: spin1, spin2
302 :
303 : CHARACTER(LEN=*), PARAMETER :: routineN = "eri_type_eri_foreach"
304 :
305 : INTEGER :: i1, i12, i12l, i2, i3, i34, i34l, i4, m1, m2, m3, m4, &
306 : irptr, nspin, nindex, nmo, proc, nonzero_elements_local, handle, dummy_int(1)
307 192 : INTEGER, ALLOCATABLE, DIMENSION(:) :: colind, offsets, nonzero_elements_global
308 192 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: erival
309 : REAL(KIND=dp) :: erint, dummy_real(1)
310 : TYPE(mp_comm_type) :: mp_group
311 :
312 192 : CALL timeset(routineN, handle)
313 :
314 192 : IF (.NOT. PRESENT(spin1)) THEN
315 0 : spin1 = nspin
316 : END IF
317 192 : IF (.NOT. PRESENT(spin2)) THEN
318 0 : spin2 = nspin
319 : END IF
320 :
321 192 : dummy_int = 0
322 192 : dummy_real = 0.0_dp
323 :
324 : ASSOCIATE (eri => this%eri(nspin)%csr_mat, norb => this%norb)
325 192 : nindex = (norb*(norb + 1))/2
326 192 : CALL mp_group%set_handle(eri%mp_group%get_handle())
327 192 : nmo = SIZE(active_orbitals, 1)
328 : ! Irrelevant in case of half-transformed integrals
329 960 : ALLOCATE (erival(nindex), colind(nindex))
330 : ALLOCATE (offsets(0:mp_group%num_pe - 1), &
331 768 : nonzero_elements_global(0:mp_group%num_pe - 1))
332 :
333 896 : DO m1 = 1, nmo
334 512 : i1 = active_orbitals(m1, spin1)
335 1732 : DO m2 = m1, nmo
336 1028 : i2 = active_orbitals(m2, spin1)
337 1028 : i12 = csr_idx_to_combined(i1, i2, norb)
338 1028 : i12l = (i12 - 1)/this%comm_exchange%num_pe + 1
339 :
340 : ! In case of half-transformed integrals, every process might carry integrals of a row
341 : ! The number of integrals varies between processes and rows (related to the randomized
342 : ! distribution of matrix blocks)
343 :
344 : ! 1) Collect the amount of local data from each process
345 1028 : nonzero_elements_local = 0
346 1028 : IF (MOD(i12 - 1, this%comm_exchange%num_pe) == this%comm_exchange%mepos) &
347 1010 : nonzero_elements_local = eri%nzerow_local(i12l)
348 1028 : CALL mp_group%allgather(nonzero_elements_local, nonzero_elements_global)
349 :
350 : ! 2) Prepare arrays for communication (calculate the offsets and the total number of elements)
351 1028 : offsets(0) = 0
352 2056 : DO proc = 1, mp_group%num_pe - 1
353 2056 : offsets(proc) = offsets(proc - 1) + nonzero_elements_global(proc - 1)
354 : END DO
355 1028 : nindex = offsets(mp_group%num_pe - 1) + nonzero_elements_global(mp_group%num_pe - 1)
356 1028 : irptr = 1
357 1028 : IF (MOD(i12 - 1, this%comm_exchange%num_pe) == this%comm_exchange%mepos) THEN
358 1010 : irptr = eri%rowptr_local(i12l)
359 :
360 : ! Exchange actual data
361 : CALL mp_group%allgatherv(eri%colind_local(irptr:irptr + nonzero_elements_local - 1), &
362 3112 : colind(1:nindex), nonzero_elements_global, offsets)
363 : CALL mp_group%allgatherv(eri%nzval_local%r_dp(irptr:irptr + nonzero_elements_local - 1), &
364 3112 : erival(1:nindex), nonzero_elements_global, offsets)
365 : ELSE
366 18 : CALL mp_group%allgatherv(dummy_int(1:0), colind(1:nindex), nonzero_elements_global, offsets)
367 18 : CALL mp_group%allgatherv(dummy_real(1:0), erival(1:nindex), nonzero_elements_global, offsets)
368 : END IF
369 :
370 5744 : DO i34l = 1, nindex
371 4204 : i34 = colind(i34l)
372 4204 : erint = erival(i34l)
373 4204 : CALL csr_idx_from_combined(i34, norb, i3, i4)
374 :
375 10972 : DO m3 = 1, nmo
376 10972 : IF (active_orbitals(m3, spin2) == i3) THEN
377 : EXIT
378 : END IF
379 : END DO
380 :
381 15024 : DO m4 = 1, nmo
382 15024 : IF (active_orbitals(m4, spin2) == i4) THEN
383 : EXIT
384 : END IF
385 : END DO
386 :
387 : ! terminate the loop prematurely if the function returns false
388 9436 : IF (.NOT. fobj%func(m1, m2, m3, m4, erint)) RETURN
389 : END DO
390 :
391 : END DO
392 : END DO
393 : END ASSOCIATE
394 :
395 192 : CALL timestop(handle)
396 384 : END SUBROUTINE eri_type_eri_foreach
397 :
398 0 : END MODULE qs_active_space_types
|
