From 866bbfbfee7b3ba863eb6972a6b9c0733f7302e8 Mon Sep 17 00:00:00 2001 From: Juerg Hutter Date: Fri, 15 Jul 2022 17:34:29 +0200 Subject: [PATCH] TDDFT/Linear Response : Add GAPW/GAPW_XC and ADMM/GAPW options (#2200) * GAPW_XC bug fix for linear response * GAPW/GAPW_XC linear response refactoring * Prettify and update regtests * Remove filter and adjust regtest value * Pretty * Linear response for ADMM/GAPW * TDDFT/ADMM/GAPW excitation energies. More regtests. Some refactoring. * Merge with upstream * Fix memory leak and non-assocoated pointer argument --- src/admm_methods.F | 8 +- src/qs_collocate_density.F | 10 -- src/qs_gapw_densities.F | 9 +- src/qs_kpp1_env_methods.F | 24 ++- src/qs_linres_kernel.F | 74 +++++++-- src/qs_linres_methods.F | 6 + src/qs_p_env_methods.F | 76 +++++++--- src/qs_p_env_types.F | 36 +++-- src/qs_tddfpt2_densities.F | 45 +++++- src/qs_tddfpt2_fhxc.F | 49 +++++- src/qs_tddfpt2_methods.F | 1 + src/qs_tddfpt2_operators.F | 3 + src/qs_tddfpt2_subgroups.F | 33 +++- src/qs_tddfpt2_types.F | 14 ++ src/qs_vxc_atom.F | 37 ++--- tests/QS/regtest-debug-1/TEST_FILES | 4 +- tests/QS/regtest-debug-1/h2o_admm_gapw.inp | 93 ++++++++++++ tests/QS/regtest-debug-1/h2o_gapw_xc.inp | 83 +++++++++++ tests/QS/regtest-tddfpt-4/TEST_FILES | 45 +++--- tests/QS/regtest-tddfpt-4/test09.inp | 7 +- tests/QS/regtest-tddfpt-4/test10.inp | 5 +- tests/QS/regtest-tddfpt-4/test13.inp | 13 +- tests/QS/regtest-tddfpt-4/test14.inp | 13 +- tests/QS/regtest-tddfpt-4/test23.inp | 166 +++++++++++++++++++++ tests/QS/regtest-tddfpt/H2O_GAPW_3.inp | 3 +- 25 files changed, 704 insertions(+), 153 deletions(-) create mode 100644 tests/QS/regtest-debug-1/h2o_admm_gapw.inp create mode 100644 tests/QS/regtest-debug-1/h2o_gapw_xc.inp create mode 100644 tests/QS/regtest-tddfpt-4/test23.inp diff --git a/src/admm_methods.F b/src/admm_methods.F index 38007eb55f..85eaf1f000 100644 --- a/src/admm_methods.F +++ b/src/admm_methods.F @@ -101,12 +101,12 @@ MODULE admm_methods PUBLIC :: admm_mo_calc_rho_aux, & admm_mo_merge_ks_matrix, & admm_mo_merge_derivs, & + admm_aux_response_density, & calc_mixed_overlap_force, & calc_aux_mo_derivs_none, & scale_dm, & admm_fit_mo_coeffs, & admm_update_ks_atom, & - admm_aux_reponse_density, & calc_admm_mo_derivatives, & calc_admm_ovlp_forces, & admm_projection_derivative @@ -2496,12 +2496,12 @@ CONTAINS !> \param dm ... !> \param dm_admm ... ! ************************************************************************************************** - SUBROUTINE admm_aux_reponse_density(qs_env, dm, dm_admm) + SUBROUTINE admm_aux_response_density(qs_env, dm, dm_admm) TYPE(qs_environment_type), INTENT(IN), POINTER :: qs_env TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN) :: dm TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT) :: dm_admm - CHARACTER(LEN=*), PARAMETER :: routineN = 'admm_aux_reponse_density' + CHARACTER(LEN=*), PARAMETER :: routineN = 'admm_aux_response_density' INTEGER :: handle, ispin, nao, nao_aux, ncol, nspins TYPE(admm_type), POINTER :: admm_env @@ -2532,6 +2532,6 @@ CONTAINS CALL timestop(handle) - END SUBROUTINE admm_aux_reponse_density + END SUBROUTINE admm_aux_response_density END MODULE admm_methods diff --git a/src/qs_collocate_density.F b/src/qs_collocate_density.F index 31051f7077..f33cab33ff 100644 --- a/src/qs_collocate_density.F +++ b/src/qs_collocate_density.F @@ -1617,15 +1617,6 @@ CONTAINS ! Figure out which task_list to use. my_soft_valid = .FALSE. IF (PRESENT(soft_valid)) my_soft_valid = soft_valid -! IF (my_soft_valid) THEN -! !TODO: Why does soft_valid overrule task_list_external? -! CALL get_ks_env(ks_env, task_list_soft=task_list) -! ELSE IF (PRESENT(task_list_external)) THEN -! task_list => task_list_external -! ELSE -! CALL get_ks_env(ks_env, task_list=task_list) -! END IF -!deb IF (PRESENT(task_list_external)) THEN task_list => task_list_external ELSEIF (my_soft_valid) THEN @@ -1633,7 +1624,6 @@ CONTAINS ELSE CALL get_ks_env(ks_env, task_list=task_list) END IF -!deb CPASSERT(ASSOCIATED(task_list)) ! Figure out which pw_env to use. diff --git a/src/qs_gapw_densities.F b/src/qs_gapw_densities.F index 3040c6a23d..b067afb9c6 100644 --- a/src/qs_gapw_densities.F +++ b/src/qs_gapw_densities.F @@ -128,12 +128,13 @@ CONTAINS CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=natom) CALL get_qs_kind(my_kind_set(ikind), paw_atom=paw_atom) -! Calculate rho1_h and rho1_s on the radial grids centered on the atomic position - IF (paw_atom) & + !Calculate rho1_h and rho1_s on the radial grids centered on the atomic position + IF (paw_atom) THEN CALL calculate_rho_atom(para_env, rho_atom_set, my_kind_set(ikind), & atom_list, natom, nspins, rho1_h_tot, rho1_s_tot) + END IF -! Calculate rho0_h and rho0_s on the radial grids centered on the atomic position + !Calculate rho0_h and rho0_s on the radial grids centered on the atomic position IF (my_do_rho0) & CALL calculate_rho0_atom(gapw_control, rho_atom_set, rho0_atom_set, rho0_mpole, & atom_list, natom, ikind, my_kind_set(ikind), rho0_h_tot) @@ -151,7 +152,7 @@ CONTAINS IF (my_do_rho0) THEN rho0_mpole%total_rho0_h = -rho0_h_tot -! Put the rho0_soft on the global grid + !Put the rho0_soft on the global grid CALL put_rho0_on_grid(qs_env, rho0_mpole, tot_rs_int) IF (ABS(rho0_h_tot) .GE. 1.0E-5_dp) THEN IF (ABS(1.0_dp - ABS(tot_rs_int/rho0_h_tot)) .GT. 1.0E-3_dp) THEN diff --git a/src/qs_kpp1_env_methods.F b/src/qs_kpp1_env_methods.F index d8ed93a6c0..049267910b 100644 --- a/src/qs_kpp1_env_methods.F +++ b/src/qs_kpp1_env_methods.F @@ -83,6 +83,7 @@ MODULE qs_kpp1_env_methods USE qs_ks_methods, ONLY: qs_ks_build_kohn_sham_matrix USE qs_p_env_types, ONLY: qs_p_env_type USE qs_rho0_ggrid, ONLY: integrate_vhg0_rspace + USE qs_rho_atom_types, ONLY: rho_atom_type USE qs_rho_types, ONLY: qs_rho_get,& qs_rho_type USE qs_vxc_atom, ONLY: calculate_xc_2nd_deriv_atom @@ -233,6 +234,7 @@ CONTAINS TYPE(pw_poisson_type), POINTER :: poisson_env TYPE(pw_pool_type), POINTER :: auxbas_pw_pool TYPE(qs_rho_type), POINTER :: rho + TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set TYPE(section_vals_type), POINTER :: input, scf_section CALL timeset(routineN, handle) @@ -276,7 +278,7 @@ CONTAINS CALL qs_rho_get(rho, rho_ao=rho_ao) CALL qs_rho_get(rho1, rho_g=rho1_g) -! gets the tmp grids + ! gets the tmp grids CPASSERT(ASSOCIATED(pw_env)) CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, & poisson_env=poisson_env) @@ -288,8 +290,7 @@ CONTAINS IF (gapw .OR. gapw_xc) & CALL prepare_gapw_den(qs_env, p_env%local_rho_set, do_rho0=(.NOT. gapw_xc)) -! *** calculate the hartree potential on the total density *** - + ! *** calculate the hartree potential on the total density *** CALL pw_pool_create_pw(auxbas_pw_pool, rho1_tot_gspace%pw, & use_data=COMPLEXDATA1D, & in_space=RECIPROCALSPACE) @@ -325,7 +326,7 @@ CONTAINS CALL pw_pool_give_back_pw(auxbas_pw_pool, rho1_tot_gspace%pw) -! *** calculate the xc potential *** + ! *** calculate the xc potential *** IF (gapw_xc) THEN CALL qs_rho_get(rho1_xc, rho_r=rho1_r, tau_r=tau1_r) ELSE @@ -390,8 +391,12 @@ CONTAINS IF (SIZE(v_xc_tau) /= nspins) CALL pw_pool_give_back_pw(auxbas_pw_pool, v_xc_tau(2)%pw) END IF - IF (gapw) CALL calculate_xc_2nd_deriv_atom(p_env%local_rho_set, qs_env, xc_section, para_env, & - do_tddft=do_tddft, do_triplet=do_triplet) + IF (gapw .OR. gapw_xc) THEN + CALL get_qs_env(qs_env, rho_atom_set=rho_atom_set) + rho1_atom_set => p_env%local_rho_set%rho_atom_set + CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, & + do_tddft=do_tddft, do_triplet=do_triplet) + END IF DO ispin = 1, SIZE(rho1_r_pw) CALL pw_release(rho1_r_pw(ispin)%pw) @@ -541,6 +546,13 @@ CONTAINS psmat(1:ns, 1:1) => rho_ao(1:ns) CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.TRUE., & rho_atom_external=p_env%local_rho_set%rho_atom_set) + ELSEIF (gapw_xc) THEN + ns = SIZE(p_env%kpp1) + ksmat(1:ns, 1:1) => p_env%kpp1(1:ns) + ns = SIZE(rho_ao) + psmat(1:ns, 1:1) => rho_ao(1:ns) + CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.TRUE., & + rho_atom_external=p_env%local_rho_set%rho_atom_set) END IF CALL pw_pool_give_back_pw(auxbas_pw_pool, v_hartree_rspace%pw) diff --git a/src/qs_linres_kernel.F b/src/qs_linres_kernel.F index e248d730d8..4760825bc8 100644 --- a/src/qs_linres_kernel.F +++ b/src/qs_linres_kernel.F @@ -50,7 +50,8 @@ MODULE qs_linres_kernel section_vals_val_get USE kg_correction, ONLY: kg_ekin_subset USE kg_environment_types, ONLY: kg_environment_type - USE kinds, ONLY: dp + USE kinds, ONLY: default_string_length,& + dp USE lri_environment_types, ONLY: lri_density_type,& lri_environment_type,& lri_kind_type @@ -89,9 +90,11 @@ MODULE qs_linres_kernel USE qs_ks_atom, ONLY: update_ks_atom USE qs_ks_types, ONLY: qs_ks_env_type USE qs_linres_types, ONLY: linres_control_type + USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type USE qs_p_env_methods, ONLY: p_env_finish_kpp1 USE qs_p_env_types, ONLY: qs_p_env_type USE qs_rho0_ggrid, ONLY: integrate_vhg0_rspace + USE qs_rho_atom_types, ONLY: rho_atom_type USE qs_rho_types, ONLY: qs_rho_get,& qs_rho_type USE qs_vxc_atom, ONLY: calculate_xc_2nd_deriv_atom @@ -200,6 +203,7 @@ CONTAINS TYPE(qs_kpp1_env_type), POINTER :: kpp1_env TYPE(qs_ks_env_type), POINTER :: ks_env TYPE(qs_rho_type), POINTER :: rho, rho1, rho1_xc, rho1a, rho_aux + TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set TYPE(section_vals_type), POINTER :: input, xc_section, xc_section_aux CALL timeset(routineN, handle) @@ -315,11 +319,15 @@ CONTAINS IF (deriv2_analytic) THEN CALL qs_rho_get(rho1a, rho_r=rho1_r, tau_r=tau1_r) CALL qs_fxc_analytic(rho, rho1_r, tau1_r, xc_section, auxbas_pw_pool, lr_triplet, v_xc, v_xc_tau) - IF (gapw) CALL calculate_xc_2nd_deriv_atom(p_env%local_rho_set, qs_env, xc_section, para_env, & - do_tddft=.FALSE., do_triplet=lr_triplet) + IF (gapw .OR. gapw_xc) THEN + CALL get_qs_env(qs_env, rho_atom_set=rho_atom_set) + rho1_atom_set => p_env%local_rho_set%rho_atom_set + CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, & + do_tddft=.FALSE., do_triplet=lr_triplet) + END IF ELSE CALL qs_fxc_fdiff(ks_env, rho, rho1a, xc_section, 6, lr_triplet, v_xc, v_xc_tau) - CPASSERT(.NOT. gapw) + CPASSERT((.NOT. gapw) .AND. (.NOT. gapw_xc)) END IF DO ispin = 1, nspins @@ -357,7 +365,7 @@ CONTAINS IF (gapw_xc) THEN ! XC and Hartree are integrated separatedly - ! XC uses the sofft basis set only + ! XC uses the soft basis set only IF (nspins == 1) THEN @@ -480,7 +488,6 @@ CONTAINS END DO IF (gapw) THEN - IF (.NOT. ((nspins == 1 .AND. lr_triplet))) THEN CALL Vh_1c_gg_integrals(qs_env, energy_hartree_1c, & p_env%hartree_local%ecoul_1c, & @@ -491,17 +498,22 @@ CONTAINS calculate_forces=.FALSE., & local_rho_set=p_env%local_rho_set) END IF - ! *** Add single atom contributions to the KS matrix *** ! remap pointer ns = SIZE(p_env%kpp1) ksmat(1:ns, 1:1) => p_env%kpp1(1:ns) ns = SIZE(rho_ao) psmat(1:ns, 1:1) => rho_ao(1:ns) - CALL update_ks_atom(qs_env, ksmat, psmat, forces=.FALSE., tddft=.TRUE., & rho_atom_external=p_env%local_rho_set%rho_atom_set) + ELSEIF (gapw_xc) THEN + ns = SIZE(p_env%kpp1) + ksmat(1:ns, 1:1) => p_env%kpp1(1:ns) + ns = SIZE(rho_ao) + psmat(1:ns, 1:1) => rho_ao(1:ns) + CALL update_ks_atom(qs_env, ksmat, psmat, forces=.FALSE., tddft=.TRUE., & + rho_atom_external=p_env%local_rho_set%rho_atom_set) END IF ! KG embedding, contribution of kinetic energy functional to kernel @@ -813,21 +825,27 @@ CONTAINS CHARACTER(len=*), PARAMETER :: routineN = 'apply_xc_admm' - INTEGER :: handle, ispin, nspins + CHARACTER(LEN=default_string_length) :: basis_type + INTEGER :: handle, ispin, ns, nspins INTEGER, DIMENSION(2, 3) :: bo - LOGICAL :: lsd + LOGICAL :: gapw, lsd REAL(KIND=dp) :: alpha TYPE(admm_type), POINTER :: admm_env + TYPE(cp_para_env_type), POINTER :: para_env TYPE(dbcsr_p_type) :: xcmat TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s + TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ksmat, psmat TYPE(dft_control_type), POINTER :: dft_control TYPE(linres_control_type), POINTER :: linres_control + TYPE(neighbor_list_set_p_type), DIMENSION(:), & + POINTER :: sab_aux_fit TYPE(pw_env_type), POINTER :: pw_env TYPE(pw_p_type), DIMENSION(:), POINTER :: rho1_aux_g, rho1_aux_r, tau_pw, v_xc, & v_xc_tau TYPE(pw_pool_type), POINTER :: auxbas_pw_pool + TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set TYPE(section_vals_type), POINTER :: xc_fun_section, xc_section - TYPE(task_list_type), POINTER :: task_list_aux_fit + TYPE(task_list_type), POINTER :: task_list TYPE(xc_rho_cflags_type) :: needs TYPE(xc_rho_set_type) :: rho1_set @@ -840,8 +858,6 @@ CONTAINS ! nothing to do ELSE CALL get_qs_env(qs_env=qs_env, linres_control=linres_control) - CPASSERT(.NOT. dft_control%qs_control%gapw) - CPASSERT(.NOT. dft_control%qs_control%gapw_xc) CPASSERT(.NOT. dft_control%qs_control%lrigpw) CPASSERT(.NOT. linres_control%lr_triplet) @@ -859,7 +875,8 @@ CONTAINS CALL dbcsr_create(xcmat%matrix, template=matrix_s(1)%matrix) CALL get_qs_env(qs_env, admm_env=admm_env) - CALL get_admm_env(admm_env, task_list_aux_fit=task_list_aux_fit) + gapw = admm_env%do_gapw + CALL qs_rho_get(p_env%rho1_admm, rho_r=rho1_aux_r, rho_g=rho1_aux_g) xc_section => admm_env%xc_section_aux bo = rho1_aux_r(1)%pw%pw_grid%bounds_local @@ -885,6 +902,21 @@ CONTAINS END IF CALL xc_rho_set_release(rho1_set) + basis_type = "AUX_FIT" + CALL get_qs_env(qs_env, para_env=para_env) + CALL get_admm_env(admm_env, task_list_aux_fit=task_list) + IF (admm_env%do_gapw) THEN + CALL prepare_gapw_den(qs_env, local_rho_set=p_env%local_rho_set_admm, & + do_rho0=.FALSE., kind_set_external=admm_env%admm_gapw_env%admm_kind_set) + rho_atom_set => admm_env%admm_gapw_env%local_rho_set%rho_atom_set + rho1_atom_set => p_env%local_rho_set%rho_atom_set + CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, & + do_tddft=.FALSE., do_triplet=.FALSE., & + kind_set_external=admm_env%admm_gapw_env%admm_kind_set) + basis_type = "AUX_FIT_SOFT" + task_list => admm_env%admm_gapw_env%task_list + END IF + alpha = 1.0_dp IF (nspins == 1) alpha = 2.0_dp @@ -894,10 +926,22 @@ CONTAINS CALL dbcsr_set(xcmat%matrix, 0.0_dp) CALL integrate_v_rspace(v_rspace=v_xc(ispin), hmat=xcmat, qs_env=qs_env, & calculate_forces=.FALSE., basis_type="AUX_FIT", & - task_list_external=task_list_aux_fit) + task_list_external=task_list) CALL dbcsr_add(p_env%kpp1_admm(ispin)%matrix, xcmat%matrix, 1.0_dp, alpha) END DO + IF (admm_env%do_gapw) THEN + CALL get_admm_env(admm_env, sab_aux_fit=sab_aux_fit) + ns = SIZE(p_env%kpp1_admm) + ksmat(1:ns, 1:1) => p_env%kpp1_admm(1:ns) + psmat(1:ns, 1:1) => p_env%p1_admm(1:ns) + CALL update_ks_atom(qs_env, ksmat, psmat, forces=.FALSE., tddft=.TRUE., & + rho_atom_external=p_env%local_rho_set_admm%rho_atom_set, & + kind_set_external=admm_env%admm_gapw_env%admm_kind_set, & + oce_external=admm_env%admm_gapw_env%oce, & + sab_external=sab_aux_fit) + END IF + DO ispin = 1, nspins CALL pw_pool_give_back_pw(auxbas_pw_pool, v_xc(ispin)%pw) END DO diff --git a/src/qs_linres_methods.F b/src/qs_linres_methods.F index e719d66959..c4ae8d91cc 100644 --- a/src/qs_linres_methods.F +++ b/src/qs_linres_methods.F @@ -66,6 +66,7 @@ MODULE qs_linres_methods USE qs_2nd_kernel_ao, ONLY: build_dm_response USE qs_environment_types, ONLY: get_qs_env,& qs_environment_type + USE qs_gapw_densities, ONLY: prepare_gapw_den USE qs_linres_kernel, ONLY: apply_op_2 USE qs_linres_types, ONLY: linres_control_type USE qs_loc_methods, ONLY: qs_loc_driver @@ -615,6 +616,11 @@ CONTAINS CALL get_qs_env(qs_env, rho=rho) ! that could be called before CALL qs_rho_update_rho(rho, qs_env=qs_env) ! that could be called before + IF (dft_control%qs_control%gapw) THEN + CALL prepare_gapw_den(qs_env) + ELSEIF (dft_control%qs_control%gapw_xc) THEN + CALL prepare_gapw_den(qs_env, do_rho0=.FALSE.) + END IF DO ispin = 1, nspins CALL dbcsr_set(p_env%kpp1(ispin)%matrix, 0.0_dp) diff --git a/src/qs_p_env_methods.F b/src/qs_p_env_methods.F index c6f2713896..b21bb558b5 100644 --- a/src/qs_p_env_methods.F +++ b/src/qs_p_env_methods.F @@ -14,8 +14,9 @@ !> 22-08-2002, TCH, started development ! ************************************************************************************************** MODULE qs_p_env_methods - USE admm_methods, ONLY: admm_aux_reponse_density - USE admm_types, ONLY: admm_type,& + USE admm_methods, ONLY: admm_aux_response_density + USE admm_types, ONLY: admm_gapw_type,& + admm_type,& get_admm_env USE atomic_kind_types, ONLY: atomic_kind_type USE cp_blacs_env, ONLY: cp_blacs_env_type @@ -68,7 +69,8 @@ MODULE qs_p_env_methods USE input_section_types, ONLY: section_vals_get,& section_vals_get_subs_vals,& section_vals_type - USE kinds, ONLY: dp + USE kinds, ONLY: default_string_length,& + dp USE preconditioner_types, ONLY: init_preconditioner USE pw_env_types, ONLY: pw_env_type USE pw_types, ONLY: pw_p_type @@ -89,10 +91,12 @@ MODULE qs_p_env_methods USE qs_matrix_pools, ONLY: mpools_get USE qs_mo_types, ONLY: get_mo_set,& mo_set_p_type + USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type USE qs_p_env_types, ONLY: qs_p_env_type USE qs_rho0_ggrid, ONLY: rho0_s_grid_create USE qs_rho0_methods, ONLY: init_rho0 - USE qs_rho_atom_methods, ONLY: allocate_rho_atom_internals + USE qs_rho_atom_methods, ONLY: allocate_rho_atom_internals,& + calculate_rho_atom_coeff USE qs_rho_methods, ONLY: qs_rho_rebuild,& qs_rho_update_rho USE qs_rho_types, ONLY: qs_rho_create,& @@ -143,6 +147,8 @@ CONTAINS CHARACTER(len=*), PARAMETER :: routineN = 'p_env_create' INTEGER :: handle, n_ao, n_mo, n_spins, natom, spin + TYPE(admm_gapw_type), POINTER :: admm_gapw_env + TYPE(admm_type), POINTER :: admm_env TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set TYPE(cp_blacs_env_type), POINTER :: blacs_env TYPE(cp_fm_pool_p_type), DIMENSION(:), POINTER :: ao_mo_fm_pools, mo_mo_fm_pools @@ -200,6 +206,14 @@ CONTAINS CALL dbcsr_set(p_env%p1_admm(spin)%matrix, 0.0_dp) END DO END IF + CALL get_qs_env(qs_env, admm_env=admm_env) + IF (admm_env%do_gapw) THEN + CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set) + admm_gapw_env => admm_env%admm_gapw_env + CALL local_rho_set_create(p_env%local_rho_set_admm) + CALL allocate_rho_atom_internals(p_env%local_rho_set_admm%rho_atom_set, atomic_kind_set, & + admm_gapw_env%admm_kind_set, dft_control, para_env) + END IF END IF CALL mpools_get(qs_env%mpools, ao_mo_fm_pools=ao_mo_fm_pools, & @@ -247,9 +261,9 @@ CONTAINS name="p_env%psi0d") END IF - !----------------------! - ! GAPW initializations ! - !----------------------! + !------------------------------! + ! GAPW/GAPW_XC initializations ! + !------------------------------! IF (dft_control%qs_control%gapw) THEN CALL get_qs_env(qs_env, & atomic_kind_set=atomic_kind_set, & @@ -265,6 +279,13 @@ CONTAINS CALL rho0_s_grid_create(pw_env, p_env%local_rho_set%rho0_mpole) CALL hartree_local_create(p_env%hartree_local) CALL init_coulomb_local(p_env%hartree_local, natom) + ELSEIF (dft_control%qs_control%gapw_xc) THEN + CALL get_qs_env(qs_env, & + atomic_kind_set=atomic_kind_set, & + qs_kind_set=qs_kind_set) + CALL local_rho_set_create(p_env%local_rho_set) + CALL allocate_rho_atom_internals(p_env%local_rho_set%rho_atom_set, atomic_kind_set, & + qs_kind_set, dft_control, para_env) END IF !------------------------! @@ -390,19 +411,23 @@ CONTAINS CHARACTER(LEN=*), PARAMETER :: routineN = 'p_env_update_rho' + CHARACTER(LEN=default_string_length) :: basis_type INTEGER :: handle, ispin + TYPE(admm_type), POINTER :: admm_env + TYPE(cp_para_env_type), POINTER :: para_env TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho1_ao TYPE(dft_control_type), POINTER :: dft_control + TYPE(neighbor_list_set_p_type), DIMENSION(:), & + POINTER :: sab_aux_fit TYPE(pw_p_type), DIMENSION(:), POINTER :: rho_g_aux, rho_r_aux TYPE(qs_ks_env_type), POINTER :: ks_env - TYPE(qs_rho_type), POINTER :: rho TYPE(task_list_type), POINTER :: task_list CALL timeset(routineN, handle) - CALL get_qs_env(qs_env, rho=rho, dft_control=dft_control) + CALL get_qs_env(qs_env, dft_control=dft_control) - IF (dft_control%do_admm) CALL admm_aux_reponse_density(qs_env, p_env%p1, p_env%p1_admm) + IF (dft_control%do_admm) CALL admm_aux_response_density(qs_env, p_env%p1, p_env%p1_admm) CALL qs_rho_get(p_env%rho1, rho_ao=rho1_ao) DO ispin = 1, SIZE(rho1_ao) @@ -410,24 +435,25 @@ CONTAINS END DO CALL qs_rho_update_rho(rho_struct=p_env%rho1, & + rho_xc_external=p_env%rho1_xc, & local_rho_set=p_env%local_rho_set, & qs_env=qs_env) IF (dft_control%do_admm) THEN IF (dft_control%admm_control%aux_exch_func /= do_admm_aux_exch_func_none) THEN - NULLIFY (ks_env, rho1_ao, rho_g_aux, rho_r_aux, task_list) - CALL get_qs_env(qs_env, & - ks_env=ks_env, & - dft_control=dft_control) + CALL get_qs_env(qs_env, ks_env=ks_env, admm_env=admm_env) + basis_type = "AUX_FIT" CALL get_admm_env(qs_env%admm_env, task_list_aux_fit=task_list) - + IF (admm_env%do_gapw) THEN + basis_type = "AUX_FIT_SOFT" + task_list => admm_env%admm_gapw_env%task_list + END IF CALL qs_rho_get(p_env%rho1_admm, & rho_ao=rho1_ao, & rho_g=rho_g_aux, & rho_r=rho_r_aux) - DO ispin = 1, SIZE(rho1_ao) CALL dbcsr_copy(rho1_ao(ispin)%matrix, p_env%p1_admm(ispin)%matrix) CALL calculate_rho_elec(ks_env=ks_env, & @@ -435,9 +461,17 @@ CONTAINS rho=rho_r_aux(ispin), & rho_gspace=rho_g_aux(ispin), & soft_valid=.FALSE., & - basis_type="AUX_FIT", & + basis_type=basis_type, & task_list_external=task_list) END DO + IF (admm_env%do_gapw) THEN + CALL get_qs_env(qs_env, para_env=para_env) + CALL get_admm_env(admm_env, sab_aux_fit=sab_aux_fit) + CALL calculate_rho_atom_coeff(qs_env, rho1_ao, & + rho_atom_set=p_env%local_rho_set_admm%rho_atom_set, & + qs_kind_set=admm_env%admm_gapw_env%admm_kind_set, & + oce=admm_env%admm_gapw_env%oce, sab=sab_aux_fit, para_env=para_env) + END IF END IF END IF @@ -798,18 +832,22 @@ CONTAINS END IF END DO - IF (ASSOCIATED(p_env%local_rho_set)) THEN + IF (gapw) THEN CALL qs_rho_update_rho(rho_struct=p_env%rho1, qs_env=qs_env, local_rho_set=p_env%local_rho_set) + ELSEIF (gapw_xc) THEN + CALL qs_rho_update_rho(rho_struct=p_env%rho1, qs_env=qs_env, rho_xc_external=p_env%rho1_xc, & + local_rho_set=p_env%local_rho_set) ELSE CALL qs_rho_update_rho(rho_struct=p_env%rho1, qs_env=qs_env) END IF IF (fdiff) THEN + CPASSERT(.NOT. (gapw .OR. gapw_xc)) CALL kpp1_calc_k_p_p1_fdiff(qs_env=qs_env, & k_p_p1=p_env%kpp1, rho=rho, rho1=p_env%rho1) ELSE CALL kpp1_calc_k_p_p1(p_env=p_env, qs_env=qs_env, & - rho1=p_env%rho1, rho1_xc=p_env%rho1) + rho1=p_env%rho1, rho1_xc=p_env%rho1_xc) END IF DO ispin = 1, n_spins diff --git a/src/qs_p_env_types.F b/src/qs_p_env_types.F index 51502eacad..dd241d2b65 100644 --- a/src/qs_p_env_types.F +++ b/src/qs_p_env_types.F @@ -58,28 +58,33 @@ MODULE qs_p_env_types !> for the moment no smearing of the orbitals. ! ************************************************************************************************** TYPE qs_p_env_type - - LOGICAL :: orthogonal_orbitals - TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: kpp1 => NULL(), kpp1_admm => NULL(), p1 => NULL(), w1 => NULL() - TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: p1_admm => NULL() - TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: m_epsilon => NULL(), & - psi0d => NULL(), S_psi0 => NULL(), Smo_inv => NULL() - TYPE(qs_kpp1_env_type), POINTER :: kpp1_env => NULL() + LOGICAL :: orthogonal_orbitals + TYPE(qs_kpp1_env_type), POINTER :: kpp1_env => NULL() + TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: kpp1 => NULL() + TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: kpp1_admm => NULL() + TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: p1 => NULL() + TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: p1_admm => NULL() + TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: w1 => NULL() + TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: m_epsilon => NULL() + TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: psi0d => NULL() + TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: S_psi0 => NULL() + TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: Smo_inv => NULL() TYPE(qs_rho_type), POINTER :: rho1 => NULL() TYPE(qs_rho_type), POINTER :: rho1_xc => NULL() TYPE(qs_rho_type), POINTER :: rho1_admm => NULL() - INTEGER, DIMENSION(2) :: n_mo, & ! no of molecular orbitals - n_ao ! no of basis functions + INTEGER, DIMENSION(2) :: n_mo, & ! no of molecular orbitals + n_ao ! no of basis functions ! GAPW stuff - TYPE(hartree_local_type), POINTER :: hartree_local => NULL() - TYPE(local_rho_type), POINTER :: local_rho_set => NULL() + TYPE(hartree_local_type), POINTER :: hartree_local => NULL() + TYPE(local_rho_type), POINTER :: local_rho_set => NULL() + TYPE(local_rho_type), POINTER :: local_rho_set_admm => NULL() ! Linear Response Modules - TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: PS_psi0 => NULL() + TYPE(cp_fm_p_type), DIMENSION(:), POINTER :: PS_psi0 => NULL() ! preconditioner matrix should be symmetric and positive definite - TYPE(preconditioner_type), DIMENSION(:), POINTER :: preconditioner => NULL() - LOGICAL :: new_preconditioner + LOGICAL :: new_preconditioner + TYPE(preconditioner_type), DIMENSION(:), POINTER :: preconditioner => NULL() END TYPE qs_p_env_type @@ -123,6 +128,9 @@ CONTAINS IF (ASSOCIATED(p_env%hartree_local)) THEN CALL hartree_local_release(p_env%hartree_local) END IF + IF (ASSOCIATED(p_env%local_rho_set_admm)) THEN + CALL local_rho_set_release(p_env%local_rho_set_admm) + END IF IF (ASSOCIATED(p_env%PS_psi0)) THEN CALL cp_fm_vect_dealloc(p_env%PS_psi0) END IF diff --git a/src/qs_tddfpt2_densities.F b/src/qs_tddfpt2_densities.F index 44d60b99b3..ab8a4fa551 100644 --- a/src/qs_tddfpt2_densities.F +++ b/src/qs_tddfpt2_densities.F @@ -6,6 +6,8 @@ !--------------------------------------------------------------------------------------------------! MODULE qs_tddfpt2_densities + USE admm_types, ONLY: admm_type,& + get_admm_env USE cp_control_types, ONLY: dft_control_type USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,& copy_fm_to_dbcsr @@ -14,19 +16,25 @@ MODULE qs_tddfpt2_densities USE cp_gemm_interface, ONLY: cp_gemm USE dbcsr_api, ONLY: dbcsr_p_type,& dbcsr_scale - USE kinds, ONLY: dp + USE kinds, ONLY: default_string_length,& + dp USE pw_env_types, ONLY: pw_env_get USE pw_pool_types, ONLY: pw_pool_type USE pw_types, ONLY: pw_p_type USE qs_collocate_density, ONLY: calculate_rho_elec USE qs_environment_types, ONLY: get_qs_env,& qs_environment_type + USE qs_gapw_densities, ONLY: prepare_gapw_den USE qs_ks_types, ONLY: qs_ks_env_type + USE qs_local_rho_types, ONLY: local_rho_type + USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type + USE qs_rho_atom_methods, ONLY: calculate_rho_atom_coeff USE qs_rho_methods, ONLY: qs_rho_copy,& qs_rho_update_rho USE qs_rho_types, ONLY: qs_rho_get,& qs_rho_type USE qs_tddfpt2_subgroups, ONLY: tddfpt_subgroup_env_type + USE task_list_types, ONLY: task_list_type #include "./base/base_uses.f90" IMPLICIT NONE @@ -102,6 +110,7 @@ CONTAINS pw_env_external=sub_env%pw_env, & task_list_external=sub_env%task_list_orb_soft, & para_env_external=sub_env%para_env) + CALL prepare_gapw_den(qs_env, local_rho_set=sub_env%local_rho_set) ELSEIF (dft_control%qs_control%gapw_xc) THEN CALL qs_rho_update_rho(rho_orb_struct, qs_env, & rho_xc_external=rho_xc_struct, & @@ -112,6 +121,7 @@ CONTAINS para_env_external=sub_env%para_env) CALL pw_env_get(sub_env%pw_env, auxbas_pw_pool=auxbas_pw_pool) CALL qs_rho_copy(rho_xc_struct, rho_orb_struct, auxbas_pw_pool, nspins) + CALL prepare_gapw_den(qs_env, local_rho_set=sub_env%local_rho_set, do_rho0=.FALSE.) ELSE CALL qs_rho_update_rho(rho_orb_struct, qs_env, & pw_env_external=sub_env%pw_env, & @@ -127,6 +137,7 @@ CONTAINS !> \brief Project a charge density expressed in primary basis set into the auxiliary basis set. !> \param rho_orb_struct response density in primary basis set !> \param rho_aux_fit_struct response density in auxiliary basis set (modified on exit) +!> \param local_rho_set GAPW density of auxiliary basis set density !> \param qs_env Quickstep environment !> \param sub_env parallel (sub)group environment !> \param wfm_rho_orb work dense matrix with shape [nao x nao] distributed among @@ -142,32 +153,47 @@ CONTAINS !> tddfpt_construct_ground_state_orb_density() and tddfpt_construct_aux_fit_density() !> in order to avoid code duplication [Sergey Chulkov] ! ************************************************************************************************** - SUBROUTINE tddfpt_construct_aux_fit_density(rho_orb_struct, rho_aux_fit_struct, qs_env, sub_env, & + SUBROUTINE tddfpt_construct_aux_fit_density(rho_orb_struct, rho_aux_fit_struct, local_rho_set, & + qs_env, sub_env, & wfm_rho_orb, wfm_rho_aux_fit, wfm_aux_orb) TYPE(qs_rho_type), POINTER :: rho_orb_struct, rho_aux_fit_struct + TYPE(local_rho_type), POINTER :: local_rho_set TYPE(qs_environment_type), POINTER :: qs_env TYPE(tddfpt_subgroup_env_type), INTENT(in) :: sub_env TYPE(cp_fm_type), POINTER :: wfm_rho_orb, wfm_rho_aux_fit, wfm_aux_orb CHARACTER(LEN=*), PARAMETER :: routineN = 'tddfpt_construct_aux_fit_density' + CHARACTER(LEN=default_string_length) :: basis_type INTEGER :: handle, ispin, nao, nao_aux, nspins REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho_aux_fit_r + TYPE(admm_type), POINTER :: admm_env TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao_aux_fit, rho_ao_orb + TYPE(neighbor_list_set_p_type), DIMENSION(:), & + POINTER :: sab_aux_fit TYPE(pw_p_type), DIMENSION(:), POINTER :: rho_aux_fit_g, rho_aux_fit_r TYPE(qs_ks_env_type), POINTER :: ks_env + TYPE(task_list_type), POINTER :: task_list CALL timeset(routineN, handle) CPASSERT(ASSOCIATED(sub_env%admm_A)) - CALL get_qs_env(qs_env, ks_env=ks_env) + CALL get_qs_env(qs_env, ks_env=ks_env, admm_env=admm_env) CALL qs_rho_get(rho_orb_struct, rho_ao=rho_ao_orb) CALL qs_rho_get(rho_aux_fit_struct, rho_ao=rho_ao_aux_fit, rho_g=rho_aux_fit_g, & rho_r=rho_aux_fit_r, tot_rho_r=tot_rho_aux_fit_r) nspins = SIZE(rho_ao_orb) + IF (admm_env%do_gapw) THEN + basis_type = "AUX_FIT_SOFT" + task_list => sub_env%task_list_aux_fit_soft + ELSE + basis_type = "AUX_FIT" + task_list => sub_env%task_list_aux_fit + END IF + CALL cp_fm_get_info(sub_env%admm_A, nrow_global=nao_aux, ncol_global=nao) DO ispin = 1, nspins ! TO DO: consider sub_env%admm_A to be a DBCSR matrix @@ -181,9 +207,18 @@ CONTAINS CALL calculate_rho_elec(matrix_p=rho_ao_aux_fit(ispin)%matrix, & rho=rho_aux_fit_r(ispin), rho_gspace=rho_aux_fit_g(ispin), & total_rho=tot_rho_aux_fit_r(ispin), ks_env=ks_env, & - soft_valid=.FALSE., basis_type="AUX_FIT", & - pw_env_external=sub_env%pw_env, task_list_external=sub_env%task_list_aux_fit) + soft_valid=.FALSE., basis_type=basis_type, & + pw_env_external=sub_env%pw_env, task_list_external=task_list) END DO + IF (admm_env%do_gapw) THEN + CALL get_admm_env(qs_env%admm_env, sab_aux_fit=sab_aux_fit) + CALL calculate_rho_atom_coeff(qs_env, rho_ao_aux_fit, & + rho_atom_set=local_rho_set%rho_atom_set, & + qs_kind_set=admm_env%admm_gapw_env%admm_kind_set, & + oce=admm_env%admm_gapw_env%oce, sab=sab_aux_fit, para_env=sub_env%para_env) + CALL prepare_gapw_den(qs_env, local_rho_set=local_rho_set, & + do_rho0=.FALSE., kind_set_external=admm_env%admm_gapw_env%admm_kind_set) + END IF CALL timestop(handle) diff --git a/src/qs_tddfpt2_fhxc.F b/src/qs_tddfpt2_fhxc.F index 8f2af196a7..853be8b7f7 100644 --- a/src/qs_tddfpt2_fhxc.F +++ b/src/qs_tddfpt2_fhxc.F @@ -6,6 +6,7 @@ !--------------------------------------------------------------------------------------------------! MODULE qs_tddfpt2_fhxc + USE admm_types, ONLY: admm_type USE cp_control_types, ONLY: dft_control_type,& stda_control_type USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl @@ -22,7 +23,9 @@ MODULE qs_tddfpt2_fhxc USE dbcsr_api, ONLY: & dbcsr_add, dbcsr_copy, dbcsr_create, dbcsr_deallocate_matrix, dbcsr_get_info, & dbcsr_p_type, dbcsr_release, dbcsr_set, dbcsr_type, dbcsr_type_symmetric - USE kinds, ONLY: dp + USE input_constants, ONLY: do_admm_aux_exch_func_none + USE kinds, ONLY: default_string_length,& + dp USE lri_environment_types, ONLY: lri_kind_type USE message_passing, ONLY: mp_sum USE pw_env_types, ONLY: pw_env_get @@ -42,6 +45,7 @@ MODULE qs_tddfpt2_fhxc integrate_v_rspace_one_center USE qs_kernel_types, ONLY: full_kernel_env_type USE qs_ks_atom, ONLY: update_ks_atom + USE qs_rho_atom_types, ONLY: rho_atom_type USE qs_rho_methods, ONLY: qs_rho_update_rho,& qs_rho_update_tddfpt USE qs_rho_types, ONLY: qs_rho_get @@ -54,6 +58,7 @@ MODULE qs_tddfpt2_fhxc USE qs_tddfpt2_subgroups, ONLY: tddfpt_subgroup_env_type USE qs_tddfpt2_types, ONLY: tddfpt_work_matrices USE qs_vxc_atom, ONLY: calculate_xc_2nd_deriv_atom + USE task_list_types, ONLY: task_list_type #include "./base/base_uses.f90" IMPLICIT NONE @@ -103,11 +108,13 @@ CONTAINS CHARACTER(LEN=*), PARAMETER :: routineN = 'fhxc_kernel' + CHARACTER(LEN=default_string_length) :: basis_type INTEGER :: handle, ikind, ispin, ivect, nao, & nao_aux, nkind, nspins, nvects INTEGER, DIMENSION(:), POINTER :: blk_sizes INTEGER, DIMENSION(maxspins) :: nactive LOGICAL :: gapw, gapw_xc + TYPE(admm_type), POINTER :: admm_env TYPE(cp_fm_type), POINTER :: work_aux_orb, work_orb_orb TYPE(cp_para_env_type), POINTER :: para_env TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: A_xc_munu_sub, rho_ia_ao, & @@ -119,6 +126,8 @@ CONTAINS rho_ia_r_aux_fit, tau_ia_r, & tau_ia_r_aux_fit, V_rspace_sub TYPE(pw_pool_type), POINTER :: auxbas_pw_pool + TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set + TYPE(task_list_type), POINTER :: task_list CALL timeset(routineN, handle) @@ -141,6 +150,7 @@ CONTAINS CALL qs_rho_get(work_matrices%rho_orb_struct_sub, rho_ao=rho_ia_ao, & rho_g=rho_ia_g, rho_r=rho_ia_r, tau_r=tau_ia_r) IF (do_hfx .AND. do_admm) THEN + CALL get_qs_env(qs_env, admm_env=admm_env) CALL qs_rho_get(work_matrices%rho_aux_fit_struct_sub, & rho_ao=rho_ia_ao_aux_fit, rho_g=rho_ia_g_aux_fit, & rho_r=rho_ia_r_aux_fit, tau_r=tau_ia_r_aux_fit) @@ -251,19 +261,21 @@ CONTAINS work_v_xc_tau=work_matrices%wpw_tau_rspace_sub) END IF IF (gapw .OR. gapw_xc) THEN - CALL calculate_xc_2nd_deriv_atom(work_matrices%local_rho_set, qs_env, kernel_env%xc_section, & + rho_atom_set => sub_env%local_rho_set%rho_atom_set + rho1_atom_set => work_matrices%local_rho_set%rho_atom_set + CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, kernel_env%xc_section, & sub_env%para_env, do_tddfpt2=.TRUE., do_triplet=is_rks_triplets) END IF ! ADMM correction - IF (do_admm) THEN + IF (do_admm .AND. dft_control%admm_control%aux_exch_func /= do_admm_aux_exch_func_none) THEN CALL tddfpt_construct_aux_fit_density(rho_orb_struct=work_matrices%rho_orb_struct_sub, & rho_aux_fit_struct=work_matrices%rho_aux_fit_struct_sub, & + local_rho_set=work_matrices%local_rho_set_admm, & qs_env=qs_env, sub_env=sub_env, & wfm_rho_orb=work_matrices%rho_ao_orb_fm_sub, & wfm_rho_aux_fit=work_matrices%rho_ao_aux_fit_fm_sub, & wfm_aux_orb=work_matrices%wfm_aux_orb_sub) - ! - C_{HF} d^{2}E_{x, ADMM}^{DFT}[\hat{\rho}] / d\hat{\rho}^2 IF (admm_symm) THEN CALL dbcsr_get_info(rho_ia_ao_aux_fit(1)%matrix, row_blk_size=blk_sizes) @@ -286,6 +298,14 @@ CONTAINS CALL pw_zero(V_rspace_sub(ispin)%pw) END DO + IF (admm_env%do_gapw) THEN + basis_type = "AUX_FIT_SOFT" + task_list => sub_env%task_list_aux_fit_soft + ELSE + basis_type = "AUX_FIT" + task_list => sub_env%task_list_aux_fit + END IF + CALL tddfpt_apply_xc(A_ia_rspace=V_rspace_sub, & kernel_env=kernel_env_admm_aux, & rho_ia_struct=work_matrices%rho_aux_fit_struct_sub, & @@ -298,9 +318,22 @@ CONTAINS hmat=A_xc_munu_sub(ispin), & qs_env=qs_env, calculate_forces=.FALSE., & pw_env_external=sub_env%pw_env, & - basis_type="AUX_FIT", & - task_list_external=sub_env%task_list_aux_fit) + basis_type=basis_type, & + task_list_external=task_list) END DO + IF (admm_env%do_gapw) THEN + rho_atom_set => sub_env%local_rho_set_admm%rho_atom_set + rho1_atom_set => work_matrices%local_rho_set_admm%rho_atom_set + CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, & + kernel_env_admm_aux%xc_section, & + sub_env%para_env, do_tddfpt2=.TRUE., do_triplet=.FALSE., & + kind_set_external=admm_env%admm_gapw_env%admm_kind_set) + CALL update_ks_atom(qs_env, A_xc_munu_sub, rho_ia_ao_aux_fit, forces=.FALSE., tddft=.TRUE., & + rho_atom_external=rho1_atom_set, & + kind_set_external=admm_env%admm_gapw_env%admm_kind_set, & + oce_external=admm_env%admm_gapw_env%oce, & + sab_external=sub_env%sab_aux_fit) + END IF ALLOCATE (dbwork) CALL dbcsr_create(dbwork, template=work_matrices%A_ia_munu_sub(1)%matrix) CALL cp_fm_create(work_aux_orb, & @@ -337,6 +370,10 @@ CONTAINS is_rks_triplets=is_rks_triplets, pw_env=sub_env%pw_env, & work_v_xc=work_matrices%wpw_rspace_sub, & work_v_xc_tau=work_matrices%wpw_tau_rspace_sub) + IF (admm_env%do_gapw) THEN + CPWARN("GAPW/ADMM needs symmetric ADMM kernel") + CPABORT("GAPW/ADMM@TDDFT") + END IF END IF END IF diff --git a/src/qs_tddfpt2_methods.F b/src/qs_tddfpt2_methods.F index d407fe6ee3..ae8d9ac822 100644 --- a/src/qs_tddfpt2_methods.F +++ b/src/qs_tddfpt2_methods.F @@ -276,6 +276,7 @@ CONTAINS CALL tddfpt_construct_aux_fit_density(rho_orb_struct=work_matrices%rho_orb_struct_sub, & rho_aux_fit_struct=work_matrices%rho_aux_fit_struct_sub, & + local_rho_set=sub_env%local_rho_set_admm, & qs_env=qs_env, sub_env=sub_env, & wfm_rho_orb=work_matrices%rho_ao_orb_fm_sub, & wfm_rho_aux_fit=work_matrices%rho_ao_aux_fit_fm_sub, & diff --git a/src/qs_tddfpt2_operators.F b/src/qs_tddfpt2_operators.F index 04b5b34b71..e8b2a43e3d 100644 --- a/src/qs_tddfpt2_operators.F +++ b/src/qs_tddfpt2_operators.F @@ -351,6 +351,9 @@ CONTAINS xc_section=kernel_env%xc_section, gapw=.FALSE., tddfpt_fac=kernel_env%beta) DEALLOCATE (rho_ia_g2, rho_ia_r2) + IF (ASSOCIATED(tau_ia_r2)) THEN + DEALLOCATE (tau_ia_r2) + END IF CALL timestop(handle) diff --git a/src/qs_tddfpt2_subgroups.F b/src/qs_tddfpt2_subgroups.F index e26acbf4f0..e8eeb63a62 100644 --- a/src/qs_tddfpt2_subgroups.F +++ b/src/qs_tddfpt2_subgroups.F @@ -146,6 +146,7 @@ MODULE qs_tddfpt2_subgroups !> GAPW local atomic grids TYPE(hartree_local_type), POINTER :: hartree_local => NULL() TYPE(local_rho_type), POINTER :: local_rho_set => NULL() + TYPE(local_rho_type), POINTER :: local_rho_set_admm => NULL() END TYPE tddfpt_subgroup_env_type ! ************************************************************************************************** @@ -311,12 +312,20 @@ CONTAINS reorder_grid_ranks=.TRUE.) END IF - IF (dft_control%do_admm) & + IF (dft_control%do_admm) THEN CALL tddfpt_build_tasklist(task_list=sub_env%task_list_aux_fit, sab=sub_env%sab_aux_fit, & basis_type="AUX_FIT", distribution_2d=sub_env%dist_2d, & pw_env=sub_env%pw_env, qs_env=qs_env, soft_valid=.FALSE., & skip_load_balance=qs_control%skip_load_balance_distributed, & reorder_grid_ranks=.FALSE.) + IF (qs_control%gapw .OR. qs_control%gapw_xc) THEN + CALL tddfpt_build_tasklist(task_list=sub_env%task_list_aux_fit_soft, sab=sub_env%sab_aux_fit, & + basis_type="AUX_FIT", distribution_2d=sub_env%dist_2d, & + pw_env=sub_env%pw_env, qs_env=qs_env, soft_valid=.TRUE., & + skip_load_balance=qs_control%skip_load_balance_distributed, & + reorder_grid_ranks=.FALSE.) + END IF + END IF IF (tddfpt_control%mgrid_is_explicit) & CALL restore_qs_mgrid(qs_control, mgrid_saved) @@ -326,9 +335,13 @@ CONTAINS CALL get_qs_env(qs_env, dbcsr_dist=sub_env%dbcsr_dist, & sab_orb=sub_env%sab_orb, task_list=sub_env%task_list_orb) - IF (dft_control%do_admm) & - CALL get_admm_env(qs_env%admm_env, sab_aux_fit=sub_env%sab_aux_fit, & + IF (dft_control%do_admm) THEN + CALL get_admm_env(admm_env, sab_aux_fit=sub_env%sab_aux_fit, & task_list_aux_fit=sub_env%task_list_aux_fit) + IF (qs_control%gapw .OR. qs_control%gapw_xc) THEN + sub_env%task_list_aux_fit_soft => admm_env%admm_gapw_env%task_list + END IF + END IF IF (qs_control%gapw .OR. qs_control%gapw_xc) THEN CALL get_qs_env(qs_env, task_list_soft=sub_env%task_list_orb_soft) END IF @@ -358,6 +371,17 @@ CONTAINS qs_kind_set, dft_control, sub_env%para_env) END IF + ! ADMM/GAPW + IF (dft_control%do_admm) THEN + IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) THEN + CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set) + CALL local_rho_set_create(sub_env%local_rho_set_admm) + CALL allocate_rho_atom_internals(sub_env%local_rho_set_admm%rho_atom_set, atomic_kind_set, & + admm_env%admm_gapw_env%admm_kind_set, & + dft_control, sub_env%para_env) + END IF + END IF + ELSE IF (kernel == tddfpt_kernel_stda) THEN sub_env%is_mgrid = .FALSE. NULLIFY (sub_env%dbcsr_dist, sub_env%dist_2d) @@ -432,6 +456,9 @@ CONTAINS IF (ASSOCIATED(sub_env%hartree_local)) THEN CALL hartree_local_release(sub_env%hartree_local) END IF + IF (ASSOCIATED(sub_env%local_rho_set_admm)) THEN + CALL local_rho_set_release(sub_env%local_rho_set_admm) + END IF ! if TDDFPT-specific plane-wave environment has not been requested, ! the pointers sub_env%dbcsr_dist, sub_env%sab_*, and sub_env%task_list_* diff --git a/src/qs_tddfpt2_types.F b/src/qs_tddfpt2_types.F index ae344b63ba..1efbaa122f 100644 --- a/src/qs_tddfpt2_types.F +++ b/src/qs_tddfpt2_types.F @@ -197,6 +197,7 @@ MODULE qs_tddfpt2_types !> GAPW local atomic grids TYPE(hartree_local_type), POINTER :: hartree_local TYPE(local_rho_type), POINTER :: local_rho_set + TYPE(local_rho_type), POINTER :: local_rho_set_admm END TYPE tddfpt_work_matrices CONTAINS @@ -229,6 +230,7 @@ CONTAINS INTEGER :: handle, igroup, ispin, istate, nao, & nao_aux, natom, ngroups, nspins INTEGER, DIMENSION(maxspins) :: nmo_occ, nmo_virt + TYPE(admm_type), POINTER :: admm_env TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set TYPE(cp_blacs_env_type), POINTER :: blacs_env TYPE(cp_fm_struct_p_type), DIMENSION(maxspins) :: fm_struct_evects @@ -258,6 +260,7 @@ CONTAINS NULLIFY (work_matrices%hartree_local) NULLIFY (work_matrices%local_rho_set) + NULLIFY (work_matrices%local_rho_set_admm) NULLIFY (work_matrices%rho_xc_struct_sub) nspins = SIZE(gs_mos) @@ -456,6 +459,13 @@ CONTAINS CALL get_qs_env(qs_env, dbcsr_dist=dbcsr_dist) CALL get_admm_env(qs_env%admm_env, sab_aux_fit=sab_hfx) dbcsr_template_hfx => matrix_s_aux_fit(1)%matrix + IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) THEN + CALL get_qs_env(qs_env, admm_env=admm_env, atomic_kind_set=atomic_kind_set) + CALL local_rho_set_create(work_matrices%local_rho_set_admm) + CALL allocate_rho_atom_internals(work_matrices%local_rho_set_admm%rho_atom_set, & + atomic_kind_set, admm_env%admm_gapw_env%admm_kind_set, & + dft_control, sub_env%para_env) + END IF ELSE CALL get_qs_env(qs_env, dbcsr_dist=dbcsr_dist, sab_orb=sab_hfx) dbcsr_template_hfx => matrix_s(1)%matrix @@ -631,6 +641,7 @@ CONTAINS NULLIFY (work_matrices%hartree_local) NULLIFY (work_matrices%local_rho_set) + NULLIFY (work_matrices%local_rho_set_admm) NULLIFY (work_matrices%rho_xc_struct_sub) CALL timestop(handle) @@ -762,6 +773,9 @@ CONTAINS IF (ASSOCIATED(work_matrices%local_rho_set)) THEN CALL local_rho_set_release(work_matrices%local_rho_set) END IF + IF (ASSOCIATED(work_matrices%local_rho_set_admm)) THEN + CALL local_rho_set_release(work_matrices%local_rho_set_admm) + END IF IF (ASSOCIATED(work_matrices%hartree_local)) THEN CALL hartree_local_release(work_matrices%hartree_local) END IF diff --git a/src/qs_vxc_atom.F b/src/qs_vxc_atom.F index 57148ed4b1..e2c4dd24d0 100644 --- a/src/qs_vxc_atom.F +++ b/src/qs_vxc_atom.F @@ -29,7 +29,6 @@ MODULE qs_vxc_atom nsoset USE paw_proj_set_types, ONLY: get_paw_proj_set,& paw_proj_set_type - USE qs_energy_types, ONLY: qs_energy_type USE qs_environment_types, ONLY: get_qs_env,& qs_environment_type USE qs_grid_atom, ONLY: grid_atom_type @@ -39,7 +38,6 @@ MODULE qs_vxc_atom has_nlcc,& qs_kind_type USE qs_linres_types, ONLY: nablavks_atom_type - USE qs_local_rho_types, ONLY: local_rho_type USE qs_rho_atom_types, ONLY: get_rho_atom,& rho_atom_coeff,& rho_atom_type @@ -462,7 +460,8 @@ CONTAINS ! ************************************************************************************************** !> \brief ... -!> \param local_rho_set ... +!> \param rho_atom_set ... +!> \param rho1_atom_set ... !> \param qs_env ... !> \param xc_section ... !> \param para_env ... @@ -470,15 +469,18 @@ CONTAINS !> 'DFT' input section !> \param do_tddfpt2 New implementation of TDDFT. !> \param do_triplet ... +!> \param kind_set_external ... ! ************************************************************************************************** - SUBROUTINE calculate_xc_2nd_deriv_atom(local_rho_set, qs_env, xc_section, para_env, & - do_tddft, do_tddfpt2, do_triplet) + SUBROUTINE calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, & + do_tddft, do_tddfpt2, do_triplet, kind_set_external) - TYPE(local_rho_type), POINTER :: local_rho_set + TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set, rho1_atom_set TYPE(qs_environment_type), POINTER :: qs_env TYPE(section_vals_type), POINTER :: xc_section TYPE(cp_para_env_type), POINTER :: para_env LOGICAL, INTENT(IN), OPTIONAL :: do_tddft, do_tddfpt2, do_triplet + TYPE(qs_kind_type), DIMENSION(:), OPTIONAL, & + POINTER :: kind_set_external CHARACTER(LEN=*), PARAMETER :: routineN = 'calculate_xc_2nd_deriv_atom' @@ -500,12 +502,10 @@ CONTAINS TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set TYPE(grid_atom_type), POINTER :: grid_atom TYPE(harmonics_atom_type), POINTER :: harmonics - TYPE(qs_energy_type), POINTER :: energy - TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set + TYPE(qs_kind_type), DIMENSION(:), POINTER :: my_kind_set, qs_kind_set TYPE(rho_atom_coeff), DIMENSION(:), POINTER :: dr1_h, dr1_s, dr_h, dr_s, r1_h, r1_s, & r_h, r_s TYPE(rho_atom_coeff), DIMENSION(:, :), POINTER :: r1_h_d, r1_s_d, r_h_d, r_s_d - TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set TYPE(rho_atom_type), POINTER :: rho1_atom, rho_atom TYPE(section_vals_type), POINTER :: input, xc_fun_section TYPE(xc_derivative_set_type) :: deriv_set @@ -517,21 +517,22 @@ CONTAINS CALL timeset(routineN, handle) - NULLIFY (qs_kind_set, energy) + NULLIFY (qs_kind_set) NULLIFY (rho_h, rho_s, drho_h, drho_s, weight) NULLIFY (rho1_h, rho1_s, drho1_h, drho1_s) NULLIFY (vxc_h, vxc_s, vxg_h, vxg_s) NULLIFY (tau_h, tau_s, tau1_h, tau1_s) CALL get_qs_env(qs_env=qs_env, & - para_env=para_env, & - energy=energy, & input=input, & qs_kind_set=qs_kind_set, & - atomic_kind_set=atomic_kind_set, & - rho_atom_set=rho_atom_set) + atomic_kind_set=atomic_kind_set) - rho1_atom_set => local_rho_set%rho_atom_set + IF (PRESENT(kind_set_external)) THEN + my_kind_set => kind_set_external + ELSE + my_kind_set => qs_kind_set + END IF my_tddft = .FALSE. IF (PRESENT(do_tddft)) my_tddft = do_tddft @@ -584,7 +585,7 @@ CONTAINS NULLIFY (atom_list, harmonics, grid_atom) CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=natom) - CALL get_qs_kind(qs_kind_set(ikind), paw_atom=paw_atom, & + CALL get_qs_kind(my_kind_set(ikind), paw_atom=paw_atom, & harmonics=harmonics, grid_atom=grid_atom) IF (.NOT. paw_atom) CYCLE @@ -716,10 +717,10 @@ CONTAINS w=weight, vxc=vxc_s, vxg=vxg_s) IF (gradient_functional) THEN - CALL gaVxcgb_GC(vxc_h, vxc_s, vxg_h, vxg_s, qs_kind_set(ikind), & + CALL gaVxcgb_GC(vxc_h, vxc_s, vxg_h, vxg_s, my_kind_set(ikind), & rho1_atom, nspins) ELSE - CALL gaVxcgb_noGC(vxc_h, vxc_s, qs_kind_set(ikind), & + CALL gaVxcgb_noGC(vxc_h, vxc_s, my_kind_set(ikind), & rho1_atom, nspins) END IF diff --git a/tests/QS/regtest-debug-1/TEST_FILES b/tests/QS/regtest-debug-1/TEST_FILES index 55ad201999..c1092e4bfb 100644 --- a/tests/QS/regtest-debug-1/TEST_FILES +++ b/tests/QS/regtest-debug-1/TEST_FILES @@ -8,5 +8,7 @@ h2o_polar.inp 87 1e-05 h2o_pdip.inp 86 1e-05 0.961018561809E+00 h2o_periodic.inp 87 1e-05 0.139741440657E+02 h2o_gga.inp 87 1e-05 0.165143521103E+02 -h2o_gapw.inp 87 4e-04 0.170220498303E+02 +h2o_gapw.inp 87 1e-05 0.170288288429E+02 +h2o_gapw_xc.inp 87 1e-05 0.170142638188E+02 +h2o_admm_gapw.inp 87 1e-05 0.165349668089E+02 #EOF diff --git a/tests/QS/regtest-debug-1/h2o_admm_gapw.inp b/tests/QS/regtest-debug-1/h2o_admm_gapw.inp new file mode 100644 index 0000000000..038f91be41 --- /dev/null +++ b/tests/QS/regtest-debug-1/h2o_admm_gapw.inp @@ -0,0 +1,93 @@ +&FORCE_EVAL + METHOD Quickstep + &PROPERTIES + &LINRES + PRECONDITIONER FULL_ALL + EPS 1.e-10 + &POLAR + DO_RAMAN T + PERIODIC_DIPOLE_OPERATOR F + &END + &END + &END + &DFT + &QS + METHOD GAPW + EPS_DEFAULT 1.e-14 + &END QS + BASIS_SET_FILE_NAME BASIS_SET + BASIS_SET_FILE_NAME BASIS_ADMM + &EFIELD + &END + &SCF + SCF_GUESS ATOMIC + &OT + PRECONDITIONER FULL_SINGLE_INVERSE + MINIMIZER DIIS + &END + &OUTER_SCF + MAX_SCF 10 + EPS_SCF 1.0E-6 + &END + MAX_SCF 10 + EPS_SCF 1.0E-6 + &END SCF + &AUXILIARY_DENSITY_MATRIX_METHOD + METHOD BASIS_PROJECTION + ADMM_PURIFICATION_METHOD NONE + EXCH_CORRECTION_FUNC PBEX + EXCH_SCALING_MODEL NONE + &END + &XC + &XC_FUNCTIONAL PBE0 + &END XC_FUNCTIONAL + &END XC + &PRINT + &MOMENTS ON + PERIODIC .FALSE. + REFERENCE COM + &END + &END + &END DFT + &SUBSYS + &CELL + ABC [angstrom] 4.0 4.0 4.0 + PERIODIC NONE + &END + &COORD + O 0.000000 0.000000 -0.065587 + H 0.000000 -0.757136 0.520545 + H 0.000000 0.757136 0.520545 + &END COORD + &TOPOLOGY + &CENTER_COORDINATES + &END + &END + &KIND H + BASIS_SET ORB DZV-GTH-PADE + BASIS_SET AUX_FIT fit3 + POTENTIAL GTH-PADE-q1 + &END KIND + &KIND O + BASIS_SET ORB DZVP-GTH-PADE + BASIS_SET AUX_FIT fit3 + POTENTIAL GTH-PADE-q6 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PRINT_LEVEL LOW + PROJECT dip + RUN_TYPE DEBUG +&END GLOBAL + +&DEBUG + DEBUG_FORCES .FALSE. + DEBUG_STRESS_TENSOR .FALSE. + DEBUG_DIPOLE .FALSE. + DEBUG_POLARIZABILITY .TRUE. + DE 0.002 + EPS_NO_ERROR_CHECK 5.e-5 +&END + + diff --git a/tests/QS/regtest-debug-1/h2o_gapw_xc.inp b/tests/QS/regtest-debug-1/h2o_gapw_xc.inp new file mode 100644 index 0000000000..e6f767de43 --- /dev/null +++ b/tests/QS/regtest-debug-1/h2o_gapw_xc.inp @@ -0,0 +1,83 @@ +&FORCE_EVAL + METHOD Quickstep + &PROPERTIES + &LINRES + PRECONDITIONER FULL_ALL + EPS 1.e-10 + &POLAR + DO_RAMAN T + PERIODIC_DIPOLE_OPERATOR F + &END + &END + &END + &DFT + &QS + METHOD GAPW_XC + EPS_DEFAULT 1.e-10 + &END QS + &EFIELD + &END + &SCF + SCF_GUESS ATOMIC + &OT + PRECONDITIONER FULL_SINGLE_INVERSE + MINIMIZER DIIS + &END + &OUTER_SCF + MAX_SCF 10 + EPS_SCF 1.0E-4 + &END + MAX_SCF 10 + EPS_SCF 1.0E-4 + &END SCF + &XC + &XC_FUNCTIONAL PBE + &END XC_FUNCTIONAL + &END XC + &PRINT + &MOMENTS ON + PERIODIC .FALSE. + REFERENCE COM + &END + &END + &END DFT + &SUBSYS + &CELL + ABC [angstrom] 4.0 4.0 4.0 + PERIODIC NONE + &END + &COORD + O 0.000000 0.000000 -0.065587 + H 0.000000 -0.757136 0.520545 + H 0.000000 0.757136 0.520545 + &END COORD + &TOPOLOGY + &CENTER_COORDINATES + &END + &END + &KIND H + BASIS_SET DZV-GTH-PADE + POTENTIAL GTH-PADE-q1 + &END KIND + &KIND O + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q6 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PRINT_LEVEL LOW + PROJECT dipole + RUN_TYPE DEBUG +&END GLOBAL + +&DEBUG + DEBUG_FORCES .FALSE. + DEBUG_STRESS_TENSOR .FALSE. + DEBUG_DIPOLE .FALSE. + DEBUG_POLARIZABILITY .TRUE. + DE 0.002 + EPS_NO_ERROR_CHECK 5.e-5 +&END + + diff --git a/tests/QS/regtest-tddfpt-4/TEST_FILES b/tests/QS/regtest-tddfpt-4/TEST_FILES index dc155dc7f4..0b0f23d44f 100644 --- a/tests/QS/regtest-tddfpt-4/TEST_FILES +++ b/tests/QS/regtest-tddfpt-4/TEST_FILES @@ -4,26 +4,27 @@ # 1 compares the last total energy in the file # for details see cp2k/tools/do_regtest # -#test01.inp 1 1.0E-11 -17.18917899734168 -#test02.inp 1 1.0E-11 -17.18917899734168 -#test03.inp 1 1.0E-11 -17.18917899734168 -#test04.inp 1 1.0E-11 -17.18917899734168 -#test05.inp 1 1.0E-11 -17.18917899734168 -#test06.inp 1 1.0E-11 -17.18917899734168 -#test07.inp 1 1.0E-11 -17.18917899734168 -#test08.inp 1 1.0E-11 -17.18917899734168 -#test09.inp 1 1.0E-11 -17.18917899734168 -#test10.inp 1 1.0E-11 -17.18917899734168 -#test11.inp 1 1.0E-11 -17.18917899734168 -#test12.inp 1 1.0E-11 -17.18917899734168 -#test13.inp 1 1.0E-11 -17.18917899734168 -#test14.inp 1 1.0E-11 -17.18917899734168 -#test15.inp 1 1.0E-11 -17.18917899734168 -#test16.inp 1 1.0E-11 -17.18917899734168 -#test17.inp 1 1.0E-11 -17.18917899734168 -#test18.inp 1 1.0E-11 -17.18917899734168 -#test19.inp 1 1.0E-11 -17.18917899734168 -#test20.inp 1 1.0E-11 -17.18917899734168 -#test21.inp 1 1.0E-11 -17.18917899734168 -#test22.inp 1 1.0E-11 -17.18917899734168 +test01.inp 37 5.0E-06 0.976371E+00 +test02.inp 37 5.0E-06 0.102132E+01 +test03.inp 37 5.0E-06 0.102132E+01 +test04.inp 37 5.0E-06 0.102641E+01 +test05.inp 37 5.0E-06 0.101974E+01 +test06.inp 37 5.0E-06 0.102171E+01 +test07.inp 37 5.0E-06 0.101792E+01 +test08.inp 37 5.0E-06 0.103120E+01 +test09.inp 37 5.0E-06 0.101371E+01 +test10.inp 37 5.0E-06 0.997064E+00 +test11.inp 37 5.0E-06 0.105541E+01 +test12.inp 37 5.0E-06 0.105058E+01 +test13.inp 37 5.0E-06 0.105541E+01 +test14.inp 37 5.0E-06 0.105058E+01 +test15.inp 37 5.0E-06 0.109231E+01 +test16.inp 37 5.0E-06 0.106367E+01 +test17.inp 37 5.0E-06 0.170669E+01 +test18.inp 37 5.0E-06 0.169670E+01 +test19.inp 37 5.0E-06 0.121034E+01 +test20.inp 37 5.0E-06 0.120983E+01 +test21.inp 37 5.0E-06 0.105522E+01 +test22.inp 37 5.0E-06 0.109138E+01 +test23.inp 37 5.0E-06 0.109395E+01 #EOF diff --git a/tests/QS/regtest-tddfpt-4/test09.inp b/tests/QS/regtest-tddfpt-4/test09.inp index 39021184f9..4561a19ac4 100644 --- a/tests/QS/regtest-tddfpt-4/test09.inp +++ b/tests/QS/regtest-tddfpt-4/test09.inp @@ -18,13 +18,10 @@ BASIS_SET_FILE_NAME BASIS_ADMM &AUXILIARY_DENSITY_MATRIX_METHOD ADMM_PURIFICATION_METHOD NONE - EXCH_CORRECTION_FUNC DEFAULT + EXCH_CORRECTION_FUNC PBEX EXCH_SCALING_MODEL NONE METHOD BASIS_PROJECTION &END - &EXCITED_STATES T - STATE 1 - &END EXCITED_STATES &SCF SCF_GUESS ATOMIC &OT @@ -113,5 +110,5 @@ &GLOBAL PRINT_LEVEL LOW PROJECT ftest - RUN_TYPE ENERGY_FORCE + RUN_TYPE ENERGY &END GLOBAL diff --git a/tests/QS/regtest-tddfpt-4/test10.inp b/tests/QS/regtest-tddfpt-4/test10.inp index 72e2785931..72c665c6e4 100644 --- a/tests/QS/regtest-tddfpt-4/test10.inp +++ b/tests/QS/regtest-tddfpt-4/test10.inp @@ -13,9 +13,6 @@ &QS METHOD GPW &END QS - &EXCITED_STATES T - STATE 1 - &END EXCITED_STATES &SCF SCF_GUESS ATOMIC &OT @@ -102,5 +99,5 @@ &GLOBAL PRINT_LEVEL LOW PROJECT ftest - RUN_TYPE ENERGY_FORCE + RUN_TYPE ENERGY &END GLOBAL diff --git a/tests/QS/regtest-tddfpt-4/test13.inp b/tests/QS/regtest-tddfpt-4/test13.inp index f472887501..52f9c34af9 100644 --- a/tests/QS/regtest-tddfpt-4/test13.inp +++ b/tests/QS/regtest-tddfpt-4/test13.inp @@ -22,9 +22,6 @@ EXCH_SCALING_MODEL NONE METHOD BASIS_PROJECTION &END - &EXCITED_STATES T - STATE 1 - &END EXCITED_STATES &SCF SCF_GUESS ATOMIC &OT @@ -45,10 +42,10 @@ SCALE_X 0.75 SCALE_C 1.0 &END - &PBE_HOLE_T_C_LR - CUTOFF_RADIUS 2.0 - SCALE_X 0.25 - &END + #&PBE_HOLE_T_C_LR + # CUTOFF_RADIUS 2.0 + # SCALE_X 0.25 + #&END &END XC_FUNCTIONAL &HF &SCREENING @@ -110,5 +107,5 @@ &GLOBAL PRINT_LEVEL LOW PROJECT ftest - RUN_TYPE ENERGY_FORCE + RUN_TYPE ENERGY &END GLOBAL diff --git a/tests/QS/regtest-tddfpt-4/test14.inp b/tests/QS/regtest-tddfpt-4/test14.inp index 3f1c95144c..bea4c49e30 100644 --- a/tests/QS/regtest-tddfpt-4/test14.inp +++ b/tests/QS/regtest-tddfpt-4/test14.inp @@ -13,9 +13,6 @@ &QS METHOD GPW &END QS - &EXCITED_STATES T - STATE 1 - &END EXCITED_STATES &SCF SCF_GUESS ATOMIC &OT @@ -36,10 +33,10 @@ SCALE_X 0.75 SCALE_C 1.0 &END - &PBE_HOLE_T_C_LR - CUTOFF_RADIUS 2.0 - SCALE_X 0.25 - &END + #&PBE_HOLE_T_C_LR + # CUTOFF_RADIUS 2.0 + # SCALE_X 0.25 + #&END &END XC_FUNCTIONAL &HF &SCREENING @@ -99,5 +96,5 @@ &GLOBAL PRINT_LEVEL LOW PROJECT ftest - RUN_TYPE ENERGY_FORCE + RUN_TYPE ENERGY &END GLOBAL diff --git a/tests/QS/regtest-tddfpt-4/test23.inp b/tests/QS/regtest-tddfpt-4/test23.inp new file mode 100644 index 0000000000..7f7803a1e7 --- /dev/null +++ b/tests/QS/regtest-tddfpt-4/test23.inp @@ -0,0 +1,166 @@ +# GPW State Excitation Transition dipole (a.u.) Oscillator +# number energy (eV) x y z strength (a.u.) +# ------------------------------------------------------------------------ +# TDDFPT| 1 9.60596 3.6419E-01 7.4416E-08 3.0317E-07 3.12150E-02 +# TDDFPT| 2 11.60274 7.7372E-08 -2.9399E-07 -1.5821E-07 3.33862E-14 +# TDDFPT| 3 12.01143 -1.4513E-07 7.2795E-08 5.7194E-01 9.62634E-02 +# TDDFPT| 4 14.53057 -3.0133E-08 3.3025E-01 -3.5556E-08 3.88271E-02 +# TDDFPT| 5 15.98254 3.9952E-08 -8.2281E-01 -1.0598E-09 2.65093E-01 + +# GAPW State Excitation Transition dipole (a.u.) Oscillator +# number energy (eV) x y z strength (a.u.) +# ------------------------------------------------------------------------ +# TDDFPT| 1 9.64373 -3.6482E-01 -7.6782E-08 -2.5542E-07 3.14463E-02 +# TDDFPT| 2 11.69324 -7.8670E-08 2.9597E-07 1.6949E-07 3.50973E-14 +# TDDFPT| 3 12.06353 4.5523E-08 -7.1777E-08 -5.7378E-01 9.73024E-02 +# TDDFPT| 4 14.58554 2.9865E-08 -3.0556E-01 8.4459E-08 3.33638E-02 +# TDDFPT| 5 16.03568 -4.1063E-08 8.3313E-01 2.7119E-07 2.72694E-01 + +# GPW State Excitation Transition dipole (a.u.) Oscillator +# ADMM number energy (eV) x y z strength (a.u.) +# ------------------------------------------------------------------------ +# TDDFPT| 1 9.96970 -3.5079E-01 -4.6816E-08 -7.6533E-08 3.00560E-02 +# TDDFPT| 2 11.99625 -1.4887E-07 1.2551E-07 1.6072E-07 1.87344E-14 +# TDDFPT| 3 12.30662 -4.5048E-08 1.3854E-08 5.6519E-01 9.63130E-02 +# TDDFPT| 4 14.85514 -5.8706E-09 3.1351E-01 -3.8501E-07 3.57708E-02 +# TDDFPT| 5 16.32599 -1.5250E-08 8.3713E-01 1.9368E-07 2.80301E-01 + +# GAPW State Excitation Transition dipole (a.u.) Oscillator +# ADMM number energy (eV) x y z strength (a.u.) +# ------------------------------------------------------------------------ +# TDDFPT| 1 10.00775 -3.5145E-01 -4.7990E-08 -8.5368E-08 3.02854E-02 +# TDDFPT| 2 12.08890 1.4932E-07 -1.2246E-07 -2.6974E-07 3.25944E-14 +# TDDFPT| 3 12.35959 -4.2690E-08 1.5301E-08 5.6713E-01 9.73934E-02 +# TDDFPT| 4 14.90739 4.9352E-09 -2.8806E-01 -2.7248E-07 3.03068E-02 +# TDDFPT| 5 16.38355 1.5262E-08 -8.4725E-01 5.9801E-07 2.88128E-01 + +# GPW State Excitation Transition dipole (a.u.) Oscillator +# ADMM number energy (eV) x y z strength (a.u.) +# PBEX ------------------------------------------------------------------------ +# TDDFPT| 1 9.95033 3.6560E-01 7.0444E-08 1.5169E-07 3.25850E-02 +# TDDFPT| 2 12.01404 1.9510E-07 -1.8384E-07 -2.3564E-07 3.74951E-14 +# TDDFPT| 3 12.29508 -6.7890E-08 9.4202E-08 5.5881E-01 9.40633E-02 +# TDDFPT| 4 14.87922 9.4347E-09 -3.1354E-01 3.0371E-07 3.58374E-02 +# TDDFPT| 5 16.35766 -2.2981E-08 8.2175E-01 1.5810E-07 2.70618E-01 + +# GAPW State Excitation Transition dipole (a.u.) Oscillator +# ADMM number energy (eV) x y z strength (a.u.) +# PBEX ------------------------------------------------------------------------ +# TDDFPT| 1 9.94584 -3.6580E-01 -6.3862E-08 6.1737E-07 3.26060E-02 +# TDDFPT| 2 12.04251 2.0416E-07 -1.6553E-07 -1.1940E-07 2.45874E-14 +# TDDFPT| 3 12.32143 1.5420E-07 -1.0626E-07 -5.5894E-01 9.43090E-02 +# TDDFPT| 4 14.89843 -9.7621E-09 3.0046E-01 -6.0159E-07 3.29514E-02 +# TDDFPT| 5 16.38306 3.1444E-08 -8.2802E-01 1.1774E-08 2.75190E-01 + +&FORCE_EVAL + METHOD Quickstep + &PROPERTIES + &TDDFPT + KERNEL FULL + ADMM_KERNEL_CORRECTION_SYMMETRIC + NSTATES 5 + MAX_ITER 50 + CONVERGENCE [eV] 1.0e-7 + RKS_TRIPLETS F + &END TDDFPT + &END PROPERTIES + &DFT + &QS + METHOD GAPW + &END QS + BASIS_SET_FILE_NAME BASIS_SET + BASIS_SET_FILE_NAME BASIS_ADMM + &AUXILIARY_DENSITY_MATRIX_METHOD + ADMM_PURIFICATION_METHOD NONE + EXCH_CORRECTION_FUNC PBEX + EXCH_SCALING_MODEL NONE + METHOD BASIS_PROJECTION + &END + &SCF + SCF_GUESS ATOMIC + &OT + PRECONDITIONER FULL_SINGLE_INVERSE + MINIMIZER DIIS + &END + &OUTER_SCF + MAX_SCF 10 + EPS_SCF 1.0E-6 + &END + MAX_SCF 50 + EPS_SCF 1.0E-6 + &END SCF + + &XC + &XC_FUNCTIONAL + &BECKE88 + SCALE_X 0.95238 + &END + &BECKE88_LR + OMEGA 0.33 + SCALE_X -0.94979 + &END + &LYP + SCALE_C 1.0 + &END + &XALPHA + SCALE_X -0.13590 + &END + &END XC_FUNCTIONAL + &HF + &SCREENING + EPS_SCHWARZ 1.0E-7 + &END + &MEMORY + MAX_MEMORY 100 + &END + &INTERACTION_POTENTIAL + POTENTIAL_TYPE MIX_CL_TRUNC + OMEGA 0.33 + SCALE_LONGRANGE 0.94979 + SCALE_COULOMB 0.18352 + CUTOFF_RADIUS 2.5 + T_C_G_DATA t_c_g.dat + &END + &END + &END XC + + &MGRID + CUTOFF 200 + REL_CUTOFF 40 + &END + &POISSON + PERIODIC NONE + POISSON_SOLVER MT + &END + &END DFT + &SUBSYS + &CELL + ABC [angstrom] 4.0 4.0 4.0 + PERIODIC NONE + &END + &COORD + O 0.000000 0.000000 0.000000 + H 0.000000 -0.757136 0.580545 + H 0.000000 0.757136 0.580545 + &END COORD + &TOPOLOGY + &CENTER_COORDINATES + &END + &END + &KIND H + BASIS_SET DZV-GTH-PADE + BASIS_SET AUX_FIT FIT3 + POTENTIAL GTH-PADE-q1 + &END KIND + &KIND O + BASIS_SET DZVP-GTH-PADE + BASIS_SET AUX_FIT FIT3 + POTENTIAL GTH-PADE-q6 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PRINT_LEVEL LOW + PROJECT ftest + RUN_TYPE ENERGY +&END GLOBAL diff --git a/tests/QS/regtest-tddfpt/H2O_GAPW_3.inp b/tests/QS/regtest-tddfpt/H2O_GAPW_3.inp index 4a34691b3c..b4f04ece76 100644 --- a/tests/QS/regtest-tddfpt/H2O_GAPW_3.inp +++ b/tests/QS/regtest-tddfpt/H2O_GAPW_3.inp @@ -29,6 +29,7 @@ METHOD Quickstep &PROPERTIES &TDDFPT + ##ADMM_KERNEL_CORRECTION_SYMMETRIC NSTATES 3 MAX_ITER 10 MAX_KV 10 @@ -47,7 +48,7 @@ &END QS &AUXILIARY_DENSITY_MATRIX_METHOD METHOD basis_projection - ADMM_PURIFICATION_METHOD none + ADMM_PURIFICATION_METHOD NONE EXCH_CORRECTION_FUNC NONE &END &SCF