Quantized versions of Jacobi and Identity preconditioners

Common/include/code_config.hpp

@@ -50,6 +50,13 @@
 #define NEVERINLINE inline
 #endif
 
+/*--- Marks a function callable from both host and device code, a no-op outside nvcc. ---*/
+#ifdef __CUDACC__
+#define SU2_CUDA_HOST_DEVICE __host__ __device__
+#else
+#define SU2_CUDA_HOST_DEVICE
+#endif
+
 #if defined(__INTEL_COMPILER)
 /*--- Disable warnings related to inline attributes. ---*/
 #pragma warning disable 2196

Common/include/linear_algebra/CMatrixInverse.hpp

@@ -30,11 +30,7 @@
 
 #include <cmath>
 
-#ifdef __CUDACC__
-#define SU2_CUDA_HOST_DEVICE __host__ __device__
-#else
-#define SU2_CUDA_HOST_DEVICE
-#endif
+#include "../code_config.hpp"
 
 namespace SU2_LinAlg {
 
@@ -95,5 +91,3 @@ SU2_CUDA_HOST_DEVICE inline void MatrixInverse(unsigned long nVar, ScalarType* m
 }
 
 }  // namespace SU2_LinAlg
-
-#undef SU2_CUDA_HOST_DEVICE

Common/include/linear_algebra/CPreconditioner.hpp

@@ -107,13 +107,25 @@ CPreconditioner<ScalarType>::~CPreconditioner() {}
 /*!
  * \class CIdentityPreconditioner
  * \brief No-op preconditioner used when Krylov solvers run without preconditioning.
+ * \note Also serves Q_IDENTITY: Build() requests quantization of the diagonal blocks, needed by
+ * the matrix-vector product shared with the Krylov solver even though this preconditioner's own
+ * operation is a no-op. CSysMatrix::QuantizeDiagonalBlocks() when quantization is off.
  */
 template <class ScalarType>
 class CIdentityPreconditioner final : public CPreconditioner<ScalarType> {
+ private:
+  CSysMatrix<ScalarType>& sparse_matrix;
+
  public:
+  inline explicit CIdentityPreconditioner(CSysMatrix<ScalarType>& matrix_ref) : sparse_matrix(matrix_ref) {}
+
+  CIdentityPreconditioner() = delete;
+
   inline void operator()(const CSysVector<ScalarType>& u, CSysVector<ScalarType>& v) const override { v = u; }
 
   inline bool IsIdentity() const override { return true; }
+
+  inline void Build() override { sparse_matrix.QuantizeDiagonalBlocks(); }
 };
 
 /*!
@@ -157,7 +169,9 @@ class CJacobiPreconditioner final : public CPreconditioner<ScalarType> {
   }
 
   /*!
-   * \note Request the associated matrix to build the preconditioner.
+   * \note Request the associated matrix to build the preconditioner. Also serves Q_JACOBI:
+   *       BuildJacobiPreconditioner() quantizes the diagonal blocks itself when the matrix was
+   *       set up for it.
    */
   inline void Build() override { sparse_matrix.BuildJacobiPreconditioner(); }
 };
@@ -248,36 +262,10 @@ class CLU_SGSPreconditioner final : public CPreconditioner<ScalarType> {
   inline void operator()(const CSysVector<ScalarType>& u, CSysVector<ScalarType>& v) const override {
     ApplyPreconditionerOnHost(u, v, [&] { sparse_matrix.ComputeLU_SGSPreconditioner(u, v, geometry, config); });
   }
-};
-
-/*!
- * \class CQuantizedLUSGSPreconditioner
- * \brief Specialization of preconditioner that uses CSysMatrix class.
- */
-template <class ScalarType>
-class CQuantizedLUSGSPreconditioner final : public CPreconditioner<ScalarType> {
- private:
-  CSysMatrix<ScalarType>& sparse_matrix;
-  CGeometry* geometry;
-  const CConfig* config;
-
- public:
-  inline CQuantizedLUSGSPreconditioner(CSysMatrix<ScalarType>& matrix_ref, CGeometry* geometry_ref,
-                                       const CConfig* config_ref)
-      : sparse_matrix(matrix_ref) {
-    if ((geometry_ref == nullptr) || (config_ref == nullptr))
-      SU2_MPI::Error("Preconditioner needs to be built with valid references.", CURRENT_FUNCTION);
-    geometry = geometry_ref;
-    config = config_ref;
-  }
-
-  CQuantizedLUSGSPreconditioner() = delete;
-
-  inline void operator()(const CSysVector<ScalarType>& u, CSysVector<ScalarType>& v) const override {
-    ApplyPreconditionerOnHost(u, v, [&] { sparse_matrix.ComputeLU_SGSPreconditioner(u, v, geometry, config); });
-  }
 
-  /*! \brief Quantize the diagonal blocks (off diagonals are quantized on the fly). */
+  /*!
+   * \note Also serves Q_LU_SGS: quantizes the diagonal blocks, no-op for plain LU_SGS.
+   */
   inline void Build() override { sparse_matrix.QuantizeDiagonalBlocks(); }
 };
 
@@ -404,19 +392,19 @@ CPreconditioner<ScalarType>* CPreconditioner<ScalarType>::Create(ENUM_LINEAR_SOL
 
   switch (kind) {
     case IDENTITY:
-      prec = new CIdentityPreconditioner<ScalarType>();
+    case Q_IDENTITY:
+      prec = new CIdentityPreconditioner<ScalarType>(jacobian);
       break;
     case JACOBI:
+    case Q_JACOBI:
       prec = new CJacobiPreconditioner<ScalarType>(jacobian, geometry, config);
       break;
     case LINELET:
       prec = new CLineletPreconditioner<ScalarType>(jacobian, geometry, config);
       break;
     case LU_SGS:
-      prec = new CLU_SGSPreconditioner<ScalarType>(jacobian, geometry, config);
-      break;
     case Q_LU_SGS:
-      prec = new CQuantizedLUSGSPreconditioner<ScalarType>(jacobian, geometry, config);
+      prec = new CLU_SGSPreconditioner<ScalarType>(jacobian, geometry, config);
       break;
     case ILU:
       prec = new CILUPreconditioner<ScalarType>(jacobian, geometry, config);