]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Enable SUNDIALS::KINSOL for LinearAlgebra::distributed::Vector 12162/head
authorPeter Munch <peterrmuench@gmail.com>
Mon, 10 May 2021 08:46:53 +0000 (10:46 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Tue, 11 May 2021 04:49:45 +0000 (06:49 +0200)
include/deal.II/sundials/copy.h
source/sundials/copy.cc
source/sundials/kinsol.cc

index 4a14ae9a5bc76bc64363ffc4987fdabea4a7149a..afc2ae3b52938d4e9c91cec972e7ce267de4a285 100644 (file)
@@ -24,6 +24,8 @@
 #    include <nvector/nvector_parallel.h>
 #  endif
 #  include <deal.II/lac/block_vector.h>
+#  include <deal.II/lac/la_parallel_block_vector.h>
+#  include <deal.II/lac/la_parallel_vector.h>
 #  include <deal.II/lac/vector.h>
 
 #  include <nvector/nvector_serial.h>
@@ -82,6 +84,18 @@ namespace SUNDIALS
     copy(Vector<double> &dst, const N_Vector &src);
     void
     copy(N_Vector &dst, const Vector<double> &src);
+
+    void
+    copy(LinearAlgebra::distributed::BlockVector<double> &dst,
+         const N_Vector &                                 src);
+    void
+    copy(N_Vector &                                             dst,
+         const LinearAlgebra::distributed::BlockVector<double> &src);
+
+    void
+    copy(LinearAlgebra::distributed::Vector<double> &dst, const N_Vector &src);
+    void
+    copy(N_Vector &dst, const LinearAlgebra::distributed::Vector<double> &src);
   } // namespace internal
 } // namespace SUNDIALS
 DEAL_II_NAMESPACE_CLOSE
index 1c18eef17f6a6cf7baa068dd5d23daead1f43efc..b11b279cc0ad0cab166d26e27ab952f2012e1096 100644 (file)
@@ -56,6 +56,77 @@ namespace SUNDIALS
       return static_cast<std::size_t>(length);
     }
 
+
+    void
+    copy(LinearAlgebra::distributed::Vector<double> &dst, const N_Vector &src)
+    {
+      const IndexSet    is = dst.locally_owned_elements();
+      const std::size_t N  = is.n_elements();
+      AssertDimension(N, N_Vector_length(src));
+      dst.zero_out_ghost_values();
+      for (std::size_t i = 0; i < N; ++i)
+        {
+#  ifdef DEAL_II_WITH_MPI
+          dst.local_element(i) = NV_Ith_P(src, i);
+#  else
+          dst.local_element(i)        = NV_Ith_S(src, i);
+#  endif
+        }
+      dst.compress(VectorOperation::insert);
+    }
+
+    void
+    copy(N_Vector &dst, const LinearAlgebra::distributed::Vector<double> &src)
+    {
+      const IndexSet    is = src.locally_owned_elements();
+      const std::size_t N  = is.n_elements();
+      AssertDimension(N, N_Vector_length(dst));
+      for (std::size_t i = 0; i < N; ++i)
+        {
+#  ifdef DEAL_II_WITH_MPI
+          NV_Ith_P(dst, i) = src.local_element(i);
+#  else
+          NV_Ith_S(dst, i)            = src.local_element(i);
+#  endif
+        }
+    }
+
+    void
+    copy(LinearAlgebra::distributed::BlockVector<double> &dst,
+         const N_Vector &                                 src)
+    {
+      const IndexSet    is = dst.locally_owned_elements();
+      const std::size_t N  = is.n_elements();
+      AssertDimension(N, N_Vector_length(src));
+      dst.zero_out_ghost_values();
+      for (std::size_t i = 0; i < N; ++i)
+        {
+#  ifdef DEAL_II_WITH_MPI
+          dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
+#  else
+          dst[is.nth_index_in_set(i)] = NV_Ith_S(src, i);
+#  endif
+        }
+      dst.compress(VectorOperation::insert);
+    }
+
+    void
+    copy(N_Vector &                                             dst,
+         const LinearAlgebra::distributed::BlockVector<double> &src)
+    {
+      IndexSet          is = src.locally_owned_elements();
+      const std::size_t N  = is.n_elements();
+      AssertDimension(N, N_Vector_length(dst));
+      for (std::size_t i = 0; i < N; ++i)
+        {
+#  ifdef DEAL_II_WITH_MPI
+          NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
+#  else
+          NV_Ith_S(dst, i)            = src[is.nth_index_in_set(i)];
+#  endif
+        }
+    }
+
 #  ifdef DEAL_II_WITH_MPI
 
 #    ifdef DEAL_II_WITH_TRILINOS
index 83327e55f27a627a4a967b4bc13e7db015dd3265..ce3cb440778308ea1b1b910b86a0c804300f7069 100644 (file)
@@ -23,6 +23,8 @@
 #  include <deal.II/base/utilities.h>
 
 #  include <deal.II/lac/block_vector.h>
+#  include <deal.II/lac/la_parallel_block_vector.h>
+#  include <deal.II/lac/la_parallel_vector.h>
 #  ifdef DEAL_II_WITH_TRILINOS
 #    include <deal.II/lac/trilinos_parallel_block_vector.h>
 #    include <deal.II/lac/trilinos_vector.h>
@@ -612,6 +614,9 @@ namespace SUNDIALS
   template class KINSOL<Vector<double>>;
   template class KINSOL<BlockVector<double>>;
 
+  template class KINSOL<LinearAlgebra::distributed::Vector<double>>;
+  template class KINSOL<LinearAlgebra::distributed::BlockVector<double>>;
+
 #  ifdef DEAL_II_WITH_MPI
 
 #    ifdef DEAL_II_WITH_TRILINOS

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.