37#include "./base/base_uses.f90"
43 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'mp2_types'
72 INTEGER,
DIMENSION(:),
ALLOCATABLE :: array
76 INTEGER,
DIMENSION(:, :),
ALLOCATABLE :: array
80 REAL(kind=
dp),
DIMENSION(:),
ALLOCATABLE :: array
84 REAL(kind=
dp),
DIMENSION(:, :),
ALLOCATABLE :: array
88 REAL(kind=
dp),
DIMENSION(:, :, :),
ALLOCATABLE :: array
92 INTEGER,
DIMENSION(:, :),
ALLOCATABLE :: index_table
96 INTEGER :: n_quadrature = -1, &
98 END TYPE mp2_laplace_type
101 LOGICAL :: big_send = .false.
105 REAL(kind=
dp) :: eps_grid = 0.0_dp, &
106 eps_filter = 0.0_dp, &
107 eps_pgf_orb_s = 0.0_dp
108 INTEGER :: print_level = 0
109 REAL(kind=
dp) :: cutoff = 0.0_dp, &
110 relative_cutoff = 0.0_dp
111 INTEGER :: size_lattice_sum = 0
115 INTEGER :: block_size = 0, &
116 number_integration_groups = 0
117 LOGICAL :: print_dgemm_info = .false.
121 INTEGER :: rpa_num_quad_points = -1, &
122 rpa_num_integ_groups = -1, &
124 TYPE(
hfx_type),
DIMENSION(:, :),
POINTER :: x_data => null()
126 xc_section_aux => null()
127 LOGICAL :: reuse_hfx = .false., &
128 minimax_quad = .false., &
129 do_ri_g0w0 = .false., &
132 print_dgemm_info = .false.
135 TYPE(
dbcsr_type),
DIMENSION(:),
POINTER :: mo_coeff_o => null(), &
138 exchange_block_size = -1
139 LOGICAL :: use_hfx_implementation = .false.
140 REAL(kind=
dp) :: ener_exchange = 0.0_dp, &
141 rse_corr_diag = 0.0_dp, &
146 INTEGER :: sigma_param = 0.0_dp
147 REAL(kind=
dp) :: e_sigma_corr = 0.0_dp
151 TYPE ri_rpa_im_time_type
152 INTEGER :: cut_memory = 0
153 LOGICAL :: memory_info = .false., &
154 make_chi_pos_definite = .false., &
155 make_overlap_mat_ao_pos_definite = .false., &
156 trunc_coulomb_ri_x = .false., &
157 keep_quad = .false., &
158 do_kpoints_from_gamma = .false., &
159 do_extrapolate_kpoints = .false.
160 REAL(kind=
dp) :: eps_filter = 0.0_dp, &
161 eps_filter_factor = 0.0_dp, &
162 eps_compress = 0.0_dp, &
163 exp_tailored_weights = 0.0_dp, &
164 regularization_ri = 0.0_dp, &
165 eps_eigval_s = 0.0_dp, &
166 eps_eigval_s_gamma = 0.0_dp, &
167 rel_cutoff_trunc_coulomb_ri_x = 0.0_dp
168 REAL(kind=
dp),
DIMENSION(:),
POINTER :: tau_tj => null(), &
172 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: weights_cos_tf_t_to_w => null(), &
173 weights_cos_tf_w_to_t => null()
174 REAL(kind=
dp),
DIMENSION(:),
ALLOCATABLE :: eigenval_gamma, &
176 INTEGER :: group_size_p = 0, &
180 INTEGER,
DIMENSION(:),
POINTER :: kp_grid => null()
181 INTEGER,
DIMENSION(3) :: kp_grid_extra = -1
182 LOGICAL :: do_im_time_kpoints = .false.
183 INTEGER :: min_bsize = 0, &
187 TYPE(
kpoint_type),
POINTER :: kpoints_g => null(), &
188 kpoints_sigma => null(), &
189 kpoints_sigma_no_xc => null()
190 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: starts_array_mc_ri, ends_array_mc_ri, &
191 starts_array_mc_block_ri, &
192 ends_array_mc_block_ri, &
193 starts_array_mc, ends_array_mc, &
194 starts_array_mc_block, &
197 END TYPE ri_rpa_im_time_type
200 INTEGER :: corr_mos_occ = 0, &
202 corr_mos_occ_beta = 0, &
203 corr_mos_virt_beta = 0, &
207 REAL(kind=
dp) :: omega_max_fit = 0.0_dp
209 REAL(kind=
dp) :: fermi_level_offset = 0.0_dp
210 INTEGER :: iter_evgw = 0, &
212 REAL(kind=
dp) :: eps_iter = 0.0_dp
213 LOGICAL :: do_hedin_shift = .false., &
214 do_ri_sigma_x = .false., &
215 do_periodic = .false., &
216 print_self_energy = .false.
217 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: vec_sigma_x_minus_vxc_gw
218 INTEGER,
DIMENSION(:),
POINTER :: kp_grid => null(), &
219 kp_grid_sigma => null()
220 INTEGER :: num_kp_grids = 0
221 REAL(kind=
dp) :: eps_kpoint = 0.0_dp
222 LOGICAL :: do_mo_coeff_gamma = .false., &
223 do_average_deg_levels = .false.
224 REAL(kind=
dp) :: eps_eigenval = 0.0_dp
225 LOGICAL :: do_extra_kpoints = .false., &
226 do_aux_bas_gw = .false.
227 REAL(kind=
dp) :: frac_aux_mos = 0.0_dp
228 INTEGER :: num_omega_points = 0
229 LOGICAL :: do_ic_model = .false., &
230 print_ic_values = .false.
231 REAL(kind=
dp) :: eps_dist = 0.0_dp
233 INTEGER :: print_exx = 0
234 LOGICAL :: do_gamma_only_sigma = .false.
235 LOGICAL :: update_xc_energy = .false., &
236 do_kpoints_sigma = .false., &
237 print_local_bandgap = .false.
238 INTEGER :: n_kp_in_kp_line = 0, &
240 nkp_self_energy = 0, &
241 nkp_self_energy_special_kp = 0, &
242 nkp_self_energy_monkh_pack = 0, &
244 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: xkp_special_kp
245 TYPE(
dbcsr_p_type),
DIMENSION(:),
ALLOCATABLE :: matrix_sigma_x_minus_vxc, &
247 REAL(kind=
dp) :: broadening_print_loc_bandgap = 0.0_dp, &
248 energy_window_print_loc_bandgap = 0.0_dp, &
249 ldos_thresh_print_loc_bandgap = 0.0_dp, &
250 energy_spacing_print_loc_bandgap = 0.0_dp, &
251 regularization_minimax = 0.0_dp, &
252 soc_energy_window = 0.0_dp
253 INTEGER,
DIMENSION(:),
POINTER :: stride_loc_bandgap => null()
256 REAL(kind=
dp) :: dos_upper = 0.0_dp, &
257 dos_lower = 0.0_dp, &
260 INTEGER :: max_level_self_energy = 0, &
261 min_level_self_energy = 0, &
264 END TYPE ri_g0w0_type
267 REAL(kind=
dp) :: di_rel = 0.0_dp, &
270 INTEGER :: max_num_iter = 0, &
272 INTEGER,
DIMENSION(:),
ALLOCATABLE :: ri_nset_per_l
273 END TYPE ri_basis_opt
279 REAL(kind=
dp),
DIMENSION(:, :),
ALLOCATABLE :: operator_half, &
283 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
ALLOCATABLE :: g_p_ia
284 TYPE(
dbcsr_p_type),
DIMENSION(:),
ALLOCATABLE :: mo_coeff_o, &
286 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:) :: p_mo, w_mo, l_jb
287 REAL(kind=
dp) :: cphf_eps_conv = 0.0_dp, &
288 scale_step_size = 0.0_dp
289 INTEGER :: cphf_max_num_iter = 0, &
290 z_solver_method = 0, &
292 LOGICAL :: enforce_decrease = .false., &
293 recalc_residual = .false., &
294 polak_ribiere = .false.
296 TYPE(
qs_force_type),
DIMENSION(:),
POINTER :: mp2_force => null()
297 REAL(kind=
dp),
DIMENSION(3, 3) :: mp2_virial = 0.0_dp
298 REAL(
dp) :: eps_canonical = 0.0_dp
299 LOGICAL :: free_hfx_buffer = .false.
300 INTEGER :: dot_blksize = 0
301 INTEGER :: max_parallel_comm = 0
305 INTEGER :: bse_spin_config = 0, &
306 bse_diag_method = 0, &
310 num_print_exc_descr = 0, &
311 screening_method = 0, &
312 num_add_start_z_space = 0, &
313 fac_max_z_space = 0, &
315 num_davidson_iter = 0, &
316 davidson_abort_cond = 0
317 REAL(kind=
dp) :: eps_res = 0.0_dp, &
318 eps_exc_en = 0.0_dp, &
320 screening_factor = 0.0_dp, &
321 bse_cutoff_occ = 0.0_dp, &
322 bse_cutoff_empty = 0.0_dp, &
323 z_space_energy_cutoff = 0.0_dp
324 LOGICAL :: do_bse = .false., &
325 bse_debug_print = .false., &
326 print_directional_exc_descr = .false., &
327 use_ks_energies = .false.
329 REAL(kind=
dp) :: bse_spectrum_freq_step_size = 0.0_dp, &
330 bse_spectrum_freq_start = 0.0_dp, &
331 bse_spectrum_freq_end = 0.0_dp
332 LOGICAL :: bse_print_spectrum = .false.
333 REAL(kind=
dp),
DIMENSION(:),
POINTER :: bse_eta_spectrum_list => null()
336 LOGICAL :: do_nto_analysis = .false., &
337 explicit_nto_list = .false.
338 REAL(kind=
dp) :: eps_nto_eigval = 0.0_dp, &
339 eps_nto_osc_str = 0.0_dp
340 INTEGER :: num_print_exc_ntos = 0
341 INTEGER,
DIMENSION(:),
POINTER :: bse_nto_state_list => null()
342 INTEGER,
DIMENSION(:),
ALLOCATABLE :: bse_nto_state_list_final
347 TYPE(mp2_laplace_type) :: ri_laplace = mp2_laplace_type()
348 TYPE(mp2_direct_type) :: direct_canonical = mp2_direct_type()
350 TYPE(mp2_gpw_type) ::
mp2_gpw = mp2_gpw_type()
351 TYPE(ri_mp2_type) :: ri_mp2 = ri_mp2_type()
352 TYPE(ri_rpa_type) :: ri_rpa = ri_rpa_type()
354#if defined(FTN_NO_DEFAULT_INIT)
355 TYPE(ri_rpa_im_time_type) :: ri_rpa_im_time = ri_rpa_im_time_type(eigenval_gamma=null(), &
357 starts_array_mc_ri=null(), &
358 ends_array_mc_ri=null(), &
359 starts_array_mc_block_ri=null(), &
360 ends_array_mc_block_ri=null(), &
361 starts_array_mc=null(), ends_array_mc=null(), &
362 starts_array_mc_block=null(), &
363 ends_array_mc_block=null())
364 TYPE(ri_g0w0_type) :: ri_g0w0 = ri_g0w0_type(vec_sigma_x_minus_vxc_gw=null(), &
365 xkp_special_kp=null(), &
366 matrix_sigma_x_minus_vxc=null(), &
368 TYPE(ri_basis_opt) :: ri_opt_param = ri_basis_opt(ri_nset_per_l=null())
369 TYPE(grad_util) :: ri_grad = grad_util(operator_half=null(), &
376 p_mo=null(), w_mo=null(), l_jb=null())
377 TYPE(bse_type) :: bse = bse_type(bse_nto_state_list_final=null())
379 TYPE(ri_rpa_im_time_type) :: ri_rpa_im_time = ri_rpa_im_time_type()
380 TYPE(ri_g0w0_type) :: ri_g0w0 = ri_g0w0_type()
381 TYPE(ri_basis_opt) :: ri_opt_param = ri_basis_opt()
382 TYPE(grad_util) :: ri_grad = grad_util()
383 TYPE(bse_type) :: bse = bse_type()
385 REAL(kind=
dp) :: mp2_memory = 0.0_dp, &
388 INTEGER :: mp2_num_proc = 0
389 INTEGER :: block_size_row = 0
390 INTEGER :: block_size_col = 0
391 LOGICAL :: calc_pq_cond_num = .false.
392 LOGICAL :: hf_fail = .false.
393 LOGICAL :: p_screen = .false.
394 LOGICAL :: not_last_hfx = .false.
395 LOGICAL :: do_im_time = .false.
398 INTEGER,
DIMENSION(:),
POINTER :: eri_blksize => null()
399 LOGICAL :: do_svd = .false.
400 REAL(kind=
dp) :: eps_svd = -1.0_dp
401 REAL(kind=
dp) :: eps_range = 0.0_dp
404 REAL(
dp) :: e_gap = 0.0_dp, &
406 LOGICAL :: ri_aux_auto_generated = .false.
410 REAL(kind=
dp),
DIMENSION(:),
ALLOCATABLE :: msg
411 INTEGER,
DIMENSION(:),
ALLOCATABLE :: sizes
412 INTEGER,
DIMENSION(:, :),
ALLOCATABLE :: indx
418 REAL(kind=
dp),
DIMENSION(:, :),
ALLOCATABLE :: msg
425 ALLOCATABLE :: elements
426 INTEGER :: n_element = 0
438 CHARACTER(LEN=*),
PARAMETER :: routinen =
'mp2_env_release'
442 CALL timeset(routinen, handle)
445 IF (.NOT. mp2_env%ri_rpa%reuse_hfx)
THEN
446 IF (
ASSOCIATED(mp2_env%ri_rpa%x_data))
CALL hfx_release(mp2_env%ri_rpa%x_data)
448 IF (
ASSOCIATED(mp2_env%ri_rpa%xc_section_aux))
CALL section_vals_release(mp2_env%ri_rpa%xc_section_aux)
449 IF (
ASSOCIATED(mp2_env%ri_rpa%xc_section_primary))
CALL section_vals_release(mp2_env%ri_rpa%xc_section_primary)
452 IF (
ASSOCIATED(mp2_env%eri_mme_param))
DEALLOCATE (mp2_env%eri_mme_param)
453 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%tau_tj))
DEALLOCATE (mp2_env%ri_rpa_im_time%tau_tj)
454 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%tau_wj))
DEALLOCATE (mp2_env%ri_rpa_im_time%tau_wj)
455 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%tj))
DEALLOCATE (mp2_env%ri_rpa_im_time%tj)
456 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%wj))
DEALLOCATE (mp2_env%ri_rpa_im_time%wj)
457 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%weights_cos_tf_t_to_w))
DEALLOCATE (mp2_env%ri_rpa_im_time%weights_cos_tf_t_to_w)
458 IF (
ASSOCIATED(mp2_env%ri_rpa_im_time%weights_cos_tf_w_to_t))
DEALLOCATE (mp2_env%ri_rpa_im_time%weights_cos_tf_w_to_t)
460 CALL mp2_env%local_gemm_ctx%destroy()
462 CALL timestop(handle)
473 CHARACTER(LEN=*),
PARAMETER :: routinen =
'mp2_env_create'
477 CALL timeset(routinen, handle)
479 cpassert(.NOT.
ASSOCIATED(mp2_env))
483 NULLIFY (mp2_env%ri_rpa%x_data)
485 CALL timestop(handle)
Interface to Minimax-Ewald method for periodic ERI's to be used in CP2K.
subroutine, public cp_eri_mme_finalize(param)
Release eri mme data. Prints some statistics on summation methods chosen.
represent a full matrix distributed on many processors
Types and set/get functions for HFX.
subroutine, public hfx_release(x_data)
This routine deallocates all data structures
Defines the basic variable types.
integer, parameter, public dp
Types and basic routines needed for a kpoint calculation.
2- and 3-center electron repulsion integral routines based on libint2 Currently available operators: ...
Interface to the message passing library MPI.
Calls routines to get RI integrals and calculate total energies.
Types needed for MP2 calculations.
subroutine, public mp2_env_create(mp2_env)
...
integer, save, public init_tshpsc_lmax
subroutine, public mp2_env_release(mp2_env)
...
basis types for the calculation of the perturbation of density theory.
stores some data used in construction of Kohn-Sham matrix
Contains information about kpoints.
Represent a qs system that is perturbed. Can calculate the linear operator and the rhs of the system ...