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 Main routines for GW + Bethe-Salpeter for computing electronic excitations
10 : !> \par History
11 : !> 04.2024 created [Maximilian Graml]
12 : ! **************************************************************************************************
13 :
14 : MODULE bse_main
15 :
16 : USE bse_full_diag, ONLY: create_A,&
17 : create_B,&
18 : create_hermitian_form_of_ABBA,&
19 : diagonalize_A,&
20 : diagonalize_C
21 : USE bse_iterative, ONLY: do_subspace_iterations,&
22 : fill_local_3c_arrays
23 : USE bse_util, ONLY: deallocate_matrices_bse,&
24 : estimate_BSE_resources,&
25 : mult_B_with_W,&
26 : print_BSE_start_flag,&
27 : truncate_BSE_matrices
28 : USE cp_fm_types, ONLY: cp_fm_release,&
29 : cp_fm_type
30 : USE cp_log_handling, ONLY: cp_get_default_logger,&
31 : cp_logger_type
32 : USE cp_output_handling, ONLY: debug_print_level
33 : USE group_dist_types, ONLY: group_dist_d1_type
34 : USE input_constants, ONLY: bse_abba,&
35 : bse_both,&
36 : bse_fulldiag,&
37 : bse_iterdiag,&
38 : bse_tda
39 : USE kinds, ONLY: dp
40 : USE message_passing, ONLY: mp_para_env_type
41 : USE mp2_types, ONLY: mp2_type
42 : USE qs_environment_types, ONLY: qs_environment_type
43 : #include "./base/base_uses.f90"
44 :
45 : IMPLICIT NONE
46 :
47 : PRIVATE
48 :
49 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'bse_main'
50 :
51 : PUBLIC :: start_bse_calculation
52 :
53 : CONTAINS
54 :
55 : ! **************************************************************************************************
56 : !> \brief Main subroutine managing BSE calculations
57 : !> \param fm_mat_S_ia_bse ...
58 : !> \param fm_mat_S_ij_bse ...
59 : !> \param fm_mat_S_ab_bse ...
60 : !> \param fm_mat_Q_static_bse_gemm ...
61 : !> \param Eigenval ...
62 : !> \param Eigenval_scf ...
63 : !> \param homo ...
64 : !> \param virtual ...
65 : !> \param dimen_RI ...
66 : !> \param dimen_RI_red ...
67 : !> \param bse_lev_virt ...
68 : !> \param gd_array ...
69 : !> \param color_sub ...
70 : !> \param mp2_env ...
71 : !> \param qs_env ...
72 : !> \param mo_coeff ...
73 : !> \param unit_nr ...
74 : ! **************************************************************************************************
75 18 : SUBROUTINE start_bse_calculation(fm_mat_S_ia_bse, fm_mat_S_ij_bse, fm_mat_S_ab_bse, &
76 : fm_mat_Q_static_bse_gemm, &
77 : Eigenval, Eigenval_scf, &
78 18 : homo, virtual, dimen_RI, dimen_RI_red, bse_lev_virt, &
79 18 : gd_array, color_sub, mp2_env, qs_env, mo_coeff, unit_nr)
80 :
81 : TYPE(cp_fm_type), INTENT(IN) :: fm_mat_S_ia_bse, fm_mat_S_ij_bse, &
82 : fm_mat_S_ab_bse
83 : TYPE(cp_fm_type), INTENT(INOUT) :: fm_mat_Q_static_bse_gemm
84 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :), &
85 : INTENT(IN) :: Eigenval, Eigenval_scf
86 : INTEGER, DIMENSION(:), INTENT(IN) :: homo, virtual
87 : INTEGER, INTENT(IN) :: dimen_RI, dimen_RI_red, bse_lev_virt
88 : TYPE(group_dist_d1_type), INTENT(IN) :: gd_array
89 : INTEGER, INTENT(IN) :: color_sub
90 : TYPE(mp2_type) :: mp2_env
91 : TYPE(qs_environment_type), POINTER :: qs_env
92 : TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: mo_coeff
93 : INTEGER, INTENT(IN) :: unit_nr
94 :
95 : CHARACTER(LEN=*), PARAMETER :: routineN = 'start_bse_calculation'
96 :
97 : INTEGER :: handle, homo_red, virtual_red
98 : LOGICAL :: my_do_abba, my_do_fulldiag, &
99 : my_do_iterat_diag, my_do_tda
100 : REAL(KIND=dp) :: diag_runtime_est
101 18 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: Eigenval_reduced
102 18 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: B_abQ_bse_local, B_bar_iaQ_bse_local, &
103 18 : B_bar_ijQ_bse_local, B_iaQ_bse_local
104 : TYPE(cp_fm_type) :: fm_A_BSE, fm_B_BSE, fm_C_BSE, fm_inv_sqrt_A_minus_B, fm_mat_S_ab_trunc, &
105 : fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, fm_mat_S_ia_trunc, fm_mat_S_ij_trunc, &
106 : fm_sqrt_A_minus_B
107 : TYPE(cp_logger_type), POINTER :: logger
108 : TYPE(mp_para_env_type), POINTER :: para_env
109 :
110 18 : CALL timeset(routineN, handle)
111 :
112 18 : para_env => fm_mat_S_ia_bse%matrix_struct%para_env
113 :
114 18 : my_do_fulldiag = .FALSE.
115 18 : my_do_iterat_diag = .FALSE.
116 18 : my_do_tda = .FALSE.
117 18 : my_do_abba = .FALSE.
118 : !Method: Iterative or full diagonalization
119 18 : SELECT CASE (mp2_env%ri_g0w0%bse_diag_method)
120 : CASE (bse_iterdiag)
121 0 : my_do_iterat_diag = .TRUE.
122 : !MG: Basics of the Davidson solver are implemented, but not rigorously checked.
123 0 : CPABORT("Iterative BSE not yet implemented")
124 : CASE (bse_fulldiag)
125 18 : my_do_fulldiag = .TRUE.
126 : END SELECT
127 : !Approximation: TDA and/or full ABBA matrix
128 30 : SELECT CASE (mp2_env%ri_g0w0%flag_tda)
129 : CASE (bse_tda)
130 12 : my_do_tda = .TRUE.
131 : CASE (bse_abba)
132 6 : my_do_abba = .TRUE.
133 : CASE (bse_both)
134 0 : my_do_tda = .TRUE.
135 18 : my_do_abba = .TRUE.
136 : END SELECT
137 :
138 18 : CALL print_BSE_start_flag(my_do_tda, my_do_abba, unit_nr)
139 :
140 : ! Link BSE debug flag against debug print level
141 18 : logger => cp_get_default_logger()
142 18 : IF (logger%iter_info%print_level == debug_print_level) THEN
143 0 : mp2_env%ri_g0w0%bse_debug_print = .TRUE.
144 : END IF
145 :
146 18 : CALL fm_mat_S_ia_bse%matrix_struct%para_env%sync()
147 : ! We apply the BSE cutoffs using the DFT Eigenenergies
148 : ! Reduce matrices in case of energy cutoff for occupied and unoccupied in A/B-BSE-matrices
149 : CALL truncate_BSE_matrices(fm_mat_S_ia_bse, fm_mat_S_ij_bse, fm_mat_S_ab_bse, &
150 : fm_mat_S_ia_trunc, fm_mat_S_ij_trunc, fm_mat_S_ab_trunc, &
151 : Eigenval_scf(:, 1, 1), Eigenval(:, 1, 1), Eigenval_reduced, &
152 : homo(1), virtual(1), dimen_RI, unit_nr, &
153 : bse_lev_virt, &
154 : homo_red, virtual_red, &
155 18 : mp2_env)
156 : ! \bar{B}^P_rs = \sum_R W_PR B^R_rs where B^R_rs = \sum_T [1/sqrt(v)]_RT (T|rs)
157 : ! r,s: MO-index, P,R,T: RI-index
158 : ! B: fm_mat_S_..., W: fm_mat_Q_...
159 : CALL mult_B_with_W(fm_mat_S_ij_trunc, fm_mat_S_ia_trunc, fm_mat_S_bar_ia_bse, &
160 : fm_mat_S_bar_ij_bse, fm_mat_Q_static_bse_gemm, &
161 18 : dimen_RI_red, homo_red, virtual_red)
162 :
163 18 : IF (my_do_iterat_diag) THEN
164 : CALL fill_local_3c_arrays(fm_mat_S_ab_trunc, fm_mat_S_ia_trunc, &
165 : fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, &
166 : B_bar_ijQ_bse_local, B_abQ_bse_local, B_bar_iaQ_bse_local, &
167 : B_iaQ_bse_local, dimen_RI_red, homo_red, virtual_red, &
168 : gd_array, color_sub, para_env)
169 : END IF
170 :
171 18 : IF (my_do_fulldiag) THEN
172 : ! Quick estimate of memory consumption and runtime of diagonalizations
173 : CALL estimate_BSE_resources(homo_red, virtual_red, unit_nr, my_do_abba, &
174 18 : para_env, diag_runtime_est)
175 : ! Matrix A constructed from GW energies and 3c-B-matrices (cf. subroutine mult_B_with_W)
176 : ! A_ia,jb = (ε_a-ε_i) δ_ij δ_ab + α * v_ia,jb - W_ij,ab
177 : ! ε_a, ε_i are GW singleparticle energies from Eigenval_reduced
178 : ! α is a spin-dependent factor
179 : ! v_ia,jb = \sum_P B^P_ia B^P_jb (unscreened Coulomb interaction)
180 : ! W_ij,ab = \sum_P \bar{B}^P_ij B^P_ab (screened Coulomb interaction)
181 : CALL create_A(fm_mat_S_ia_trunc, fm_mat_S_bar_ij_bse, fm_mat_S_ab_trunc, &
182 : fm_A_BSE, Eigenval_reduced, unit_nr, &
183 : homo_red, virtual_red, dimen_RI, mp2_env, &
184 18 : para_env)
185 18 : IF (my_do_abba) THEN
186 : ! Matrix B constructed from 3c-B-matrices (cf. subroutine mult_B_with_W)
187 : ! B_ia,jb = α * v_ia,jb - W_ib,aj
188 : ! α is a spin-dependent factor
189 : ! v_ia,jb = \sum_P B^P_ia B^P_jb (unscreened Coulomb interaction)
190 : ! W_ib,aj = \sum_P \bar{B}^P_ib B^P_aj (screened Coulomb interaction)
191 : CALL create_B(fm_mat_S_ia_trunc, fm_mat_S_bar_ia_bse, fm_B_BSE, &
192 6 : homo_red, virtual_red, dimen_RI, unit_nr, mp2_env)
193 : ! Construct Matrix C=(A-B)^0.5 (A+B) (A-B)^0.5 to solve full BSE matrix as a hermitian problem
194 : ! (cf. Eq. (A7) in F. Furche J. Chem. Phys., Vol. 114, No. 14, (2001)).
195 : ! We keep fm_sqrt_A_minus_B and fm_inv_sqrt_A_minus_B for print of singleparticle transitions
196 : ! of ABBA as described in Eq. (A10) in F. Furche J. Chem. Phys., Vol. 114, No. 14, (2001).
197 : CALL create_hermitian_form_of_ABBA(fm_A_BSE, fm_B_BSE, fm_C_BSE, &
198 : fm_sqrt_A_minus_B, fm_inv_sqrt_A_minus_B, &
199 6 : homo_red, virtual_red, unit_nr, mp2_env, diag_runtime_est)
200 : END IF
201 18 : CALL cp_fm_release(fm_B_BSE)
202 18 : IF (my_do_tda) THEN
203 : ! Solving the hermitian eigenvalue equation A X^n = Ω^n X^n
204 : CALL diagonalize_A(fm_A_BSE, homo_red, virtual_red, homo(1), &
205 12 : unit_nr, diag_runtime_est, mp2_env, qs_env, mo_coeff)
206 : END IF
207 : ! Release to avoid faulty use of changed A matrix
208 18 : CALL cp_fm_release(fm_A_BSE)
209 18 : IF (my_do_abba) THEN
210 : ! Solving eigenvalue equation C Z^n = (Ω^n)^2 Z^n .
211 : ! Here, the eigenvectors Z^n relate to X^n via
212 : ! Eq. (A10) in F. Furche J. Chem. Phys., Vol. 114, No. 14, (2001).
213 : CALL diagonalize_C(fm_C_BSE, homo_red, virtual_red, homo(1), &
214 : fm_sqrt_A_minus_B, fm_inv_sqrt_A_minus_B, &
215 6 : unit_nr, diag_runtime_est, mp2_env, qs_env, mo_coeff)
216 : END IF
217 : ! Release to avoid faulty use of changed C matrix
218 18 : CALL cp_fm_release(fm_C_BSE)
219 : END IF
220 :
221 : CALL deallocate_matrices_bse(fm_mat_S_bar_ia_bse, fm_mat_S_bar_ij_bse, &
222 : fm_mat_S_ia_trunc, fm_mat_S_ij_trunc, fm_mat_S_ab_trunc, &
223 18 : fm_mat_Q_static_bse_gemm)
224 18 : DEALLOCATE (Eigenval_reduced)
225 18 : IF (my_do_iterat_diag) THEN
226 : ! Contains untested Block-Davidson algorithm
227 : CALL do_subspace_iterations(B_bar_ijQ_bse_local, B_abQ_bse_local, B_bar_iaQ_bse_local, &
228 : B_iaQ_bse_local, homo(1), virtual(1), mp2_env%ri_g0w0%bse_spin_config, unit_nr, &
229 0 : Eigenval(:, 1, 1), para_env, mp2_env)
230 : ! Deallocate local 3c-B-matrices
231 0 : DEALLOCATE (B_bar_ijQ_bse_local, B_abQ_bse_local, B_bar_iaQ_bse_local, B_iaQ_bse_local)
232 : END IF
233 :
234 18 : IF (unit_nr > 0) THEN
235 9 : WRITE (unit_nr, '(T2,A4,T7,A53)') 'BSE|', 'The BSE was successfully calculated. Have a nice day!'
236 : END IF
237 :
238 18 : CALL timestop(handle)
239 :
240 36 : END SUBROUTINE start_bse_calculation
241 :
242 : END MODULE bse_main
|
