From: David Wells Date: Sun, 17 Sep 2017 20:55:58 +0000 (-0400) Subject: Add a test for copying SUNDIALS vectors. X-Git-Tag: v9.0.0-rc1~1062^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bd1c495aa1ff95c39bf3c8dab184be718c92a544;p=dealii.git Add a test for copying SUNDIALS vectors. --- diff --git a/tests/sundials/copy_01.cc b/tests/sundials/copy_01.cc new file mode 100644 index 0000000000..11cb7a68ca --- /dev/null +++ b/tests/sundials/copy_01.cc @@ -0,0 +1,73 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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. +// +// --------------------------------------------------------------------- + +// Verify that we can correctly copy a PETSc vector to a parallel (non-PETSc) +// SUNDIALS vector. + +#include "../tests.h" + +#include +#include +#include + +#include +#include + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + MPILogInitAll log_all; + + const unsigned int current_process = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + const unsigned int n_processes = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + // give the zeroth process 1 dof, the first 2, etc + const types::global_dof_index n_local_dofs = 1 + current_process; + const types::global_dof_index begin_dof = (current_process*(current_process + 1))/2; + const types::global_dof_index end_dof = begin_dof + n_local_dofs; + + const types::global_dof_index n_global_dofs = ((n_processes - 1)*(n_processes))/2 + + n_processes; + IndexSet local_dofs(n_global_dofs); + local_dofs.add_range(begin_dof, end_dof); + local_dofs.compress(); + AssertDimension(local_dofs.n_elements(), n_local_dofs); + + + PETScWrappers::MPI::Vector vec(local_dofs, MPI_COMM_WORLD); + for (const types::global_dof_index local_dof : local_dofs) + { + vec[local_dof] = 2*local_dof - current_process; + } + vec.compress(VectorOperation::insert); + const std::size_t n_dofs = vec.size(); + + N_Vector sundials_vector = N_VNew_Parallel(MPI_COMM_WORLD, + n_local_dofs, + n_dofs); + SUNDIALS::internal::copy(sundials_vector, vec); + + PETScWrappers::MPI::Vector vec2(local_dofs, MPI_COMM_WORLD); + SUNDIALS::internal::copy(vec2, sundials_vector); + + AssertThrow(vec2 == vec, + ExcMessage("The two PETSc vectors should be equal since they are " + "copies (by means of an intermediate SUNDIALs vector)")); + + deallog << "n_local_dofs: " + << n_local_dofs + << std::endl; + deallog << "OK" << std::endl; +} diff --git a/tests/sundials/copy_01.with_mpi=true.with_petsc=true.mpirun=4.output b/tests/sundials/copy_01.with_mpi=true.with_petsc=true.mpirun=4.output new file mode 100644 index 0000000000..ea6a78c882 --- /dev/null +++ b/tests/sundials/copy_01.with_mpi=true.with_petsc=true.mpirun=4.output @@ -0,0 +1,15 @@ + +DEAL:0::n_local_dofs: 1 +DEAL:0::OK + +DEAL:1::n_local_dofs: 2 +DEAL:1::OK + + +DEAL:2::n_local_dofs: 3 +DEAL:2::OK + + +DEAL:3::n_local_dofs: 4 +DEAL:3::OK +