]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make the SUNDIALS vector copy more generic.
authorDavid Wells <wellsd2@rpi.edu>
Sat, 16 Sep 2017 18:00:30 +0000 (14:00 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Sat, 16 Sep 2017 19:07:00 +0000 (15:07 -0400)
source/sundials/copy.cc

index bf4e3c7ab235e73ecebe520c13cb90cb788425e5..27c554d063959263436e408400ba2cdbd7a1bcb5 100644 (file)
@@ -1,15 +1,15 @@
 //-----------------------------------------------------------
 //
-//    Copyright (C) 2015 by the deal.II authors
+// Copyright (C) 2017 by the deal.II authors
 //
-//    This file is part of the deal.II library.
+// 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.
+// 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.
 //
 //-----------------------------------------------------------
 
@@ -22,6 +22,41 @@ namespace SUNDIALS
 {
   namespace internal
   {
+    namespace
+    {
+      /**
+       * SUNDIALS provides different macros for getting the local length of a
+       * vector for serial and parallel vectors (as well as various parallel
+       * vectors that are not yet supported by deal.II). This function provides
+       * a generic interface to both and does a (checked) conversion from long
+       * int (the type SUNDIALS uses for lengths) to std::size_t.
+       */
+      inline
+      std::size_t
+      N_Vector_length(const N_Vector &vec)
+      {
+        const N_Vector_ID id = N_VGetVectorID(vec);
+        long int length = -1;
+        switch (id)
+          {
+          case SUNDIALS_NVEC_SERIAL:
+          {
+            length = NV_LENGTH_S(vec);
+            break;
+          }
+          case SUNDIALS_NVEC_PARALLEL:
+          {
+            length = NV_LOCLENGTH_P(vec);
+            break;
+          }
+          default:
+            Assert(false, ExcNotImplemented());
+          }
+
+        Assert(length >= 0, ExcInternalError());
+        return static_cast<std::size_t>(length);
+      }
+    }
 
 #ifdef DEAL_II_WITH_MPI
 
@@ -31,7 +66,7 @@ namespace SUNDIALS
     void copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
     {
       IndexSet is = dst.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(src));
+      AssertDimension(is.n_elements(), N_Vector_length(src));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
@@ -42,7 +77,7 @@ namespace SUNDIALS
     void copy(N_Vector &dst, const TrilinosWrappers::MPI::Vector &src)
     {
       IndexSet is = src.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(dst));
+      AssertDimension(is.n_elements(), N_Vector_length(dst));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
@@ -52,7 +87,7 @@ namespace SUNDIALS
     void copy(TrilinosWrappers::MPI::BlockVector &dst, const N_Vector &src)
     {
       IndexSet is = dst.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(src));
+      AssertDimension(is.n_elements(), N_Vector_length(src));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
@@ -63,7 +98,7 @@ namespace SUNDIALS
     void copy(N_Vector &dst, const TrilinosWrappers::MPI::BlockVector &src)
     {
       IndexSet is = src.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(dst));
+      AssertDimension(is.n_elements(), N_Vector_length(dst));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
@@ -77,7 +112,7 @@ namespace SUNDIALS
     void copy(PETScWrappers::MPI::Vector &dst, const N_Vector &src)
     {
       IndexSet is = dst.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(src));
+      AssertDimension(is.n_elements(), N_Vector_length(src));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
@@ -88,7 +123,7 @@ namespace SUNDIALS
     void copy(N_Vector &dst, const PETScWrappers::MPI::Vector &src)
     {
       IndexSet is = src.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(dst));
+      AssertDimension(is.n_elements(), N_Vector_length(dst));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
@@ -98,7 +133,7 @@ namespace SUNDIALS
     void copy(PETScWrappers::MPI::BlockVector &dst, const N_Vector &src)
     {
       IndexSet is = dst.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(src));
+      AssertDimension(is.n_elements(), N_Vector_length(src));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
@@ -109,7 +144,7 @@ namespace SUNDIALS
     void copy(N_Vector &dst, const PETScWrappers::MPI::BlockVector &src)
     {
       IndexSet is = src.locally_owned_elements();
-      AssertDimension(is.n_elements(), NV_LOCLENGTH_P(dst));
+      AssertDimension(is.n_elements(), N_Vector_length(dst));
       for (unsigned int i=0; i<is.n_elements(); ++i)
         {
           NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
@@ -122,7 +157,7 @@ namespace SUNDIALS
 
     void copy(BlockVector<double> &dst, const N_Vector &src)
     {
-      AssertDimension((unsigned int)NV_LENGTH_S(src), dst.size());
+      AssertDimension(N_Vector_length(src), dst.size());
       for (unsigned int i=0; i<dst.size(); ++i)
         {
           dst[i] = NV_Ith_S(src, i);
@@ -131,7 +166,7 @@ namespace SUNDIALS
 
     void copy(N_Vector &dst, const BlockVector<double> &src)
     {
-      AssertDimension((unsigned int)NV_LENGTH_S(dst), src.size());
+      AssertDimension(N_Vector_length(dst), src.size());
       for (unsigned int i=0; i<src.size(); ++i)
         {
           NV_Ith_S(dst, i) = src[i];
@@ -140,7 +175,7 @@ namespace SUNDIALS
 
     void copy(Vector<double> &dst, const N_Vector &src)
     {
-      AssertDimension((unsigned int)NV_LENGTH_S(src), dst.size());
+      AssertDimension(N_Vector_length(src), dst.size());
       for (unsigned int i=0; i<dst.size(); ++i)
         {
           dst[i] = NV_Ith_S(src, i);
@@ -149,7 +184,7 @@ namespace SUNDIALS
 
     void copy(N_Vector &dst, const Vector<double> &src)
     {
-      AssertDimension((unsigned int)NV_LENGTH_S(dst), src.size());
+      AssertDimension(N_Vector_length(dst), src.size());
       for (unsigned int i=0; i<src.size(); ++i)
         {
           NV_Ith_S(dst, i) = src[i];
@@ -160,4 +195,3 @@ namespace SUNDIALS
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
-

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.