From 9081568142b3ad6fa66841c2b42b700a58240177 Mon Sep 17 00:00:00 2001 From: "Thomas D. Kuehne" Date: Sun, 24 May 2026 19:36:24 +0200 Subject: [PATCH] Expose spin-summed RKS setting through C APIs --- include/gauxc/c/xc_integrator.h | 71 ++++++++++- src/c-api/c_xc_integrator.cxx | 192 ++++++++++++++++++++++++++++++ src/fortran-api/xc_integrator.f90 | 156 ++++++++++++++++++++++-- tests/c_api_test.cxx | 28 +++++ 4 files changed, 434 insertions(+), 13 deletions(-) diff --git a/include/gauxc/c/xc_integrator.h b/include/gauxc/c/xc_integrator.h index e0a5ecd3..893456f5 100644 --- a/include/gauxc/c/xc_integrator.h +++ b/include/gauxc/c/xc_integrator.h @@ -36,6 +36,14 @@ typedef struct GauXCIntegrator { void* ptr; ///< Pointer to the XCIntegrator instance. } GauXCIntegrator; +/** + * @brief Settings for Kohn-Sham XC integrations. + */ +typedef struct GauXCKSSettings { + double gks_dtol; ///< Density tolerance used by GKS integrations. + bool rks_density_matrix_is_spin_summed; ///< True if an RKS matrix contains the spin-summed density. +} GauXCKSSettings; + /** * @brief Delete an XCIntegrator instance. * @param status Status object to capture any errors. @@ -109,6 +117,20 @@ extern void gauxc_integrator_eval_exc_rks( double* exc ); +/** + * @brief Evaluate the exchange-correlation energy for RKS with explicit settings. + */ +extern void gauxc_integrator_eval_exc_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + const int64_t m, + const int64_t n, + const double* density_matrix, + const int64_t ldp, + const GauXCKSSettings* settings, + double* exc +); + /** * @brief Evaluate the exchange-correlation energy for UKS. * @param status Status object to capture any errors. @@ -189,6 +211,22 @@ extern void gauxc_integrator_eval_exc_vxc_rks( const int64_t vxc_ld ); +/** + * @brief Evaluate the exchange-correlation energy and potential for RKS with explicit settings. + */ +extern void gauxc_integrator_eval_exc_vxc_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + const int64_t m, + const int64_t n, + const double* density_matrix, + const int64_t ldp, + const GauXCKSSettings* settings, + double* exc, + double* vxc_matrix, + const int64_t vxc_ld +); + /** * @brief Evaluate the exchange-correlation energy and potential for UKS. * @param status Status object to capture any errors. @@ -323,6 +361,20 @@ extern void gauxc_integrator_eval_exc_grad_rks( double* exc_grad ); +/** + * @brief Evaluate the exchange-correlation nuclear gradient for RKS with explicit settings. + */ +extern void gauxc_integrator_eval_exc_grad_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + const int64_t m, + const int64_t n, + const double* density_matrix, + const int64_t ldp, + const GauXCKSSettings* settings, + double* exc_grad +); + /** * @brief Evaluate the exchange-correlation energy gradient for UKS. * @param status Status object to capture any errors. @@ -421,6 +473,23 @@ extern void gauxc_integrator_eval_fxc_contraction_rks( const int64_t ldfxc ); +/** + * @brief Evaluate the FXC contraction for RKS with explicit settings. + */ +extern void gauxc_integrator_eval_fxc_contraction_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + const int64_t m, + const int64_t n, + const double* density_matrix, + const int64_t ldp, + const double* t_density_matrix, + const int64_t ltdp, + const GauXCKSSettings* settings, + double* fxc, + const int64_t ldfxc +); + /** * @brief Evaluate the FXC contraction for UKS. * @param status Status object to capture any errors. @@ -462,4 +531,4 @@ extern void gauxc_integrator_eval_fxc_contraction_uks( #ifdef __cplusplus } // namespace GauXC::C } // extern "C" -#endif \ No newline at end of file +#endif diff --git a/src/c-api/c_xc_integrator.cxx b/src/c-api/c_xc_integrator.cxx index 534eb3d0..b091711b 100644 --- a/src/c-api/c_xc_integrator.cxx +++ b/src/c-api/c_xc_integrator.cxx @@ -69,6 +69,31 @@ get_replicated_integrator( const char* integrator_kernel_name, namespace GauXC::C { + +namespace { + +IntegratorSettingsKS convert_ks_settings( const GauXCKSSettings* c_settings ) { + IntegratorSettingsKS settings; + if( c_settings != nullptr ) { + settings.gks_dtol = c_settings->gks_dtol; + settings.rks_density_matrix_is_spin_summed = + c_settings->rks_density_matrix_is_spin_summed; + } + return settings; +} + +IntegratorSettingsEXC_GRAD convert_exc_grad_settings( const GauXCKSSettings* c_settings ) { + IntegratorSettingsEXC_GRAD settings; + if( c_settings != nullptr ) { + settings.gks_dtol = c_settings->gks_dtol; + settings.rks_density_matrix_is_spin_summed = + c_settings->rks_density_matrix_is_spin_summed; + } + return settings; +} + +} // namespace + extern "C" { void gauxc_integrator_delete( @@ -211,6 +236,44 @@ void gauxc_integrator_eval_exc_rks( } } +void gauxc_integrator_eval_exc_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + int64_t m, + int64_t n, + const double* density_matrix, + int64_t ldp, + const GauXCKSSettings* settings, + double* exc +) { + detail::gauxc_status_init(status); + if (integrator.ptr == nullptr || integrator.hdr.type != GauXC_Type_Integrator) { + detail::gauxc_status_handle(status, 1, "Invalid Integrator handle"); + return; + } + if (density_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "Density matrix pointer cannot be null"); + return; + } + if (settings == nullptr) { + detail::gauxc_status_handle(status, 1, "KS settings pointer cannot be null"); + return; + } + if (exc == nullptr) { + detail::gauxc_status_handle(status, 1, "Exc output pointer cannot be null"); + return; + } + try { + detail::get_xc_integrator_ptr(integrator)->eval_exc( + m, n, + density_matrix, ldp, + exc, + convert_ks_settings(settings) ); + } catch (std::exception& e) { + detail::gauxc_status_handle(status, 1, e.what()); + } +} + void gauxc_integrator_eval_exc_uks( GauXCStatus* status, const GauXCIntegrator integrator, @@ -345,6 +408,51 @@ void gauxc_integrator_eval_exc_vxc_rks( } } +void gauxc_integrator_eval_exc_vxc_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + int64_t m, + int64_t n, + const double* density_matrix, + int64_t ldp, + const GauXCKSSettings* settings, + double* exc, + double* vxc_matrix, + int64_t vxc_ld +) { + detail::gauxc_status_init(status); + if (integrator.ptr == nullptr || integrator.hdr.type != GauXC_Type_Integrator) { + detail::gauxc_status_handle(status, 1, "Invalid Integrator handle"); + return; + } + if (density_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "Density matrix pointer cannot be null"); + return; + } + if (settings == nullptr) { + detail::gauxc_status_handle(status, 1, "KS settings pointer cannot be null"); + return; + } + if (exc == nullptr) { + detail::gauxc_status_handle(status, 1, "Exc output pointer cannot be null"); + return; + } + if (vxc_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "VXC matrix pointer cannot be null"); + return; + } + try { + detail::get_xc_integrator_ptr(integrator)->eval_exc_vxc( + m, n, + density_matrix, ldp, + vxc_matrix, vxc_ld, + exc, + convert_ks_settings(settings) ); + } catch (std::exception& e) { + detail::gauxc_status_handle(status, 1, e.what()); + } +} + void gauxc_integrator_eval_exc_vxc_uks( GauXCStatus* status, const GauXCIntegrator integrator, @@ -576,6 +684,44 @@ void gauxc_integrator_eval_exc_grad_rks( } } +void gauxc_integrator_eval_exc_grad_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + int64_t m, + int64_t n, + const double* density_matrix, + int64_t ldp, + const GauXCKSSettings* settings, + double* exc_grad +) { + detail::gauxc_status_init(status); + if (integrator.ptr == nullptr || integrator.hdr.type != GauXC_Type_Integrator) { + detail::gauxc_status_handle(status, 1, "Invalid Integrator handle"); + return; + } + if (density_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "Density matrix pointer cannot be null"); + return; + } + if (settings == nullptr) { + detail::gauxc_status_handle(status, 1, "KS settings pointer cannot be null"); + return; + } + if (exc_grad == nullptr) { + detail::gauxc_status_handle(status, 1, "Exc gradient output pointer cannot be null"); + return; + } + try { + detail::get_xc_integrator_ptr(integrator)->eval_exc_grad( + m, n, + density_matrix, ldp, + exc_grad, + convert_exc_grad_settings(settings) ); + } catch (std::exception& e) { + detail::gauxc_status_handle(status, 1, e.what()); + } +} + void gauxc_integrator_eval_exc_grad_uks( GauXCStatus* status, const GauXCIntegrator integrator, @@ -738,6 +884,52 @@ void gauxc_integrator_eval_fxc_contraction_rks( } } +void gauxc_integrator_eval_fxc_contraction_rks_with_settings( + GauXCStatus* status, + const GauXCIntegrator integrator, + int64_t m, + int64_t n, + const double* density_matrix, + int64_t ldp, + const double* t_density_matrix, + int64_t ldtp, + const GauXCKSSettings* settings, + double* fxc, + int64_t ldfxc +) { + detail::gauxc_status_init(status); + if (integrator.ptr == nullptr || integrator.hdr.type != GauXC_Type_Integrator) { + detail::gauxc_status_handle(status, 1, "Invalid Integrator handle"); + return; + } + if (density_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "Density matrix pointer cannot be null"); + return; + } + if (t_density_matrix == nullptr) { + detail::gauxc_status_handle(status, 1, "Density matrix pointer cannot be null"); + return; + } + if (settings == nullptr) { + detail::gauxc_status_handle(status, 1, "KS settings pointer cannot be null"); + return; + } + if (fxc == nullptr) { + detail::gauxc_status_handle(status, 1, "FXC output pointer cannot be null"); + return; + } + try { + detail::get_xc_integrator_ptr(integrator)->eval_fxc_contraction( + m, n, + density_matrix, ldp, + t_density_matrix, ldtp, + fxc, ldfxc, + convert_ks_settings(settings) ); + } catch (std::exception& e) { + detail::gauxc_status_handle(status, 1, e.what()); + } +} + void gauxc_integrator_eval_fxc_contraction_uks( GauXCStatus* status, const GauXCIntegrator integrator, diff --git a/src/fortran-api/xc_integrator.f90 b/src/fortran-api/xc_integrator.f90 index 9d186580..7fbcd9a0 100644 --- a/src/fortran-api/xc_integrator.f90 +++ b/src/fortran-api/xc_integrator.f90 @@ -50,6 +50,12 @@ module gauxc_integrator type(c_ptr) :: ptr = c_null_ptr end type gauxc_integrator_type + !> @brief Settings for Kohn-Sham XC integrations + type, bind(c), public :: gauxc_ks_settings_type + real(c_double) :: gks_dtol + logical(c_bool) :: rks_density_matrix_is_spin_summed + end type gauxc_ks_settings_type + interface gauxc_delete !> @brief Delete an XCIntegrator instance. subroutine gauxc_integrator_delete(status, integrator) bind(c) @@ -136,6 +142,28 @@ subroutine gauxc_integrator_eval_exc_rks_c(status, integrator, m, n, density_mat !> @param exc Pointer to store the exchange-correlation energy. real(c_double), intent(out) :: exc end subroutine gauxc_integrator_eval_exc_rks_c + subroutine gauxc_integrator_eval_exc_rks_with_settings_c(status, integrator, m, n, density_matrix, & + & ldp, settings, exc) & + & bind(c, name="gauxc_integrator_eval_exc_rks_with_settings") + import :: gauxc_status_type, gauxc_integrator_type, gauxc_ks_settings_type, c_double, c_int64_t + implicit none + !> @param status Status object to capture any errors. + type(gauxc_status_type), intent(inout) :: status + !> @param integrator Handle to the XCIntegrator. + type(gauxc_integrator_type), value :: integrator + !> @param m Number of rows in the density matrix. + integer(c_int64_t), value :: m + !> @param n Number of columns in the density matrix. + integer(c_int64_t), value :: n + !> @param density_matrix Pointer to the density matrix data. + real(c_double), intent(in) :: density_matrix(ldp, *) + !> @param ldp Leading dimension of the density matrix. + integer(c_int64_t), value :: ldp + !> @param settings Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in) :: settings + !> @param exc Pointer to store the exchange-correlation energy. + real(c_double), intent(out) :: exc + end subroutine gauxc_integrator_eval_exc_rks_with_settings_c module procedure gauxc_integrator_eval_exc_rks !> @brief Evaluate the exchange-correlation energy for UKS. @@ -228,6 +256,32 @@ subroutine gauxc_integrator_eval_exc_vxc_rks_c(status, integrator, m, n, & !> @param ldvxc Leading dimension of the potential matrix. integer(c_int64_t), value :: ldvxc end subroutine gauxc_integrator_eval_exc_vxc_rks_c + subroutine gauxc_integrator_eval_exc_vxc_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, settings, exc, vxc_matrix, ldvxc) & + & bind(c, name="gauxc_integrator_eval_exc_vxc_rks_with_settings") + import :: gauxc_status_type, gauxc_integrator_type, gauxc_ks_settings_type, c_double, c_int64_t + implicit none + !> @param status Status object to capture any errors. + type(gauxc_status_type), intent(inout) :: status + !> @param integrator Handle to the XCIntegrator. + type(gauxc_integrator_type), value :: integrator + !> @param m Number of rows in the density matrix. + integer(c_int64_t), value :: m + !> @param n Number of columns in the density matrix. + integer(c_int64_t), value :: n + !> @param density_matrix Pointer to the density matrix data. + real(c_double), intent(in) :: density_matrix(ldp, *) + !> @param ldp Leading dimension of the density matrix. + integer(c_int64_t), value :: ldp + !> @param settings Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in) :: settings + !> @param exc Pointer to store the exchange-correlation energy. + real(c_double), intent(out) :: exc + !> @param vxc_matrix Pointer to the potential matrix data. + real(c_double), intent(out) :: vxc_matrix(ldvxc, *) + !> @param ldvxc Leading dimension of the potential matrix. + integer(c_int64_t), value :: ldvxc + end subroutine gauxc_integrator_eval_exc_vxc_rks_with_settings_c module procedure gauxc_integrator_eval_exc_vxc_rks !> @brief Evaluate the exchange-correlation energy and potential for UKS. @@ -378,6 +432,28 @@ subroutine gauxc_integrator_eval_exc_grad_rks_c(status, integrator, m, n, & !> @param exc_grad Pointer to the gradient. real(c_double), intent(out) :: exc_grad(*) end subroutine gauxc_integrator_eval_exc_grad_rks_c + subroutine gauxc_integrator_eval_exc_grad_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, settings, exc_grad) & + & bind(c, name="gauxc_integrator_eval_exc_grad_rks_with_settings") + import :: gauxc_status_type, gauxc_integrator_type, gauxc_ks_settings_type, c_double, c_int64_t + implicit none + !> @param status Status object to capture any errors. + type(gauxc_status_type), intent(inout) :: status + !> @param integrator Handle to the XCIntegrator. + type(gauxc_integrator_type), value :: integrator + !> @param m Number of rows in the density matrix. + integer(c_int64_t), value :: m + !> @param n Number of columns in the density matrix. + integer(c_int64_t), value :: n + !> @param density_matrix Pointer to the density matrix data. + real(c_double), intent(in) :: density_matrix(ldp, *) + !> @param ldp Leading dimension of the density matrix. + integer(c_int64_t), value :: ldp + !> @param settings Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in) :: settings + !> @param exc_grad Pointer to the gradient. + real(c_double), intent(out) :: exc_grad(*) + end subroutine gauxc_integrator_eval_exc_grad_rks_with_settings_c module procedure gauxc_integrator_eval_exc_grad_rks subroutine gauxc_integrator_eval_exc_grad_uks_c(status, integrator, m, n, & @@ -487,6 +563,34 @@ subroutine gauxc_integrator_eval_fxc_contraction_rks_c(status, integrator, m, n, !> @param ldfxc Leading dimension of the fxc contraction matrix. integer(c_int64_t), value :: ldfxc end subroutine gauxc_integrator_eval_fxc_contraction_rks_c + subroutine gauxc_integrator_eval_fxc_contraction_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, t_density_matrix, ldtp, settings, fxc_contraction_matrix, ldfxc) & + & bind(c, name="gauxc_integrator_eval_fxc_contraction_rks_with_settings") + import :: gauxc_status_type, gauxc_integrator_type, gauxc_ks_settings_type, c_double, c_int64_t + implicit none + !> @param status Status object to capture any errors. + type(gauxc_status_type), intent(inout) :: status + !> @param integrator Handle to the XCIntegrator. + type(gauxc_integrator_type), value :: integrator + !> @param m Number of rows in the density matrix. + integer(c_int64_t), value :: m + !> @param n Number of columns in the density matrix. + integer(c_int64_t), value :: n + !> @param density_matrix Pointer to the density matrix data. + real(c_double), intent(in) :: density_matrix(ldp, *) + !> @param ldp Leading dimension of the density matrix. + integer(c_int64_t), value :: ldp + !> @param t_density_matrix Pointer to the density matrix data. + real(c_double), intent(in) :: t_density_matrix(ldtp, *) + !> @param ldtp Leading dimension of the density matrix. + integer(c_int64_t), value :: ldtp + !> @param settings Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in) :: settings + !> @param fxc_contraction_matrix Pointer to the fxc contraction matrix data. + real(c_double), intent(out) :: fxc_contraction_matrix(ldfxc, *) + !> @param ldfxc Leading dimension of the fxc contraction matrix. + integer(c_int64_t), value :: ldfxc + end subroutine gauxc_integrator_eval_fxc_contraction_rks_with_settings_c module procedure gauxc_integrator_eval_fxc_contraction_rks subroutine gauxc_integrator_eval_fxc_contraction_uks_c(status, integrator, m, n, & @@ -621,7 +725,7 @@ subroutine gauxc_integrator_integrate_den(status, integrator, density_matrix, de end subroutine gauxc_integrator_integrate_den !> @brief Evaluate the exchange-correlation energy for RKS. - subroutine gauxc_integrator_eval_exc_rks(status, integrator, density_matrix, exc) + subroutine gauxc_integrator_eval_exc_rks(status, integrator, density_matrix, exc, settings) !> @param status Status object to capture any errors. type(gauxc_status_type), intent(inout) :: status !> @param integrator Handle to the XCIntegrator. @@ -630,13 +734,20 @@ subroutine gauxc_integrator_eval_exc_rks(status, integrator, density_matrix, exc real(c_double), contiguous, intent(in) :: density_matrix(:, :) !> @param exc Pointer to store the exchange-correlation energy. real(c_double), intent(out) :: exc + !> @param settings Optional Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in), optional :: settings integer(c_int64_t) :: m, n, ldp m = size(density_matrix, 1, kind=c_int64_t) n = size(density_matrix, 2, kind=c_int64_t) ldp = size(density_matrix, 1, kind=c_int64_t) - call gauxc_integrator_eval_exc_rks_c(status, integrator, m, n, density_matrix, ldp, exc) + if (present(settings)) then + call gauxc_integrator_eval_exc_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, settings, exc) + else + call gauxc_integrator_eval_exc_rks_c(status, integrator, m, n, density_matrix, ldp, exc) + end if end subroutine gauxc_integrator_eval_exc_rks !> @brief Evaluate the exchange-correlation energy for UKS. @@ -696,7 +807,7 @@ end subroutine gauxc_integrator_eval_exc_gks !> @brief Evaluate the exchange-correlation energy and potential for RKS. subroutine gauxc_integrator_eval_exc_vxc_rks(status, integrator, & - & density_matrix, exc, vxc_matrix) + & density_matrix, exc, vxc_matrix, settings) !> @param status Status object to capture any errors. type(gauxc_status_type), intent(inout) :: status !> @param integrator Handle to the XCIntegrator. @@ -707,6 +818,8 @@ subroutine gauxc_integrator_eval_exc_vxc_rks(status, integrator, & real(c_double), intent(out) :: exc !> @param vxc_matrix Pointer to the potential matrix data. real(c_double), contiguous, intent(out) :: vxc_matrix(:, :) + !> @param settings Optional Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in), optional :: settings integer(c_int64_t) :: m, n, ldp, ldvxc @@ -715,8 +828,13 @@ subroutine gauxc_integrator_eval_exc_vxc_rks(status, integrator, & ldp = size(density_matrix, 1, kind=c_int64_t) ldvxc = size(vxc_matrix, 1, kind=c_int64_t) - call gauxc_integrator_eval_exc_vxc_rks_c(status, integrator, m, n, & - & density_matrix, ldp, exc, vxc_matrix, ldvxc) + if (present(settings)) then + call gauxc_integrator_eval_exc_vxc_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, settings, exc, vxc_matrix, ldvxc) + else + call gauxc_integrator_eval_exc_vxc_rks_c(status, integrator, m, n, & + & density_matrix, ldp, exc, vxc_matrix, ldvxc) + end if end subroutine gauxc_integrator_eval_exc_vxc_rks !> @brief Evaluate the exchange-correlation energy and potential for UKS. @@ -834,7 +952,7 @@ subroutine gauxc_integrator_eval_exc_vxc_gks(status, integrator, & end subroutine gauxc_integrator_eval_exc_vxc_gks subroutine gauxc_integrator_eval_exc_grad_rks(status, integrator, & - & density_matrix, exc_grad) + & density_matrix, exc_grad, settings) !> @param status Status object to capture any errors. type(gauxc_status_type), intent(inout) :: status !> @param integrator Handle to the XCIntegrator. @@ -843,6 +961,8 @@ subroutine gauxc_integrator_eval_exc_grad_rks(status, integrator, & real(c_double), contiguous, intent(in) :: density_matrix(:, :) !> @param exc_grad Pointer to the gradient. real(c_double), contiguous, intent(out) :: exc_grad(:) + !> @param settings Optional Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in), optional :: settings integer(c_int64_t) :: m, n, ldp @@ -850,8 +970,13 @@ subroutine gauxc_integrator_eval_exc_grad_rks(status, integrator, & n = size(density_matrix, 2, kind=c_int64_t) ldp = size(density_matrix, 1, kind=c_int64_t) - call gauxc_integrator_eval_exc_grad_rks_c(status, integrator, m, n, & - & density_matrix, ldp, exc_grad) + if (present(settings)) then + call gauxc_integrator_eval_exc_grad_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, settings, exc_grad) + else + call gauxc_integrator_eval_exc_grad_rks_c(status, integrator, m, n, & + & density_matrix, ldp, exc_grad) + end if end subroutine gauxc_integrator_eval_exc_grad_rks subroutine gauxc_integrator_eval_exc_grad_uks(status, integrator, & @@ -927,7 +1052,7 @@ subroutine gauxc_integrator_eval_exx_rks(status, integrator, & end subroutine gauxc_integrator_eval_exx_rks subroutine gauxc_integrator_eval_fxc_contraction_rks(status, integrator, & - & density_matrix, t_density_matrix, fxc_contraction_matrix) + & density_matrix, t_density_matrix, fxc_contraction_matrix, settings) !> @param status Status object to capture any errors. type(gauxc_status_type), intent(inout) :: status !> @param integrator Handle to the XCIntegrator. @@ -938,6 +1063,8 @@ subroutine gauxc_integrator_eval_fxc_contraction_rks(status, integrator, & real(c_double), contiguous, intent(in) :: t_density_matrix(:, :) !> @param fxc_contraction_matrix Pointer to the fxc contraction matrix data. real(c_double), contiguous, intent(out) :: fxc_contraction_matrix(:, :) + !> @param settings Optional Kohn-Sham integration settings. + type(gauxc_ks_settings_type), intent(in), optional :: settings integer(c_int64_t) :: m, n, ldp, ldtp, ldfxc @@ -947,8 +1074,13 @@ subroutine gauxc_integrator_eval_fxc_contraction_rks(status, integrator, & ldtp = size(t_density_matrix, 1, kind=c_int64_t) ldfxc = size(fxc_contraction_matrix, 1, kind=c_int64_t) - call gauxc_integrator_eval_fxc_contraction_rks_c(status, integrator, m, n, & - & density_matrix, ldp, t_density_matrix, ldtp, fxc_contraction_matrix, ldfxc) + if (present(settings)) then + call gauxc_integrator_eval_fxc_contraction_rks_with_settings_c(status, integrator, m, n, & + & density_matrix, ldp, t_density_matrix, ldtp, settings, fxc_contraction_matrix, ldfxc) + else + call gauxc_integrator_eval_fxc_contraction_rks_c(status, integrator, m, n, & + & density_matrix, ldp, t_density_matrix, ldtp, fxc_contraction_matrix, ldfxc) + end if end subroutine gauxc_integrator_eval_fxc_contraction_rks subroutine gauxc_integrator_eval_fxc_contraction_uks(status, integrator, & @@ -990,4 +1122,4 @@ subroutine gauxc_integrator_eval_fxc_contraction_uks(status, integrator, & & fxc_contraction_matrix_s, ldfxc_s, fxc_contraction_matrix_z, ldfxc_z) end subroutine gauxc_integrator_eval_fxc_contraction_uks -end module gauxc_integrator \ No newline at end of file +end module gauxc_integrator diff --git a/tests/c_api_test.cxx b/tests/c_api_test.cxx index e91941a4..a158c18c 100644 --- a/tests/c_api_test.cxx +++ b/tests/c_api_test.cxx @@ -1169,6 +1169,34 @@ TEST_CASE("C-API XCIntegrator", "[c-api]") { p.cleanup(); } + SECTION("eval_exc_rks with spin-summed density (Benzene / SVWN5 / cc-pVDZ)") { + std::string ref_file = GAUXC_REF_DATA_PATH "/benzene_svwn5_cc-pvdz_ufg_ssf.hdf5"; + auto p = make_c_api_pipeline(ref_file, "SVWN5"); + + HighFive::File file(ref_file, HighFive::File::ReadOnly); + auto dset = file.getDataSet("/DENSITY"); + auto dims = dset.getDimensions(); + int64_t m = dims[0], n = dims[1]; + std::vector P(m * n); + dset.read(P.data()); + for(auto& x : P) x *= 2.0; + + double EXC_ref; + file.getDataSet("/EXC").read(&EXC_ref); + + C::GauXCKSSettings settings{}; + settings.gks_dtol = 1e-12; + settings.rks_density_matrix_is_spin_summed = true; + + double exc = 0.0; + C::gauxc_integrator_eval_exc_rks_with_settings( + &p.status, p.integrator, m, n, P.data(), m, &settings, &exc); + CHECK(p.status.code == 0); + CHECK(exc == Approx(EXC_ref)); + + p.cleanup(); + } + SECTION("eval_exc_vxc_rks (Benzene / SVWN5 / cc-pVDZ)") { std::string ref_file = GAUXC_REF_DATA_PATH "/benzene_svwn5_cc-pvdz_ufg_ssf.hdf5"; auto p = make_c_api_pipeline(ref_file, "SVWN5");