*/
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
bool
SparsityPattern::is_compressed () const
{
- return compressed;
+ return graph->Filled();
}
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)
// 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),
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,
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);
SparsityPattern::SparsityPattern (const IndexSet ¶llel_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);
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);
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);
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);
"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
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,
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."));
graph->FillComplete();
nonlocal_graph.reset();
-
- compressed = true;
}
ierr = graph->OptimizeStorage ();
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-
- compressed = true;
}
--- /dev/null
+// ---------------------------------------------------------------------
+// $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;
+ };
+}
--- /dev/null
+
+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