Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Routines to compute the Coulomb integral V_(alpha beta)(k) for a k-point k using lattice
10 : !> summation in real space. These integrals are e.g. needed in periodic RI for RPA, GW
11 : !> \par History
12 : !> 2018.05 created [Jan Wilhelm]
13 : !> \author Jan Wilhelm
14 : ! **************************************************************************************************
15 : MODULE kpoint_coulomb_2c
16 : USE atomic_kind_types, ONLY: atomic_kind_type,&
17 : get_atomic_kind_set
18 : USE basis_set_types, ONLY: gto_basis_set_type
19 : USE cell_types, ONLY: cell_type,&
20 : get_cell,&
21 : pbc
22 : USE dbcsr_api, ONLY: &
23 : dbcsr_create, dbcsr_init_p, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
24 : dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_p_type, &
25 : dbcsr_release_p, dbcsr_reserve_all_blocks, dbcsr_set, dbcsr_type, dbcsr_type_no_symmetry
26 : USE generic_shg_integrals, ONLY: int_operators_r12_ab_shg
27 : USE generic_shg_integrals_init, ONLY: contraction_matrix_shg
28 : USE kinds, ONLY: dp
29 : USE kpoint_types, ONLY: get_kpoint_info,&
30 : kpoint_type
31 : USE mathconstants, ONLY: twopi
32 : USE particle_types, ONLY: particle_type
33 : USE qs_kind_types, ONLY: get_qs_kind,&
34 : qs_kind_type
35 : #include "./base/base_uses.f90"
36 :
37 : IMPLICIT NONE
38 :
39 : PRIVATE
40 :
41 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'kpoint_coulomb_2c'
42 :
43 : PUBLIC :: build_2c_coulomb_matrix_kp
44 :
45 : ! **************************************************************************************************
46 :
47 : TYPE two_d_util_type
48 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: block
49 : END TYPE
50 :
51 : CONTAINS
52 :
53 : ! **************************************************************************************************
54 : !> \brief ...
55 : !> \param matrix_v_kp ...
56 : !> \param kpoints ...
57 : !> \param basis_type ...
58 : !> \param cell ...
59 : !> \param particle_set ...
60 : !> \param qs_kind_set ...
61 : !> \param atomic_kind_set ...
62 : !> \param size_lattice_sum ...
63 : !> \param operator_type ...
64 : !> \param ikp_start ...
65 : !> \param ikp_end ...
66 : ! **************************************************************************************************
67 138 : SUBROUTINE build_2c_coulomb_matrix_kp(matrix_v_kp, kpoints, basis_type, cell, particle_set, qs_kind_set, &
68 : atomic_kind_set, size_lattice_sum, operator_type, ikp_start, ikp_end)
69 :
70 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_kp
71 : TYPE(kpoint_type), POINTER :: kpoints
72 : CHARACTER(LEN=*), INTENT(IN) :: basis_type
73 : TYPE(cell_type), POINTER :: cell
74 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
75 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
76 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
77 : INTEGER :: size_lattice_sum, operator_type, &
78 : ikp_start, ikp_end
79 :
80 : CHARACTER(LEN=*), PARAMETER :: routineN = 'build_2c_coulomb_matrix_kp'
81 :
82 : INTEGER :: handle, total_periodicity
83 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
84 :
85 138 : CALL timeset(routineN, handle)
86 :
87 138 : CALL check_periodicity(cell, kpoints, total_periodicity)
88 :
89 138 : CALL allocate_tmp(matrix_v_L_tmp, matrix_v_kp, ikp_start)
90 :
91 : CALL lattice_sum(matrix_v_kp, kpoints, basis_type, cell, particle_set, &
92 : qs_kind_set, atomic_kind_set, size_lattice_sum, matrix_v_L_tmp, &
93 138 : operator_type, ikp_start, ikp_end)
94 :
95 138 : CALL deallocate_tmp(matrix_v_L_tmp)
96 :
97 138 : CALL timestop(handle)
98 :
99 138 : END SUBROUTINE build_2c_coulomb_matrix_kp
100 :
101 : ! **************************************************************************************************
102 : !> \brief ...
103 : !> \param matrix_v_kp ...
104 : !> \param kpoints ...
105 : !> \param basis_type ...
106 : !> \param cell ...
107 : !> \param particle_set ...
108 : !> \param qs_kind_set ...
109 : !> \param atomic_kind_set ...
110 : !> \param size_lattice_sum ...
111 : !> \param matrix_v_L_tmp ...
112 : !> \param operator_type ...
113 : !> \param ikp_start ...
114 : !> \param ikp_end ...
115 : ! **************************************************************************************************
116 138 : SUBROUTINE lattice_sum(matrix_v_kp, kpoints, basis_type, cell, particle_set, &
117 : qs_kind_set, atomic_kind_set, size_lattice_sum, matrix_v_L_tmp, &
118 : operator_type, ikp_start, ikp_end)
119 :
120 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_kp
121 : TYPE(kpoint_type), POINTER :: kpoints
122 : CHARACTER(LEN=*), INTENT(IN) :: basis_type
123 : TYPE(cell_type), POINTER :: cell
124 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
125 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
126 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
127 : INTEGER :: size_lattice_sum
128 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
129 : INTEGER :: operator_type, ikp_start, ikp_end
130 :
131 : CHARACTER(LEN=*), PARAMETER :: routineN = 'lattice_sum'
132 :
133 : INTEGER :: factor, handle, i_block, i_x, i_x_inner, i_x_outer, ik, j_y, j_y_inner, &
134 : j_y_outer, k_z, k_z_inner, k_z_outer, nkp, x_max, x_min, y_max, y_min, z_max, z_min
135 : INTEGER, DIMENSION(3) :: nkp_grid, periodic
136 : REAL(KIND=dp) :: coskl, sinkl
137 : REAL(KIND=dp), DIMENSION(3) :: vec_L, vec_s
138 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat
139 138 : TYPE(two_d_util_type), ALLOCATABLE, DIMENSION(:) :: blocks_v_L, blocks_v_L_store
140 : TYPE(two_d_util_type), ALLOCATABLE, &
141 138 : DIMENSION(:, :, :) :: blocks_v_kp
142 :
143 138 : CALL timeset(routineN, handle)
144 :
145 138 : CALL get_kpoint_info(kpoints, nkp_grid=nkp_grid, nkp=nkp)
146 138 : CALL get_cell(cell=cell, h=hmat, periodic=periodic)
147 :
148 138 : IF (MODULO(nkp_grid(1), 2) == 1) THEN
149 74 : factor = 3**(size_lattice_sum - 1)
150 64 : ELSE IF (MODULO(nkp_grid(1), 2) == 0) THEN
151 64 : factor = 2**(size_lattice_sum - 1)
152 : END IF
153 :
154 138 : IF (MODULO(nkp_grid(1), 2) == 1) THEN
155 74 : x_min = -(factor*nkp_grid(1) - 1)/2
156 74 : x_max = (factor*nkp_grid(1) - 1)/2
157 64 : ELSE IF (MODULO(nkp_grid(1), 2) == 0) THEN
158 64 : x_min = -factor*nkp_grid(1)/2
159 64 : x_max = factor*nkp_grid(1)/2 - 1
160 : END IF
161 138 : IF (periodic(1) == 0) THEN
162 74 : x_min = 0
163 74 : x_max = 0
164 : END IF
165 :
166 138 : IF (MODULO(nkp_grid(2), 2) == 1) THEN
167 12 : y_min = -(factor*nkp_grid(2) - 1)/2
168 12 : y_max = (factor*nkp_grid(2) - 1)/2
169 126 : ELSE IF (MODULO(nkp_grid(2), 2) == 0) THEN
170 126 : y_min = -factor*nkp_grid(2)/2
171 126 : y_max = factor*nkp_grid(2)/2 - 1
172 : END IF
173 138 : IF (periodic(2) == 0) THEN
174 12 : y_min = 0
175 12 : y_max = 0
176 : END IF
177 :
178 138 : IF (MODULO(nkp_grid(3), 2) == 1) THEN
179 56 : z_min = -(factor*nkp_grid(3) - 1)/2
180 56 : z_max = (factor*nkp_grid(3) - 1)/2
181 82 : ELSE IF (MODULO(nkp_grid(3), 2) == 0) THEN
182 82 : z_min = -factor*nkp_grid(3)/2
183 82 : z_max = factor*nkp_grid(3)/2 - 1
184 : END IF
185 138 : IF (periodic(3) == 0) THEN
186 56 : z_min = 0
187 56 : z_max = 0
188 : END IF
189 :
190 138 : CALL allocate_blocks_v_kp(blocks_v_kp, matrix_v_kp, ikp_start, ikp_end)
191 138 : CALL allocate_blocks_v_L(blocks_v_L, matrix_v_L_tmp)
192 138 : CALL allocate_blocks_v_L(blocks_v_L_store, matrix_v_L_tmp)
193 :
194 974 : DO i_x_inner = 0, 2*nkp_grid(1) - 1
195 9406 : DO j_y_inner = 0, 2*nkp_grid(2) - 1
196 40628 : DO k_z_inner = 0, 2*nkp_grid(3) - 1
197 :
198 93184 : DO i_x_outer = x_min, x_max + nkp_grid(1), 2*nkp_grid(1)
199 295424 : DO j_y_outer = y_min, y_max + nkp_grid(2), 2*nkp_grid(2)
200 766848 : DO k_z_outer = z_min, z_max + nkp_grid(3), 2*nkp_grid(3)
201 :
202 502784 : i_x = i_x_inner + i_x_outer
203 502784 : j_y = j_y_inner + j_y_outer
204 502784 : k_z = k_z_inner + k_z_outer
205 :
206 : IF (i_x > x_max .OR. i_x < x_min .OR. &
207 : j_y > y_max .OR. j_y < y_min .OR. &
208 502784 : k_z > z_max .OR. k_z < z_min) CYCLE
209 :
210 710528 : vec_s = [REAL(i_x, dp), REAL(j_y, dp), REAL(k_z, dp)]
211 :
212 2309216 : vec_L = MATMUL(hmat, vec_s)
213 :
214 : ! Compute (P 0 | Q vec_L) and store it in matrix_v_L_tmp
215 : CALL compute_v_transl(matrix_v_L_tmp, blocks_v_L, vec_L, particle_set, &
216 : qs_kind_set, atomic_kind_set, basis_type, cell, &
217 177632 : operator_type)
218 :
219 1284932 : DO i_block = 1, SIZE(blocks_v_L)
220 : blocks_v_L_store(i_block)%block(:, :) = blocks_v_L_store(i_block)%block(:, :) &
221 143739292 : + blocks_v_L(i_block)%block(:, :)
222 : END DO
223 :
224 : END DO
225 : END DO
226 : END DO
227 :
228 : ! add exp(iq*vec_L) * (P 0 | Q vec_L) to V_PQ(q)
229 230272 : DO ik = ikp_start, ikp_end
230 :
231 : ! coskl and sinkl are identical for all i_x_outer, j_y_outer, k_z_outer
232 795648 : coskl = COS(twopi*DOT_PRODUCT(vec_s(1:3), kpoints%xkp(1:3, ik)))
233 795648 : sinkl = SIN(twopi*DOT_PRODUCT(vec_s(1:3), kpoints%xkp(1:3, ik)))
234 :
235 1578560 : DO i_block = 1, SIZE(blocks_v_L)
236 :
237 : blocks_v_kp(ik, 1, i_block)%block(:, :) = blocks_v_kp(ik, 1, i_block)%block(:, :) &
238 301555840 : + coskl*blocks_v_L_store(i_block)%block(:, :)
239 : blocks_v_kp(ik, 2, i_block)%block(:, :) = blocks_v_kp(ik, 2, i_block)%block(:, :) &
240 301754752 : + sinkl*blocks_v_L_store(i_block)%block(:, :)
241 :
242 : END DO
243 :
244 : END DO
245 :
246 311120 : DO i_block = 1, SIZE(blocks_v_L)
247 :
248 22406336 : blocks_v_L_store(i_block)%block(:, :) = 0.0_dp
249 :
250 : END DO
251 :
252 : END DO
253 : END DO
254 : END DO
255 :
256 138 : CALL set_blocks_to_matrix_v_kp(matrix_v_kp, blocks_v_kp, ikp_start, ikp_end)
257 :
258 138 : CALL deallocate_blocks_v_kp(blocks_v_kp)
259 138 : CALL deallocate_blocks_v_L(blocks_v_L)
260 138 : CALL deallocate_blocks_v_L(blocks_v_L_store)
261 :
262 138 : CALL timestop(handle)
263 :
264 138 : END SUBROUTINE
265 :
266 : ! **************************************************************************************************
267 : !> \brief ...
268 : !> \param matrix_v_kp ...
269 : !> \param blocks_v_kp ...
270 : !> \param ikp_start ...
271 : !> \param ikp_end ...
272 : ! **************************************************************************************************
273 138 : SUBROUTINE set_blocks_to_matrix_v_kp(matrix_v_kp, blocks_v_kp, ikp_start, ikp_end)
274 :
275 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_kp
276 : TYPE(two_d_util_type), ALLOCATABLE, &
277 : DIMENSION(:, :, :) :: blocks_v_kp
278 : INTEGER :: ikp_start, ikp_end
279 :
280 : CHARACTER(LEN=*), PARAMETER :: routineN = 'set_blocks_to_matrix_v_kp'
281 :
282 : INTEGER :: col, handle, i_block, i_real_im, ik, row
283 138 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: data_block
284 : TYPE(dbcsr_iterator_type) :: iter
285 :
286 138 : CALL timeset(routineN, handle)
287 :
288 978 : DO ik = ikp_start, ikp_end
289 :
290 2658 : DO i_real_im = 1, 2
291 :
292 1680 : i_block = 1
293 :
294 1680 : CALL dbcsr_iterator_start(iter, matrix_v_kp(ik, i_real_im)%matrix)
295 :
296 13036 : DO WHILE (dbcsr_iterator_blocks_left(iter))
297 :
298 11356 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
299 :
300 2575432 : data_block(:, :) = blocks_v_kp(ik, i_real_im, i_block)%block(:, :)
301 :
302 11356 : i_block = i_block + 1
303 :
304 : END DO
305 :
306 2520 : CALL dbcsr_iterator_stop(iter)
307 :
308 : END DO
309 :
310 : END DO
311 :
312 138 : CALL timestop(handle)
313 :
314 138 : END SUBROUTINE
315 :
316 : ! **************************************************************************************************
317 : !> \brief ...
318 : !> \param matrix_v_L_tmp ...
319 : !> \param blocks_v_L ...
320 : !> \param vec_L ...
321 : !> \param particle_set ...
322 : !> \param qs_kind_set ...
323 : !> \param atomic_kind_set ...
324 : !> \param basis_type ...
325 : !> \param cell ...
326 : !> \param operator_type ...
327 : ! **************************************************************************************************
328 177632 : SUBROUTINE compute_v_transl(matrix_v_L_tmp, blocks_v_L, vec_L, particle_set, &
329 : qs_kind_set, atomic_kind_set, basis_type, cell, operator_type)
330 :
331 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
332 : TYPE(two_d_util_type), ALLOCATABLE, DIMENSION(:) :: blocks_v_L
333 : REAL(KIND=dp), DIMENSION(3) :: vec_L
334 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
335 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
336 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
337 : CHARACTER(LEN=*), INTENT(IN) :: basis_type
338 : TYPE(cell_type), POINTER :: cell
339 : INTEGER :: operator_type
340 :
341 : CHARACTER(LEN=*), PARAMETER :: routineN = 'compute_v_transl'
342 :
343 : INTEGER :: col, handle, i_block, kind_a, kind_b, row
344 177632 : INTEGER, ALLOCATABLE, DIMENSION(:) :: kind_of
345 : REAL(dp), DIMENSION(3) :: ra, rab_L, rb
346 177632 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: data_block
347 177632 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: contr_a, contr_b
348 : TYPE(dbcsr_iterator_type) :: iter
349 : TYPE(gto_basis_set_type), POINTER :: basis_set_a, basis_set_b
350 :
351 177632 : CALL timeset(routineN, handle)
352 :
353 177632 : NULLIFY (basis_set_a, basis_set_b, data_block)
354 :
355 177632 : CALL get_atomic_kind_set(atomic_kind_set, kind_of=kind_of)
356 :
357 177632 : CALL dbcsr_set(matrix_v_L_tmp, 0.0_dp)
358 :
359 177632 : i_block = 1
360 :
361 177632 : CALL dbcsr_iterator_start(iter, matrix_v_L_tmp)
362 :
363 1082692 : DO WHILE (dbcsr_iterator_blocks_left(iter))
364 :
365 905060 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
366 :
367 905060 : kind_a = kind_of(row)
368 905060 : kind_b = kind_of(col)
369 :
370 905060 : CALL get_qs_kind(qs_kind=qs_kind_set(kind_a), basis_set=basis_set_a, basis_type=basis_type)
371 905060 : CALL get_qs_kind(qs_kind=qs_kind_set(kind_b), basis_set=basis_set_b, basis_type=basis_type)
372 :
373 905060 : ra(1:3) = pbc(particle_set(row)%r(1:3), cell)
374 905060 : rb(1:3) = pbc(particle_set(col)%r(1:3), cell)
375 :
376 3620240 : rab_L(1:3) = rb(1:3) - ra(1:3) + vec_L(1:3)
377 :
378 905060 : CALL contraction_matrix_shg(basis_set_a, contr_a)
379 905060 : CALL contraction_matrix_shg(basis_set_b, contr_b)
380 :
381 143236508 : blocks_v_L(i_block)%block = 0.0_dp
382 :
383 : CALL int_operators_r12_ab_shg(operator_type, blocks_v_L(i_block)%block, rab=rab_L, &
384 : fba=basis_set_a, fbb=basis_set_b, scona_shg=contr_a, sconb_shg=contr_b, &
385 905060 : calculate_forces=.FALSE.)
386 :
387 905060 : i_block = i_block + 1
388 :
389 905060 : DEALLOCATE (contr_a, contr_b)
390 :
391 : END DO
392 :
393 177632 : CALL dbcsr_iterator_stop(iter)
394 :
395 177632 : DEALLOCATE (kind_of)
396 :
397 177632 : CALL timestop(handle)
398 :
399 177632 : END SUBROUTINE
400 :
401 : ! **************************************************************************************************
402 : !> \brief ...
403 : !> \param blocks_v_kp ...
404 : ! **************************************************************************************************
405 138 : SUBROUTINE deallocate_blocks_v_kp(blocks_v_kp)
406 : TYPE(two_d_util_type), ALLOCATABLE, &
407 : DIMENSION(:, :, :) :: blocks_v_kp
408 :
409 : CHARACTER(LEN=*), PARAMETER :: routineN = 'deallocate_blocks_v_kp'
410 :
411 : INTEGER :: handle, i_block, i_real_img, ik
412 :
413 138 : CALL timeset(routineN, handle)
414 :
415 1254 : DO ik = LBOUND(blocks_v_kp, 1), UBOUND(blocks_v_kp, 1)
416 2658 : DO i_real_img = 1, SIZE(blocks_v_kp, 2)
417 13876 : DO i_block = 1, SIZE(blocks_v_kp, 3)
418 13036 : DEALLOCATE (blocks_v_kp(ik, i_real_img, i_block)%block)
419 : END DO
420 : END DO
421 : END DO
422 :
423 11494 : DEALLOCATE (blocks_v_kp)
424 :
425 138 : CALL timestop(handle)
426 :
427 138 : END SUBROUTINE
428 :
429 : ! **************************************************************************************************
430 : !> \brief ...
431 : !> \param blocks_v_L ...
432 : ! **************************************************************************************************
433 276 : SUBROUTINE deallocate_blocks_v_L(blocks_v_L)
434 : TYPE(two_d_util_type), ALLOCATABLE, DIMENSION(:) :: blocks_v_L
435 :
436 : CHARACTER(LEN=*), PARAMETER :: routineN = 'deallocate_blocks_v_L'
437 :
438 : INTEGER :: handle, i_block
439 :
440 276 : CALL timeset(routineN, handle)
441 :
442 2590 : DO i_block = 1, SIZE(blocks_v_L, 1)
443 2590 : DEALLOCATE (blocks_v_L(i_block)%block)
444 : END DO
445 :
446 2590 : DEALLOCATE (blocks_v_L)
447 :
448 276 : CALL timestop(handle)
449 :
450 276 : END SUBROUTINE
451 :
452 : ! **************************************************************************************************
453 : !> \brief ...
454 : !> \param blocks_v_L ...
455 : !> \param matrix_v_L_tmp ...
456 : ! **************************************************************************************************
457 276 : SUBROUTINE allocate_blocks_v_L(blocks_v_L, matrix_v_L_tmp)
458 : TYPE(two_d_util_type), ALLOCATABLE, DIMENSION(:) :: blocks_v_L
459 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
460 :
461 : CHARACTER(LEN=*), PARAMETER :: routineN = 'allocate_blocks_v_L'
462 :
463 : INTEGER :: col, handle, i_block, nblocks, row
464 276 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: data_block
465 : TYPE(dbcsr_iterator_type) :: iter
466 :
467 276 : CALL timeset(routineN, handle)
468 :
469 276 : nblocks = 0
470 :
471 276 : CALL dbcsr_iterator_start(iter, matrix_v_L_tmp)
472 :
473 2590 : DO WHILE (dbcsr_iterator_blocks_left(iter))
474 :
475 2314 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
476 :
477 2314 : nblocks = nblocks + 1
478 :
479 : END DO
480 :
481 276 : CALL dbcsr_iterator_stop(iter)
482 :
483 3142 : ALLOCATE (blocks_v_L(nblocks))
484 :
485 276 : i_block = 1
486 :
487 276 : CALL dbcsr_iterator_start(iter, matrix_v_L_tmp)
488 :
489 2590 : DO WHILE (dbcsr_iterator_blocks_left(iter))
490 :
491 2314 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
492 :
493 9256 : ALLOCATE (blocks_v_L(i_block)%block(SIZE(data_block, 1), SIZE(data_block, 2)))
494 236140 : blocks_v_L(i_block)%block = 0.0_dp
495 :
496 4904 : i_block = i_block + 1
497 :
498 : END DO
499 :
500 276 : CALL dbcsr_iterator_stop(iter)
501 :
502 276 : CALL timestop(handle)
503 :
504 276 : END SUBROUTINE
505 :
506 : ! **************************************************************************************************
507 : !> \brief ...
508 : !> \param blocks_v_kp ...
509 : !> \param matrix_v_kp ...
510 : !> \param ikp_start ...
511 : !> \param ikp_end ...
512 : ! **************************************************************************************************
513 138 : SUBROUTINE allocate_blocks_v_kp(blocks_v_kp, matrix_v_kp, ikp_start, ikp_end)
514 : TYPE(two_d_util_type), ALLOCATABLE, &
515 : DIMENSION(:, :, :) :: blocks_v_kp
516 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_kp
517 : INTEGER :: ikp_start, ikp_end
518 :
519 : CHARACTER(LEN=*), PARAMETER :: routineN = 'allocate_blocks_v_kp'
520 :
521 : INTEGER :: col, handle, i_block, i_real_img, ik, &
522 : nblocks, row
523 138 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: data_block
524 : TYPE(dbcsr_iterator_type) :: iter
525 :
526 138 : CALL timeset(routineN, handle)
527 :
528 138 : nblocks = 0
529 :
530 138 : CALL dbcsr_iterator_start(iter, matrix_v_kp(ikp_start, 1)%matrix)
531 :
532 1295 : DO WHILE (dbcsr_iterator_blocks_left(iter))
533 :
534 1157 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
535 :
536 1157 : nblocks = nblocks + 1
537 :
538 : END DO
539 :
540 138 : CALL dbcsr_iterator_stop(iter)
541 :
542 15517 : ALLOCATE (blocks_v_kp(ikp_start:ikp_end, SIZE(matrix_v_kp, 2), nblocks))
543 :
544 978 : DO ik = ikp_start, ikp_end
545 :
546 2658 : DO i_real_img = 1, SIZE(matrix_v_kp, 2)
547 :
548 1680 : i_block = 1
549 :
550 1680 : CALL dbcsr_iterator_start(iter, matrix_v_kp(ik, i_real_img)%matrix)
551 :
552 13036 : DO WHILE (dbcsr_iterator_blocks_left(iter))
553 :
554 11356 : CALL dbcsr_iterator_next_block(iter, row, col, data_block)
555 :
556 0 : ALLOCATE (blocks_v_kp(ik, i_real_img, i_block)%block(SIZE(data_block, 1), &
557 45424 : SIZE(data_block, 2)))
558 2575432 : blocks_v_kp(ik, i_real_img, i_block)%block = 0.0_dp
559 :
560 24392 : i_block = i_block + 1
561 :
562 : END DO
563 :
564 2520 : CALL dbcsr_iterator_stop(iter)
565 :
566 : END DO
567 :
568 : END DO
569 :
570 138 : CALL timestop(handle)
571 :
572 138 : END SUBROUTINE
573 :
574 : ! **************************************************************************************************
575 : !> \brief ...
576 : !> \param cell ...
577 : !> \param kpoints ...
578 : !> \param total_periodicity ...
579 : ! **************************************************************************************************
580 138 : SUBROUTINE check_periodicity(cell, kpoints, total_periodicity)
581 : TYPE(cell_type), POINTER :: cell
582 : TYPE(kpoint_type), POINTER :: kpoints
583 : INTEGER :: total_periodicity
584 :
585 : CHARACTER(LEN=*), PARAMETER :: routineN = 'check_periodicity'
586 :
587 : INTEGER :: handle
588 : INTEGER, DIMENSION(3) :: nkp_grid, periodic
589 :
590 138 : CALL timeset(routineN, handle)
591 :
592 138 : CALL get_cell(cell=cell, periodic=periodic)
593 138 : CALL get_kpoint_info(kpoints, nkp_grid=nkp_grid)
594 :
595 138 : IF (periodic(1) == 0) THEN
596 74 : CPASSERT(nkp_grid(1) == 1)
597 : END IF
598 138 : IF (periodic(2) == 0) THEN
599 12 : CPASSERT(nkp_grid(2) == 1)
600 : END IF
601 138 : IF (periodic(3) == 0) THEN
602 56 : CPASSERT(nkp_grid(3) == 1)
603 : END IF
604 :
605 552 : total_periodicity = SUM(periodic)
606 :
607 138 : CALL timestop(handle)
608 :
609 138 : END SUBROUTINE
610 :
611 : ! **************************************************************************************************
612 : !> \brief ...
613 : !> \param matrix_v_L_tmp ...
614 : !> \param matrix_v_kp ...
615 : !> \param ikp_start ...
616 : ! **************************************************************************************************
617 138 : SUBROUTINE allocate_tmp(matrix_v_L_tmp, matrix_v_kp, ikp_start)
618 :
619 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
620 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_v_kp
621 : INTEGER :: ikp_start
622 :
623 : CHARACTER(LEN=*), PARAMETER :: routineN = 'allocate_tmp'
624 :
625 : INTEGER :: handle
626 :
627 138 : CALL timeset(routineN, handle)
628 :
629 138 : NULLIFY (matrix_v_L_tmp)
630 138 : CALL dbcsr_init_p(matrix_v_L_tmp)
631 : CALL dbcsr_create(matrix=matrix_v_L_tmp, &
632 : template=matrix_v_kp(ikp_start, 1)%matrix, &
633 138 : matrix_type=dbcsr_type_no_symmetry)
634 138 : CALL dbcsr_reserve_all_blocks(matrix_v_L_tmp)
635 138 : CALL dbcsr_set(matrix_v_L_tmp, 0.0_dp)
636 :
637 138 : CALL timestop(handle)
638 :
639 138 : END SUBROUTINE
640 :
641 : ! **************************************************************************************************
642 : !> \brief ...
643 : !> \param matrix_v_L_tmp ...
644 : ! **************************************************************************************************
645 138 : SUBROUTINE deallocate_tmp(matrix_v_L_tmp)
646 :
647 : TYPE(dbcsr_type), POINTER :: matrix_v_L_tmp
648 :
649 : CHARACTER(LEN=*), PARAMETER :: routineN = 'deallocate_tmp'
650 :
651 : INTEGER :: handle
652 :
653 138 : CALL timeset(routineN, handle)
654 :
655 138 : CALL dbcsr_release_p(matrix_v_L_tmp)
656 :
657 138 : CALL timestop(handle)
658 :
659 138 : END SUBROUTINE
660 :
661 0 : END MODULE kpoint_coulomb_2c
|
