]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use SUNContext as needed in SUNDIALS tests, part 1.
authorDavid Wells <drwells@email.unc.edu>
Sun, 15 May 2022 22:02:56 +0000 (18:02 -0400)
committerDavid Wells <drwells@email.unc.edu>
Tue, 17 May 2022 13:01:02 +0000 (09:01 -0400)
This new context object is mostly for logging (i.e., its API is for profiling
and logging only at the moment) so we just need to make sure it exists.

tests/sundials/n_vector.cc

index 14eb7a5a363964d5906b6742063a24ab5c18cbdf..21802a148d47a13fde569102dc33423e202d6d03 100644 (file)
 
 using namespace SUNDIALS::internal;
 
+// For better compatibility with older versions of SUNDIALS just define this as
+// a static member so we can use it whenever we need it
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+static SUNContext global_nvector_context;
+#endif
+
 // anonymous namespace groups helper functions for testing
 namespace
 {
@@ -199,10 +205,20 @@ template <typename VectorType>
 void
 test_nvector_view_unwrap()
 {
-  auto       vector         = create_test_vector<VectorType>();
-  const auto const_vector   = create_test_vector<VectorType>();
-  auto       n_vector       = make_nvector_view(vector);
-  auto       const_n_vector = make_nvector_view(const_vector);
+  auto       vector       = create_test_vector<VectorType>();
+  const auto const_vector = create_test_vector<VectorType>();
+  auto       n_vector     = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
+  auto const_n_vector = make_nvector_view(const_vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                          ,
+                                          global_nvector_context
+#endif
+  );
 
   Assert(n_vector != nullptr, NVectorTestError());
   Assert((n_vector)->content != nullptr, NVectorTestError());
@@ -248,8 +264,14 @@ void
 test_get_vector_id()
 {
   auto vector   = create_test_vector<VectorType>();
-  auto n_vector = make_nvector_view(vector);
-  auto id       = N_VGetVectorID(n_vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+
+  );
+  auto id = N_VGetVectorID(n_vector);
   Assert(id == SUNDIALS_NVEC_CUSTOM, NVectorTestError());
   deallog << "test_get_vector_id OK" << std::endl;
 }
@@ -261,8 +283,13 @@ void
 test_clone()
 {
   auto vector   = create_test_vector<VectorType>();
-  auto n_vector = make_nvector_view(vector);
-  auto cloned   = N_VClone(n_vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
+  auto cloned = N_VClone(n_vector);
 
   Assert(cloned != nullptr, NVectorTestError());
   AssertDimension(unwrap_nvector<VectorType>(cloned)->size(), vector.size());
@@ -279,7 +306,12 @@ test_destroy()
 {
   GrowingVectorMemory<VectorType>                   mem;
   typename GrowingVectorMemory<VectorType>::Pointer vector(mem);
-  auto n_vector = make_nvector_view(*vector);
+  auto n_vector = make_nvector_view(*vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
 
   Assert(n_vector != nullptr, NVectorTestError());
   auto cloned = N_VClone(n_vector);
@@ -313,7 +345,12 @@ void
 test_get_communicator()
 {
   auto vector   = create_test_vector<VectorType>();
-  auto n_vector = make_nvector_view(vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
   // required by SUNDIALS: MPI-unaware vectors should return the nullptr
   Assert(N_VGetCommunicator(n_vector) == nullptr, NVectorTestError());
 
@@ -328,7 +365,12 @@ void
 test_get_communicator()
 {
   auto vector   = create_test_vector<VectorType>();
-  auto n_vector = make_nvector_view(vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
   Assert(*static_cast<MPI_Comm *>(N_VGetCommunicator(n_vector)) ==
            MPI_COMM_WORLD,
          NVectorTestError());
@@ -343,7 +385,12 @@ void
 test_length()
 {
   auto vector   = create_test_vector<VectorType>();
-  auto n_vector = make_nvector_view(vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
   Assert(N_VGetLength(n_vector) == static_cast<int>(vector.size()),
          NVectorTestError());
 
@@ -361,9 +408,24 @@ test_linear_sum()
   auto vc       = create_test_vector<VectorType>();
   auto expected = create_test_vector<VectorType>();
 
-  auto nv_a = make_nvector_view(va);
-  auto nv_b = make_nvector_view(vb);
-  auto nv_c = make_nvector_view(vc);
+  auto nv_a = make_nvector_view(va
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vb
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_c = make_nvector_view(vc
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   va = 1.0;
   vb = 2.0;
@@ -398,8 +460,18 @@ test_dot_product()
   const auto vb   = create_test_vector<VectorType>(3.0);
   const auto size = va.size();
 
-  auto nv_a = make_nvector_view(va);
-  auto nv_b = make_nvector_view(vb);
+  auto nv_a = make_nvector_view(va
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vb
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   auto result   = N_VDotProd(nv_a, nv_b);
   auto expected = 6.0 * size;
@@ -423,7 +495,12 @@ test_set_constant()
   auto expected = create_test_vector<VectorType>();
   expected      = 1.0;
 
-  auto n_vector = make_nvector_view(vector);
+  auto n_vector = make_nvector_view(vector
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                    ,
+                                    global_nvector_context
+#endif
+  );
 
   N_VConst(1.0, n_vector);
   Assert(vector_equal(vector, expected), NVectorTestError());
@@ -445,8 +522,18 @@ test_add_constant()
   auto       result   = create_test_vector<VectorType>(-1.0);
   auto       expected = create_test_vector<VectorType>(3.0);
 
-  auto nv_a      = make_nvector_view(vector_a);
-  auto nv_result = make_nvector_view(result);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_result = make_nvector_view(result
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     ,
+                                     global_nvector_context
+#endif
+  );
 
   N_VAddConst(nv_a, 1.0, nv_result);
   Assert(vector_equal(result, expected), NVectorTestError());
@@ -470,9 +557,24 @@ test_elementwise_product()
   auto       vector_result = create_test_vector<VectorType>();
   auto       expected      = create_test_vector<VectorType>(6.0);
 
-  auto nv_a      = make_nvector_view(vector_a);
-  auto nv_b      = make_nvector_view(vector_b);
-  auto nv_result = make_nvector_view(vector_result);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_result = make_nvector_view(vector_result
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     ,
+                                     global_nvector_context
+#endif
+  );
 
   // result = a.*b
   N_VProd(nv_a, nv_b, nv_result);
@@ -480,7 +582,12 @@ test_elementwise_product()
 
   auto vector_c = create_test_vector<VectorType>(2.0);
 
-  auto nv_c = make_nvector_view(vector_c);
+  auto nv_c = make_nvector_view(vector_c
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   // c .*= b
   N_VProd(nv_c, nv_b, nv_c);
@@ -509,9 +616,24 @@ test_elementwise_div()
   auto       vector_result = create_test_vector<VectorType>();
   auto       expected      = create_test_vector<VectorType>(2.0);
 
-  auto nv_a      = make_nvector_view(vector_a);
-  auto nv_b      = make_nvector_view(vector_b);
-  auto nv_result = make_nvector_view(vector_result);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_result = make_nvector_view(vector_result
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     ,
+                                     global_nvector_context
+#endif
+  );
 
   N_VDiv(nv_a, nv_b, nv_result);
   Assert(vector_equal(vector_result, expected), NVectorTestError());
@@ -529,8 +651,18 @@ test_elementwise_inv()
   auto       vector_result = create_test_vector<VectorType>();
   auto       expected      = create_test_vector<VectorType>(1.0 / 6.0);
 
-  auto nv_a      = make_nvector_view(vector_a);
-  auto nv_result = make_nvector_view(vector_result);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_result = make_nvector_view(vector_result
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     ,
+                                     global_nvector_context
+#endif
+  );
 
   N_VInv(nv_a, nv_result);
   Assert(vector_equal(vector_result, expected), NVectorTestError());
@@ -553,8 +685,18 @@ test_elementwise_abs()
   auto       vector_result = create_test_vector<VectorType>();
   auto       expected      = create_test_vector<VectorType>(2.0);
 
-  auto nv_a      = make_nvector_view(vector_a);
-  auto nv_result = make_nvector_view(vector_result);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_result = make_nvector_view(vector_result
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     ,
+                                     global_nvector_context
+#endif
+  );
 
   N_VAbs(nv_a, nv_result);
   Assert(vector_equal(vector_result, expected), NVectorTestError());
@@ -576,8 +718,18 @@ test_weighted_rms_norm()
   const auto vector_a = create_test_vector<VectorType>(2.0);
   const auto vector_b = create_test_vector<VectorType>(3.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   const auto result = N_VWrmsNorm(nv_a, nv_b);
   Assert(std::fabs(result - 6.0) < 1e-12, NVectorTestError());
@@ -596,10 +748,30 @@ test_weighted_rms_norm_mask()
   const auto vector_c = create_test_vector<VectorType>(1.0);
   const auto vector_d = create_test_vector<VectorType>(0.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
-  auto nv_c = make_nvector_view(vector_c);
-  auto nv_d = make_nvector_view(vector_d);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_c = make_nvector_view(vector_c
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_d = make_nvector_view(vector_d
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
   {
     const auto result = N_VWrmsNormMask(nv_a, nv_b, nv_c);
     Assert(std::fabs(result - 6.0) < 1e-12, NVectorTestError());
@@ -621,8 +793,18 @@ test_weighted_l2_norm()
   const auto vector_a = create_test_vector<VectorType>(2.0);
   const auto vector_b = create_test_vector<VectorType>(3.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   const auto result = N_VWL2Norm(nv_a, nv_b);
 
@@ -645,8 +827,18 @@ test_max_norm()
   const auto vector_a = create_test_vector<VectorType>(2.0);
   const auto vector_b = create_test_vector<VectorType>(-3.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   auto result = N_VMaxNorm(nv_a);
   Assert(std::fabs(result - 2.0) < 1e-12, NVectorTestError());
@@ -666,8 +858,18 @@ test_l1_norm()
   const auto vector_a = create_test_vector<VectorType>(2.0);
   const auto vector_b = create_test_vector<VectorType>(-3.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   auto result = N_VL1Norm(nv_a);
   Assert(std::fabs(result - vector_a.l1_norm()) < 1e-12, NVectorTestError());
@@ -687,8 +889,18 @@ test_min_element()
   const auto vector_a = create_test_vector<VectorType>(2.0);
   const auto vector_b = create_test_vector<VectorType>(-3.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   auto result = N_VMin(nv_a);
   Assert(std::fabs(result - 2.0) < 1e-12, NVectorTestError());
@@ -709,8 +921,18 @@ test_scale()
   auto       vector_b = create_test_vector<VectorType>(2.0);
   const auto expected = create_test_vector<VectorType>(4.0);
 
-  auto nv_a = make_nvector_view(vector_a);
-  auto nv_b = make_nvector_view(vector_b);
+  auto nv_a = make_nvector_view(vector_a
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
+  auto nv_b = make_nvector_view(vector_b
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                ,
+                                global_nvector_context
+#endif
+  );
 
   // b *= 2
   N_VScale(2.0, nv_b, nv_b);
@@ -761,6 +983,11 @@ main(int argc, char **argv)
 {
   Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
   MPILogInitAll                    log_all;
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+  MPI_Comm  communicator = MPI_COMM_WORLD;
+  const int ierr = SUNContext_Create(&communicator, &global_nvector_context);
+  AssertThrow(ierr == 0, ExcMessage("unable to create SUNContext object"));
+#endif
 
   using VectorType = Vector<double>;
 
@@ -782,4 +1009,8 @@ main(int argc, char **argv)
   // which is invoked first
   GrowingVectorMemory<PETScWrappers::MPI::Vector>::release_unused_memory();
   GrowingVectorMemory<PETScWrappers::MPI::BlockVector>::release_unused_memory();
+
+#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+  SUNContext_Free(&global_nvector_context);
+#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.