From: Sebastian Proell Date: Sat, 4 Nov 2023 11:39:18 +0000 (+0100) Subject: SUNDIALS N_Vector: implement efficient operations X-Git-Tag: relicensing~31^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F16645%2Fhead;p=dealii.git SUNDIALS N_Vector: implement efficient operations These operation are only supported by deal.II vectors. --- diff --git a/doc/news/changes/minor/20240215Proell b/doc/news/changes/minor/20240215Proell new file mode 100644 index 0000000000..2204b0ff0b --- /dev/null +++ b/doc/news/changes/minor/20240215Proell @@ -0,0 +1,4 @@ +New: Add potentially more efficient vector operation to SUNDIALS wrappers. These operations are used automatically +for deal.II's distributed vectors. +
+(Sebastian Proell, 2024/02/15) diff --git a/include/deal.II/sundials/n_vector.templates.h b/include/deal.II/sundials/n_vector.templates.h index 12455a9c22..da2d2fb57a 100644 --- a/include/deal.II/sundials/n_vector.templates.h +++ b/include/deal.II/sundials/n_vector.templates.h @@ -25,6 +25,7 @@ #ifdef DEAL_II_WITH_SUNDIALS # include +# include # include # include @@ -170,6 +171,84 @@ namespace SUNDIALS # endif ); + + /** + * Variable template that is true for distributed deal.II vectors and false + * otherwise. + */ + template + constexpr bool is_dealii_compatible_distributed_vector = + std::is_same_v< + VectorType, + LinearAlgebra::distributed::Vector> || + std::is_same_v>; + + template ::value> * = nullptr> + unsigned int + n_blocks(const VectorType &) + { + return 1; + } + + + + template ::value> * = nullptr> + unsigned int + n_blocks(const VectorType &vector) + { + return vector.n_blocks(); + } + + + + template ::value> * = nullptr> + VectorType & + block(VectorType &vector, const unsigned int b) + { + AssertDimension(b, 0); + (void)b; + return vector; + } + + + + template ::value> * = nullptr> + const VectorType & + block(const VectorType &vector, const unsigned int b) + { + AssertDimension(b, 0); + (void)b; + return vector; + } + + + + template ::value> * = nullptr> + typename VectorType::BlockType & + block(VectorType &vector, const unsigned int b) + { + return vector.block(b); + } + + + + template ::value> * = nullptr> + const typename VectorType::BlockType & + block(const VectorType &vector, const unsigned int b) + { + return vector.block(b); + } + + /** * Collection of all operations specified by SUNDIALS N_Vector * documentation. These functions are attached to the generic N_Vector @@ -203,10 +282,39 @@ namespace SUNDIALS N_Vector y, N_Vector z); + template < + typename VectorType, + std::enable_if_t> + * = nullptr> + int + linear_combination(int nv, realtype *c, N_Vector *x, N_Vector z); + template SUNDIALS::realtype dot_product(N_Vector x, N_Vector y); + template < + typename VectorType, + std::enable_if_t> + * = nullptr> + int + dot_product_multi(int nv, N_Vector x, N_Vector *Y, realtype *d); + + template < + typename VectorType, + std::enable_if_t> + * = nullptr> + int + dot_product_multi_local(int nv, N_Vector x, N_Vector *y, realtype *d); + + + template < + typename VectorType, + std::enable_if_t> + * = nullptr> + int + dot_product_multi_all_reduce(int nv, N_Vector x, realtype *d); + template SUNDIALS::realtype weighted_l2_norm(N_Vector x, N_Vector y); @@ -674,6 +782,33 @@ namespace SUNDIALS + template < + typename VectorType, + std::enable_if_t> *> + int + linear_combination(int nv, realtype *c, N_Vector *x, N_Vector z) + { + std::vector unwrapped_vectors(nv); + for (int i = 0; i < nv; ++i) + unwrapped_vectors[i] = unwrap_nvector_const(x[i]); + auto *z_dealii = unwrap_nvector(z); + + // N.B. The first pointer may alias with z. + for (unsigned int b = 0; b < n_blocks(*z_dealii); ++b) + for (unsigned int j = 0; j < block(*z_dealii, b).locally_owned_size(); + ++j) + { + double temp = 0.; + for (int i = 0; i < nv; ++i) + temp += block(*unwrapped_vectors[i], b).local_element(j) * c[i]; + block(*z_dealii, b).local_element(j) = temp; + } + + return 0; + } + + + template SUNDIALS::realtype dot_product(N_Vector x, N_Vector y) @@ -684,6 +819,62 @@ namespace SUNDIALS + template < + typename VectorType, + std::enable_if_t> *> + int + dot_product_multi(int nv, N_Vector x, N_Vector *y, realtype *d) + { + const int status = dot_product_multi_local(nv, x, y, d); + if (status != 0) + return status; + + return dot_product_multi_all_reduce(nv, x, d); + } + + + + template < + typename VectorType, + std::enable_if_t> *> + int + dot_product_multi_local(int nv, N_Vector x, N_Vector *y, realtype *d) + { + std::vector unwrapped_vectors(nv); + for (int i = 0; i < nv; ++i) + unwrapped_vectors[i] = unwrap_nvector_const(y[i]); + const VectorType *x_dealii = unwrap_nvector_const(x); + + std::fill(d, d + nv, 0.); + + for (unsigned int b = 0; b < n_blocks(*x_dealii); ++b) + for (unsigned int j = 0; j < block(*x_dealii, b).locally_owned_size(); + ++j) + { + for (int i = 0; i < nv; ++i) + d[i] += block(*x_dealii, b).local_element(j) * + block(*unwrapped_vectors[i], b).local_element(j); + } + return 0; + } + + + + template < + typename VectorType, + std::enable_if_t> *> + int + dot_product_multi_all_reduce(int nv, N_Vector x, realtype *d) + { + ArrayView products(d, nv); + Utilities::MPI::sum(products, + get_communicator(x), + products); + return 0; + } + + + template SUNDIALS::realtype weighted_l2_norm(N_Vector x, N_Vector w) @@ -1069,23 +1260,87 @@ namespace SUNDIALS v->ops->nvl1norm = &NVectorOperations::l1_norm; v->ops->nvwrmsnormmask = &NVectorOperations::weighted_rms_norm_mask; - // v->ops->nvcompare = undef; - // v->ops->nvinvtest = undef; - // v->ops->nvconstrmask = undef; - // v->ops->nvminquotient = undef; - - /* fused and vector array operations are disabled (NULL) by default */ - - /* local reduction operations */ - // v->ops->nvdotprodlocal = undef; - // v->ops->nvmaxnormlocal = undef; - // v->ops->nvminlocal = undef; - // v->ops->nvl1normlocal = undef; - // v->ops->nvinvtestlocal = undef; - // v->ops->nvconstrmasklocal = undef; - // v->ops->nvminquotientlocal = undef; - // v->ops->nvwsqrsumlocal = undef; - // v->ops->nvwsqrsummasklocal = undef; + + /* The following are declared as standard by SUNDIALS but were not + * necessary so far. + */ + v->ops->nvcompare = nullptr; + v->ops->nvinvtest = nullptr; + v->ops->nvconstrmask = nullptr; + v->ops->nvminquotient = nullptr; + + /* fused and vector array operations are disabled by default */ + + /* OPTIONAL fused vector operations */ + + // These operations are only implemented for deal.II vectors. SUNDIALS + // will automatically fall back to a less optimized version implemented + // through the standard functions above for other VectorTypes. + if constexpr (is_dealii_compatible_distributed_vector) + { + v->ops->nvlinearcombination = + &NVectorOperations::linear_combination; + v->ops->nvdotprodmulti = + &NVectorOperations::dot_product_multi; + } + else + { + v->ops->nvlinearcombination = nullptr; + v->ops->nvdotprodmulti = nullptr; + } + + v->ops->nvscaleaddmulti = nullptr; + + /* OPTIONAL vector array operations */ + v->ops->nvlinearsumvectorarray = nullptr; + v->ops->nvscalevectorarray = nullptr; + v->ops->nvconstvectorarray = nullptr; + v->ops->nvwrmsnormvectorarray = nullptr; + v->ops->nvwrmsnormmaskvectorarray = nullptr; + v->ops->nvscaleaddmultivectorarray = nullptr; + v->ops->nvlinearcombinationvectorarray = nullptr; + + /* Local reduction kernels (no parallel communication) */ + v->ops->nvdotprodlocal = nullptr; + v->ops->nvmaxnormlocal = nullptr; + v->ops->nvminlocal = nullptr; + v->ops->nvl1normlocal = nullptr; + v->ops->nvinvtestlocal = nullptr; + v->ops->nvconstrmasklocal = nullptr; + v->ops->nvminquotientlocal = nullptr; + v->ops->nvwsqrsumlocal = nullptr; + v->ops->nvwsqrsummasklocal = nullptr; + +# if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0) + /* Single buffer reduction operations */ + + // These operations are only implemented for deal.II vectors. SUNDIALS + // will automatically fall back to a less optimized version implemented + // through the standard functions above for other VectorTypes. + if constexpr (is_dealii_compatible_distributed_vector) + { + v->ops->nvdotprodmultilocal = + &NVectorOperations::dot_product_multi_local; + v->ops->nvdotprodmultiallreduce = + &NVectorOperations::dot_product_multi_all_reduce; + } + else + { + v->ops->nvdotprodmultilocal = nullptr; + v->ops->nvdotprodmultiallreduce = nullptr; + } +# endif + + /* XBraid interface operations */ + v->ops->nvbufsize = nullptr; + v->ops->nvbufpack = nullptr; + v->ops->nvbufunpack = nullptr; + + /* Debugging functions (called when SUNDIALS_DEBUG_PRINTVEC is defined). + */ + v->ops->nvprint = nullptr; + v->ops->nvprintfile = nullptr; + return v; } diff --git a/tests/sundials/n_vector.cc b/tests/sundials/n_vector.cc index 496d7d0fa9..72c18e0c13 100644 --- a/tests/sundials/n_vector.cc +++ b/tests/sundials/n_vector.cc @@ -16,8 +16,6 @@ // Test SUNDIALS' vector operations on N_Vector implementation. The N_Vectors // are created by calling NVectorView on one of the internal vector types. -#include "../../include/deal.II/sundials/n_vector.h" - #include #include @@ -26,8 +24,8 @@ #include #include -#include "../../include/deal.II/sundials/n_vector.templates.h" #include +#include #ifdef DEAL_II_WITH_MPI # include @@ -466,6 +464,66 @@ test_linear_sum() +template +void +test_linear_combination() +{ + auto va = create_test_vector(); + auto vb = create_test_vector(); + auto vc = create_test_vector(); + auto v_dst = create_test_vector(); + auto expected = create_test_vector(); + + 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 + ); + auto nv_dst = make_nvector_view(v_dst +#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0) + , + global_nvector_context +#endif + ); + + va = 1.0; + vb = 2.0; + vc = 3.0; + + expected = 1. * 1. + 2. * 2. + 3. * 3.; + + // Make a contiguous array of N_Vectors for consumption by SUNDIALS + // N.B. The array needs to contain N_Vector and not our NVectorView, so we use + // the conversion operator. + std::array n_vectors{ + {(N_Vector)nv_a, (N_Vector)nv_b, (N_Vector)nv_c}}; + + std::array weights{{1., 2., 3.}}; + // test sum three vectors + N_VLinearCombination(3, weights.data(), n_vectors.data(), nv_dst); + Assert(vector_equal(v_dst, expected), NVectorTestError()); + // repeat to test that sum overwrites initial content + N_VLinearCombination(3, weights.data(), n_vectors.data(), nv_dst); + Assert(vector_equal(v_dst, expected), NVectorTestError()); + + deallog << "test_linear_combination OK" << std::endl; +} + + + template void test_dot_product() @@ -501,6 +559,49 @@ test_dot_product() +template +void +test_dot_product_multi() +{ + const auto va = create_test_vector(1.0); + const auto vb = create_test_vector(2.0); + const auto vc = create_test_vector(3.0); + const auto size = va.size(); + + 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 + ); + + std::array nv_array{{(N_Vector)nv_b, (N_Vector)nv_c}}; + std::array dst; + + int status = + N_VDotProdMulti(nv_array.size(), nv_a, nv_array.data(), dst.data()); + Assert(status == 0, NVectorTestError()); + + Assert(std::fabs(dst[0] - 2. * size) < 1e-12, NVectorTestError()); + Assert(std::fabs(dst[1] - 3. * size) < 1e-12, NVectorTestError()); + + deallog << "test_dot_product_multi OK" << std::endl; +} + + + template void test_set_constant() @@ -976,7 +1077,9 @@ run_all_tests(const std::string &prefix) test_get_communicator(); test_length(); test_linear_sum(); + test_linear_combination(); test_dot_product(); + test_dot_product_multi(); test_set_constant(); test_add_constant(); test_elementwise_product(); diff --git a/tests/sundials/n_vector.mpirun=1.with_trilinos=off.with_petsc=on.output b/tests/sundials/n_vector.mpirun=1.with_trilinos=off.with_petsc=on.output index 0a1ccaeb3f..bc7fa901fc 100644 --- a/tests/sundials/n_vector.mpirun=1.with_trilinos=off.with_petsc=on.output +++ b/tests/sundials/n_vector.mpirun=1.with_trilinos=off.with_petsc=on.output @@ -6,7 +6,9 @@ DEAL:0:Vector::test_destroy OK DEAL:0:Vector::test_get_communicator OK DEAL:0:Vector::test_length OK DEAL:0:Vector::test_linear_sum OK +DEAL:0:Vector::test_linear_combination OK DEAL:0:Vector::test_dot_product OK +DEAL:0:Vector::test_dot_product_multi OK DEAL:0:Vector::test_set_constant OK DEAL:0:Vector::test_add_constant OK DEAL:0:Vector::test_elementwise_product OK @@ -28,7 +30,9 @@ DEAL:0:BlockVector::test_destroy OK DEAL:0:BlockVector::test_get_communicator OK DEAL:0:BlockVector::test_length OK DEAL:0:BlockVector::test_linear_sum OK +DEAL:0:BlockVector::test_linear_combination OK DEAL:0:BlockVector::test_dot_product OK +DEAL:0:BlockVector::test_dot_product_multi OK DEAL:0:BlockVector::test_set_constant OK DEAL:0:BlockVector::test_add_constant OK DEAL:0:BlockVector::test_elementwise_product OK @@ -50,7 +54,9 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_destroy OK DEAL:0:LinearAlgebra::distributed::Vector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::Vector::test_length OK DEAL:0:LinearAlgebra::distributed::Vector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::Vector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::Vector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_product OK @@ -72,7 +78,9 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_destroy OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_length OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_product OK @@ -94,7 +102,9 @@ DEAL:0:PETScWrappers::MPI::Vector::test_destroy OK DEAL:0:PETScWrappers::MPI::Vector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::Vector::test_length OK DEAL:0:PETScWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::Vector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::Vector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::Vector::test_set_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_add_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_product OK @@ -116,7 +126,9 @@ DEAL:0:PETScWrappers::MPI::BlockVector::test_destroy OK DEAL:0:PETScWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::BlockVector::test_length OK DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_product OK diff --git a/tests/sundials/n_vector.mpirun=1.with_trilinos=on.with_petsc=off.output b/tests/sundials/n_vector.mpirun=1.with_trilinos=on.with_petsc=off.output index 0f190076ae..dd000e1e8c 100644 --- a/tests/sundials/n_vector.mpirun=1.with_trilinos=on.with_petsc=off.output +++ b/tests/sundials/n_vector.mpirun=1.with_trilinos=on.with_petsc=off.output @@ -6,7 +6,9 @@ DEAL:0:Vector::test_destroy OK DEAL:0:Vector::test_get_communicator OK DEAL:0:Vector::test_length OK DEAL:0:Vector::test_linear_sum OK +DEAL:0:Vector::test_linear_combination OK DEAL:0:Vector::test_dot_product OK +DEAL:0:Vector::test_dot_product_multi OK DEAL:0:Vector::test_set_constant OK DEAL:0:Vector::test_add_constant OK DEAL:0:Vector::test_elementwise_product OK @@ -28,7 +30,9 @@ DEAL:0:BlockVector::test_destroy OK DEAL:0:BlockVector::test_get_communicator OK DEAL:0:BlockVector::test_length OK DEAL:0:BlockVector::test_linear_sum OK +DEAL:0:BlockVector::test_linear_combination OK DEAL:0:BlockVector::test_dot_product OK +DEAL:0:BlockVector::test_dot_product_multi OK DEAL:0:BlockVector::test_set_constant OK DEAL:0:BlockVector::test_add_constant OK DEAL:0:BlockVector::test_elementwise_product OK @@ -50,7 +54,9 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_destroy OK DEAL:0:LinearAlgebra::distributed::Vector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::Vector::test_length OK DEAL:0:LinearAlgebra::distributed::Vector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::Vector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::Vector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_product OK @@ -72,7 +78,9 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_destroy OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_length OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_product OK @@ -94,7 +102,9 @@ DEAL:0:TrilinosWrappers::MPI::Vector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::Vector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::Vector::test_length OK DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::Vector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_product OK @@ -116,7 +126,9 @@ DEAL:0:TrilinosWrappers::MPI::BlockVector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_length OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_product OK diff --git a/tests/sundials/n_vector.mpirun=1.with_trilinos=true.with_petsc=true.output b/tests/sundials/n_vector.mpirun=1.with_trilinos=true.with_petsc=true.output index 36ed4bc06b..a635b34129 100644 --- a/tests/sundials/n_vector.mpirun=1.with_trilinos=true.with_petsc=true.output +++ b/tests/sundials/n_vector.mpirun=1.with_trilinos=true.with_petsc=true.output @@ -6,7 +6,9 @@ DEAL:0:Vector::test_destroy OK DEAL:0:Vector::test_get_communicator OK DEAL:0:Vector::test_length OK DEAL:0:Vector::test_linear_sum OK +DEAL:0:Vector::test_linear_combination OK DEAL:0:Vector::test_dot_product OK +DEAL:0:Vector::test_dot_product_multi OK DEAL:0:Vector::test_set_constant OK DEAL:0:Vector::test_add_constant OK DEAL:0:Vector::test_elementwise_product OK @@ -28,7 +30,9 @@ DEAL:0:BlockVector::test_destroy OK DEAL:0:BlockVector::test_get_communicator OK DEAL:0:BlockVector::test_length OK DEAL:0:BlockVector::test_linear_sum OK +DEAL:0:BlockVector::test_linear_combination OK DEAL:0:BlockVector::test_dot_product OK +DEAL:0:BlockVector::test_dot_product_multi OK DEAL:0:BlockVector::test_set_constant OK DEAL:0:BlockVector::test_add_constant OK DEAL:0:BlockVector::test_elementwise_product OK @@ -50,7 +54,9 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_destroy OK DEAL:0:LinearAlgebra::distributed::Vector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::Vector::test_length OK DEAL:0:LinearAlgebra::distributed::Vector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::Vector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::Vector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_product OK @@ -72,7 +78,9 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_destroy OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_length OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_product OK @@ -94,7 +102,9 @@ DEAL:0:TrilinosWrappers::MPI::Vector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::Vector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::Vector::test_length OK DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::Vector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_product OK @@ -116,7 +126,9 @@ DEAL:0:TrilinosWrappers::MPI::BlockVector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_length OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_product OK @@ -138,7 +150,9 @@ DEAL:0:PETScWrappers::MPI::Vector::test_destroy OK DEAL:0:PETScWrappers::MPI::Vector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::Vector::test_length OK DEAL:0:PETScWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::Vector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::Vector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::Vector::test_set_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_add_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_product OK @@ -160,7 +174,9 @@ DEAL:0:PETScWrappers::MPI::BlockVector::test_destroy OK DEAL:0:PETScWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::BlockVector::test_length OK DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_product OK diff --git a/tests/sundials/n_vector.mpirun=2.with_trilinos=true.with_petsc=true.output b/tests/sundials/n_vector.mpirun=2.with_trilinos=true.with_petsc=true.output index 327f8d85a7..8011b1b7bf 100644 --- a/tests/sundials/n_vector.mpirun=2.with_trilinos=true.with_petsc=true.output +++ b/tests/sundials/n_vector.mpirun=2.with_trilinos=true.with_petsc=true.output @@ -6,7 +6,9 @@ DEAL:0:Vector::test_destroy OK DEAL:0:Vector::test_get_communicator OK DEAL:0:Vector::test_length OK DEAL:0:Vector::test_linear_sum OK +DEAL:0:Vector::test_linear_combination OK DEAL:0:Vector::test_dot_product OK +DEAL:0:Vector::test_dot_product_multi OK DEAL:0:Vector::test_set_constant OK DEAL:0:Vector::test_add_constant OK DEAL:0:Vector::test_elementwise_product OK @@ -28,7 +30,9 @@ DEAL:0:BlockVector::test_destroy OK DEAL:0:BlockVector::test_get_communicator OK DEAL:0:BlockVector::test_length OK DEAL:0:BlockVector::test_linear_sum OK +DEAL:0:BlockVector::test_linear_combination OK DEAL:0:BlockVector::test_dot_product OK +DEAL:0:BlockVector::test_dot_product_multi OK DEAL:0:BlockVector::test_set_constant OK DEAL:0:BlockVector::test_add_constant OK DEAL:0:BlockVector::test_elementwise_product OK @@ -50,7 +54,9 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_destroy OK DEAL:0:LinearAlgebra::distributed::Vector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::Vector::test_length OK DEAL:0:LinearAlgebra::distributed::Vector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::Vector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::Vector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::Vector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_product OK @@ -72,7 +78,9 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_destroy OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_get_communicator OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_length OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_sum OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_linear_combination OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_dot_product_multi OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_set_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_add_constant OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_product OK @@ -94,7 +102,9 @@ DEAL:0:TrilinosWrappers::MPI::Vector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::Vector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::Vector::test_length OK DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::Vector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_product OK @@ -116,7 +126,9 @@ DEAL:0:TrilinosWrappers::MPI::BlockVector::test_destroy OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_length OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_product OK @@ -138,7 +150,9 @@ DEAL:0:PETScWrappers::MPI::Vector::test_destroy OK DEAL:0:PETScWrappers::MPI::Vector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::Vector::test_length OK DEAL:0:PETScWrappers::MPI::Vector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::Vector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::Vector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::Vector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::Vector::test_set_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_add_constant OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_product OK @@ -160,7 +174,9 @@ DEAL:0:PETScWrappers::MPI::BlockVector::test_destroy OK DEAL:0:PETScWrappers::MPI::BlockVector::test_get_communicator OK DEAL:0:PETScWrappers::MPI::BlockVector::test_length OK DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_linear_combination OK DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:0:PETScWrappers::MPI::BlockVector::test_set_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_add_constant OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_product OK @@ -183,7 +199,9 @@ DEAL:1:Vector::test_destroy OK DEAL:1:Vector::test_get_communicator OK DEAL:1:Vector::test_length OK DEAL:1:Vector::test_linear_sum OK +DEAL:1:Vector::test_linear_combination OK DEAL:1:Vector::test_dot_product OK +DEAL:1:Vector::test_dot_product_multi OK DEAL:1:Vector::test_set_constant OK DEAL:1:Vector::test_add_constant OK DEAL:1:Vector::test_elementwise_product OK @@ -205,7 +223,9 @@ DEAL:1:BlockVector::test_destroy OK DEAL:1:BlockVector::test_get_communicator OK DEAL:1:BlockVector::test_length OK DEAL:1:BlockVector::test_linear_sum OK +DEAL:1:BlockVector::test_linear_combination OK DEAL:1:BlockVector::test_dot_product OK +DEAL:1:BlockVector::test_dot_product_multi OK DEAL:1:BlockVector::test_set_constant OK DEAL:1:BlockVector::test_add_constant OK DEAL:1:BlockVector::test_elementwise_product OK @@ -227,7 +247,9 @@ DEAL:1:LinearAlgebra::distributed::Vector::test_destroy OK DEAL:1:LinearAlgebra::distributed::Vector::test_get_communicator OK DEAL:1:LinearAlgebra::distributed::Vector::test_length OK DEAL:1:LinearAlgebra::distributed::Vector::test_linear_sum OK +DEAL:1:LinearAlgebra::distributed::Vector::test_linear_combination OK DEAL:1:LinearAlgebra::distributed::Vector::test_dot_product OK +DEAL:1:LinearAlgebra::distributed::Vector::test_dot_product_multi OK DEAL:1:LinearAlgebra::distributed::Vector::test_set_constant OK DEAL:1:LinearAlgebra::distributed::Vector::test_add_constant OK DEAL:1:LinearAlgebra::distributed::Vector::test_elementwise_product OK @@ -249,7 +271,9 @@ DEAL:1:LinearAlgebra::distributed::BlockVector::test_destroy OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_get_communicator OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_length OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_linear_sum OK +DEAL:1:LinearAlgebra::distributed::BlockVector::test_linear_combination OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_dot_product OK +DEAL:1:LinearAlgebra::distributed::BlockVector::test_dot_product_multi OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_set_constant OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_add_constant OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_elementwise_product OK @@ -271,7 +295,9 @@ DEAL:1:TrilinosWrappers::MPI::Vector::test_destroy OK DEAL:1:TrilinosWrappers::MPI::Vector::test_get_communicator OK DEAL:1:TrilinosWrappers::MPI::Vector::test_length OK DEAL:1:TrilinosWrappers::MPI::Vector::test_linear_sum OK +DEAL:1:TrilinosWrappers::MPI::Vector::test_linear_combination OK DEAL:1:TrilinosWrappers::MPI::Vector::test_dot_product OK +DEAL:1:TrilinosWrappers::MPI::Vector::test_dot_product_multi OK DEAL:1:TrilinosWrappers::MPI::Vector::test_set_constant OK DEAL:1:TrilinosWrappers::MPI::Vector::test_add_constant OK DEAL:1:TrilinosWrappers::MPI::Vector::test_elementwise_product OK @@ -293,7 +319,9 @@ DEAL:1:TrilinosWrappers::MPI::BlockVector::test_destroy OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_get_communicator OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_length OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:1:TrilinosWrappers::MPI::BlockVector::test_linear_combination OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_dot_product OK +DEAL:1:TrilinosWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_set_constant OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_add_constant OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_elementwise_product OK @@ -315,7 +343,9 @@ DEAL:1:PETScWrappers::MPI::Vector::test_destroy OK DEAL:1:PETScWrappers::MPI::Vector::test_get_communicator OK DEAL:1:PETScWrappers::MPI::Vector::test_length OK DEAL:1:PETScWrappers::MPI::Vector::test_linear_sum OK +DEAL:1:PETScWrappers::MPI::Vector::test_linear_combination OK DEAL:1:PETScWrappers::MPI::Vector::test_dot_product OK +DEAL:1:PETScWrappers::MPI::Vector::test_dot_product_multi OK DEAL:1:PETScWrappers::MPI::Vector::test_set_constant OK DEAL:1:PETScWrappers::MPI::Vector::test_add_constant OK DEAL:1:PETScWrappers::MPI::Vector::test_elementwise_product OK @@ -337,7 +367,9 @@ DEAL:1:PETScWrappers::MPI::BlockVector::test_destroy OK DEAL:1:PETScWrappers::MPI::BlockVector::test_get_communicator OK DEAL:1:PETScWrappers::MPI::BlockVector::test_length OK DEAL:1:PETScWrappers::MPI::BlockVector::test_linear_sum OK +DEAL:1:PETScWrappers::MPI::BlockVector::test_linear_combination OK DEAL:1:PETScWrappers::MPI::BlockVector::test_dot_product OK +DEAL:1:PETScWrappers::MPI::BlockVector::test_dot_product_multi OK DEAL:1:PETScWrappers::MPI::BlockVector::test_set_constant OK DEAL:1:PETScWrappers::MPI::BlockVector::test_add_constant OK DEAL:1:PETScWrappers::MPI::BlockVector::test_elementwise_product OK