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