From 2c0c8aa7cf2e70ee759f5b9bac7e3ecabca5c863 Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Mon, 15 Jan 2024 17:11:38 +0200 Subject: [PATCH 1/4] Expose partial Cholesky orthogonalization in header --- src/linalg.h | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/linalg.h b/src/linalg.h index 8f299003..c8e87eca 100644 --- a/src/linalg.h +++ b/src/linalg.h @@ -104,6 +104,9 @@ arma::mat CanonicalOrth(const arma::mat & S, double cutoff); /// Same, but use computed decomposition arma::mat CanonicalOrth(const arma::mat & Svec, const arma::vec & Sval, double cutoff); +/// Pivoted Cholesky orthogonalization +arma::mat PartialCholeskyOrth(const arma::mat & S, double cholcut, double scut); + /// Automatic orthonormalization. arma::mat BasOrth(const arma::mat & S, bool verbose); /// Orthogonalize basis From f8eb6ad4965607b4b484d635598db70263f22fc2 Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Mon, 15 Jan 2024 17:16:13 +0200 Subject: [PATCH 2/4] Expose another way to compute nuclear attraction matrices --- src/basis.cpp | 31 +++++++++++++++++++------------ src/basis.h | 2 ++ 2 files changed, 21 insertions(+), 12 deletions(-) diff --git a/src/basis.cpp b/src/basis.cpp index ba470aa6..8d4915b4 100644 --- a/src/basis.cpp +++ b/src/basis.cpp @@ -2318,7 +2318,23 @@ arma::mat BasisSet::kinetic() const { } arma::mat BasisSet::nuclear() const { - // Form nuclear attraction matrix + std::vector> nuclear_data; + for(size_t inuc=0;inuc> & nuclear_data) const { // Size of basis set size_t N=get_Nbf(); @@ -2332,18 +2348,9 @@ arma::mat BasisSet::nuclear() const { #pragma omp parallel for schedule(dynamic) #endif for(size_t ip=0;ip> & nuclei) const; /// Calculate electric potential matrix arma::mat potential(coords_t r) const; /// Calculate superposition of atomic potentials From c7436176202ef36ed556deb23f026280468cacfe Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Mon, 15 Jan 2024 17:19:03 +0200 Subject: [PATCH 3/4] Use pivoted Cholesky in determining auxiliary basis for density fitting --- src/density_fitting.cpp | 26 +++----------------------- src/density_fitting.h | 2 +- src/localization.cpp | 6 +++--- src/scf-base.cpp | 6 ++++-- src/scf.h | 2 ++ src/settings.cpp | 1 + 6 files changed, 14 insertions(+), 29 deletions(-) diff --git a/src/density_fitting.cpp b/src/density_fitting.cpp index 1efa1a97..5d128370 100644 --- a/src/density_fitting.cpp +++ b/src/density_fitting.cpp @@ -52,7 +52,7 @@ bool DensityFit::Bmat_enabled() const { return Bmat; } -size_t DensityFit::fill(const BasisSet & orbbas, const BasisSet & auxbas, bool dir, double erithr, double linthr, bool bmat) { +size_t DensityFit::fill(const BasisSet & orbbas, const BasisSet & auxbas, bool dir, double erithr, double linthr, double cholthr, bool bmat) { // Construct density fitting basis // Store amount of functions @@ -125,28 +125,8 @@ size_t DensityFit::fill(const BasisSet & orbbas, const BasisSet & auxbas, bool d } if(Bmat) { - // Form ab^-1 and ab^-1/2 - arma::mat abvec; - arma::vec abval; - eig_sym_ordered(abval,abvec,ab); - - // Count linearly independent vectors - size_t Nind=0; - for(size_t i=0;i=linthr) - Nind++; - - // and drop the linearly dependent ones - abval=abval.subvec(abval.n_elem-Nind,abval.n_elem-1); - abvec=abvec.cols(abvec.n_cols-Nind,abvec.n_cols-1); - - // Form matrices - ab_inv.zeros(abvec.n_rows,abvec.n_rows); - ab_invh.zeros(abvec.n_rows,abvec.n_rows); - for(size_t i=0;iget_unique_shellpairs().size()); diff --git a/src/scf.h b/src/scf.h index ff3689e2..3349b4ab 100644 --- a/src/scf.h +++ b/src/scf.h @@ -287,6 +287,8 @@ class SCF { size_t fitmem; /// Threshold for density fitting double fitthr; + /// Pivoted Cholesky threshold for density fitting + double fitcholthr; /// Cholesky calculation? bool cholesky; diff --git a/src/settings.cpp b/src/settings.cpp index 12a01258..c68ed51c 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -142,6 +142,7 @@ void Settings::add_scf_settings() { add_int("FittingMemory", "Amount of memory in MB to use for exchange fitting",1000); // Threshold for screening eigenvectors add_double("FittingThreshold", "Linear dependence threshold for Coulomb integrals in density fitting",1e-8); + add_double("FittingCholeskyThreshold", "Linear dependence threshold for pivoted Cholesky of Coulomb integrals in density fitting",1e-9); // SAP basis add_string("SAPBasis", "Tabulated atomic effective potential \"basis set\"","helfem_large.gbs"); From 88084170c8ea3f01de05114705e7555f637af9f3 Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Mon, 15 Jan 2024 17:02:55 +0200 Subject: [PATCH 4/4] Density matrix digestion --- src/density_fitting.cpp | 9 +++++++++ src/density_fitting.h | 2 ++ 2 files changed, 11 insertions(+) diff --git a/src/density_fitting.cpp b/src/density_fitting.cpp index 5d128370..6a42d6fd 100644 --- a/src/density_fitting.cpp +++ b/src/density_fitting.cpp @@ -783,6 +783,15 @@ arma::mat DensityFit::calcJ(const arma::mat & P) const { // Get the expansion coefficients arma::vec c=compute_expansion(P); + return digestJ(c); +} + +arma::mat DensityFit::calcJ(const arma::vec & c) const { + if(c.n_elem != Naux) { + std::ostringstream oss; + oss << "Error in DensityFit: Naux = " << Naux << ", c.n_elem = " << c.n_elem << "!\n"; + throw std::logic_error(oss.str()); + } arma::mat J(Nbf,Nbf); J.zeros(); diff --git a/src/density_fitting.h b/src/density_fitting.h index 3445a06a..5d6337d1 100644 --- a/src/density_fitting.h +++ b/src/density_fitting.h @@ -149,6 +149,8 @@ class DensityFit { arma::mat calcJ(const arma::mat & P) const; /// Get Coulomb matrix from P std::vector calcJ(const std::vector & P) const; + /// Digest J matrix from computed expansion + arma::mat digest(const arma::vec & gamma) const; /// Calculate force from P arma::vec forceJ(const arma::mat & P);