]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Generalise Scalapack scale_columns and scale_rows functions
authorJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 15 Feb 2018 08:50:16 +0000 (09:50 +0100)
committerBenjamin Brands <benjamin.brands@fau.de>
Thu, 15 Feb 2018 11:45:25 +0000 (12:45 +0100)
include/deal.II/base/array_view.h
include/deal.II/lac/scalapack.h
source/lac/CMakeLists.txt
source/lac/scalapack.cc
source/lac/scalapack.inst.in [new file with mode: 0644]

index 4321e61fd5d9dab3c4a90b082b7e2bd003291c84..84cfeb6627af33ef0d53cf02220a2452682fa5c1 100644 (file)
@@ -524,6 +524,46 @@ make_array_view (ElementType *const begin, ElementType *const end)
 
 
 
+/**
+ * Create a view from an ArrayView itself.
+ *
+ * This function is used for @p const references to objects of ArrayView type.
+ * It only exists for compatibility purposes.
+ *
+ * @param[in] array_view The ArrayView that we wish to make a copy of.
+ *
+ * @relates ArrayView
+ */
+template <typename Number>
+inline
+ArrayView<const Number>
+make_array_view (const ArrayView<Number> &array_view)
+{
+  return make_array_view (array_view.cbegin(), array_view.cend());
+}
+
+
+
+/**
+ * Create a view from an ArrayView itself.
+ *
+ * This function is used for non-@p const references to objects of ArrayView type.
+ * It only exists for compatibility purposes.
+ *
+ * @param[in] array_view The ArrayView that we wish to make a copy of.
+ *
+ * @relates ArrayView
+ */
+template <typename Number>
+inline
+ArrayView<Number>
+make_array_view (ArrayView<Number> &array_view)
+{
+  return make_array_view (array_view.begin(), array_view.end());
+}
+
+
+
 /**
  * Create a view to an entire Tensor object. This is equivalent to initializing
  * an ArrayView object with a pointer to the first element and the size of the
index 0cc03e9b87f570de4d5f5152ad0f69d2a288a919..68e673bb56a18809bd4df78d2d7a9ef97b127948 100644 (file)
@@ -26,7 +26,6 @@
 #include <deal.II/lac/lapack_support.h>
 #include <deal.II/base/mpi.h>
 #include <deal.II/base/thread_management.h>
-#include <deal.II/base/array_view.h>
 #include <mpi.h>
 
 #include <memory>
@@ -543,8 +542,12 @@ public:
    * The array @p factors must have as many entries as the matrix columns.
    *
    * Copies of @p factors have to be available on all processes of the underlying MPI communicator.
+   *
+   * @note The fundamental prerequisite for the @p InputVector is that it must be possible to
+   * create an ArrayView from it; this is satisfied by the @p std::vector and Vector classes.
    */
-  void scale_columns(const ArrayView<const NumberType> &factors);
+  template <class InputVector>
+  void scale_columns(const InputVector &factors);
 
   /**
    * Scale the rows of the distributed matrix by the scalars provided in the array @p factors.
@@ -552,8 +555,12 @@ public:
    * The array @p factors must have as many entries as the matrix rows.
    *
    * Copies of @p factors have to be available on all processes of the underlying MPI communicator.
+   *
+   * @note The fundamental prerequisite for the @p InputVector is that it must be possible to
+   * create an ArrayView from it; this is satisfied by the @p std::vector and Vector classes.
    */
-  void scale_rows(const ArrayView<const NumberType> &factors);
+  template <class InputVector>
+  void scale_rows(const InputVector &factors);
 
 private:
 
index defcbb0fb2b69f6346183214ee5920e62e673fa0..f98f40e702e9972c4b67fed9c41fc56be5f6ab80 100644 (file)
@@ -73,6 +73,7 @@ SET(_inst
   precondition_block.inst.in
   relaxation_block.inst.in
   read_write_vector.inst.in
+  scalapack.inst.in
   solver.inst.in
   sparse_matrix_ez.inst.in
   sparse_matrix.inst.in
index 864962d2df355b474123fe5a95775926d4618c33..67e068e0cc61ed49c384eb140ee0f11d5fbae655 100644 (file)
@@ -1506,43 +1506,72 @@ void ScaLAPACKMatrix<NumberType>::load_parallel(const char *filename)
 
 
 
-template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::scale_columns(const ArrayView<const NumberType> &factors)
+namespace internal
 {
-  Assert(n_columns==(int)factors.size(),ExcDimensionMismatch(n_columns,factors.size()));
+  namespace
+  {
 
-  if (this->grid->mpi_process_is_active)
-    for (int i=0; i<n_local_columns; ++i)
-      {
-        const NumberType s = factors[global_column(i)];
+    template <typename NumberType>
+    void scale_columns(ScaLAPACKMatrix<NumberType>       &matrix,
+                       const ArrayView<const NumberType> &factors,
+                       const bool                         grid_mpi_process_is_active)
+    {
+      Assert(matrix.n()==factors.size(),ExcDimensionMismatch(matrix.n(),factors.size()));
+
+      if (grid_mpi_process_is_active)
+        for (unsigned int i=0; i<matrix.local_n(); ++i)
+          {
+            const NumberType s = factors[matrix.global_column(i)];
+
+            for (unsigned int j=0; j<matrix.local_m(); ++j)
+              matrix.local_el(j,i) *= s;
+          }
+    }
+
+    template <typename NumberType>
+    void scale_rows(ScaLAPACKMatrix<NumberType>       &matrix,
+                    const ArrayView<const NumberType> &factors,
+                    const bool                         grid_mpi_process_is_active)
+    {
+      Assert(matrix.m()==factors.size(),ExcDimensionMismatch(matrix.m(),factors.size()));
 
-        for (int j=0; j<n_local_rows; ++j)
-          local_el(j,i) *= s;
-      }
+      if (grid_mpi_process_is_active)
+        for (unsigned int i=0; i<matrix.local_m(); ++i)
+          {
+            const NumberType s = factors[matrix.global_row(i)];
+
+            for (unsigned int j=0; j<matrix.local_n(); ++j)
+              matrix.local_el(i,j) *= s;
+          }
+    }
+
+  }
 }
 
 
 
 template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::scale_rows(const ArrayView<const NumberType> &factors)
+template <class InputVector>
+void ScaLAPACKMatrix<NumberType>::scale_columns(const InputVector &factors)
 {
-  Assert(n_rows==(int)factors.size(),ExcDimensionMismatch(n_rows,factors.size()));
+  internal::scale_columns(*this, make_array_view(factors),
+                          this->grid->mpi_process_is_active);
+}
+
 
-  if (this->grid->mpi_process_is_active)
-    for (int i=0; i<n_local_rows; ++i)
-      {
-        const NumberType s = factors[global_row(i)];
 
-        for (int j=0; j<n_local_columns; ++j)
-          local_el(i,j) *= s;
-      }
+template <typename NumberType>
+template <class InputVector>
+void ScaLAPACKMatrix<NumberType>::scale_rows(const InputVector &factors)
+{
+  internal::scale_rows(*this, make_array_view(factors),
+                       this->grid->mpi_process_is_active);
 }
 
 
 
 // instantiations
-template class ScaLAPACKMatrix<double>;
-template class ScaLAPACKMatrix<float>;
+#include "scalapack.inst"
 
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/source/lac/scalapack.inst.in b/source/lac/scalapack.inst.in
new file mode 100644 (file)
index 0000000..0e3ca68
--- /dev/null
@@ -0,0 +1,36 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+for (Number : REAL_SCALARS)
+{
+    template class ScaLAPACKMatrix<Number>;
+
+    template
+    void ScaLAPACKMatrix<Number>::scale_columns<Vector<Number>> (const Vector<Number>&);
+
+    template
+    void ScaLAPACKMatrix<Number>::scale_rows<Vector<Number>> (const Vector<Number>&);
+}
+
+for (VEC : GENERAL_CONTAINER_TYPES;
+        Number : REAL_SCALARS)
+{
+    template
+    void ScaLAPACKMatrix<Number>::scale_columns<VEC<Number>> (const VEC<Number>&);
+
+    template
+    void ScaLAPACKMatrix<Number>::scale_rows<VEC<Number>> (const VEC<Number>&);
+}

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.