From: bangerth Date: Mon, 2 Sep 2013 22:03:39 +0000 (+0000) Subject: Add tests for: Patch by Fahad Alrasched: Add extract_subvector_to() function to all... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d1acba1e2822d23cbac5f10129cf8acc8e99d66c;p=dealii-svn.git Add tests for: Patch by Fahad Alrasched: Add extract_subvector_to() function to all vectors. Use it in DoFAccessor::get_dof_values(). git-svn-id: https://svn.dealii.org/trunk@30556 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/gla/extract_subvector_to.cc b/tests/gla/extract_subvector_to.cc new file mode 100644 index 0000000000..e96721758c --- /dev/null +++ b/tests/gla/extract_subvector_to.cc @@ -0,0 +1,128 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2004 - 2013 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. +// +// --------------------------------------------------------------------- + + + +// test the various VECTOR::extract_subvector_to functions + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + + +template +void test (Vector &vector) +{ + const unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + + for (unsigned int i=0; i indices; + for (unsigned int j=0; j values1 (indices.size()); + vector.extract_subvector_to (indices, values1); + for (unsigned int j=0; j values2 (indices.size()); + vector.extract_subvector_to (indices.begin(), + indices.end(), + values2.begin()); + for (unsigned int j=0; j v(17); + test (v); + deallog.pop(); + } + + { + deallog.push("PETSc"); + PETScWrappers::Vector v(17); + test (v); + deallog.pop(); + } + + { + deallog.push("Trilinos"); + TrilinosWrappers::Vector v(17); + test (v); + deallog.pop(); + } + + + { + deallog.push("deal.II"); + BlockVector v(3); + v.block(0).reinit(7); + v.block(1).reinit(5); + v.block(2).reinit(3); + v.collect_sizes(); + test (v); + deallog.pop(); + } + + { + deallog.push("PETSc"); + PETScWrappers::BlockVector v(3); + v.block(0).reinit(7); + v.block(1).reinit(5); + v.block(2).reinit(3); + v.collect_sizes(); + test (v); + deallog.pop(); + } + + { + deallog.push("Trilinos"); + TrilinosWrappers::BlockVector v(3); + v.block(0).reinit(7); + v.block(1).reinit(5); + v.block(2).reinit(3); + v.collect_sizes(); + test (v); + deallog.pop(); + } + + } + +} diff --git a/tests/gla/extract_subvector_to/ncpu_1/cmp/generic b/tests/gla/extract_subvector_to/ncpu_1/cmp/generic new file mode 100644 index 0000000000..784d66bf3f --- /dev/null +++ b/tests/gla/extract_subvector_to/ncpu_1/cmp/generic @@ -0,0 +1,7 @@ + +DEAL:0:deal.II::OK +DEAL:0:PETSc::OK +DEAL:0:Trilinos::OK +DEAL:0:deal.II::OK +DEAL:0:PETSc::OK +DEAL:0:Trilinos::OK diff --git a/tests/gla/extract_subvector_to_parallel.cc b/tests/gla/extract_subvector_to_parallel.cc new file mode 100644 index 0000000000..c859a1a9b0 --- /dev/null +++ b/tests/gla/extract_subvector_to_parallel.cc @@ -0,0 +1,183 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2004 - 2013 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. +// +// --------------------------------------------------------------------- + + + +// test the various VECTOR::extract_subvector_to functions for +// parallel vectors and block vectors + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + + +template +void set (Vector &vector) +{ + for (unsigned int i=0; i +void test (Vector &vector) +{ + const unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + + // select every other element + std::vector indices; + for (unsigned int j=0; j values1 (indices.size()); + vector.extract_subvector_to (indices, values1); + for (unsigned int j=0; j values2 (indices.size()); + vector.extract_subvector_to (indices.begin(), + indices.end(), + values2.begin()); + for (unsigned int j=0; j w(local, dense_local, MPI_COMM_WORLD); + // set (w); + // parallel::distributed::Vector v(local, dense_local, MPI_COMM_WORLD); + // v = w; // get copy of vector including ghost elements + // test (v); + // deallog.pop(); + // } + + { + deallog.push("PETSc"); + PETScWrappers::MPI::Vector w(local, MPI_COMM_WORLD); + set (w); + PETScWrappers::MPI::Vector v(local, dense_local, MPI_COMM_WORLD); + v = w; // get copy of vector including ghost elements + test (v); + deallog.pop(); + } + + { + deallog.push("Trilinos"); + TrilinosWrappers::MPI::Vector w(local, MPI_COMM_WORLD); + set (w); + TrilinosWrappers::MPI::Vector v(local, dense_local, MPI_COMM_WORLD); + v = w; // get copy of vector including ghost elements + test (v); + deallog.pop(); + } + + + std::vector partitioning; + { + IndexSet block1(10); + if (myid==0) + block1.add_range(0,7); + if (myid==1) + block1.add_range(7,10); + + IndexSet block2(6); + if (myid==0) + block2.add_range(0,2); + if (myid==1) + block2.add_range(2,6); + + partitioning.push_back(block1); + partitioning.push_back(block2); + } + + std::vector dense_partitioning; + { + IndexSet block1(10); + block1.add_range(0,10); + + IndexSet block2(6); + block2.add_range(0,6); + + dense_partitioning.push_back(block1); + dense_partitioning.push_back(block2); + } + + + // { + // deallog.push("deal.II"); + // parallel::distributed::BlockVector w(partitioning, MPI_COMM_WORLD); + // set (w); + // parallel::distributed::BlockVector v(partitioning, dense_partitioning, MPI_COMM_WORLD); + // v = w; // get copy of vector including ghost elements + // test (v); + // deallog.pop(); + // } + + { + deallog.push("PETSc"); + PETScWrappers::MPI::BlockVector w(partitioning, MPI_COMM_WORLD); + set (w); + PETScWrappers::MPI::BlockVector v(partitioning, dense_partitioning, MPI_COMM_WORLD); + v = w; // get copy of vector including ghost elements + test (v); + deallog.pop(); + } + + { + deallog.push("Trilinos"); + TrilinosWrappers::MPI::BlockVector w(partitioning, MPI_COMM_WORLD); + set (w); + TrilinosWrappers::MPI::BlockVector v(partitioning, dense_partitioning, MPI_COMM_WORLD); + v = w; // get copy of vector including ghost elements + test (v); + deallog.pop(); + } + + } + +} diff --git a/tests/gla/extract_subvector_to_parallel/ncpu_2/cmp/generic b/tests/gla/extract_subvector_to_parallel/ncpu_2/cmp/generic new file mode 100644 index 0000000000..ae4bdbf19f --- /dev/null +++ b/tests/gla/extract_subvector_to_parallel/ncpu_2/cmp/generic @@ -0,0 +1,7 @@ + +DEAL:0:PETSc::OK +DEAL:0:Trilinos::OK +DEAL:0:PETSc::OK +DEAL:0:Trilinos::OK + +