Support for more UQ types via Jacobi eigensolver

src/scf/driver/scf_loop.cpp

@@ -28,23 +28,14 @@ struct Kernel {
 
     template<typename FloatType>
     auto operator()(const std::span<FloatType>& a) {
-        using tensorwrapper::types::fabs;
-        using clean_type = std::decay_t<FloatType>;
-        bool in_noise    = false;
-        bool converged   = false;
-        if constexpr(tensorwrapper::types::is_interval_v<clean_type>) {
-            auto abs_val = fabs(a[0]);
-            // If uncertainty is greater than the tolerance, return true
-            in_noise  = abs_val.radius() > m_tol;
-            converged = abs_val.upper() < m_tol;
-        } else if constexpr(tensorwrapper::types::is_uncertain_v<clean_type>) {
-            auto abs_val = fabs(a[0]);
-            in_noise     = abs_val.sd() > m_tol;
-            converged    = (abs_val.mean() + abs_val.sd()) < m_tol;
-        } else {
-            converged = fabs(a[0]) < m_tol;
-        }
-        return std::make_pair(converged, in_noise);
+        using tensorwrapper::types::uq_center;
+        /// For UQ, the center/median of the difference is our best estimate
+        /// for the uncertainty caused by incomplete convergence. It plus/minus
+        // the radius/standard deviation is the total range of uncertainty.
+        // Convergence happens when the center/median is within the tolerance,
+        // regardless of the radius/standard deviation.
+        auto abs_val = std::fabs(uq_center(a[0]));
+        return abs_val < m_tol;
     }
 };
 
@@ -290,11 +281,7 @@ MODULE_RUN(SCFLoop) {
             auto e_conv  = check_tolerance(de.buffer(), e_tol);
             auto g_conv  = check_tolerance(grad_norm.buffer(), g_tol);
             auto dp_conv = check_tolerance(dp_norm.buffer(), dp_tol);
-            if(e_conv.first && g_conv.first && dp_conv.first) converged = true;
-            if(e_conv.second) {
-                logger.log("  Energy has converged to within uncertainty");
-                converged = true;
-            }
+            if(e_conv && g_conv && dp_conv) converged = true;
 
             // If using DIIS and not converged, extrapolate new Fock matrix
             if(diis_on && !converged) { F = diis.extrapolate(F, grad); }

src/scf/eigen_solver/eigen_generalized.cpp

@@ -44,7 +44,7 @@ struct Kernel {
         using clean_t = std::decay_t<FloatType>;
         // Convert to Eigen buffers
 
-        if constexpr(tensorwrapper::types::is_interval_v<clean_t>) {
+        if constexpr(tensorwrapper::types::is_uq_type_v<clean_t>) {
             throw std::runtime_error(
               "EigenGeneralized Kernel: Interval types not supported");
         } else {