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) :: alpha = 0.0_dp
105 : INTEGER :: as_mixing_method = direct_mixing_nr
106 : INTEGER :: as_mixing_iter = 0
107 : INTEGER :: as_mixing_dim = 0
108 : LOGICAL :: do_scf_embedding = .FALSE.
109 : LOGICAL :: qcschema = .FALSE.
110 : LOGICAL :: fcidump = .FALSE.
111 : CHARACTER(LEN=default_path_length) :: qcschema_filename = ''
112 : TYPE(eri_type) :: eri = eri_type()
113 : REAL(KIND=dp), DIMENSION(:), POINTER :: as_mix_r_old => NULL()
114 : REAL(KIND=dp), DIMENSION(:), POINTER :: as_mix_weight => NULL()
115 : REAL(KIND=dp), DIMENSION(:), POINTER :: as_mix_x_old => NULL()
116 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: as_mix_r_buffer => NULL()
117 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: as_mix_x_buffer => NULL()
118 : TYPE(mixing_storage_type), POINTER :: as_mixing_store => NULL()
119 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_active => NULL()
120 : TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_inactive => NULL()
121 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: p_active => NULL()
122 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: ks_sub => NULL()
123 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: vxc_sub => NULL()
124 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: h_sub => NULL()
125 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: fock_sub => NULL()
126 : TYPE(cp_fm_type), DIMENSION(:), POINTER :: sab_sub => NULL()
127 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: pmat_inactive => NULL()
128 : TYPE(section_vals_type), POINTER :: xc_section => NULL()
129 : END TYPE active_space_type
130 :
131 : ABSTRACT INTERFACE
132 : ! **************************************************************************************************
133 : !> \brief The function signature to be implemented by a child of `eri_type_eri_element_func`
134 : !> \param this object reference
135 : !> \param i i-index
136 : !> \param j j-index
137 : !> \param k k-index
138 : !> \param l l-index
139 : !> \param val value of the integral at (i,j,k.l)
140 : !> \return True if the ERI foreach loop should continue, false, if not
141 : ! **************************************************************************************************
142 : LOGICAL FUNCTION eri_type_eri_element_func_interface(this, i, j, k, l, val)
143 : IMPORT :: eri_type_eri_element_func, dp
144 : CLASS(eri_type_eri_element_func), INTENT(inout) :: this
145 : INTEGER, INTENT(in) :: i, j, k, l
146 : REAL(KIND=dp), INTENT(in) :: val
147 : END FUNCTION eri_type_eri_element_func_interface
148 : END INTERFACE
149 :
150 : ! **************************************************************************************************
151 :
152 : CONTAINS
153 :
154 : ! **************************************************************************************************
155 : !> \brief Creates an active space environment type, nullifying all quantities.
156 : !> \param active_space_env the active space environment to be initialized
157 : ! **************************************************************************************************
158 82 : SUBROUTINE create_active_space_type(active_space_env)
159 : TYPE(active_space_type), POINTER :: active_space_env
160 :
161 82 : IF (ASSOCIATED(active_space_env)) THEN
162 0 : CALL release_active_space_type(active_space_env)
163 : END IF
164 :
165 984 : ALLOCATE (active_space_env)
166 : NULLIFY (active_space_env%active_orbitals, active_space_env%inactive_orbitals)
167 : NULLIFY (active_space_env%mos_active, active_space_env%mos_inactive)
168 : NULLIFY (active_space_env%ks_sub, active_space_env%p_active)
169 : NULLIFY (active_space_env%vxc_sub, active_space_env%h_sub)
170 : NULLIFY (active_space_env%fock_sub, active_space_env%pmat_inactive)
171 : NULLIFY (active_space_env%as_mixing_store)
172 : NULLIFY (active_space_env%as_mix_r_old, active_space_env%as_mix_weight)
173 : NULLIFY (active_space_env%as_mix_x_old)
174 : NULLIFY (active_space_env%as_mix_r_buffer, active_space_env%as_mix_x_buffer)
175 :
176 82 : END SUBROUTINE create_active_space_type
177 :
178 : ! **************************************************************************************************
179 : !> \brief Releases all quantities in the active space environment.
180 : !> \param active_space_env the active space environment to be released
181 : ! **************************************************************************************************
182 82 : SUBROUTINE release_active_space_type(active_space_env)
183 : TYPE(active_space_type), POINTER :: active_space_env
184 :
185 : INTEGER :: imo
186 :
187 82 : IF (ASSOCIATED(active_space_env)) THEN
188 :
189 82 : IF (ASSOCIATED(active_space_env%active_orbitals)) THEN
190 82 : DEALLOCATE (active_space_env%active_orbitals)
191 : END IF
192 :
193 82 : IF (ASSOCIATED(active_space_env%inactive_orbitals)) THEN
194 82 : DEALLOCATE (active_space_env%inactive_orbitals)
195 : END IF
196 :
197 82 : IF (ASSOCIATED(active_space_env%mos_active)) THEN
198 182 : DO imo = 1, SIZE(active_space_env%mos_active)
199 182 : CALL deallocate_mo_set(active_space_env%mos_active(imo))
200 : END DO
201 82 : DEALLOCATE (active_space_env%mos_active)
202 : END IF
203 :
204 82 : IF (ASSOCIATED(active_space_env%mos_inactive)) THEN
205 182 : DO imo = 1, SIZE(active_space_env%mos_inactive)
206 182 : CALL deallocate_mo_set(active_space_env%mos_inactive(imo))
207 : END DO
208 82 : DEALLOCATE (active_space_env%mos_inactive)
209 : END IF
210 :
211 82 : CALL release_eri_type(active_space_env%eri)
212 :
213 82 : CALL cp_fm_release(active_space_env%p_active)
214 82 : CALL cp_fm_release(active_space_env%ks_sub)
215 82 : CALL cp_fm_release(active_space_env%vxc_sub)
216 82 : CALL cp_fm_release(active_space_env%h_sub)
217 82 : CALL cp_fm_release(active_space_env%fock_sub)
218 82 : CALL cp_fm_release(active_space_env%sab_sub)
219 :
220 82 : IF (ASSOCIATED(active_space_env%as_mixing_store)) THEN
221 82 : CALL mixing_storage_release(active_space_env%as_mixing_store)
222 82 : DEALLOCATE (active_space_env%as_mixing_store)
223 : END IF
224 82 : IF (ASSOCIATED(active_space_env%as_mix_r_old)) THEN
225 0 : DEALLOCATE (active_space_env%as_mix_r_old)
226 : END IF
227 82 : IF (ASSOCIATED(active_space_env%as_mix_weight)) THEN
228 0 : DEALLOCATE (active_space_env%as_mix_weight)
229 : END IF
230 82 : IF (ASSOCIATED(active_space_env%as_mix_x_old)) THEN
231 0 : DEALLOCATE (active_space_env%as_mix_x_old)
232 : END IF
233 82 : IF (ASSOCIATED(active_space_env%as_mix_r_buffer)) THEN
234 0 : DEALLOCATE (active_space_env%as_mix_r_buffer)
235 : END IF
236 82 : IF (ASSOCIATED(active_space_env%as_mix_x_buffer)) THEN
237 0 : DEALLOCATE (active_space_env%as_mix_x_buffer)
238 : END IF
239 :
240 82 : IF (ASSOCIATED(active_space_env%pmat_inactive)) &
241 82 : CALL dbcsr_deallocate_matrix_set(active_space_env%pmat_inactive)
242 :
243 82 : DEALLOCATE (active_space_env)
244 : END IF
245 :
246 82 : END SUBROUTINE release_active_space_type
247 :
248 : ! **************************************************************************************************
249 : !> \brief Releases the ERI environment type.
250 : !> \param eri_env the ERI environment to be released
251 : ! **************************************************************************************************
252 82 : SUBROUTINE release_eri_type(eri_env)
253 : TYPE(eri_type) :: eri_env
254 :
255 : INTEGER :: i
256 :
257 82 : IF (ASSOCIATED(eri_env%eri)) THEN
258 :
259 192 : DO i = 1, SIZE(eri_env%eri)
260 110 : CALL dbcsr_csr_destroy(eri_env%eri(i)%csr_mat)
261 192 : DEALLOCATE (eri_env%eri(i)%csr_mat)
262 : END DO
263 82 : CALL mp_para_env_release(eri_env%para_env_sub)
264 82 : CALL eri_env%comm_exchange%free()
265 82 : DEALLOCATE (eri_env%eri)
266 :
267 : END IF
268 :
269 82 : END SUBROUTINE release_eri_type
270 :
271 : ! **************************************************************************************************
272 : !> \brief calculates combined index (ij)
273 : !> \param i Index j
274 : !> \param j Index i
275 : !> \param n Dimension in i or j direction
276 : !> \returns The combined index
277 : !> \par History
278 : !> 04.2016 created [JGH]
279 : ! **************************************************************************************************
280 5504 : INTEGER FUNCTION csr_idx_to_combined(i, j, n) RESULT(ij)
281 : INTEGER, INTENT(IN) :: i, j, n
282 :
283 5504 : CPASSERT(i <= j)
284 5504 : CPASSERT(i <= n)
285 5504 : CPASSERT(j <= n)
286 :
287 5504 : ij = (i - 1)*n - ((i - 1)*(i - 2))/2 + (j - i + 1)
288 :
289 5504 : CPASSERT(ij <= (n*(n + 1))/2 .AND. 0 <= ij)
290 :
291 5504 : END FUNCTION csr_idx_to_combined
292 :
293 : ! **************************************************************************************************
294 : !> \brief extracts indices i and j from combined index ij
295 : !> \param ij The combined index
296 : !> \param n Dimension in i or j direction
297 : !> \param i Resulting i index
298 : !> \param j Resulting j index
299 : !> \par History
300 : !> 04.2016 created [JGH]
301 : ! **************************************************************************************************
302 6949 : SUBROUTINE csr_idx_from_combined(ij, n, i, j)
303 : INTEGER, INTENT(IN) :: ij, n
304 : INTEGER, INTENT(OUT) :: i, j
305 :
306 : INTEGER :: m, m0
307 :
308 6949 : m = MAX(ij/n, 1)
309 18793 : DO i = m, n
310 18793 : m0 = (i - 1)*n - ((i - 1)*(i - 2))/2
311 18793 : j = ij - m0 + i - 1
312 18793 : IF (j <= n) EXIT
313 : END DO
314 :
315 6949 : CPASSERT(i > 0 .AND. i <= n)
316 6949 : CPASSERT(j > 0 .AND. j <= n)
317 6949 : CPASSERT(i <= j)
318 :
319 6949 : END SUBROUTINE csr_idx_from_combined
320 :
321 : ! **************************************************************************************************
322 : !> \brief Calls the provided function for each element in the ERI
323 : !> \param this object reference
324 : !> \param nspin The spin number
325 : !> \param active_orbitals the active orbital indices
326 : !> \param fobj The function object from which to call `func(i, j, k, l, val)`
327 : !> \param spin1 the first spin value
328 : !> \param spin2 the second spin value
329 : !> \par History
330 : !> 04.2016 created [JHU]
331 : !> 06.2016 factored out from qs_a_s_methods:fcidump [TMU]
332 : !> \note Calls MPI, must be executed on all ranks.
333 : ! **************************************************************************************************
334 200 : SUBROUTINE eri_type_eri_foreach(this, nspin, active_orbitals, fobj, spin1, spin2)
335 : CLASS(eri_type), INTENT(in) :: this
336 : CLASS(eri_type_eri_element_func) :: fobj
337 : INTEGER, DIMENSION(:, :), INTENT(IN) :: active_orbitals
338 : INTEGER, OPTIONAL :: spin1, spin2
339 :
340 : CHARACTER(LEN=*), PARAMETER :: routineN = "eri_type_eri_foreach"
341 :
342 : INTEGER :: i1, i12, i12l, i2, i3, i34, i34l, i4, m1, m2, m3, m4, &
343 : irptr, nspin, nindex, nmo, proc, nonzero_elements_local, handle, dummy_int(1)
344 200 : INTEGER, ALLOCATABLE, DIMENSION(:) :: colind, offsets, nonzero_elements_global
345 200 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: erival
346 : REAL(KIND=dp) :: erint, dummy_real(1)
347 : TYPE(mp_comm_type) :: mp_group
348 :
349 200 : CALL timeset(routineN, handle)
350 :
351 200 : IF (.NOT. PRESENT(spin1)) THEN
352 0 : spin1 = nspin
353 : END IF
354 200 : IF (.NOT. PRESENT(spin2)) THEN
355 0 : spin2 = nspin
356 : END IF
357 :
358 200 : dummy_int = 0
359 200 : dummy_real = 0.0_dp
360 :
361 : ASSOCIATE (eri => this%eri(nspin)%csr_mat, norb => this%norb)
362 200 : nindex = (norb*(norb + 1))/2
363 200 : CALL mp_group%set_handle(eri%mp_group%get_handle())
364 200 : nmo = SIZE(active_orbitals, 1)
365 : ! Irrelevant in case of half-transformed integrals
366 1000 : ALLOCATE (erival(nindex), colind(nindex))
367 : ALLOCATE (offsets(0:mp_group%num_pe - 1), &
368 800 : nonzero_elements_global(0:mp_group%num_pe - 1))
369 :
370 944 : DO m1 = 1, nmo
371 544 : i1 = active_orbitals(m1, spin1)
372 1868 : DO m2 = m1, nmo
373 1124 : i2 = active_orbitals(m2, spin1)
374 1124 : i12 = csr_idx_to_combined(i1, i2, norb)
375 1124 : i12l = (i12 - 1)/this%comm_exchange%num_pe + 1
376 :
377 : ! In case of half-transformed integrals, every process might carry integrals of a row
378 : ! The number of integrals varies between processes and rows (related to the randomized
379 : ! distribution of matrix blocks)
380 :
381 : ! 1) Collect the amount of local data from each process
382 1124 : nonzero_elements_local = 0
383 1124 : IF (MOD(i12 - 1, this%comm_exchange%num_pe) == this%comm_exchange%mepos) &
384 1106 : nonzero_elements_local = eri%nzerow_local(i12l)
385 1124 : CALL mp_group%allgather(nonzero_elements_local, nonzero_elements_global)
386 :
387 : ! 2) Prepare arrays for communication (calculate the offsets and the total number of elements)
388 1124 : offsets(0) = 0
389 2248 : DO proc = 1, mp_group%num_pe - 1
390 2248 : offsets(proc) = offsets(proc - 1) + nonzero_elements_global(proc - 1)
391 : END DO
392 1124 : nindex = offsets(mp_group%num_pe - 1) + nonzero_elements_global(mp_group%num_pe - 1)
393 1124 : irptr = 1
394 1124 : IF (MOD(i12 - 1, this%comm_exchange%num_pe) == this%comm_exchange%mepos) THEN
395 1106 : irptr = eri%rowptr_local(i12l)
396 :
397 : ! Exchange actual data
398 : CALL mp_group%allgatherv(eri%colind_local(irptr:irptr + nonzero_elements_local - 1), &
399 3676 : colind(1:nindex), nonzero_elements_global, offsets)
400 : CALL mp_group%allgatherv(eri%nzval_local%r_dp(irptr:irptr + nonzero_elements_local - 1), &
401 3676 : erival(1:nindex), nonzero_elements_global, offsets)
402 : ELSE
403 18 : CALL mp_group%allgatherv(dummy_int(1:0), colind(1:nindex), nonzero_elements_global, offsets)
404 18 : CALL mp_group%allgatherv(dummy_real(1:0), erival(1:nindex), nonzero_elements_global, offsets)
405 : END IF
406 :
407 6808 : DO i34l = 1, nindex
408 5140 : i34 = colind(i34l)
409 5140 : erint = erival(i34l)
410 5140 : CALL csr_idx_from_combined(i34, norb, i3, i4)
411 :
412 14180 : DO m3 = 1, nmo
413 14180 : IF (active_orbitals(m3, spin2) == i3) THEN
414 : EXIT
415 : END IF
416 : END DO
417 :
418 19540 : DO m4 = 1, nmo
419 19540 : IF (active_orbitals(m4, spin2) == i4) THEN
420 : EXIT
421 : END IF
422 : END DO
423 :
424 : ! terminate the loop prematurely if the function returns false
425 11404 : IF (.NOT. fobj%func(m1, m2, m3, m4, erint)) RETURN
426 : END DO
427 :
428 : END DO
429 : END DO
430 : END ASSOCIATE
431 :
432 200 : CALL timestop(handle)
433 400 : END SUBROUTINE eri_type_eri_foreach
434 :
435 0 : END MODULE qs_active_space_types
|
