]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add new wrappers to initialize Paralution and modify the current tests to use them.
authorturcksin <turcksin@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Jun 2014 21:13:28 +0000 (21:13 +0000)
committerturcksin <turcksin@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Jun 2014 21:13:28 +0000 (21:13 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_paralution@33028 0785d39b-7218-0410-832d-ea1e28bc413d

13 files changed:
deal.II/include/deal.II/base/paralution.h [new file with mode: 0644]
deal.II/include/deal.II/lac/paralution_sparse_matrix.h
deal.II/source/lac/paralution_sparse_matrix.cc
tests/lac/paralution_solver_01.cc
tests/lac/paralution_solver_02.cc
tests/lac/paralution_sparse_matrix_01.cc
tests/lac/paralution_vector_01.cc
tests/lac/paralution_vector_02.cc
tests/lac/paralution_vector_03.cc
tests/lac/paralution_vector_04.cc
tests/lac/paralution_vector_05.cc
tests/lac/paralution_vector_06.cc
tests/tests.h

diff --git a/deal.II/include/deal.II/base/paralution.h b/deal.II/include/deal.II/base/paralution.h
new file mode 100644 (file)
index 0000000..315df97
--- /dev/null
@@ -0,0 +1,175 @@
+// ---------------------------------------------------------------------
+// $Id: paralution.h 32926 2014-05-16 17:12:16Z turcksin $
+//
+// Copyright (C) 2014 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#ifndef __deal2__paralution_h
+#define __deal2__paralution_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_PARALUTION
+
+#include <deal.II/base/multithread_info.h>
+
+#include <paralution.hpp>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+namespace Utilities
+{
+  /**
+   * A namespace for utility functions that abstract certain operations using
+   * the Paralution library.
+   */
+  namespace Paralution
+  {
+    /**
+     * A class that is used to initialize and to stop Paralution. If a program
+     * uses Paralution one would typically just create an object of this type at
+     * the beginninh of <code>main()</code>. The constructor of this class then
+     * runs <code>paralution::init_paralution()</code>. At the end of the
+     * program, the compiler will invoke the destructor of this object which in
+     * turns calls <code>paralution::stop_paralution()</code>.
+     */
+    class Paralution_InitFinalize
+    {
+    public :
+      /**
+       * Constructor. Initialize Paralution by calling <tt>paralution::init_paralution</tt>
+       * and set the number of threads used by Paralution to the given parameters.
+       */
+      Paralution_InitFinalize(const unsigned int max_num_threads);
+
+      /**
+       * Destructor. Calls <tt>paralution::stop_paralution</tt>.
+       */
+      ~Paralution_InitFinalize();
+
+      /**
+       * Set a device.
+       */
+      void set_device(const unsigned int device) const;
+
+      /**
+       * Set a specific gpu.
+       */
+      void set_gpu_cuda(const unsigned int gpu) const;
+
+      /**
+       * Set OpenCL compute units.
+       */
+      void set_ocl_compute_units(const unsigned int compute_unit) const;
+
+      /**
+       * Set a specific OpenCL platform.
+       */
+      void set_ocl_platform(const unsigned int platform) const;
+
+      /**
+       * Set a specific OpenCL platform and device.
+       */
+      void set_ocl(const unsigned int platform, const unsigned int device) const;
+
+      /**
+       * Set the number of OpenMP thread.
+       */
+      void set_omp_threads(const unsigned int n_threads) const;
+
+      /**
+       * Print information about the platform.
+       */
+      void info() const;
+    };
+
+
+
+// ------------------- inline functions --------------
+    inline Paralution_InitFinalize::Paralution_InitFinalize(const unsigned int max_num_threads)
+    {
+      // Initialize paralution
+      paralution::init_paralution();
+
+      // Set the number of OpenMP threads
+      set_omp_threads(max_num_threads);
+
+      // Set the number of TBB threads
+      multithread_info.set_thread_limit(max_num_threads);
+    }
+
+
+    inline Paralution_InitFinalize::~Paralution_InitFinalize()
+    {
+      paralution::stop_paralution();
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_device(const unsigned int device) const
+    {
+      paralution::set_device_paralution(device);
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_gpu_cuda(const unsigned int gpu) const
+    {
+      paralution::set_gpu_cuda_paralution(gpu);
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_ocl_compute_units(const unsigned int compute_unit) const
+    {
+      paralution::set_ocl_compute_units_paralution(compute_unit);
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_ocl_platform(const unsigned int platform) const
+    {
+      paralution::set_ocl_platform_paralution(platform);
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_ocl(const unsigned int platform,
+                                                 const unsigned int device) const
+    {
+      paralution::set_ocl_paralution(platform, device);
+    }
+
+
+
+    inline void Paralution_InitFinalize::set_omp_threads(const unsigned int n_threads) const
+    {
+      paralution::set_omp_threads_paralution(n_threads);
+    }
+
+
+
+    inline void Paralution_InitFinalize::info() const
+    {
+      paralution::info_paralution();
+    }
+  }
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_PARALUTION
+
+#endif
index c98e9d84f18324f0ca267cab814b5e7c37c8aff2..40f365597bdaedf281954244450db859a171b4bb 100644 (file)
@@ -24,7 +24,6 @@
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/subscriptor.h>
 #include <deal.II/lac/sparse_matrix.h>
-#include <algorithm>
 
 #include <paralution.hpp>
 
@@ -371,6 +370,12 @@ namespace ParalutionWrappers
      */
     void sync();
 
+    /**
+     * Convert the paralution matrix to a different format. This function should
+     * only be called after convert_to_paralution_csr.
+     */
+    void convert_format(matrix_format format);
+
     /**
      * Return a constant reference to the underlying dealii::SparseMatrix.
      */
index 7f20e32150ff268e5948e491170085df3b55e800..530f34b1a829789c176eee4c76796a9abaee30da 100644 (file)
@@ -137,6 +137,56 @@ namespace ParalutionWrappers
     else
       local_matrix.CopyFromAsync(sparse_matrix.paralution_matrix());
   }
+
+  template <typename Number>
+  void SparseMatrix<Number>::convert_format(matrix_format format)
+  {
+    if (is_local_matrix==false)
+    {
+      switch (format)
+      {
+        case DENSE :
+          {
+            local_matrix.ConvertToDense();
+            break;
+          }
+        case CSR :
+          {
+            local_matrix.ConvertToCSR();
+            break;
+          }
+        case MCSR :
+          {
+            local_matrix.ConvertToMCSR();
+            break;
+          }
+        case BCSR :
+          {
+            local_matrix.ConvertToBCSR();
+            break;
+          }
+        case COO :
+          {
+            local_matrix.ConvertToDIA();
+            break;
+          }
+        case ELL :
+          {
+            local_matrix.ConvertToELL();
+            break;
+          }
+        case HYB :
+          {
+            local_matrix.ConvertToHYB();
+            break;
+          }
+        default :
+          {
+            AssertThrow(false,ExcMessage("Wrong format of Paralution matrix."));
+          }
+      }
+    }
+  }
 }
 
 // Explicit instantiations
index 720f95e07416c5f004d6825af8b66fec038e0235..8057394373cad5b5d399dc97eca8c64b5007e168 100644 (file)
@@ -23,6 +23,7 @@
 #include <fstream>
 #include <iomanip>
 #include <deal.II/base/logstream.h>
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_sparse_matrix.h>
 #include <deal.II/lac/paralution_vector.h>
 #include <deal.II/lac/paralution_solver.h>
@@ -30,7 +31,7 @@
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -69,6 +70,4 @@ int main()
     }
     solver.solve(A,u,f,p);
   }
-
-  paralution::stop_paralution();
 }
index 7bc614a68e104fa1a7c05ba5456d99e762f2746c..34b3cca2cde13e5b7cb04d5b4a05b733a309b22b 100644 (file)
@@ -23,6 +23,7 @@
 #include <fstream>
 #include <iomanip>
 #include <deal.II/base/logstream.h>
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_sparse_matrix.h>
 #include <deal.II/lac/paralution_vector.h>
 #include <deal.II/lac/paralution_solver.h>
@@ -30,7 +31,7 @@
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -69,6 +70,4 @@ int main()
     }
     solver.solve(A,u,f,p);
   }
-
-  paralution::stop_paralution();
 }
index 44335a3de8c33def6b6fba2133ea4f96a779b94f..ca8b7661ea6138f594406e1e270d061f02153986 100644 (file)
@@ -18,9 +18,9 @@
 // Test ParalutionWrappers::SparseMatrix
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/sparsity_pattern.h>
 #include <deal.II/lac/paralution_sparse_matrix.h>
-#include "paralution.hpp"
 
 template <typename Number>
 void check()
@@ -55,7 +55,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -66,7 +66,5 @@ int main()
 
   check<float>();
   check<double>();
-
-  paralution::stop_paralution();
 }
 
index 73311620d68b27f4caa7ed08257fc80815de008b..e7997c38d98e0f1f1e886d6053726c76ceb6561b 100644 (file)
@@ -18,8 +18,8 @@
 // Test ParalutionWrappers::Vector constructor
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 template <typename Number>
 void check()
@@ -35,7 +35,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -47,6 +47,4 @@ int main()
   check<int>();
   check<float>();
   check<double>();
-
-  paralution::stop_paralution();
 }
index d3d30684652601328822ef19ee3516c09363ea78..05adabc943ba4667c819127cc271d74db948ca20 100644 (file)
@@ -18,8 +18,8 @@
 // Test copy constructor and access elements of ParalutionWrappers::Vector
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 void check()
 {
@@ -39,7 +39,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -49,6 +49,4 @@ int main()
   deallog.threshold_double(1.e-10);
 
   check();
-
-  paralution::stop_paralution();
 }
index d4f298679dbeb282918c3d9d6c79f61116f58ca2..51ff86d8029ea798cdcc4a74d16aceeae008bd06 100644 (file)
@@ -18,8 +18,8 @@
 // Test ParalutionWrappers::Vector iterator
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 void check()
 {
@@ -40,7 +40,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -50,6 +50,4 @@ int main()
   deallog.threshold_double(1.e-10);
 
   check();
-
-  paralution::stop_paralution();
 }
index c21c6fd4183fe4e12f375d7935f2b756d172239f..6a77a6db57a8c5f0d217e57b1df1d85d04822b06 100644 (file)
@@ -18,8 +18,8 @@
 // Test reinit for ParalutionWrappers::Vector
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 void check()
 {
@@ -39,7 +39,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -49,6 +49,4 @@ int main()
   deallog.threshold_double(1.e-10);
 
   check();
-
-  paralution::stop_paralution();
 }
index 6e55260b08be90b7097e1752d416a6b167a6a8b8..79bd2df3645cd8b70145ae68f873d59ecbf54fbb 100644 (file)
@@ -19,8 +19,8 @@
 // operator-=
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 void check_1()
 {
@@ -58,7 +58,7 @@ void check_2()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -69,6 +69,4 @@ int main()
 
   check_1();
   check_2();
-
-  paralution::stop_paralution();
 }
index fe00d96ed600ba1ca366d8c153fc90aec42b8477..6ac90a3ab0dfba0109f4c17e65d7c0c6ba397367 100644 (file)
@@ -19,8 +19,8 @@
 
 
 #include "../tests.h"
+#include <deal.II/base/paralution.h>
 #include <deal.II/lac/paralution_vector.h>
-#include "paralution.hpp"
 
 void check()
 {
@@ -37,7 +37,7 @@ void check()
 
 int main()
 {
-  paralution::init_paralution();
+  Utilities::Paralution::Paralution_InitFinalize paralution(1);
 
   std::ofstream logfile("output");
   deallog << std::fixed;
@@ -47,6 +47,4 @@ int main()
   deallog.threshold_double(1.e-10);
 
   check();
-
-  paralution::stop_paralution();
 }
index ef37ca2046a016c9758b4900e5b45f7f451029db..9d8680ed62e7b0fc3419d53da7b5d3226d6c0600 100644 (file)
@@ -157,9 +157,10 @@ void unify_pretty_function (const std::string &filename)
  * Note that we can't do this if we run in MPI mode because then
  * MPI_InitFinalize already calls this function. Since every test
  * calls MPI_InitFinalize itself, we can't adjust the thread count
- * for this here.
+ * for this here. The same is true for Paralution, Paralution_InitFinalize
+ * calls set_thread_limit.
  */
-#ifndef DEAL_II_WITH_MPI
+#if !(defined DEAL_II_WITH_MPI) && !(defined DEAL_II_WITH_PARALUTION)
 struct LimitConcurrency
 {
   LimitConcurrency ()

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.