From a119fabe599c399c8d076e51ee67875ca909a71e Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Mon, 12 Jan 2015 11:30:29 -0500 Subject: [PATCH] add test for reinit() and copy_from() --- tests/gla/mat_05.cc | 138 +++++++++++++++++++++++++++++++ tests/gla/mat_05.mpirun=3.output | 49 +++++++++++ 2 files changed, 187 insertions(+) create mode 100644 tests/gla/mat_05.cc create mode 100644 tests/gla/mat_05.mpirun=3.output diff --git a/tests/gla/mat_05.cc b/tests/gla/mat_05.cc new file mode 100644 index 0000000000..6ab75b11b1 --- /dev/null +++ b/tests/gla/mat_05.cc @@ -0,0 +1,138 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// reinit and copy_from for LA::MPI::SparseMatrix + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + +#include "gla.h" + +template +void test () +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + if (myid==0) + deallog << "numproc=" << numproc << std::endl; + + // each processor owns 2 indices and all + // are ghosting Element 1 (the second) + + IndexSet local_active(numproc*2); + local_active.add_range(myid*2,myid*2+2); + IndexSet local_relevant= local_active; + local_relevant.add_range(0,1); + + CompressedSimpleSparsityPattern csp (local_relevant); + + for (unsigned int i=0; i<2*numproc; ++i) + if (local_relevant.is_element(i)) + csp.add(i,i); + + csp.add(0,1); + + + typename LA::MPI::SparseMatrix mat; + mat.reinit (local_active, local_active, csp, MPI_COMM_WORLD); + + Assert(mat.n()==numproc*2, ExcInternalError()); + Assert(mat.m()==numproc*2, ExcInternalError()); + + // set local values + mat.set(myid*2,myid*2, myid*2.0); + mat.set(myid*2+1,myid*2+1, myid*2.0+1.0); + + mat.compress(VectorOperation::insert); + + mat.add(0,1,1.0); + + mat.compress(VectorOperation::add); + + deallog << "l1-norm: " << mat.l1_norm() << std::endl; + if (myid==0) + deallog << "mat(0,1): " << mat(0,1) << std::endl; + + { + // disabled: + //typename LA::MPI::SparseMatrix mat2(mat); + } + { + // disabled: + //typename LA::MPI::SparseMatrix mat2; + //mat2 = mat; + } + + { + deallog << "* reinit(other):" << std::endl; + typename LA::MPI::SparseMatrix mat2; + mat2.reinit(mat); + deallog << "l1-norm: " << mat2.l1_norm() << std::endl; + if (myid==0) + deallog << "mat(0,1): " << mat2(0,1) << std::endl; + mat2.add(0,1,0.3); + mat2.compress(VectorOperation::add); + if (myid==0) + deallog << "after adding mat(0,1): " << mat2(0,1) << std::endl; + } + { + deallog << "* copy_from():" << std::endl; + typename LA::MPI::SparseMatrix mat2; + // PETSc needs a reinit() befor copy_from() + mat2.reinit (local_active, local_active, csp, MPI_COMM_WORLD); + mat2.copy_from(mat); + deallog << "l1-norm: " << mat2.l1_norm() << std::endl; + if (myid==0) + deallog << "mat(0,1): " << mat2(0,1) << std::endl; + mat2.add(0,1,0.3); + mat2.compress(VectorOperation::add); + if (myid==0) + deallog << "after adding mat(0,1): " << mat2(0,1) << std::endl; + } + + // done + if (myid==0) + deallog << "OK" << std::endl; +} + + + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + { + deallog.push("PETSc"); + test(); + deallog.pop(); + deallog.push("Trilinos"); + test(); + deallog.pop(); + } + + // compile, don't run + //if (myid==9999) + // test(); + + +} diff --git a/tests/gla/mat_05.mpirun=3.output b/tests/gla/mat_05.mpirun=3.output new file mode 100644 index 0000000000..b32c68e1e6 --- /dev/null +++ b/tests/gla/mat_05.mpirun=3.output @@ -0,0 +1,49 @@ + +DEAL:0:PETSc::numproc=3 +DEAL:0:PETSc::l1-norm: 5.00000 +DEAL:0:PETSc::mat(0,1): 3.00000 +DEAL:0:PETSc::* reinit(other): +DEAL:0:PETSc::l1-norm: 0.00000 +DEAL:0:PETSc::mat(0,1): 0.00000 +DEAL:0:PETSc::after adding mat(0,1): 0.900000 +DEAL:0:PETSc::* copy_from(): +DEAL:0:PETSc::l1-norm: 5.00000 +DEAL:0:PETSc::mat(0,1): 3.00000 +DEAL:0:PETSc::after adding mat(0,1): 3.90000 +DEAL:0:PETSc::OK +DEAL:0:Trilinos::numproc=3 +DEAL:0:Trilinos::l1-norm: 5.00000 +DEAL:0:Trilinos::mat(0,1): 3.00000 +DEAL:0:Trilinos::* reinit(other): +DEAL:0:Trilinos::l1-norm: 0.00000 +DEAL:0:Trilinos::mat(0,1): 0.00000 +DEAL:0:Trilinos::after adding mat(0,1): 0.900000 +DEAL:0:Trilinos::* copy_from(): +DEAL:0:Trilinos::l1-norm: 5.00000 +DEAL:0:Trilinos::mat(0,1): 3.00000 +DEAL:0:Trilinos::after adding mat(0,1): 3.90000 +DEAL:0:Trilinos::OK + +DEAL:1:PETSc::l1-norm: 5.00000 +DEAL:1:PETSc::* reinit(other): +DEAL:1:PETSc::l1-norm: 0.00000 +DEAL:1:PETSc::* copy_from(): +DEAL:1:PETSc::l1-norm: 5.00000 +DEAL:1:Trilinos::l1-norm: 5.00000 +DEAL:1:Trilinos::* reinit(other): +DEAL:1:Trilinos::l1-norm: 0.00000 +DEAL:1:Trilinos::* copy_from(): +DEAL:1:Trilinos::l1-norm: 5.00000 + + +DEAL:2:PETSc::l1-norm: 5.00000 +DEAL:2:PETSc::* reinit(other): +DEAL:2:PETSc::l1-norm: 0.00000 +DEAL:2:PETSc::* copy_from(): +DEAL:2:PETSc::l1-norm: 5.00000 +DEAL:2:Trilinos::l1-norm: 5.00000 +DEAL:2:Trilinos::* reinit(other): +DEAL:2:Trilinos::l1-norm: 0.00000 +DEAL:2:Trilinos::* copy_from(): +DEAL:2:Trilinos::l1-norm: 5.00000 + -- 2.39.5