[GSoC] Feature: Jacobi preconditioner on GPU and CUDA Unified preconditiong matrix

Common/include/linear_algebra/CPreconditioner.hpp

@@ -114,13 +114,39 @@ class CJacobiPreconditioner final : public CPreconditioner<ScalarType> {
    * \param[out] v - CSysVector that is the result of the preconditioning
    */
   inline void operator()(const CSysVector<ScalarType>& u, CSysVector<ScalarType>& v) const override {
-    sparse_matrix.ComputeJacobiPreconditioner(u, v, geometry, config);
+    if (config->GetCUDA()) {
+#ifdef HAVE_CUDA
+      sparse_matrix.GPUComputeJacobiPreconditioner(u, v, geometry, config);
+#else
+      SU2_MPI::Error(
+          "\nError in launching sparse matrix Preconditioner Function\nENABLE_CUDA is set to YES\nPlease compile with "
+          "CUDA "
+          "options enabled in Meson to access GPU Functions",
+          CURRENT_FUNCTION);
+#endif
+    } else {
+      sparse_matrix.ComputeJacobiPreconditioner(u, v, geometry, config);
+    }
   }
 
   /*!
    * \note Request the associated matrix to build the preconditioner.
    */
-  inline void Build() override { sparse_matrix.BuildJacobiPreconditioner(); }
+  inline void Build() override {
+    if (config->GetCUDA()) {
+#ifdef HAVE_CUDA
+      sparse_matrix.GPUBuildJacobiPreconditioner();
+#else
+      SU2_MPI::Error(
+          "\nError in building sparse matrix Preconditioner Function\nENABLE_CUDA is set to YES\nPlease compile with "
+          "CUDA "
+          "options enabled in Meson to access GPU Functions",
+          CURRENT_FUNCTION);
+#endif
+    } else {
+      sparse_matrix.BuildJacobiPreconditioner();
+    }
+  }
 };
 
 /*!

Common/include/linear_algebra/CSysMatrix.hpp

@@ -176,7 +176,9 @@ class CSysMatrix {
   /*!< \brief Level structure for alternative shared memory parallelization of ILU. */
   CCompressedSparsePatternUL levels_ilu;
 
-  ScalarType* invM; /*!< \brief Inverse of (Jacobi) preconditioner. */
+  ScalarType* invM;             /*!< \brief Inverse of (Jacobi) preconditioner. */
+  ScalarType* d_invM;           /*!< \brief Inverse of (Jacobi) preconditioner on device. */
+  bool invM_is_managed = false; /*!< \brief Boolean that indicates whether GPU supports Unified Memory or not */
 
   /*--- Temporary (hence mutable) working memory used in the Linelet preconditioner, outer vector is for threads ---*/
   mutable vector<vector<const ScalarType*> >
@@ -924,6 +926,11 @@ class CSysMatrix {
    */
   void BuildJacobiPreconditioner();
 
+  /*!
+   * \brief Build the Jacobi preconditioner on GPU.
+   */
+  void GPUBuildJacobiPreconditioner();
+
   /*!
    * \brief Multiply CSysVector by the preconditioner
    * \param[in] vec - CSysVector to be multiplied by the preconditioner.
@@ -934,6 +941,16 @@ class CSysMatrix {
   void ComputeJacobiPreconditioner(const CSysVector<ScalarType>& vec, CSysVector<ScalarType>& prod, CGeometry* geometry,
                                    const CConfig* config) const;
 
+  /*!
+   * \brief Multiply CSysVector by the preconditioner
+   * \param[in] vec - CSysVector to be multiplied by the preconditioner.
+   * \param[out] prod - Result of the product A*vec.
+   * \param[in] geometry - Geometrical definition of the problem.
+   * \param[in] config - Definition of the particular problem.
+   */
+  void GPUComputeJacobiPreconditioner(const CSysVector<ScalarType>& vec, CSysVector<ScalarType>& prod,
+                                      CGeometry* geometry, const CConfig* config) const;
+
   /*!
    * \brief Build the ILU preconditioner.
    */