From b501a6d18ce4c0fefa4f040b40e8f37869dec7e8 Mon Sep 17 00:00:00 2001 From: Thomas Kuehne Date: Sun, 5 Jul 2026 15:20:56 +0200 Subject: [PATCH] Fix GauXC GAPW pseudopotential corrections Enable GauXC MODEL NONE for pseudopotential GAPW/GAPW_XC one-center corrections through mapped LibXC evaluations, skip empty-core GTH one-center terms, and keep the conventional GauXC RKS density normalization intact. Refresh native-grid GAPW GTH/ECP regression references after the force, virial, and k-point fixes. Checks on spark: cmake --build build-host -j20; QS/regtest-gauxc-api 7/7; QS/regtest-gauxc-gapw-gth 8/8; broad pre-format smoke QS/regtest-gauxc-* 113/113. --- src/qs_ks_methods.F | 27 ++-- src/qs_vxc_atom.F | 20 ++- src/torch_c_api.cpp | 11 ++ src/xc/xc_derivatives.F | 147 +++++++++++++++++- src/xc/xc_gauxc_functional.F | 7 - src/xc/xc_gauxc_interface.F | 4 - src/xc/xc_libxc.F | 66 ++++++-- .../QS/regtest-gauxc-gapw-ecp/TEST_FILES.toml | 10 +- .../regtest-gauxc-gapw-gth-kp/TEST_FILES.toml | 6 +- .../QS/regtest-gauxc-gapw-gth/TEST_FILES.toml | 6 +- 10 files changed, 251 insertions(+), 53 deletions(-) diff --git a/src/qs_ks_methods.F b/src/qs_ks_methods.F index 652655dae1..9f6a420799 100644 --- a/src/qs_ks_methods.F +++ b/src/qs_ks_methods.F @@ -143,6 +143,7 @@ MODULE qs_ks_methods get_gauxc_section,& xc_section_uses_native_skala_grid USE smeagol_interface, ONLY: smeagol_shift_v_hartree + USE string_utilities, ONLY: uppercase USE surface_dipole, ONLY: calc_dipsurf_potential USE tblite_ks_matrix, ONLY: build_tblite_ks_matrix USE virial_types, ONLY: virial_type @@ -199,14 +200,15 @@ CONTAINS CHARACTER(LEN=*), PARAMETER :: routineN = 'qs_ks_build_kohn_sham_matrix' - CHARACTER(len=default_string_length) :: name + CHARACTER(len=default_string_length) :: gauxc_model_name, name INTEGER :: ace_rebuild_frequency, atom_a, handle, & iatom, ikind, img, ispin, natom, & nimages, nspins, output_unit INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_of_kind, kind_of LOGICAL :: ace_active, do_adiabatic_rescaling, do_ddapc, do_hfx, do_kpoints, do_ppl, dokp, & - gapw, gapw_xc, just_energy_xc, lrigpw, my_print, native_grid_diagnostics, & - native_grid_use_cuda, native_skala_restore_exc, rigpw, use_gauxc_matrix, use_virial + gapw, gapw_xc, gauxc_model_none, just_energy_xc, lrigpw, my_print, & + native_grid_diagnostics, native_grid_use_cuda, native_skala_restore_exc, rigpw, & + use_gauxc_matrix, use_virial LOGICAL, SAVE :: native_grid_cpu_kpoints_warned = .FALSE. REAL(KIND=dp) :: ecore_ppl, edisp, ee_ener, ekin_mol, & mulliken_order_p, & @@ -922,12 +924,19 @@ CONTAINS NULLIFY (rho1) IF (dft_control%use_gauxc .AND. (gapw .OR. gapw_xc) .AND. & .NOT. xc_section_uses_native_skala_grid(xc_section)) THEN - ! Molecular GauXC evaluates the XC term outside xc_derivatives. - ! The accurate-XCINT force correction would otherwise try to - ! evaluate the GAUXC section through CP2K's local functional path. - ! Native-grid SKALA can provide this correction through the - ! CP2K grid path and needs it for GAPW_XC force consistency. - CONTINUE + gauxc_model_none = .FALSE. + gauxc_section => get_gauxc_section(xc_section) + IF (ASSOCIATED(gauxc_section)) THEN + CALL section_vals_val_get(gauxc_section, "MODEL", c_val=gauxc_model_name) + gauxc_model_name = ADJUSTL(gauxc_model_name) + CALL uppercase(gauxc_model_name) + gauxc_model_none = (TRIM(gauxc_model_name) == "" .OR. & + TRIM(gauxc_model_name) == "NONE") + END IF + IF (gauxc_model_none .AND. & + (gapw_xc .OR. gauxc_gapw_has_paw_pseudopotentials(qs_kind_set))) THEN + CALL accint_weight_force(qs_env, rho_struct, rho1, 0, xc_section) + END IF ELSE CALL accint_weight_force(qs_env, rho_struct, rho1, 0, xc_section) END IF diff --git a/src/qs_vxc_atom.F b/src/qs_vxc_atom.F index 77ca466195..84b2bd4516 100644 --- a/src/qs_vxc_atom.F +++ b/src/qs_vxc_atom.F @@ -15,6 +15,8 @@ MODULE qs_vxc_atom USE basis_set_types, ONLY: get_gto_basis_set,& gto_basis_set_type USE cp_control_types, ONLY: dft_control_type + USE external_potential_types, ONLY: gth_potential_type,& + sgp_potential_type USE input_constants, ONLY: xc_none USE input_section_types, ONLY: section_vals_get_subs_vals,& section_vals_type,& @@ -120,14 +122,16 @@ CONTAINS INTEGER :: bo(2), gapw_density_partition, handle, & iat, iatom, idir, ikind, ir, jdir, & - myfun, na, natom, nr, nspins, num_pe + myfun, na, natom, nr, nspins, num_pe, & + zatom INTEGER, DIMENSION(2, 3) :: bounds INTEGER, DIMENSION(:), POINTER :: atom_list LOGICAL :: accint, donlcc, evaluate_hard, evaluate_soft, gradient_f, lsd, & my_calculate_forces, nlcc, paw_atom, skala_atom_grid, tau_f, use_virial REAL(dp) :: agr, alpha, density_cut, exc_h, exc_s, & gradient_cut, & - my_adiabatic_rescale_factor, tau_cut + my_adiabatic_rescale_factor, tau_cut, & + zeff REAL(dp), DIMENSION(1, 1, 1) :: tau_d REAL(dp), DIMENSION(1, 1, 1, 1) :: rho_d REAL(dp), DIMENSION(3) :: skala_atom_force_h, skala_atom_force_s @@ -140,6 +144,7 @@ CONTAINS TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set TYPE(dft_control_type), POINTER :: dft_control TYPE(grid_atom_type), POINTER :: grid_atom + TYPE(gth_potential_type), POINTER :: gth_potential TYPE(gto_basis_set_type), POINTER :: basis_1c TYPE(harmonics_atom_type), POINTER :: harmonics TYPE(mp_para_env_type), POINTER :: para_env @@ -151,6 +156,7 @@ CONTAINS TYPE(rho_atom_type), DIMENSION(:), POINTER :: my_rho_atom_set TYPE(rho_atom_type), POINTER :: rho_atom TYPE(section_vals_type), POINTER :: input, my_xc_section, xc_fun_section + TYPE(sgp_potential_type), POINTER :: sgp_potential TYPE(tau_basis_cache_type) :: tau_basis_cache TYPE(virial_type), POINTER :: virial TYPE(xc_derivative_set_type) :: deriv_set @@ -165,6 +171,7 @@ CONTAINS NULLIFY (my_kind_set) NULLIFY (atomic_kind_set) NULLIFY (grid_atom) + NULLIFY (gth_potential) NULLIFY (force) NULLIFY (harmonics) NULLIFY (input) @@ -173,6 +180,7 @@ CONTAINS NULLIFY (rho_atom) NULLIFY (my_rho_atom_set) NULLIFY (rho_nlcc) + NULLIFY (sgp_potential) NULLIFY (virial) my_calculate_forces = .FALSE. IF (PRESENT(calculate_forces)) my_calculate_forces = calculate_forces @@ -250,11 +258,17 @@ CONTAINS DO ikind = 1, SIZE(atomic_kind_set) CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=natom) + NULLIFY (gth_potential, sgp_potential) CALL get_qs_kind(my_kind_set(ikind), paw_atom=paw_atom, & - harmonics=harmonics, grid_atom=grid_atom) + gth_potential=gth_potential, harmonics=harmonics, & + grid_atom=grid_atom, sgp_potential=sgp_potential, & + zatom=zatom, zeff=zeff) CALL get_qs_kind(my_kind_set(ikind), basis_set=basis_1c, basis_type="GAPW_1C") IF (.NOT. paw_atom) CYCLE + IF (skala_atom_grid .AND. & + (ASSOCIATED(gth_potential) .OR. ASSOCIATED(sgp_potential)) .AND. & + ABS(zeff - REAL(zatom, dp)) <= 1.0E-10_dp) CYCLE nr = grid_atom%nr na = grid_atom%ng_sphere diff --git a/src/torch_c_api.cpp b/src/torch_c_api.cpp index 9349653ebf..4b09327cf1 100644 --- a/src/torch_c_api.cpp +++ b/src/torch_c_api.cpp @@ -7,6 +7,7 @@ #if defined(__LIBTORCH) +#include #include #include #include @@ -47,6 +48,16 @@ private: ******************************************************************************/ static bool use_cuda_if_available = true; +class TorchCpuThreadGuard { +public: + TorchCpuThreadGuard() { + at::set_num_threads(1); + at::set_num_interop_threads(1); + } +}; + +static TorchCpuThreadGuard torch_cpu_thread_guard; + static torch::Device get_device() { if (!use_cuda_if_available || !torch::cuda::is_available()) { return torch::kCPU; diff --git a/src/xc/xc_derivatives.F b/src/xc/xc_derivatives.F index 431bb25b6e..6121e69d4d 100644 --- a/src/xc/xc_derivatives.F +++ b/src/xc/xc_derivatives.F @@ -11,7 +11,9 @@ MODULE xc_derivatives USE input_section_types, ONLY: section_vals_get_subs_vals2,& section_vals_type,& section_vals_val_get - USE kinds, ONLY: dp + USE kinds, ONLY: default_string_length,& + dp + USE string_utilities, ONLY: uppercase USE xc_b97, ONLY: b97_lda_eval,& b97_lda_info,& b97_lsd_eval,& @@ -281,9 +283,14 @@ CONTAINS CALL pbe_lda_info(functional, reference, shortform, needs, max_deriv) END IF CASE ("GAUXC") - CALL skala_info(functional, lsd, reference, shortform, needs, max_deriv) - ! Note: SKALA functional routes through apply_gauxc in qs_ks_methods.F - ! when USE_GAUXC = .TRUE. (requires dft_control%use_gauxc to be set) + IF (gauxc_model_none_selected(functional)) THEN + CALL gauxc_model_none_xc_info(functional, lsd, reference, shortform, & + needs, max_deriv, print_warn) + ELSE + CALL skala_info(functional, lsd, reference, shortform, needs, max_deriv) + ! Note: SKALA functional routes through apply_gauxc in qs_ks_methods.F + ! when USE_GAUXC = .TRUE. (requires dft_control%use_gauxc to be set) + END IF CASE ("XWPBE") IF (lsd) THEN CALL xwpbe_lsd_info(reference, shortform, needs, max_deriv) @@ -507,7 +514,11 @@ CONTAINS CALL pbe_lda_eval(rho_set, deriv_set, deriv_order, functional) END IF CASE ("GAUXC") - CPABORT(abort_message_skala) + IF (gauxc_model_none_selected(functional)) THEN + CALL gauxc_model_none_xc_eval(functional, lsd, rho_set, deriv_set, deriv_order) + ELSE + CPABORT(abort_message_skala) + END IF CASE ("XWPBE") IF (lsd) THEN CALL xwpbe_lsd_eval(rho_set, deriv_set, deriv_order, functional) @@ -552,6 +563,132 @@ CONTAINS CALL timestop(handle) END SUBROUTINE xc_functional_eval +! ************************************************************************************************** +!> \brief true for GAUXC sections that wrap a conventional LibXC functional +!> \param functional the GAUXC section +!> \return whether MODEL NONE is active +! ************************************************************************************************** + FUNCTION gauxc_model_none_selected(functional) + TYPE(section_vals_type), POINTER :: functional + LOGICAL :: gauxc_model_none_selected + + CHARACTER(LEN=default_string_length) :: model_key, model_name + + CALL section_vals_val_get(functional, "MODEL", c_val=model_name) + model_key = ADJUSTL(model_name) + CALL uppercase(model_key) + gauxc_model_none_selected = (TRIM(model_key) == "" .OR. TRIM(model_key) == "NONE") + END FUNCTION gauxc_model_none_selected + +! ************************************************************************************************** +!> \brief map GAUXC MODEL NONE shorthand names to LibXC exchange/correlation components +!> \param functional the GAUXC section +!> \param xc_fun_name the GAUXC FUNCTIONAL value +!> \param libxc_names LibXC section names to evaluate and add +!> \param nfunc number of LibXC components +! ************************************************************************************************** + SUBROUTINE gauxc_model_none_libxc_names(functional, xc_fun_name, libxc_names, nfunc) + TYPE(section_vals_type), POINTER :: functional + CHARACTER(LEN=*), INTENT(OUT) :: xc_fun_name + CHARACTER(LEN=*), DIMENSION(:), INTENT(OUT) :: libxc_names + INTEGER, INTENT(OUT) :: nfunc + + CHARACTER(LEN=default_string_length) :: xc_key + + CALL section_vals_val_get(functional, "FUNCTIONAL", c_val=xc_fun_name) + xc_key = ADJUSTL(xc_fun_name) + CALL uppercase(xc_key) + libxc_names(:) = "" + SELECT CASE (TRIM(xc_key)) + CASE ("PBE") + nfunc = 2 + libxc_names(1) = "GGA_X_PBE" + libxc_names(2) = "GGA_C_PBE" + CASE ("TPSS") + nfunc = 2 + libxc_names(1) = "MGGA_X_TPSS" + libxc_names(2) = "MGGA_C_TPSS" + CASE ("R2SCAN") + nfunc = 2 + libxc_names(1) = "MGGA_X_R2SCAN" + libxc_names(2) = "MGGA_C_R2SCAN" + CASE default + CPABORT("GAUXC MODEL NONE one-center GAPW corrections currently know PBE, TPSS, and R2SCAN.") + END SELECT + END SUBROUTINE gauxc_model_none_libxc_names + +! ************************************************************************************************** +!> \brief needs information for GAUXC MODEL NONE one-center GAPW corrections +!> \param functional the GAUXC section +!> \param lsd whether spin-polarized derivatives are needed +!> \param reference reference string for the wrapped functional +!> \param shortform short name for printout +!> \param needs density ingredients needed by the wrapped functional +!> \param max_deriv maximum implemented derivative order +!> \param print_warn whether LibXC should print development warnings +! ************************************************************************************************** + SUBROUTINE gauxc_model_none_xc_info(functional, lsd, reference, shortform, & + needs, max_deriv, print_warn) + TYPE(section_vals_type), POINTER :: functional + LOGICAL, INTENT(in) :: lsd + CHARACTER(LEN=*), INTENT(OUT), OPTIONAL :: reference, shortform + TYPE(xc_rho_cflags_type), INTENT(inout), OPTIONAL :: needs + INTEGER, INTENT(out), OPTIONAL :: max_deriv + LOGICAL, INTENT(IN), OPTIONAL :: print_warn + + CHARACTER(LEN=default_string_length) :: libxc_names(2), xc_fun_name + INTEGER :: ifunc, max_deriv_i, max_deriv_min, nfunc + + CALL gauxc_model_none_libxc_names(functional, xc_fun_name, libxc_names, nfunc) + max_deriv_min = HUGE(max_deriv_min) + DO ifunc = 1, nfunc + IF (lsd) THEN + CALL libxc_lsd_info(functional, needs=needs, max_deriv=max_deriv_i, & + print_warn=print_warn, & + func_name_override=TRIM(libxc_names(ifunc))) + ELSE + CALL libxc_lda_info(functional, needs=needs, max_deriv=max_deriv_i, & + print_warn=print_warn, & + func_name_override=TRIM(libxc_names(ifunc))) + END IF + max_deriv_min = MIN(max_deriv_min, max_deriv_i) + END DO + IF (PRESENT(max_deriv)) max_deriv = max_deriv_min + IF (PRESENT(reference)) & + reference = "Functional computed by GauXC (underlying: "//TRIM(xc_fun_name)//")" + IF (PRESENT(shortform)) shortform = "GAUXC ("//TRIM(xc_fun_name)//")" + END SUBROUTINE gauxc_model_none_xc_info + +! ************************************************************************************************** +!> \brief one-center GAPW correction evaluation for GAUXC MODEL NONE +!> \param functional the GAUXC section +!> \param lsd whether spin-polarized derivatives are evaluated +!> \param rho_set density ingredients +!> \param deriv_set derivative accumulator +!> \param deriv_order derivative order to evaluate +! ************************************************************************************************** + SUBROUTINE gauxc_model_none_xc_eval(functional, lsd, rho_set, deriv_set, deriv_order) + TYPE(section_vals_type), POINTER :: functional + LOGICAL, INTENT(in) :: lsd + TYPE(xc_rho_set_type), INTENT(IN) :: rho_set + TYPE(xc_derivative_set_type), INTENT(IN) :: deriv_set + INTEGER, INTENT(IN) :: deriv_order + + CHARACTER(LEN=default_string_length) :: libxc_names(2), xc_fun_name + INTEGER :: ifunc, nfunc + + CALL gauxc_model_none_libxc_names(functional, xc_fun_name, libxc_names, nfunc) + DO ifunc = 1, nfunc + IF (lsd) THEN + CALL libxc_lsd_eval(rho_set, deriv_set, deriv_order, functional, & + func_name_override=TRIM(libxc_names(ifunc))) + ELSE + CALL libxc_lda_eval(rho_set, deriv_set, deriv_order, functional, & + func_name_override=TRIM(libxc_names(ifunc))) + END IF + END DO + END SUBROUTINE gauxc_model_none_xc_eval + ! ************************************************************************************************** !> \brief ... !> \param functionals a section containing the functional combination to be diff --git a/src/xc/xc_gauxc_functional.F b/src/xc/xc_gauxc_functional.F index 3c12631d2f..09d669e818 100644 --- a/src/xc/xc_gauxc_functional.F +++ b/src/xc/xc_gauxc_functional.F @@ -1146,13 +1146,6 @@ CONTAINS CALL uppercase(gradient_runtime_key) use_gauxc_model = (TRIM(model_key) /= "" .AND. TRIM(model_key) /= "NONE") use_skala_model = (INDEX(TRIM(model_key), "SKALA") > 0) - IF (gapw_pseudopotentials .AND. .NOT. use_gauxc_model) THEN - CALL cp_abort(__LOCATION__, & - "GauXC with METHOD GAPW/GAPW_XC and pseudopotentials is supported only for "// & - "Skala-style models that replace the molecular XC term. "// & - "Use POTENTIAL ALL for local/semi-local GauXC GAPW validation or METHOD GPW "// & - "with pseudopotentials.") - END IF IF (gapw_pseudopotentials .AND. use_gauxc_model .AND. .NOT. dft_control%qs_control%gapw_xc .AND. & .NOT. gapw_paw_pseudopotentials .AND. para_env%mepos == 0 .AND. ASSOCIATED(scf_env)) THEN IF (scf_env%iter_count == 1) THEN diff --git a/src/xc/xc_gauxc_interface.F b/src/xc/xc_gauxc_interface.F index ed1356ebde..df71596bdb 100644 --- a/src/xc/xc_gauxc_interface.F +++ b/src/xc/xc_gauxc_interface.F @@ -809,10 +809,6 @@ CONTAINS END IF IF (nspins == 1) THEN - ! xmat factor 2 is applied by both CP2K and GauXC - ! "unapply" it here to even things back out. - ! This is NOT necessary in the Skala branch. - density_scalar = 0.5_dp*density_scalar CALL gauxc_integrator_eval_exc_vxc_rks( & status%status, & integrator%integrator, & diff --git a/src/xc/xc_libxc.F b/src/xc/xc_libxc.F index 3faed929f3..4c58f994a0 100644 --- a/src/xc/xc_libxc.F +++ b/src/xc/xc_libxc.F @@ -307,9 +307,11 @@ CONTAINS !> true (does not set the unneeded components to false) !> \param max_deriv maximum implemented derivative of the xc functional !> \param print_warn whether to print warning about development status of a functional +!> \param func_name_override optional LibXC functional name overriding the section name !> \author F. Tran ! ************************************************************************************************** - SUBROUTINE libxc_lda_info(libxc_params, reference, shortform, needs, max_deriv, print_warn) + SUBROUTINE libxc_lda_info(libxc_params, reference, shortform, needs, max_deriv, print_warn, & + func_name_override) TYPE(section_vals_type), POINTER :: libxc_params CHARACTER(LEN=*), INTENT(OUT), OPTIONAL :: reference, shortform @@ -317,6 +319,7 @@ CONTAINS INTENT(inout), OPTIONAL :: needs INTEGER, INTENT(out), OPTIONAL :: max_deriv LOGICAL, INTENT(IN), OPTIONAL :: print_warn + CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: func_name_override #if defined (__LIBXC) CHARACTER(LEN=128) :: s1, s2 @@ -326,8 +329,13 @@ CONTAINS TYPE(xc_f03_func_t) :: xc_func TYPE(xc_f03_func_info_t) :: xc_info - func_name = libxc_params%section%name - CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + IF (PRESENT(func_name_override)) THEN + func_name = func_name_override + func_scale = 1.0_dp + ELSE + func_name = libxc_params%section%name + CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + END IF CALL cite_reference(Marques2012) CALL cite_reference(Lehtola2018) @@ -398,6 +406,7 @@ CONTAINS MARK_USED(needs) MARK_USED(max_deriv) MARK_USED(print_warn) + MARK_USED(func_name_override) CALL cp_abort(__LOCATION__, "Unknown functional! If you are asking "// & "for a functional of the LibXC library, "// & @@ -415,9 +424,11 @@ CONTAINS !> true (does not set the unneeded components to false) !> \param max_deriv maximum implemented derivative of the xc functional !> \param print_warn whether to print warning about development status of a functional +!> \param func_name_override optional LibXC functional name overriding the section name !> \author F. Tran ! ************************************************************************************************** - SUBROUTINE libxc_lsd_info(libxc_params, reference, shortform, needs, max_deriv, print_warn) + SUBROUTINE libxc_lsd_info(libxc_params, reference, shortform, needs, max_deriv, print_warn, & + func_name_override) TYPE(section_vals_type), POINTER :: libxc_params CHARACTER(LEN=*), INTENT(OUT), OPTIONAL :: reference, shortform @@ -425,6 +436,7 @@ CONTAINS INTENT(inout), OPTIONAL :: needs INTEGER, INTENT(out), OPTIONAL :: max_deriv LOGICAL, INTENT(IN), OPTIONAL :: print_warn + CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: func_name_override #if defined (__LIBXC) CHARACTER(LEN=128) :: s1, s2 @@ -434,8 +446,13 @@ CONTAINS TYPE(xc_f03_func_t) :: xc_func TYPE(xc_f03_func_info_t) :: xc_info - func_name = libxc_params%section%name - CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + IF (PRESENT(func_name_override)) THEN + func_name = func_name_override + func_scale = 1.0_dp + ELSE + func_name = libxc_params%section%name + CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + END IF CALL cite_reference(Marques2012) CALL cite_reference(Lehtola2018) @@ -508,6 +525,7 @@ CONTAINS MARK_USED(needs) MARK_USED(max_deriv) MARK_USED(print_warn) + MARK_USED(func_name_override) CALL cp_abort(__LOCATION__, "Unknown functional! If you are "// & "asking for a functional of the LibXC library, "// & @@ -542,14 +560,16 @@ CONTAINS !> if positive all the derivatives up to the given degree are evaluated, !> if negative only the given degree is calculated !> \param libxc_params input parameter (functional name, scaling and parameters) +!> \param func_name_override optional LibXC functional name overriding the section name !> \author F. Tran ! ************************************************************************************************** - SUBROUTINE libxc_lda_eval(rho_set, deriv_set, grad_deriv, libxc_params) + SUBROUTINE libxc_lda_eval(rho_set, deriv_set, grad_deriv, libxc_params, func_name_override) TYPE(xc_rho_set_type), INTENT(IN) :: rho_set TYPE(xc_derivative_set_type), INTENT(IN) :: deriv_set INTEGER, INTENT(in) :: grad_deriv TYPE(section_vals_type), POINTER :: libxc_params + CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: func_name_override #if defined (__LIBXC) CHARACTER(len=*), PARAMETER :: routineN = 'libxc_lda_eval' @@ -574,15 +594,22 @@ CONTAINS NULLIFY (dummy) NULLIFY (rho, norm_drho, laplace_rho, tau) - func_name = libxc_params%section%name - CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + IF (PRESENT(func_name_override)) THEN + func_name = func_name_override + func_scale = 1.0_dp + ELSE + func_name = libxc_params%section%name + CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + END IF IF (ABS(func_scale - 1.0_dp) < 1.0e-10_dp) func_scale = 1.0_dp func_id = xc_libxc_wrap_functional_get_number(func_name) CALL xc_f03_func_init(xc_func, func_id, XC_UNPOLARIZED) xc_info = xc_f03_func_get_info(xc_func) - CALL xc_libxc_wrap_functional_set_params(xc_func, xc_info, libxc_params, no_exc) + no_exc = .FALSE. + IF (.NOT. PRESENT(func_name_override)) & + CALL xc_libxc_wrap_functional_set_params(xc_func, xc_info, libxc_params, no_exc) CALL xc_rho_set_get(rho_set, can_return_null=.TRUE., & rho=rho, norm_drho=norm_drho, laplace_rho=laplace_rho, & @@ -762,6 +789,7 @@ CONTAINS MARK_USED(deriv_set) MARK_USED(grad_deriv) MARK_USED(libxc_params) + MARK_USED(func_name_override) CALL cp_abort(__LOCATION__, "Unknown functional! If you are asking "// & "for a functional of the LibXC library, "// & "you have to download and install the library!") @@ -777,14 +805,16 @@ CONTAINS !> if positive all the derivatives up to the given degree are evaluated, !> if negative only the given degree is calculated !> \param libxc_params input parameter (functional name, scaling and parameters) +!> \param func_name_override optional LibXC functional name overriding the section name !> \author F. Tran ! ************************************************************************************************** - SUBROUTINE libxc_lsd_eval(rho_set, deriv_set, grad_deriv, libxc_params) + SUBROUTINE libxc_lsd_eval(rho_set, deriv_set, grad_deriv, libxc_params, func_name_override) TYPE(xc_rho_set_type), INTENT(IN) :: rho_set TYPE(xc_derivative_set_type), INTENT(IN) :: deriv_set INTEGER, INTENT(in) :: grad_deriv TYPE(section_vals_type), POINTER :: libxc_params + CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: func_name_override #if defined (__LIBXC) CHARACTER(len=*), PARAMETER :: routineN = 'libxc_lsd_eval' @@ -823,15 +853,22 @@ CONTAINS NULLIFY (rhoa, rhob, norm_drho, norm_drhoa, norm_drhob, laplace_rhoa, & laplace_rhob, tau_a, tau_b) - func_name = libxc_params%section%name - CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + IF (PRESENT(func_name_override)) THEN + func_name = func_name_override + func_scale = 1.0_dp + ELSE + func_name = libxc_params%section%name + CALL section_vals_val_get(libxc_params, "scale", r_val=func_scale) + END IF IF (ABS(func_scale - 1.0_dp) < 1.0e-10_dp) func_scale = 1.0_dp func_id = xc_libxc_wrap_functional_get_number(func_name) CALL xc_f03_func_init(xc_func, func_id, XC_POLARIZED) xc_info = xc_f03_func_get_info(xc_func) - CALL xc_libxc_wrap_functional_set_params(xc_func, xc_info, libxc_params, no_exc) + no_exc = .FALSE. + IF (.NOT. PRESENT(func_name_override)) & + CALL xc_libxc_wrap_functional_set_params(xc_func, xc_info, libxc_params, no_exc) CALL xc_rho_set_get(rho_set, can_return_null=.TRUE., & rhoa=rhoa, rhob=rhob, norm_drho=norm_drho, & @@ -1290,6 +1327,7 @@ CONTAINS MARK_USED(deriv_set) MARK_USED(grad_deriv) MARK_USED(libxc_params) + MARK_USED(func_name_override) CALL cp_abort(__LOCATION__, "Unknown functional! If you are asking "// & "for a functional of the LibXC library, "// & diff --git a/tests/QS/regtest-gauxc-gapw-ecp/TEST_FILES.toml b/tests/QS/regtest-gauxc-gapw-ecp/TEST_FILES.toml index 274511571a..babec12d78 100644 --- a/tests/QS/regtest-gauxc-gapw-ecp/TEST_FILES.toml +++ b/tests/QS/regtest-gauxc-gapw-ecp/TEST_FILES.toml @@ -1,8 +1,8 @@ "HCl_GAPW_SKALA_ECP_ENERGY.inp" = [{matcher="E_total", tol=5e-6, ref=-15.464581508021762}] -"HCl_NATIVE_SKALA_GAPW_ECP_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.45263701}] +"HCl_NATIVE_SKALA_GAPW_ECP_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.44142596263930}] "HCl_NATIVE_SKALA_GAPW_GPWTYPE_ECP_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.36078871}] -"HCl_NATIVE_SKALA_GAPW_XC_ECP_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.36291913}] +"HCl_NATIVE_SKALA_GAPW_XC_ECP_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.35170807520943}] "H2_NATIVE_SKALA_GAPW_ECP_STRESS_DEBUG.inp" = [{matcher="DEBUG_stress_sum", tol=6e-1, ref=0.0}] -"HCl_NATIVE_SKALA_GAPW_ECP_KP_INV_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.44577629}] -"HCl_NATIVE_SKALA_GAPW_ECP_KP_SYM_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.44577629}] -"HCl_NATIVE_SKALA_GAPW_ECP_KP_SPGLIB_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.44577629}] +"HCl_NATIVE_SKALA_GAPW_ECP_KP_INV_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.43466946670853}] +"HCl_NATIVE_SKALA_GAPW_ECP_KP_SYM_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.43466946670853}] +"HCl_NATIVE_SKALA_GAPW_ECP_KP_SPGLIB_STRESS.inp" = [{matcher="E_total", tol=5e-6, ref=-15.43466946670853}] diff --git a/tests/QS/regtest-gauxc-gapw-gth-kp/TEST_FILES.toml b/tests/QS/regtest-gauxc-gapw-gth-kp/TEST_FILES.toml index 80fac5377e..dea40ed44e 100644 --- a/tests/QS/regtest-gauxc-gapw-gth-kp/TEST_FILES.toml +++ b/tests/QS/regtest-gauxc-gapw-gth-kp/TEST_FILES.toml @@ -1,3 +1,3 @@ -"H2_NATIVE_SKALA_GAPW_GTH_KP_INV_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.97063781133388}] -"H2_NATIVE_SKALA_GAPW_GTH_KP_SYM_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.97063781133388}] -"H2_NATIVE_SKALA_GAPW_GTH_KP_SPGLIB_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.97063781133388}] +"H2_NATIVE_SKALA_GAPW_GTH_KP_INV_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.94951222851689}] +"H2_NATIVE_SKALA_GAPW_GTH_KP_SYM_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.94951222851689}] +"H2_NATIVE_SKALA_GAPW_GTH_KP_SPGLIB_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.94951222851689}] diff --git a/tests/QS/regtest-gauxc-gapw-gth/TEST_FILES.toml b/tests/QS/regtest-gauxc-gapw-gth/TEST_FILES.toml index 300b8bb9d2..fbeaab6fca 100644 --- a/tests/QS/regtest-gauxc-gapw-gth/TEST_FILES.toml +++ b/tests/QS/regtest-gauxc-gapw-gth/TEST_FILES.toml @@ -1,8 +1,8 @@ -"H2_NATIVE_SKALA_GAPW_GTH_STRESS.inp" = [{matcher="E_total", tol=5e-8, ref=-0.9912833172}] -"H2_NATIVE_SKALA_GAPW_GTH_STRESS_CHUNK_REQUEST.inp" = [{matcher="E_total", tol=5e-8, ref=-0.9912833169}] +"H2_NATIVE_SKALA_GAPW_GTH_STRESS.inp" = [{matcher="E_total", tol=5e-8, ref=-0.97019845356640}] +"H2_NATIVE_SKALA_GAPW_GTH_STRESS_CHUNK_REQUEST.inp" = [{matcher="E_total", tol=5e-8, ref=-0.97019845356640}] "H2_NATIVE_SKALA_GAPW_GPWTYPE_STRESS.inp" = [{matcher="E_total", tol=1e-8, ref=-0.96567258919119}, {matcher="M072", tol=1e-4, ref=5.37252718E-05}, {matcher="M031", tol=1e-5, ref=-8.58976831594E+04}] -"H2_NATIVE_SKALA_GAPW_XC_GTH_STRESS.inp" = [{matcher="E_total", tol=5e-8, ref=-0.9698629747}] +"H2_NATIVE_SKALA_GAPW_XC_GTH_STRESS.inp" = [{matcher="E_total", tol=5e-8, ref=-0.94877811114650}] "H2_NATIVE_SKALA_GAPW_GTH_FORCE_DEBUG.inp" = [{matcher="DEBUG_force_sum", tol=2e-1, ref=0.0}] "H2_NATIVE_SKALA_GAPW_GTH_STRESS_DEBUG.inp" = [{matcher="DEBUG_stress_sum", tol=2e-1, ref=0.0}]