From eca3f638c874033865c707a4a3b7dee6405cd714 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Thu, 12 Jun 2014 08:03:40 +0000 Subject: [PATCH] Simplify TrilinosSparsityPattern::is_compressed(), add test git-svn-id: https://svn.dealii.org/trunk@33039 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/lac/trilinos_sparsity_pattern.h | 9 +- .../source/lac/trilinos_sparsity_pattern.cc | 20 ---- tests/trilinos/trilinos_sparsity_pattern.cc | 110 ++++++++++++++++++ .../trilinos/trilinos_sparsity_pattern.output | 21 ++++ 4 files changed, 132 insertions(+), 28 deletions(-) create mode 100644 tests/trilinos/trilinos_sparsity_pattern.cc create mode 100644 tests/trilinos/trilinos_sparsity_pattern.output diff --git a/deal.II/include/deal.II/lac/trilinos_sparsity_pattern.h b/deal.II/include/deal.II/lac/trilinos_sparsity_pattern.h index 01d6d16676..83644e0500 100644 --- a/deal.II/include/deal.II/lac/trilinos_sparsity_pattern.h +++ b/deal.II/include/deal.II/lac/trilinos_sparsity_pattern.h @@ -1094,12 +1094,6 @@ namespace TrilinosWrappers */ std_cxx1x::shared_ptr column_space_map; - /** - * A boolean variable to hold information on whether the vector is - * compressed or not. - */ - bool compressed; - /** * A sparsity pattern object in Trilinos to be used for finite element * based problems which allows for adding non-local elements to the @@ -1363,7 +1357,7 @@ namespace TrilinosWrappers bool SparsityPattern::is_compressed () const { - return compressed; + return graph->Filled(); } @@ -1413,7 +1407,6 @@ namespace TrilinosWrappers TrilinosWrappers::types::int_type *col_index_ptr = (TrilinosWrappers::types::int_type *)(&*begin); const int n_cols = static_cast(end - begin); - compressed = false; int ierr; if ( graph->RowMap().LID(static_cast(row)) != -1) diff --git a/deal.II/source/lac/trilinos_sparsity_pattern.cc b/deal.II/source/lac/trilinos_sparsity_pattern.cc index b406c7300f..072b1ceaa2 100644 --- a/deal.II/source/lac/trilinos_sparsity_pattern.cc +++ b/deal.II/source/lac/trilinos_sparsity_pattern.cc @@ -167,8 +167,6 @@ namespace TrilinosWrappers // MPI will still get a parallel // interface. SparsityPattern::SparsityPattern () - : - compressed (true) { column_space_map.reset(new Epetra_Map (TrilinosWrappers::types::int_type(0), TrilinosWrappers::types::int_type(0), @@ -240,7 +238,6 @@ namespace TrilinosWrappers column_space_map (new Epetra_Map(TrilinosWrappers::types::int_type(0), TrilinosWrappers::types::int_type(0), Utilities::Trilinos::comm_self())), - compressed (false), graph (new Epetra_FECrsGraph(View, *column_space_map, *column_space_map, @@ -255,8 +252,6 @@ namespace TrilinosWrappers SparsityPattern::SparsityPattern (const IndexSet ¶llel_partitioning, const MPI_Comm &communicator, const size_type n_entries_per_row) - : - compressed (false) { Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator, false); @@ -268,8 +263,6 @@ namespace TrilinosWrappers SparsityPattern::SparsityPattern (const IndexSet ¶llel_partitioning, const MPI_Comm &communicator, const std::vector &n_entries_per_row) - : - compressed (false) { Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator, false); @@ -282,8 +275,6 @@ namespace TrilinosWrappers const IndexSet &col_parallel_partitioning, const MPI_Comm &communicator, const size_type n_entries_per_row) - : - compressed (false) { Epetra_Map row_map = row_parallel_partitioning.make_trilinos_map (communicator, false); @@ -299,8 +290,6 @@ namespace TrilinosWrappers const IndexSet &col_parallel_partitioning, const MPI_Comm &communicator, const std::vector &n_entries_per_row) - : - compressed (false) { Epetra_Map row_map = row_parallel_partitioning.make_trilinos_map (communicator, false); @@ -317,8 +306,6 @@ namespace TrilinosWrappers const IndexSet &writable_rows, const MPI_Comm &communicator, const size_type n_max_entries_per_row) - : - compressed (false) { reinit (row_parallel_partitioning, col_parallel_partitioning, writable_rows, communicator, n_max_entries_per_row); @@ -366,7 +353,6 @@ namespace TrilinosWrappers "the maps of different processors.")); graph.reset (); column_space_map.reset (new Epetra_Map (input_col_map)); - compressed = false; // for more than one processor, need to // specify only row map first and let the @@ -423,7 +409,6 @@ namespace TrilinosWrappers static_cast(n_global_elements(input_row_map))); column_space_map.reset (new Epetra_Map (input_col_map)); - compressed = false; if (input_row_map.Comm().NumProc() > 1) graph.reset(new Epetra_FECrsGraph(Copy, input_row_map, @@ -609,7 +594,6 @@ namespace TrilinosWrappers static_cast(n_global_elements(input_col_map))); column_space_map.reset (new Epetra_Map (input_col_map)); - compressed = false; Assert (input_row_map.LinearMap() == true, ExcMessage ("This function is not efficient if the map is not contiguous.")); @@ -708,8 +692,6 @@ namespace TrilinosWrappers graph->FillComplete(); nonlocal_graph.reset(); - - compressed = true; } @@ -749,8 +731,6 @@ namespace TrilinosWrappers ierr = graph->OptimizeStorage (); AssertThrow (ierr == 0, ExcTrilinosError(ierr)); - - compressed = true; } diff --git a/tests/trilinos/trilinos_sparsity_pattern.cc b/tests/trilinos/trilinos_sparsity_pattern.cc new file mode 100644 index 0000000000..501a9fd810 --- /dev/null +++ b/tests/trilinos/trilinos_sparsity_pattern.cc @@ -0,0 +1,110 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2014 - 2014 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. +// +// --------------------------------------------------------------------- + + + +// Tests basic stuff of Trilinos sparsity patterns + +#include "../tests.h" +#include +#include +#include + + +void test () +{ + TrilinosWrappers::SparsityPattern sp; + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + + deallog << "Creating entries..." << std::endl; + + sp.reinit(5,7,3); + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + + for (unsigned int i=0; i<5; ++i) + for (unsigned int j=0; j<7; ++j) + if ((i+2*j+1) % 3 == 0) + sp.add (i,j); + + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + + sp.compress (); + + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + deallog << "Number of entries: " << sp.n_nonzero_elements() << std::endl; + deallog << "Number of rows: " << sp.n_rows() << std::endl; + deallog << "Number of colums: " << sp.n_cols() << std::endl; + deallog << "Local size: " << sp.local_size() << std::endl; + deallog << "Max row length: " << sp.max_entries_per_row() << std::endl; + deallog << "SP::row_length(0): " << sp.row_length(0) << std::endl; + deallog << "Bandwidth: " << sp.bandwidth() << std::endl; + deallog << "SP::empty(): " << sp.empty() << std::endl; + + sp.compress (); + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + + deallog << "Clearing..." << std::endl; + + sp.clear(); + + deallog << "SP::is_compressed(): " << sp.is_compressed() << std::endl; + deallog << "Bandwidth: " << sp.bandwidth() << std::endl; + deallog << "SP::empty(): " << sp.empty() << std::endl; + deallog << "Number of rows: " << sp.n_rows() << std::endl; + + deallog << "OK" << std::endl; +} + + + +int main (int argc, char **argv) +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv); + + try + { + test (); + } + catch (std::exception &exc) + { + deallog << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + deallog << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + + return 1; + } + catch (...) + { + deallog << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + deallog << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} diff --git a/tests/trilinos/trilinos_sparsity_pattern.output b/tests/trilinos/trilinos_sparsity_pattern.output new file mode 100644 index 0000000000..245037ab2c --- /dev/null +++ b/tests/trilinos/trilinos_sparsity_pattern.output @@ -0,0 +1,21 @@ + +DEAL::SP::is_compressed(): 1 +DEAL::Creating entries... +DEAL::SP::is_compressed(): 0 +DEAL::SP::is_compressed(): 0 +DEAL::SP::is_compressed(): 1 +DEAL::Number of entries: 11 +DEAL::Number of rows: 5 +DEAL::Number of colums: 7 +DEAL::Local size: 5 +DEAL::Max row length: 3 +DEAL::SP::row_length(0): 2 +DEAL::Bandwidth: 4 +DEAL::SP::empty(): 0 +DEAL::SP::is_compressed(): 1 +DEAL::Clearing... +DEAL::SP::is_compressed(): 1 +DEAL::Bandwidth: 0 +DEAL::SP::empty(): 1 +DEAL::Number of rows: 0 +DEAL::OK -- 2.39.5