]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify TrilinosSparsityPattern::is_compressed(), add test
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 12 Jun 2014 08:03:40 +0000 (08:03 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 12 Jun 2014 08:03:40 +0000 (08:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@33039 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/trilinos_sparsity_pattern.h
deal.II/source/lac/trilinos_sparsity_pattern.cc
tests/trilinos/trilinos_sparsity_pattern.cc [new file with mode: 0644]
tests/trilinos/trilinos_sparsity_pattern.output [new file with mode: 0644]

index 01d6d16676600c06a878ff43313b3da793c1a575..83644e05008d204e656e4e83a0745497e3d415d6 100644 (file)
@@ -1094,12 +1094,6 @@ namespace TrilinosWrappers
      */
     std_cxx1x::shared_ptr<Epetra_Map> 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<int>(end - begin);
-    compressed = false;
 
     int ierr;
     if ( graph->RowMap().LID(static_cast<TrilinosWrappers::types::int_type>(row)) != -1)
index b406c7300f20e1cc5b971e2e86a230337a80d345..072b1ceaa2d29afc3ba4a3b4c4b9b583c2e12d5a 100644 (file)
@@ -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  &parallel_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     &parallel_partitioning,
                                      const MPI_Comm     &communicator,
                                      const std::vector<size_type> &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<size_type> &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<size_type>(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<size_type>(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 (file)
index 0000000..501a9fd
--- /dev/null
@@ -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 <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <fstream>
+#include <iomanip>
+
+
+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 (file)
index 0000000..245037a
--- /dev/null
@@ -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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.