38 dbt_batched_contract_finalize, dbt_batched_contract_init, dbt_contract, dbt_copy, &
39 dbt_copy_matrix_to_tensor, dbt_copy_tensor_to_matrix, dbt_create, dbt_destroy, dbt_filter, &
40 dbt_get_info, dbt_nblks_total, dbt_nd_mp_comm, dbt_pgrid_destroy, dbt_pgrid_type, dbt_type
66#include "./base/base_uses.f90"
72 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'rpa_im_time'
97 TYPE(
dbcsr_p_type),
DIMENSION(:, :, :),
POINTER :: propagator
98 INTEGER,
INTENT(IN) :: ntime
99 TYPE(
dbcsr_type),
INTENT(IN) :: matrix_template
100 INTEGER,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell
102 INTEGER :: icell, isector, jtime, ncell
105 cpassert(
SIZE(index_to_cell, 1) == 3)
106 ncell =
SIZE(index_to_cell, 2)
114 ALLOCATE (propagator(isector, jtime, icell)%matrix)
115 CALL dbcsr_create(matrix=propagator(isector, jtime, icell)%matrix, &
116 template=matrix_template, &
117 matrix_type=dbcsr_type_no_symmetry)
164 t_3c_M, t_3c_O, t_3c_O_compressed, t_3c_O_ind, &
165 starts_array_mc, ends_array_mc, &
166 starts_array_mc_block, ends_array_mc_block, &
167 grid, e_fermi, eps_filter, &
168 alpha, eps_filter_im_time, Eigenval, nmo, &
169 cut_memory, unit_nr, &
171 qs_env, do_kpoints_from_Gamma, &
172 index_to_cell_3c, cell_to_index_3c, &
173 has_mat_P_blocks, do_ri_sos_laplace_mp2, &
174 dbcsr_time, dbcsr_nflop)
175 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
INTENT(IN) :: mat_p_omega
177 INTEGER,
INTENT(IN) :: homo
180 INTEGER,
INTENT(IN) :: ispin
181 TYPE(dbt_type),
INTENT(INOUT) :: t_3c_m
182 TYPE(dbt_type),
DIMENSION(:, :),
INTENT(INOUT) :: t_3c_o
184 INTENT(INOUT) :: t_3c_o_compressed
186 INTENT(INOUT) :: t_3c_o_ind
187 INTEGER,
DIMENSION(:),
INTENT(IN) :: starts_array_mc, ends_array_mc, &
188 starts_array_mc_block, &
191 REAL(kind=
dp),
INTENT(IN) :: e_fermi, eps_filter, alpha, &
193 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval
194 INTEGER,
INTENT(IN) :: nmo, cut_memory, unit_nr
198 LOGICAL,
INTENT(IN) :: do_kpoints_from_gamma
199 INTEGER,
ALLOCATABLE,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell_3c
200 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :), &
201 INTENT(IN) :: cell_to_index_3c
202 LOGICAL,
DIMENSION(:, :, :, :, :),
INTENT(INOUT) :: has_mat_p_blocks
203 LOGICAL,
INTENT(IN) :: do_ri_sos_laplace_mp2
204 REAL(
dp),
INTENT(INOUT) :: dbcsr_time
205 INTEGER(int_8),
INTENT(INOUT) :: dbcsr_nflop
207 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_mat_P_omega'
209 INTEGER :: comm_2d_handle, handle, handle2, handle3, i, i_cell, i_cell_r_1, &
210 i_cell_r_1_minus_s, i_cell_r_1_minus_t, i_cell_r_2, i_cell_r_2_minus_s_minus_t, i_cell_s, &
211 i_cell_t, i_mem, iquad, j, j_mem, jquad, num_3c_repl, num_cells_dm, num_integ_points, &
213 INTEGER(int_8) :: nze, nze_dm_occ, nze_dm_virt, nze_m_occ, &
215 INTEGER(KIND=int_8) :: flops_1_occ, flops_1_virt, flops_2
216 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: dist_1, dist_2, mc_ranges, size_dm, &
218 INTEGER,
DIMENSION(2) :: pdims_2d
219 INTEGER,
DIMENSION(2, 1) :: ibounds_2, jbounds_2
220 INTEGER,
DIMENSION(2, 2) :: ibounds_1, jbounds_1
221 INTEGER,
DIMENSION(3) :: bounds_3c
222 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell_dm
223 LOGICAL :: do_gamma_rpa, do_kpoints_cubic_rpa, first_cycle_im_time, first_cycle_omega_loop, &
224 memory_info, r_1_minus_s_needed, r_1_minus_t_needed, r_2_minus_s_minus_t_needed
225 REAL(
dp) :: occ, occ_dm_occ, occ_dm_virt, occ_m_occ, &
226 occ_m_virt, occ_o, t1_flop
227 REAL(kind=
dp) :: omega, omega_old, t1, t2, tau, weight, &
230 TYPE(
dbcsr_p_type),
DIMENSION(:, :, :),
POINTER :: propagator
231 TYPE(dbt_pgrid_type) :: pgrid_2d
232 TYPE(dbt_type) :: t_3c_m_occ, t_3c_m_occ_tmp, t_3c_m_virt, &
233 t_3c_m_virt_tmp, t_dm, t_dm_tmp, t_p, &
235 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:) :: t_dm_occ, t_dm_virt
236 TYPE(dbt_type),
ALLOCATABLE,
DIMENSION(:, :) :: t_3c_o_occ, t_3c_o_virt
239 CALL timeset(routinen, handle)
241 IF (do_ri_sos_laplace_mp2)
THEN
242 num_integ_points =
SIZE(grid%imaginary_time)
244 num_integ_points =
SIZE(grid%frequency)
249 memory_info = mp2_env%ri_rpa_im_time%memory_info
250 IF (memory_info)
THEN
251 unit_nr_dbcsr = unit_nr
256 do_kpoints_cubic_rpa = qs_env%mp2_env%ri_rpa_im_time%do_im_time_kpoints
257 do_gamma_rpa = .NOT. do_kpoints_cubic_rpa
258 num_3c_repl = maxval(cell_to_index_3c)
260 first_cycle_im_time = .true.
261 ALLOCATE (t_3c_o_occ(
SIZE(t_3c_o, 1),
SIZE(t_3c_o, 2)), t_3c_o_virt(
SIZE(t_3c_o, 1),
SIZE(t_3c_o, 2)))
262 DO i = 1,
SIZE(t_3c_o, 1)
263 DO j = 1,
SIZE(t_3c_o, 2)
264 CALL dbt_create(t_3c_o(i, j), t_3c_o_occ(i, j))
265 CALL dbt_create(t_3c_o(i, j), t_3c_o_virt(i, j))
269 CALL dbt_create(t_3c_m, t_3c_m_occ, name=
"M occ (RI | AO AO)")
270 CALL dbt_create(t_3c_m, t_3c_m_virt, name=
"M virt (RI | AO AO)")
272 ALLOCATE (mc_ranges(cut_memory + 1))
273 mc_ranges(:cut_memory) = starts_array_mc_block(:)
274 mc_ranges(cut_memory + 1) = ends_array_mc_block(cut_memory) + 1
276 DO jquad = 1, num_integ_points
283 eigenval, e_fermi, eps_filter, memory_info, unit_nr, &
284 jquad, do_kpoints_cubic_rpa, do_kpoints_from_gamma, qs_env, &
285 num_cells_dm, index_to_cell_dm, para_env)
287 ALLOCATE (t_dm_virt(num_cells_dm))
288 ALLOCATE (t_dm_occ(num_cells_dm))
291 CALL comm_2d%set_handle(comm_2d_handle)
293 pgrid_2d = dbt_nd_mp_comm(comm_2d, [1], [2], pdims_2d=pdims_2d)
294 ALLOCATE (size_p(dbt_nblks_total(t_3c_m, 1)))
295 CALL dbt_get_info(t_3c_m, blk_size_1=size_p)
297 ALLOCATE (size_dm(dbt_nblks_total(t_3c_o(1, 1), 3)))
298 CALL dbt_get_info(t_3c_o(1, 1), blk_size_3=size_dm)
299 CALL create_2c_tensor(t_dm, dist_1, dist_2, pgrid_2d, size_dm, size_dm, name=
"D (AO | AO)")
301 DEALLOCATE (dist_1, dist_2)
302 CALL create_2c_tensor(t_p, dist_1, dist_2, pgrid_2d, size_p, size_p, name=
"P (RI | RI)")
304 DEALLOCATE (dist_1, dist_2)
305 CALL dbt_pgrid_destroy(pgrid_2d)
307 DO i_cell = 1, num_cells_dm
308 CALL dbt_create(t_dm, t_dm_virt(i_cell), name=
"D virt (AO | AO)")
311 CALL dbt_copy(t_dm_tmp, t_dm_virt(i_cell), move_data=.true.)
314 CALL dbt_create(t_dm, t_dm_occ(i_cell), name=
"D occ (AO | AO)")
316 CALL dbt_copy(t_dm_tmp, t_dm_occ(i_cell), move_data=.true.)
317 CALL dbt_destroy(t_dm_tmp)
324 CALL dbt_destroy(t_dm)
326 CALL dbt_create(t_3c_o_occ(1, 1), t_3c_m_occ_tmp, name=
"M (RI AO | AO)")
327 CALL dbt_create(t_3c_o_virt(1, 1), t_3c_m_virt_tmp, name=
"M (RI AO | AO)")
329 CALL timeset(routinen//
"_contract", handle2)
334 DO i = 1,
SIZE(t_3c_o_occ, 1)
335 DO j = 1,
SIZE(t_3c_o_occ, 2)
336 CALL dbt_batched_contract_init(t_3c_o_occ(i, j), batch_range_2=mc_ranges, batch_range_3=mc_ranges)
339 DO i = 1,
SIZE(t_3c_o_virt, 1)
340 DO j = 1,
SIZE(t_3c_o_virt, 2)
341 CALL dbt_batched_contract_init(t_3c_o_virt(i, j), batch_range_2=mc_ranges, batch_range_3=mc_ranges)
344 CALL dbt_batched_contract_init(t_3c_m_occ_tmp, batch_range_2=mc_ranges, batch_range_3=mc_ranges)
345 CALL dbt_batched_contract_init(t_3c_m_virt_tmp, batch_range_2=mc_ranges, batch_range_3=mc_ranges)
346 CALL dbt_batched_contract_init(t_3c_m_occ, batch_range_2=mc_ranges, batch_range_3=mc_ranges)
347 CALL dbt_batched_contract_init(t_3c_m_virt, batch_range_2=mc_ranges, batch_range_3=mc_ranges)
349 DO i_cell_t = 1, num_cells_dm/2 + 1
351 IF (.NOT. any(has_mat_p_blocks(i_cell_t, :, :, :, :))) cycle
353 CALL dbt_batched_contract_init(t_p)
355 IF (do_gamma_rpa)
THEN
364 DO j_mem = 1, cut_memory
366 CALL dbt_get_info(t_3c_o_occ(1, 1), nfull_total=bounds_3c)
368 jbounds_1(:, 1) = [1, bounds_3c(1)]
369 jbounds_1(:, 2) = [starts_array_mc(j_mem), ends_array_mc(j_mem)]
371 jbounds_2(:, 1) = [starts_array_mc(j_mem), ends_array_mc(j_mem)]
373 IF (do_gamma_rpa)
CALL dbt_batched_contract_init(t_dm_virt(1))
375 DO i_mem = 1, cut_memory
377 IF (.NOT. any(has_mat_p_blocks(i_cell_t, i_mem, j_mem, :, :))) cycle
379 ibounds_1(:, 1) = [1, bounds_3c(1)]
380 ibounds_1(:, 2) = [starts_array_mc(i_mem), ends_array_mc(i_mem)]
382 ibounds_2(:, 1) = [starts_array_mc(i_mem), ends_array_mc(i_mem)]
384 IF (unit_nr_dbcsr > 0)
WRITE (unit=unit_nr_dbcsr, fmt=
"(T3,A,I3,1X,I3)") &
385 "RPA_LOW_SCALING_INFO| Memory Cut iteration", i_mem, j_mem
387 DO i_cell_r_1 = 1, num_3c_repl
389 DO i_cell_r_2 = 1, num_3c_repl
391 IF (.NOT. has_mat_p_blocks(i_cell_t, i_mem, j_mem, i_cell_r_1, i_cell_r_2)) cycle
393 CALL get_diff_index_3c(i_cell_r_1, i_cell_t, i_cell_r_1_minus_t, &
394 index_to_cell_3c, cell_to_index_3c, index_to_cell_dm, &
395 r_1_minus_t_needed, do_kpoints_cubic_rpa)
397 IF (do_gamma_rpa)
CALL dbt_batched_contract_init(t_dm_occ(1))
398 DO i_cell_s = 1, num_cells_dm
399 CALL get_diff_index_3c(i_cell_r_1, i_cell_s, i_cell_r_1_minus_s, index_to_cell_3c, &
400 cell_to_index_3c, index_to_cell_dm, r_1_minus_s_needed, &
401 do_kpoints_cubic_rpa)
402 IF (r_1_minus_s_needed)
THEN
404 CALL timeset(routinen//
"_calc_M_occ_t", handle3)
406 t_3c_o_ind(i_cell_r_1_minus_s, i_cell_r_2, j_mem)%ind, &
407 t_3c_o_compressed(i_cell_r_1_minus_s, i_cell_r_2, j_mem), &
408 qs_env%mp2_env%ri_rpa_im_time%eps_compress)
410 IF (do_gamma_rpa .AND. i_mem == 1)
THEN
416 CALL dbt_copy(t_3c_o(i_cell_r_1_minus_s, i_cell_r_2), &
417 t_3c_o_occ(i_cell_r_1_minus_s, i_cell_r_2), move_data=.true.)
419 CALL dbt_contract(alpha=1.0_dp, &
420 tensor_1=t_3c_o_occ(i_cell_r_1_minus_s, i_cell_r_2), &
421 tensor_2=t_dm_occ(i_cell_s), &
423 tensor_3=t_3c_m_occ_tmp, &
424 contract_1=[3], notcontract_1=[1, 2], &
425 contract_2=[2], notcontract_2=[1], &
426 map_1=[1, 2], map_2=[3], &
427 bounds_2=jbounds_1, bounds_3=ibounds_2, &
428 filter_eps=eps_filter, unit_nr=unit_nr_dbcsr, &
430 CALL timestop(handle3)
432 dbcsr_nflop = dbcsr_nflop + flops_1_occ
437 IF (do_gamma_rpa)
CALL dbt_batched_contract_finalize(t_dm_occ(1))
442 CALL timeset(routinen//
"_copy_M_occ_t", handle3)
443 CALL dbt_copy(t_3c_m_occ_tmp, t_3c_m_occ, order=[1, 3, 2], move_data=.true.)
444 CALL dbt_filter(t_3c_m_occ, eps_filter)
445 CALL timestop(handle3)
447 IF (do_gamma_rpa)
THEN
449 nze_m_occ = nze_m_occ + nze
450 occ_m_occ = occ_m_occ + occ
453 DO i_cell_s = 1, num_cells_dm
454 CALL get_diff_diff_index_3c(i_cell_r_2, i_cell_s, i_cell_t, i_cell_r_2_minus_s_minus_t, &
455 index_to_cell_3c, cell_to_index_3c, index_to_cell_dm, &
456 r_2_minus_s_minus_t_needed, do_kpoints_cubic_rpa)
458 IF (r_1_minus_t_needed .AND. r_2_minus_s_minus_t_needed)
THEN
460 t_3c_o_ind(i_cell_r_2_minus_s_minus_t, i_cell_r_1_minus_t, i_mem)%ind, &
461 t_3c_o_compressed(i_cell_r_2_minus_s_minus_t, i_cell_r_1_minus_t, i_mem), &
462 qs_env%mp2_env%ri_rpa_im_time%eps_compress)
464 CALL dbt_copy(t_3c_o(i_cell_r_2_minus_s_minus_t, i_cell_r_1_minus_t), &
465 t_3c_o_virt(i_cell_r_2_minus_s_minus_t, i_cell_r_1_minus_t), move_data=.true.)
467 CALL timeset(routinen//
"_calc_M_virt_t", handle3)
468 CALL dbt_contract(alpha=alpha/2.0_dp, &
469 tensor_1=t_3c_o_virt( &
470 i_cell_r_2_minus_s_minus_t, i_cell_r_1_minus_t), &
471 tensor_2=t_dm_virt(i_cell_s), &
473 tensor_3=t_3c_m_virt_tmp, &
474 contract_1=[3], notcontract_1=[1, 2], &
475 contract_2=[2], notcontract_2=[1], &
476 map_1=[1, 2], map_2=[3], &
477 bounds_2=ibounds_1, bounds_3=jbounds_2, &
478 filter_eps=eps_filter, unit_nr=unit_nr_dbcsr, &
480 CALL timestop(handle3)
482 dbcsr_nflop = dbcsr_nflop + flops_1_virt
487 CALL timeset(routinen//
"_copy_M_virt_t", handle3)
488 CALL dbt_copy(t_3c_m_virt_tmp, t_3c_m_virt, move_data=.true.)
489 CALL dbt_filter(t_3c_m_virt, eps_filter)
490 CALL timestop(handle3)
492 IF (do_gamma_rpa)
THEN
494 nze_m_virt = nze_m_virt + nze
495 occ_m_virt = occ_m_virt + occ
500 CALL timeset(routinen//
"_calc_P_t", handle3)
502 CALL dbt_contract(alpha=1.0_dp, tensor_1=t_3c_m_occ, &
503 tensor_2=t_3c_m_virt, &
506 contract_1=[2, 3], notcontract_1=[1], &
507 contract_2=[2, 3], notcontract_2=[1], &
508 map_1=[1], map_2=[2], &
509 filter_eps=eps_filter_im_time/real(cut_memory**2, kind=
dp), &
512 unit_nr=unit_nr_dbcsr)
514 CALL timestop(handle3)
516 first_cycle_im_time = .false.
518 IF (jquad == 1 .AND. flops_2 == 0)
THEN
519 has_mat_p_blocks(i_cell_t, i_mem, j_mem, i_cell_r_1, i_cell_r_2) = .false.
525 IF (do_gamma_rpa)
CALL dbt_batched_contract_finalize(t_dm_virt(1))
528 CALL dbt_batched_contract_finalize(t_p, unit_nr=unit_nr_dbcsr)
530 CALL dbt_create(mat_p_global%matrix, t_p_tmp)
531 CALL dbt_copy(t_p, t_p_tmp, move_data=.true.)
532 CALL dbt_copy_tensor_to_matrix(t_p_tmp, mat_p_global%matrix)
533 CALL dbt_destroy(t_p_tmp)
535 IF (do_ri_sos_laplace_mp2)
THEN
539 CALL dbcsr_add(mat_p_omega(jquad, i_cell_t)%matrix, mat_p_global%matrix, 1.0_dp, 1.0_dp)
541 CALL timeset(routinen//
"_Fourier_transform", handle3)
544 first_cycle_omega_loop = .true.
546 tau = grid%imaginary_time(jquad)
548 DO iquad = 1, num_integ_points
550 omega = grid%frequency(iquad)
551 weight = grid%cosine_time_to_frequency_weights(iquad, jquad)
553 IF (first_cycle_omega_loop)
THEN
556 CALL dbcsr_scale(mat_p_global%matrix, cos(omega*tau)*weight)
558 CALL dbcsr_scale(mat_p_global%matrix, cos(omega*tau)/cos(omega_old*tau)*weight/weight_old)
561 CALL dbcsr_add(mat_p_omega(iquad, i_cell_t)%matrix, mat_p_global%matrix, 1.0_dp, 1.0_dp)
563 first_cycle_omega_loop = .false.
570 CALL timestop(handle3)
575 CALL timestop(handle2)
577 CALL dbt_batched_contract_finalize(t_3c_m_occ_tmp)
578 CALL dbt_batched_contract_finalize(t_3c_m_virt_tmp)
579 CALL dbt_batched_contract_finalize(t_3c_m_occ)
580 CALL dbt_batched_contract_finalize(t_3c_m_virt)
582 DO i = 1,
SIZE(t_3c_o_occ, 1)
583 DO j = 1,
SIZE(t_3c_o_occ, 2)
584 CALL dbt_batched_contract_finalize(t_3c_o_occ(i, j))
588 DO i = 1,
SIZE(t_3c_o_virt, 1)
589 DO j = 1,
SIZE(t_3c_o_virt, 2)
590 CALL dbt_batched_contract_finalize(t_3c_o_virt(i, j))
594 CALL dbt_destroy(t_p)
595 DO i_cell = 1, num_cells_dm
596 CALL dbt_destroy(t_dm_virt(i_cell))
597 CALL dbt_destroy(t_dm_occ(i_cell))
600 CALL dbt_destroy(t_3c_m_occ_tmp)
601 CALL dbt_destroy(t_3c_m_virt_tmp)
602 DEALLOCATE (t_dm_virt)
603 DEALLOCATE (t_dm_occ)
608 dbcsr_time = dbcsr_time + t2 - t1_flop
610 IF (unit_nr > 0)
THEN
611 WRITE (unit_nr,
'(/T3,A,1X,I3)') &
612 'RPA_LOW_SCALING_INFO| Info for time point', jquad
613 WRITE (unit_nr,
'(T6,A,T56,F25.1)') &
614 'Execution time (s):', t2 - t1
615 WRITE (unit_nr,
'(T6,A,T63,ES7.1,1X,A1,1X,F7.3,A1)') &
616 'Occupancy of D occ:', real(nze_dm_occ,
dp),
'/', occ_dm_occ*100,
'%'
617 WRITE (unit_nr,
'(T6,A,T63,ES7.1,1X,A1,1X,F7.3,A1)') &
618 'Occupancy of D virt:', real(nze_dm_virt,
dp),
'/', occ_dm_virt*100,
'%'
619 IF (do_gamma_rpa)
THEN
620 WRITE (unit_nr,
'(T6,A,T63,ES7.1,1X,A1,1X,F7.3,A1)') &
621 'Occupancy of 3c ints:', real(nze_o,
dp),
'/', occ_o*100,
'%'
622 WRITE (unit_nr,
'(T6,A,T63,ES7.1,1X,A1,1X,F7.3,A1)') &
623 'Occupancy of M occ:', real(nze_m_occ,
dp),
'/', occ_m_occ*100,
'%'
624 WRITE (unit_nr,
'(T6,A,T63,ES7.1,1X,A1,1X,F7.3,A1)') &
625 'Occupancy of M virt:', real(nze_m_virt,
dp),
'/', occ_m_virt*100,
'%'
633 CALL dbt_destroy(t_3c_m_occ)
634 CALL dbt_destroy(t_3c_m_virt)
636 DO i = 1,
SIZE(t_3c_o, 1)
637 DO j = 1,
SIZE(t_3c_o, 2)
638 CALL dbt_destroy(t_3c_o_occ(i, j))
639 CALL dbt_destroy(t_3c_o_virt(i, j))
645 CALL timestop(handle)
654 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
INTENT(IN) :: mat_p_omega
656 INTEGER :: i_kp, jquad
658 DO jquad = 1,
SIZE(mat_p_omega, 1)
659 DO i_kp = 1,
SIZE(mat_p_omega, 2)
661 CALL dbcsr_set(mat_p_omega(jquad, i_kp)%matrix, 0.0_dp)
692 Eigenval, e_fermi, eps_filter, memory_info, unit_nr, &
693 jquad, do_kpoints_cubic_RPA, do_kpoints_from_Gamma, qs_env, &
694 num_cells_dm, index_to_cell_dm, para_env)
697 INTEGER,
INTENT(IN) :: nmo
699 INTEGER,
INTENT(IN) :: homo
701 INTENT(INOUT),
POINTER :: propagator
702 TYPE(
dbcsr_p_type),
DIMENSION(:),
INTENT(IN) :: matrix_s
703 INTEGER,
INTENT(IN) :: ispin
704 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval
705 REAL(kind=
dp),
INTENT(IN) :: e_fermi, eps_filter
706 LOGICAL,
INTENT(IN) :: memory_info
707 INTEGER,
INTENT(IN) :: unit_nr, jquad
708 LOGICAL,
INTENT(IN) :: do_kpoints_cubic_rpa, &
709 do_kpoints_from_gamma
711 INTEGER,
INTENT(OUT) :: num_cells_dm
712 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell_dm
715 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_mat_dm_global'
717 INTEGER :: handle, num_integ_points
718 INTEGER,
DIMENSION(3, 1) :: index_to_cell_zero
721 index_to_cell_zero = 0
722 num_integ_points =
SIZE(grid%imaginary_time)
724 CALL timeset(routinen, handle)
726 IF (memory_info .AND. unit_nr > 0)
WRITE (unit=unit_nr, fmt=
"(T3,A,T75,i6)") &
727 "RPA_LOW_SCALING_INFO| Started with time point: ", jquad
729 tau = grid%imaginary_time(jquad)
731 IF (do_kpoints_cubic_rpa)
THEN
734 ispin, num_integ_points, jquad, e_fermi, tau, &
735 eps_filter, num_cells_dm, index_to_cell_dm, &
739 ispin, num_integ_points, jquad, e_fermi, tau, &
740 eps_filter, num_cells_dm, index_to_cell_dm, &
743 ELSE IF (do_kpoints_from_gamma)
THEN
746 ispin, num_integ_points, jquad, e_fermi, tau, &
750 ispin, num_integ_points, jquad, e_fermi, tau, &
765 eps_filter, e_fermi, tau, para_env)
777 CALL timestop(handle)
795 eps_filter, e_fermi, tau, para_env)
797 INTENT(INOUT),
POINTER :: propagator
798 INTEGER,
INTENT(IN) :: jquad
800 INTEGER,
INTENT(IN) :: homo
801 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: eigenval
802 INTEGER,
INTENT(IN) :: nmo
803 REAL(kind=
dp),
INTENT(IN) :: eps_filter, e_fermi, tau
806 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_gamma_propagator'
807 REAL(kind=
dp),
PARAMETER :: stabilize_exp = 70.0_dp
809 INTEGER :: handle, i_global, iib, jjb, nao, &
810 ncol_local, nrow_local
811 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
813 cfm_mo_coeff_virt_scaled, &
814 cfm_scaled_dm_occ_tau, &
815 cfm_scaled_dm_virt_tau
817 fm_scaled_dm_virt_tau
818 TYPE(
dbcsr_type),
POINTER :: mat_occ, mat_virt
820 CALL timeset(routinen, handle)
828 nrow_local=nrow_local, &
829 ncol_local=ncol_local, &
830 row_indices=row_indices, &
831 col_indices=col_indices)
833 CALL cp_cfm_create(cfm_scaled_dm_occ_tau, cfm_mo_coeff%matrix_struct, nrow=nao, ncol=nao)
834 CALL cp_cfm_create(cfm_scaled_dm_virt_tau, cfm_mo_coeff%matrix_struct, nrow=nao, ncol=nao)
835 CALL cp_cfm_create(cfm_mo_coeff_occ_scaled, cfm_mo_coeff%matrix_struct)
836 CALL cp_cfm_create(cfm_mo_coeff_virt_scaled, cfm_mo_coeff%matrix_struct)
838 DO jjb = 1, nrow_local
839 DO iib = 1, ncol_local
840 i_global = col_indices(iib)
841 IF (i_global <= homo .AND. abs(tau*0.5_dp*(eigenval(i_global) - e_fermi)) < stabilize_exp)
THEN
842 cfm_mo_coeff_occ_scaled%local_data(jjb, iib) = &
843 cfm_mo_coeff%local_data(jjb, iib)*exp(tau*0.5_dp*(eigenval(i_global) - e_fermi))
845 cfm_mo_coeff_occ_scaled%local_data(jjb, iib) = (0.0_dp, 0.0_dp)
850 DO jjb = 1, nrow_local
851 DO iib = 1, ncol_local
852 i_global = col_indices(iib)
853 IF (i_global > homo .AND. abs(tau*0.5_dp*(eigenval(i_global) - e_fermi)) < stabilize_exp)
THEN
854 cfm_mo_coeff_virt_scaled%local_data(jjb, iib) = &
855 cfm_mo_coeff%local_data(jjb, iib)*exp(-tau*0.5_dp*(eigenval(i_global) - e_fermi))
857 cfm_mo_coeff_virt_scaled%local_data(jjb, iib) = (0.0_dp, 0.0_dp)
864 CALL parallel_gemm(transa=
"N", transb=
"C", m=nao, n=nao, k=nmo, alpha=(1.0_dp, 0.0_dp), &
865 matrix_a=cfm_mo_coeff_occ_scaled, matrix_b=cfm_mo_coeff_occ_scaled, beta=(0.0_dp, 0.0_dp), &
866 matrix_c=cfm_scaled_dm_occ_tau)
867 CALL parallel_gemm(transa=
"N", transb=
"C", m=nao, n=nao, k=nmo, alpha=(1.0_dp, 0.0_dp), &
868 matrix_a=cfm_mo_coeff_virt_scaled, matrix_b=cfm_mo_coeff_virt_scaled, beta=(0.0_dp, 0.0_dp), &
869 matrix_c=cfm_scaled_dm_virt_tau)
871 CALL cp_fm_create(fm_scaled_dm_occ_tau, cfm_scaled_dm_occ_tau%matrix_struct)
872 CALL cp_fm_create(fm_scaled_dm_virt_tau, cfm_scaled_dm_virt_tau%matrix_struct)
873 CALL cp_cfm_to_fm(cfm_scaled_dm_occ_tau, fm_scaled_dm_occ_tau)
874 CALL cp_cfm_to_fm(cfm_scaled_dm_virt_tau, fm_scaled_dm_virt_tau)
880 CALL copy_fm_to_dbcsr(fm_scaled_dm_virt_tau, mat_virt, keep_sparsity=.false.)
889 CALL timestop(handle)
903 SUBROUTINE build_kpoint_density_matrix_rpa(kpoint, ikpgr, ispin, tau, e_fermi, sector, density)
906 INTEGER,
INTENT(IN) :: ikpgr, ispin
907 REAL(kind=
dp),
INTENT(IN) :: tau, e_fermi
908 INTEGER,
INTENT(IN) :: sector
911 CHARACTER(LEN=*),
PARAMETER :: routinen =
'build_kpoint_density_matrix_rpa'
912 REAL(kind=
dp),
PARAMETER :: stabilize_exp = 70.0_dp
914 INTEGER :: handle, i_mo, nao, nmo, nspin
915 REAL(kind=
dp),
DIMENSION(:),
POINTER :: eigenvalues, exp_scaling, occupation
920 CALL timeset(routinen, handle)
923 cpassert(kpoint%use_real_wfn .EQV. .false.)
925 kp => kpoint%kp_env(ikpgr)%kpoint_env
926 nspin =
SIZE(kp%mos, 2)
927 cpassert(ispin >= 1 .AND. ispin <= nspin)
928 mo_set => kp%mos(1, ispin)
929 mo_set_im => kp%mos(2, ispin)
930 CALL get_mo_set(mo_set, nao=nao, nmo=nmo, eigenvalues=eigenvalues, occupation_numbers=occupation)
936 ALLOCATE (exp_scaling(nmo))
937 CALL cp_cfm_create(cfm_mo_coeff, mo_set%mo_coeff%matrix_struct)
939 CALL cp_fm_to_cfm(msourcer=mo_set%mo_coeff, msourcei=mo_set_im%mo_coeff, &
940 mtarget=cfm_mo_coeff)
947 cmplx(2.0_dp/real(nspin, kind=
dp) - occupation, 0.0_dp, kind=
dp))
955 IF (abs(tau*0.5_dp*(eigenvalues(i_mo) - e_fermi)) < stabilize_exp)
THEN
956 exp_scaling(i_mo) = exp(-abs(tau*(eigenvalues(i_mo) - e_fermi)))
958 exp_scaling(i_mo) = 0.0_dp
963 CALL parallel_gemm(
"N",
"C", nao, nao, nmo, (1.0_dp, 0.0_dp), cfm_mo_coeff, cfm_work, &
964 (0.0_dp, 0.0_dp), density)
968 DEALLOCATE (exp_scaling)
970 CALL timestop(handle)
972 END SUBROUTINE build_kpoint_density_matrix_rpa
988 SUBROUTINE compute_transl_dm(propagator, qs_env, ispin, num_integ_points, jquad, e_fermi, tau, &
989 eps_filter, num_cells_dm, index_to_cell_dm, &
992 INTENT(INOUT),
POINTER :: propagator
994 INTEGER,
INTENT(IN) :: ispin, num_integ_points, jquad
995 REAL(kind=
dp),
INTENT(IN) :: e_fermi, tau, eps_filter
996 INTEGER,
INTENT(OUT) :: num_cells_dm
997 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell_dm
998 INTEGER,
INTENT(IN) :: sector
1000 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_transl_dm'
1002 INTEGER :: handle, i_dim, i_img, nao
1003 INTEGER,
DIMENSION(3) :: cell_grid_dm
1006 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: mat_dm_global_work
1007 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_s_kp
1010 CALL timeset(routinen, handle)
1013 matrix_s_kp=matrix_s_kp, &
1020 cell_grid_dm(i_dim) = (kpoints%nkp_grid(i_dim)/2)*2 - 1
1023 num_cells_dm = cell_grid_dm(1)*cell_grid_dm(2)*cell_grid_dm(3)
1025 NULLIFY (mat_dm_global_work)
1029 DO i_img = 1, num_cells_dm
1031 ALLOCATE (mat_dm_global_work(i_img)%matrix)
1032 CALL dbcsr_create(matrix=mat_dm_global_work(i_img)%matrix, &
1033 template=matrix_s_kp(1, 1)%matrix, &
1035 matrix_type=dbcsr_type_no_symmetry)
1042 CALL get_mo_set(kpoints%kp_env(1)%kpoint_env%mos(1, 1), nao=nao)
1043 CALL cp_cfm_create(cfm_density, kpoints%kp_env(1)%kpoint_env%mos(1, 1)%mo_coeff%matrix_struct, &
1047 CALL init_cell_index_rpa(cell_grid_dm, kpoints%cell_to_index, kpoints%index_to_cell, cell)
1051 CALL density_matrix_from_kp_to_transl(kpoints, cfm_density, mat_dm_global_work, kpoints%index_to_cell, &
1052 ispin, tau, e_fermi, sector)
1057 index_to_cell_dm => kpoints%index_to_cell
1060 IF (.NOT.
ASSOCIATED(propagator))
THEN
1062 kpoints%index_to_cell)
1065 DO i_img = 1, num_cells_dm
1068 CALL dbcsr_filter(mat_dm_global_work(i_img)%matrix, eps_filter)
1070 CALL dbcsr_copy(propagator(sector, jquad, i_img)%matrix, &
1071 mat_dm_global_work(i_img)%matrix)
1077 CALL timestop(handle)
1095 INTENT(INOUT),
POINTER :: propagator
1097 INTEGER,
INTENT(IN) :: ispin, num_integ_points, jquad
1098 REAL(kind=
dp),
INTENT(IN) :: e_fermi, tau
1099 INTEGER,
INTENT(IN) :: sector
1101 CHARACTER(LEN=*),
PARAMETER :: routinen =
'compute_periodic_dm'
1103 INTEGER :: handle, nao, num_cells_dm
1104 INTEGER,
DIMENSION(3, 1) :: index_to_cell_zero
1106 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: mat_dm_global_work
1107 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_s_kp
1110 CALL timeset(routinen, handle)
1112 NULLIFY (matrix_s_kp)
1115 matrix_s_kp=matrix_s_kp)
1117 kpoints_g => qs_env%mp2_env%ri_rpa_im_time%kpoints_G
1120 index_to_cell_zero = 0
1122 NULLIFY (mat_dm_global_work)
1127 IF (.NOT.
ASSOCIATED(propagator))
THEN
1132 ALLOCATE (mat_dm_global_work(1)%matrix)
1133 CALL dbcsr_create(matrix=mat_dm_global_work(1)%matrix, &
1134 template=matrix_s_kp(1, 1)%matrix, &
1135 matrix_type=dbcsr_type_no_symmetry)
1139 CALL dbcsr_set(mat_dm_global_work(1)%matrix, 0.0_dp)
1142 CALL get_mo_set(kpoints_g%kp_env(1)%kpoint_env%mos(1, 1), nao=nao)
1143 CALL cp_cfm_create(cfm_density, kpoints_g%kp_env(1)%kpoint_env%mos(1, 1)%mo_coeff%matrix_struct, &
1146 CALL density_matrix_from_kp_to_mic(kpoints_g, cfm_density, mat_dm_global_work, qs_env, ispin, tau, e_fermi, sector)
1150 CALL dbcsr_copy(propagator(sector, jquad, 1)%matrix, &
1151 mat_dm_global_work(1)%matrix)
1155 CALL timestop(handle)
1170 SUBROUTINE density_matrix_from_kp_to_mic(kpoints_G, density, mat_dm_global_work, qs_env, ispin, tau, e_fermi, sector)
1174 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: mat_dm_global_work
1176 INTEGER,
INTENT(IN) :: ispin
1177 REAL(kind=
dp),
INTENT(IN) :: tau, e_fermi
1178 INTEGER,
INTENT(IN) :: sector
1180 CHARACTER(LEN=*),
PARAMETER :: routinen =
'density_matrix_from_kp_to_mic'
1182 COMPLEX(KIND=dp),
CONTIGUOUS,
DIMENSION(:, :), &
1184 INTEGER :: handle, iatom, iatom_old, ik, ik_global, &
1185 ik_old, irow, jatom, jatom_old, jcol, &
1186 nao, ncol_local, nrow_local, num_cells
1187 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: atom_from_ao_index
1188 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
1189 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell
1190 REAL(kind=
dp) :: contribution, weight_im, weight_re
1191 REAL(kind=
dp),
DIMENSION(3, 3) :: hmat
1192 REAL(kind=
dp),
DIMENSION(:),
POINTER :: wkp
1193 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: xkp
1199 CALL timeset(routinen, handle)
1204 index_to_cell => kpoints_g%index_to_cell
1205 num_cells =
SIZE(index_to_cell, 2)
1207 CALL cp_cfm_get_info(density, nrow_global=nao, nrow_local=nrow_local, ncol_local=ncol_local, &
1208 row_indices=row_indices, col_indices=col_indices)
1209 CALL cp_fm_create(fm_mat_work, density%matrix_struct, set_zero=.true.)
1211 ALLOCATE (atom_from_ao_index(nao))
1215 NULLIFY (cell, particle_set)
1216 CALL get_qs_env(qs_env, cell=cell, particle_set=particle_set)
1223 CALL dbcsr_set(mat_dm_global_work(1)%matrix, 0.0_dp)
1225 DO ik = 1,
SIZE(kpoints_g%kp_env)
1227 kp => kpoints_g%kp_env(ik)%kpoint_env
1228 ik_global = kp%nkpoint
1229 CALL build_kpoint_density_matrix_rpa(kpoints_g, ik, ispin, tau, e_fermi, sector, density)
1232 DO irow = 1, nrow_local
1233 DO jcol = 1, ncol_local
1235 iatom = atom_from_ao_index(row_indices(irow))
1236 jatom = atom_from_ao_index(col_indices(jcol))
1238 IF (ik_global /= ik_old .OR. iatom /= iatom_old .OR. jatom /= jatom_old)
THEN
1241 num_cells, iatom, jatom, xkp(1:3, ik_global), wkp(ik_global), &
1242 cell, index_to_cell, hmat, particle_set)
1250 contribution = real(cmplx(weight_re, -weight_im, kind=
dp)*cfm_data(irow, jcol), kind=
dp)
1252 fm_mat_work%local_data(irow, jcol) = fm_mat_work%local_data(irow, jcol) + contribution
1259 CALL copy_fm_to_dbcsr(fm_mat_work, mat_dm_global_work(1)%matrix, keep_sparsity=.false.)
1264 DEALLOCATE (atom_from_ao_index)
1266 CALL timestop(handle)
1268 END SUBROUTINE density_matrix_from_kp_to_mic
1281 SUBROUTINE density_matrix_from_kp_to_transl(kpoints, density, mat_dm_global_work, index_to_cell, &
1282 ispin, tau, e_fermi, sector)
1286 TYPE(
dbcsr_p_type),
DIMENSION(:),
INTENT(IN) :: mat_dm_global_work
1287 INTEGER,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell
1288 INTEGER,
INTENT(IN) :: ispin
1289 REAL(kind=
dp),
INTENT(IN) :: tau, e_fermi
1290 INTEGER,
INTENT(IN) :: sector
1292 CHARACTER(LEN=*),
PARAMETER :: routinen =
'density_matrix_from_kp_to_transl'
1294 INTEGER :: handle, icell, ik, ik_global, xcell, &
1296 REAL(kind=
dp) :: arg, coskl, sinkl
1297 REAL(kind=
dp),
DIMENSION(:),
POINTER :: wkp
1298 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: xkp
1300 TYPE(
dbcsr_type),
POINTER :: mat_work_im, mat_work_re
1303 CALL timeset(routinen, handle)
1307 NULLIFY (mat_work_re)
1310 template=mat_dm_global_work(1)%matrix, &
1311 matrix_type=dbcsr_type_no_symmetry)
1313 NULLIFY (mat_work_im)
1316 template=mat_dm_global_work(1)%matrix, &
1317 matrix_type=dbcsr_type_no_symmetry)
1324 cpassert(
SIZE(mat_dm_global_work) ==
SIZE(index_to_cell, 2))
1326 DO icell = 1,
SIZE(mat_dm_global_work)
1328 CALL dbcsr_set(mat_dm_global_work(icell)%matrix, 0.0_dp)
1332 DO ik = 1,
SIZE(kpoints%kp_env)
1334 kp => kpoints%kp_env(ik)%kpoint_env
1335 ik_global = kp%nkpoint
1336 CALL build_kpoint_density_matrix_rpa(kpoints, ik, ispin, tau, e_fermi, sector, density)
1342 DO icell = 1,
SIZE(mat_dm_global_work)
1344 xcell = index_to_cell(1, icell)
1345 ycell = index_to_cell(2, icell)
1346 zcell = index_to_cell(3, icell)
1348 arg = real(xcell,
dp)*xkp(1, ik_global) + real(ycell,
dp)*xkp(2, ik_global) + &
1349 REAL(zcell,
dp)*xkp(3, ik_global)
1350 coskl = wkp(ik_global)*cos(
twopi*arg)
1351 sinkl = wkp(ik_global)*sin(
twopi*arg)
1353 CALL dbcsr_add(mat_dm_global_work(icell)%matrix, mat_work_re, 1.0_dp, coskl)
1355 CALL dbcsr_add(mat_dm_global_work(icell)%matrix, mat_work_im, 1.0_dp, -sinkl)
1366 CALL timestop(handle)
1368 END SUBROUTINE density_matrix_from_kp_to_transl
1378 INTEGER,
DIMENSION(3),
INTENT(IN) :: cell_grid
1379 INTEGER,
DIMENSION(:, :, :),
POINTER :: cell_to_index
1380 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell
1381 TYPE(
cell_type),
INTENT(IN),
POINTER :: cell
1383 CHARACTER(LEN=*),
PARAMETER :: routinen =
'init_cell_index_rpa'
1385 INTEGER :: cell_counter, handle, i_cell, &
1386 index_min_dist, num_cells, xcell, &
1388 INTEGER,
DIMENSION(3) :: itm
1389 INTEGER,
DIMENSION(:, :),
POINTER :: index_to_cell_unsorted
1390 INTEGER,
DIMENSION(:, :, :),
POINTER :: cell_to_index_unsorted
1391 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: abs_cell_vectors
1392 REAL(kind=
dp),
DIMENSION(3) :: cell_vector
1393 REAL(kind=
dp),
DIMENSION(3, 3) :: hmat
1395 CALL timeset(routinen, handle)
1399 num_cells = cell_grid(1)*cell_grid(2)*cell_grid(3)
1400 itm(:) = cell_grid(:)/2
1404 cpassert(cell_grid(1) /= itm(1)*2)
1405 cpassert(cell_grid(2) /= itm(2)*2)
1406 cpassert(cell_grid(3) /= itm(3)*2)
1408 IF (
ASSOCIATED(cell_to_index))
DEALLOCATE (cell_to_index)
1409 IF (
ASSOCIATED(index_to_cell))
DEALLOCATE (index_to_cell)
1411 ALLOCATE (cell_to_index_unsorted(-itm(1):itm(1), -itm(2):itm(2), -itm(3):itm(3)))
1412 cell_to_index_unsorted(:, :, :) = 0
1414 ALLOCATE (index_to_cell_unsorted(3, num_cells))
1415 index_to_cell_unsorted(:, :) = 0
1417 ALLOCATE (cell_to_index(-itm(1):itm(1), -itm(2):itm(2), -itm(3):itm(3)))
1418 cell_to_index(:, :, :) = 0
1420 ALLOCATE (index_to_cell(3, num_cells))
1421 index_to_cell(:, :) = 0
1423 ALLOCATE (abs_cell_vectors(1:num_cells))
1427 DO xcell = -itm(1), itm(1)
1428 DO ycell = -itm(2), itm(2)
1429 DO zcell = -itm(3), itm(3)
1431 cell_counter = cell_counter + 1
1432 cell_to_index_unsorted(xcell, ycell, zcell) = cell_counter
1434 index_to_cell_unsorted(1, cell_counter) = xcell
1435 index_to_cell_unsorted(2, cell_counter) = ycell
1436 index_to_cell_unsorted(3, cell_counter) = zcell
1438 cell_vector(1:3) = matmul(hmat, real(index_to_cell_unsorted(1:3, cell_counter),
dp))
1440 abs_cell_vectors(cell_counter) = sqrt(cell_vector(1)**2 + cell_vector(2)**2 + cell_vector(3)**2)
1448 DO i_cell = 1, num_cells/2 + 1
1450 index_min_dist = minloc(abs_cell_vectors(1:num_cells/2 + 1), dim=1)
1452 xcell = index_to_cell_unsorted(1, index_min_dist)
1453 ycell = index_to_cell_unsorted(2, index_min_dist)
1454 zcell = index_to_cell_unsorted(3, index_min_dist)
1456 index_to_cell(1, i_cell) = xcell
1457 index_to_cell(2, i_cell) = ycell
1458 index_to_cell(3, i_cell) = zcell
1460 cell_to_index(xcell, ycell, zcell) = i_cell
1462 abs_cell_vectors(index_min_dist) = 10000000000.0_dp
1467 DO i_cell = num_cells/2 + 2, num_cells
1469 index_min_dist = minloc(abs_cell_vectors(1:num_cells), dim=1)
1471 xcell = index_to_cell_unsorted(1, index_min_dist)
1472 ycell = index_to_cell_unsorted(2, index_min_dist)
1473 zcell = index_to_cell_unsorted(3, index_min_dist)
1475 index_to_cell(1, i_cell) = xcell
1476 index_to_cell(2, i_cell) = ycell
1477 index_to_cell(3, i_cell) = zcell
1479 cell_to_index(xcell, ycell, zcell) = i_cell
1481 abs_cell_vectors(index_min_dist) = 10000000000.0_dp
1485 DEALLOCATE (index_to_cell_unsorted, cell_to_index_unsorted, abs_cell_vectors)
1487 CALL timestop(handle)
1502 SUBROUTINE get_diff_index_3c(i_cell_R, i_cell_S, i_cell_R_minus_S, index_to_cell_3c, &
1503 cell_to_index_3c, index_to_cell_dm, R_minus_S_needed, &
1504 do_kpoints_cubic_RPA)
1506 INTEGER,
INTENT(IN) :: i_cell_r, i_cell_s
1507 INTEGER,
INTENT(OUT) :: i_cell_r_minus_s
1508 INTEGER,
ALLOCATABLE,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell_3c
1509 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :), &
1510 INTENT(IN) :: cell_to_index_3c
1511 INTEGER,
DIMENSION(:, :),
INTENT(IN),
POINTER :: index_to_cell_dm
1512 LOGICAL,
INTENT(OUT) :: r_minus_s_needed
1513 LOGICAL,
INTENT(IN) :: do_kpoints_cubic_rpa
1515 CHARACTER(LEN=*),
PARAMETER :: routinen =
'get_diff_index_3c'
1517 INTEGER :: handle, x_cell_r, x_cell_r_minus_s, x_cell_s, y_cell_r, y_cell_r_minus_s, &
1518 y_cell_s, z_cell_r, z_cell_r_minus_s, z_cell_s
1520 CALL timeset(routinen, handle)
1522 IF (do_kpoints_cubic_rpa)
THEN
1524 x_cell_r = index_to_cell_3c(1, i_cell_r)
1525 y_cell_r = index_to_cell_3c(2, i_cell_r)
1526 z_cell_r = index_to_cell_3c(3, i_cell_r)
1528 x_cell_s = index_to_cell_dm(1, i_cell_s)
1529 y_cell_s = index_to_cell_dm(2, i_cell_s)
1530 z_cell_s = index_to_cell_dm(3, i_cell_s)
1532 x_cell_r_minus_s = x_cell_r - x_cell_s
1533 y_cell_r_minus_s = y_cell_r - y_cell_s
1534 z_cell_r_minus_s = z_cell_r - z_cell_s
1536 IF (x_cell_r_minus_s >= lbound(cell_to_index_3c, 1) .AND. &
1537 x_cell_r_minus_s <= ubound(cell_to_index_3c, 1) .AND. &
1538 y_cell_r_minus_s >= lbound(cell_to_index_3c, 2) .AND. &
1539 y_cell_r_minus_s <= ubound(cell_to_index_3c, 2) .AND. &
1540 z_cell_r_minus_s >= lbound(cell_to_index_3c, 3) .AND. &
1541 z_cell_r_minus_s <= ubound(cell_to_index_3c, 3))
THEN
1543 i_cell_r_minus_s = cell_to_index_3c(x_cell_r_minus_s, y_cell_r_minus_s, z_cell_r_minus_s)
1546 IF (i_cell_r_minus_s == 0)
THEN
1548 r_minus_s_needed = .false.
1549 i_cell_r_minus_s = 0
1553 r_minus_s_needed = .true.
1559 i_cell_r_minus_s = 0
1560 r_minus_s_needed = .false.
1566 r_minus_s_needed = .true.
1567 i_cell_r_minus_s = 1
1571 CALL timestop(handle)
1573 END SUBROUTINE get_diff_index_3c
1587 SUBROUTINE get_diff_diff_index_3c(i_cell_R, i_cell_S, i_cell_T, i_cell_R_minus_S_minus_T, &
1588 index_to_cell_3c, cell_to_index_3c, index_to_cell_dm, &
1589 R_minus_S_minus_T_needed, &
1590 do_kpoints_cubic_RPA)
1592 INTEGER,
INTENT(IN) :: i_cell_r, i_cell_s, i_cell_t
1593 INTEGER,
INTENT(OUT) :: i_cell_r_minus_s_minus_t
1594 INTEGER,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell_3c
1595 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :), &
1596 INTENT(IN) :: cell_to_index_3c
1597 INTEGER,
DIMENSION(:, :),
INTENT(IN) :: index_to_cell_dm
1598 LOGICAL,
INTENT(OUT) :: r_minus_s_minus_t_needed
1599 LOGICAL,
INTENT(IN) :: do_kpoints_cubic_rpa
1601 CHARACTER(LEN=*),
PARAMETER :: routinen =
'get_diff_diff_index_3c'
1603 INTEGER :: handle, x_cell_r, x_cell_r_minus_s_minus_t, x_cell_s, x_cell_t, y_cell_r, &
1604 y_cell_r_minus_s_minus_t, y_cell_s, y_cell_t, z_cell_r, z_cell_r_minus_s_minus_t, &
1607 CALL timeset(routinen, handle)
1609 IF (do_kpoints_cubic_rpa)
THEN
1611 x_cell_r = index_to_cell_3c(1, i_cell_r)
1612 y_cell_r = index_to_cell_3c(2, i_cell_r)
1613 z_cell_r = index_to_cell_3c(3, i_cell_r)
1615 x_cell_s = index_to_cell_dm(1, i_cell_s)
1616 y_cell_s = index_to_cell_dm(2, i_cell_s)
1617 z_cell_s = index_to_cell_dm(3, i_cell_s)
1619 x_cell_t = index_to_cell_dm(1, i_cell_t)
1620 y_cell_t = index_to_cell_dm(2, i_cell_t)
1621 z_cell_t = index_to_cell_dm(3, i_cell_t)
1623 x_cell_r_minus_s_minus_t = x_cell_r - x_cell_s - x_cell_t
1624 y_cell_r_minus_s_minus_t = y_cell_r - y_cell_s - y_cell_t
1625 z_cell_r_minus_s_minus_t = z_cell_r - z_cell_s - z_cell_t
1627 IF (x_cell_r_minus_s_minus_t >= lbound(cell_to_index_3c, 1) .AND. &
1628 x_cell_r_minus_s_minus_t <= ubound(cell_to_index_3c, 1) .AND. &
1629 y_cell_r_minus_s_minus_t >= lbound(cell_to_index_3c, 2) .AND. &
1630 y_cell_r_minus_s_minus_t <= ubound(cell_to_index_3c, 2) .AND. &
1631 z_cell_r_minus_s_minus_t >= lbound(cell_to_index_3c, 3) .AND. &
1632 z_cell_r_minus_s_minus_t <= ubound(cell_to_index_3c, 3))
THEN
1634 i_cell_r_minus_s_minus_t = cell_to_index_3c(x_cell_r_minus_s_minus_t, &
1635 y_cell_r_minus_s_minus_t, &
1636 z_cell_r_minus_s_minus_t)
1639 IF (i_cell_r_minus_s_minus_t == 0)
THEN
1641 r_minus_s_minus_t_needed = .false.
1645 r_minus_s_minus_t_needed = .true.
1651 i_cell_r_minus_s_minus_t = 0
1652 r_minus_s_minus_t_needed = .false.
1659 r_minus_s_minus_t_needed = .true.
1660 i_cell_r_minus_s_minus_t = 1
1664 CALL timestop(handle)
1666 END SUBROUTINE get_diff_diff_index_3c
Handles all functions related to the CELL.
subroutine, public get_cell(cell, alpha, beta, gamma, deth, orthorhombic, abc, periodic, h, h_inv, symmetry_id, tag)
Get informations about a simulation cell.
Basic linear algebra operations for complex full matrices.
subroutine, public cp_cfm_column_scale(matrix_a, scaling)
Scales columns of the full matrix by corresponding factors.
Represents a complex full matrix distributed on many processors.
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
subroutine, public cp_fm_to_cfm(msourcer, msourcei, mtarget)
Construct a complex full matrix by taking its real and imaginary parts from two separate real-value f...
subroutine, public cp_cfm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
Creates a new full matrix with the given structure.
subroutine, public cp_cfm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, matrix_struct, para_env)
Returns information about a full matrix.
subroutine, public cp_cfm_to_fm(msource, mtargetr, mtargeti)
Copy real and imaginary parts of a complex full matrix into separate real-value full matrices.
subroutine, public dbcsr_release_p(matrix)
...
subroutine, public dbcsr_scale(matrix, alpha_scalar)
...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_init_p(matrix)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_clear(matrix)
...
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
subroutine, public dbcsr_distribution_get(dist, row_dist, col_dist, nrows, ncols, has_threads, group, mynode, numnodes, nprows, npcols, myprow, mypcol, pgrid, subgroups_defined, prow_group, pcol_group)
...
subroutine, public dbcsr_reserve_all_blocks(matrix)
Reserves all blocks.
DBCSR operations in CP2K.
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
represent a full matrix distributed on many processors
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
This is the start of a dbt_api, all publically needed functions are exported here....
Types and set/get functions for HFX.
Defines the basic variable types.
integer, parameter, public int_8
integer, parameter, public dp
Types and basic routines needed for a kpoint calculation.
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym, inversion_symmetry_only, symmetry_backend, symmetry_reduction_method, gamma_centered, lattice_fft)
Retrieve information from a kpoint environment.
Machine interface based on Fortran 2003 and POSIX.
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
real(kind=dp) function, public m_walltime()
returns time from a real-time clock, protected against rolling early/easily
Definition of mathematical constants and functions.
real(kind=dp), parameter, public twopi
Interface to the message passing library MPI.
Types needed for MP2 calculations.
basic linear algebra operations for full matrixes
Define the data structure for the particle information.
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, mimic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_vhxc, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Definition and initialisation of the mo data type.
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Utility methods to build 3-center integral tensors of various types.
subroutine, public create_2c_tensor(t2c, dist_1, dist_2, pgrid, sizes_1, sizes_2, order, name)
...
Utility methods to build 3-center integral tensors of various types.
subroutine, public decompress_tensor(tensor, blk_indices, compressed, eps)
...
subroutine, public get_tensor_occupancy(tensor, nze, occ)
...
Utility routines for GW with imaginary time.
subroutine, public compute_weight_re_im(weight_re, weight_im, num_cells, iatom, jatom, xkp, wkp_w, cell, index_to_cell, hmat, particle_set)
...
subroutine, public get_atom_index_from_basis_function_index(qs_env, atom_from_basis_index, basis_size, basis_type, first_bf_from_atom)
...
Routines for low-scaling RPA/GW with imaginary time.
subroutine, public zero_mat_p_omega(mat_p_omega)
...
subroutine, public create_propagator_matrix_set(propagator, ntime, matrix_template, index_to_cell)
Creates the sector-, time-, and cell-resolved propagator matrix set.
subroutine, public compute_mat_dm_global(grid, nmo, cfm_mo_coeff, homo, propagator, matrix_s, ispin, eigenval, e_fermi, eps_filter, memory_info, unit_nr, jquad, do_kpoints_cubic_rpa, do_kpoints_from_gamma, qs_env, num_cells_dm, index_to_cell_dm, para_env)
...
integer, parameter, public propagator_sector_virtual
subroutine, public compute_periodic_dm(propagator, qs_env, ispin, num_integ_points, jquad, e_fermi, tau, sector)
...
subroutine, public init_cell_index_rpa(cell_grid, cell_to_index, index_to_cell, cell)
...
integer, parameter, public propagator_sector_occupied
subroutine, public compute_transl_dm(propagator, qs_env, ispin, num_integ_points, jquad, e_fermi, tau, eps_filter, num_cells_dm, index_to_cell_dm, sector)
...
subroutine, public compute_gamma_propagator(propagator, jquad, cfm_mo_coeff, homo, eigenval, nmo, eps_filter, e_fermi, tau, para_env)
Builds the Gamma-point occupied and virtual propagators.
subroutine, public compute_mat_p_omega(mat_p_omega, cfm_mo_coeff, homo, mat_p_global, matrix_s, ispin, t_3c_m, t_3c_o, t_3c_o_compressed, t_3c_o_ind, starts_array_mc, ends_array_mc, starts_array_mc_block, ends_array_mc_block, grid, e_fermi, eps_filter, alpha, eps_filter_im_time, eigenval, nmo, cut_memory, unit_nr, mp2_env, para_env, qs_env, do_kpoints_from_gamma, index_to_cell_3c, cell_to_index_3c, has_mat_p_blocks, do_ri_sos_laplace_mp2, dbcsr_time, dbcsr_nflop)
...
integer, parameter, public propagator_number_of_sectors
Definition and construction of time/frequency grids for correlation methods.
Type defining parameters related to the simulation cell.
Represent a complex full matrix.
Keeps information about a specific k-point.
Contains information about kpoints.
stores all the informations relevant to an mpi environment