From 95babf0f543aa4fc91b4a917613accd6b5b7579a Mon Sep 17 00:00:00 2001 From: kronbichler Date: Fri, 30 Jan 2009 08:15:18 +0000 Subject: [PATCH] Optimized constructor of Trilinos matrices on one processor. This makes assembly and some other initialization operations a bit faster. git-svn-id: https://svn.dealii.org/trunk@18322 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/source/trilinos_sparse_matrix.cc | 79 ++++++--- .../lac/source/trilinos_sparsity_pattern.cc | 150 +++++++++++++----- 2 files changed, 168 insertions(+), 61 deletions(-) diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc index 5af9a75a70..580bcb552b 100755 --- a/deal.II/lac/source/trilinos_sparse_matrix.cc +++ b/deal.II/lac/source/trilinos_sparse_matrix.cc @@ -93,7 +93,7 @@ namespace TrilinosWrappers last_action (Zero), compressed (true), matrix (std::auto_ptr - (new Epetra_FECrsMatrix(View, row_map, 0))) + (new Epetra_FECrsMatrix(View, row_map, col_map, 0))) { matrix->FillComplete(); } @@ -153,7 +153,7 @@ namespace TrilinosWrappers last_action (Zero), compressed (true), matrix (std::auto_ptr - (new Epetra_FECrsMatrix(Copy, row_map, + (new Epetra_FECrsMatrix(Copy, row_map, (int*)const_cast(&(n_entries_per_row[0])), false))) {} @@ -173,9 +173,18 @@ namespace TrilinosWrappers #endif last_action (Zero), compressed (true), + // on one processor only, we know how the + // columns of the matrix will be + // distributed (everything on one + // processor), so we can hand in this + // information to the constructor. we + // can't do so in parallel, where the + // information from columns is only + // available when entries have been added matrix (std::auto_ptr - (new Epetra_FECrsMatrix(Copy, row_map, - int(n_max_entries_per_row), false))) + (new Epetra_FECrsMatrix(Copy, row_map, col_map, + int(n_max_entries_per_row), + false))) {} @@ -194,9 +203,9 @@ namespace TrilinosWrappers last_action (Zero), compressed (true), matrix (std::auto_ptr - (new Epetra_FECrsMatrix(Copy, row_map, - (int*)const_cast(&(n_entries_per_row[0])), - false))) + (new Epetra_FECrsMatrix(Copy, row_map, col_map, + (int*)const_cast(&(n_entries_per_row[0])), + false))) {} @@ -232,8 +241,7 @@ namespace TrilinosWrappers SparseMatrix::~SparseMatrix () - { - } + {} @@ -331,9 +339,26 @@ namespace TrilinosWrappers // it were square. The call to // GlobalAssemble later will set the // correct values. - Epetra_CrsGraph graph (Copy, row_map, - &n_entries_per_row[input_row_map.MinMyGID()], true); + // for more than one processor, need to + // specify only row map first and let the + // matrix entries decide about the column + // map (which says which columns are + // present in the matrix, not to be + // confused with the col_map that tells + // how the domain dofs of the matrix will + // be distributed). for only one + // processor, we can directly assign the + // columns as well. + std::auto_ptr graph; + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_CrsGraph (Copy, row_map, + &n_entries_per_row[input_row_map.MinMyGID()], true)); + else + graph = std::auto_ptr + (new Epetra_CrsGraph (Copy, row_map, col_map, + &n_entries_per_row[input_row_map.MinMyGID()], true)); // TODO: As of now, assume that the // sparsity pattern sits at the all @@ -341,8 +366,8 @@ namespace TrilinosWrappers // each processor set its rows. Since // this is wasteful, a better solution // should be found in the future. - Assert (graph.NumGlobalRows() == (int)sparsity_pattern.n_rows(), - ExcDimensionMismatch (graph.NumGlobalRows(), + Assert (graph->NumGlobalRows() == (int)sparsity_pattern.n_rows(), + ExcDimensionMismatch (graph->NumGlobalRows(), sparsity_pattern.n_rows())); // Trilinos has a bug for rectangular @@ -368,16 +393,16 @@ namespace TrilinosWrappers for (int col=0; col < row_length; ++col) row_indices[col] = sparsity_pattern.column_number (row, col); - graph.InsertGlobalIndices (row, row_length, &row_indices[0]); + graph->InsertGlobalIndices (row, row_length, &row_indices[0]); } // Now, fill the graph (sort indices, make // memory contiguous, etc). - graph.FillComplete(col_map, row_map); + graph->FillComplete(col_map, row_map); // And now finally generate the matrix. matrix = std::auto_ptr - (new Epetra_FECrsMatrix(Copy, graph, false)); + (new Epetra_FECrsMatrix(Copy, *graph, false)); last_action = Zero; @@ -443,9 +468,15 @@ namespace TrilinosWrappers // it were square. The call to // GlobalAssemble later will set the // correct values. - Epetra_CrsGraph graph (Copy, row_map, - &n_entries_per_row[input_row_map.MinMyGID()], true); - + std::auto_ptr graph; + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_CrsGraph (Copy, row_map, + &n_entries_per_row[input_row_map.MinMyGID()], true)); + else + graph = std::auto_ptr + (new Epetra_CrsGraph (Copy, row_map, col_map, + &n_entries_per_row[input_row_map.MinMyGID()], true)); // TODO: As of now, assume that the // sparsity pattern sits at the all @@ -453,8 +484,8 @@ namespace TrilinosWrappers // each processor set its rows. Since // this is wasteful, a better solution // should be found in the future. - Assert (graph.NumGlobalRows() == (int)sparsity_pattern.n_rows(), - ExcDimensionMismatch (graph.NumGlobalRows(), + Assert (graph->NumGlobalRows() == (int)sparsity_pattern.n_rows(), + ExcDimensionMismatch (graph->NumGlobalRows(), sparsity_pattern.n_rows())); // Trilinos has a bug for rectangular @@ -485,16 +516,16 @@ namespace TrilinosWrappers ++col_num, ++col) row_indices[col] = *col_num; - graph.InsertGlobalIndices (row, row_length, &row_indices[0]); + graph->InsertGlobalIndices (row, row_length, &row_indices[0]); } // Now, fill the graph (sort indices, make // memory contiguous, etc). - graph.FillComplete(col_map, row_map); + graph->FillComplete(col_map, row_map); // And now finally generate the matrix. matrix = std::auto_ptr - (new Epetra_FECrsMatrix(Copy, graph, false)); + (new Epetra_FECrsMatrix(Copy, *graph, false)); last_action = Zero; diff --git a/deal.II/lac/source/trilinos_sparsity_pattern.cc b/deal.II/lac/source/trilinos_sparsity_pattern.cc index 7b005997cf..a4de59a8f6 100755 --- a/deal.II/lac/source/trilinos_sparsity_pattern.cc +++ b/deal.II/lac/source/trilinos_sparsity_pattern.cc @@ -90,7 +90,7 @@ namespace TrilinosWrappers col_map (row_map), compressed (true), graph (std::auto_ptr - (new Epetra_FECrsGraph(View, row_map, 0))), + (new Epetra_FECrsGraph(View, row_map, col_map, 0))), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -104,24 +104,51 @@ namespace TrilinosWrappers row_map (InputMap), col_map (row_map), compressed (false), - graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - int(n_entries_per_row), false))), + // for more than one processor, need to + // specify only row map first and let the + // matrix entries decide about the column + // map (which says which columns are + // present in the matrix, not to be + // confused with the col_map that tells + // how the domain dofs of the matrix will + // be distributed). for only one + // processor, we can directly assign the + // columns as well. + graph (row_map.Comm().NumProc() > 1 ? + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + (int)n_entries_per_row, + false))) + : + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + (int)n_entries_per_row, + false))) + ), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) {} - + SparsityPattern::SparsityPattern (const Epetra_Map &InputMap, const std::vector &n_entries_per_row) : row_map (InputMap), col_map (row_map), compressed (false), - graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - (int*)const_cast(&(n_entries_per_row[0])), - false))), + graph (row_map.Comm().NumProc() > 1 ? + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + (int*)const_cast + (&(n_entries_per_row[row_map.MinMyGID()])), + false))) + : + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + (int*)const_cast + (&(n_entries_per_row[row_map.MinMyGID()])), + false))) + ), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -134,9 +161,17 @@ namespace TrilinosWrappers row_map (InputRowMap), col_map (InputColMap), compressed (false), - graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - int(n_entries_per_row), false))), + graph (row_map.Comm().NumProc() > 1 ? + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + (int)n_entries_per_row, + false))) + : + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + (int)n_entries_per_row, + false))) + ), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -149,10 +184,19 @@ namespace TrilinosWrappers row_map (InputRowMap), col_map (InputColMap), compressed (false), - graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - (int*)const_cast(&(n_entries_per_row[0])), - false))), + graph (row_map.Comm().NumProc() > 1 ? + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + (int*)const_cast + (&(n_entries_per_row[row_map.MinMyGID()])), + false))) + : + (std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + (int*)const_cast + (&(n_entries_per_row[row_map.MinMyGID()])), + false))) + ), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -171,8 +215,8 @@ namespace TrilinosWrappers #endif compressed (false), graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - int(n_entries_per_row), false))), + (new Epetra_FECrsGraph(Copy, row_map, col_map, + int(n_entries_per_row), false))), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -191,9 +235,9 @@ namespace TrilinosWrappers #endif compressed (false), graph (std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - (int*)const_cast(&(n_entries_per_row[0])), - false))), + (new Epetra_FECrsGraph(Copy, row_map, col_map, + (int*)const_cast(&(n_entries_per_row[0])), + false))), cached_row_indices (1), n_cached_elements (0), row_in_cache (0) @@ -201,8 +245,8 @@ namespace TrilinosWrappers // Copy function is currently not working // because the Trilinos Epetra_FECrsGraph - // does not implement a reinit function - // from another graph. + // does not implement a constructor from + // another graph. /* SparsityPattern::SparsityPattern (const SparsityPattern &InputSP) : @@ -261,8 +305,22 @@ namespace TrilinosWrappers row_map = input_row_map; col_map = input_col_map; - graph = std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, n_entries_per_row, false)); + // for more than one processor, need to + // specify only row map first and let the + // matrix entries decide about the column + // map (which says which columns are + // present in the matrix, not to be + // confused with the col_map that tells + // how the domain dofs of the matrix will + // be distributed). for only one + // processor, we can directly assign the + // columns as well. + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, n_entries_per_row, false)); + else + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, n_entries_per_row, false)); n_cached_elements = 0; row_in_cache = 0; @@ -312,10 +370,16 @@ namespace TrilinosWrappers row_map = input_row_map; col_map = input_col_map; - graph = std::auto_ptr - (new Epetra_FECrsGraph(Copy, row_map, - n_entries_per_row[input_row_map.MinMyGID()], - false)); + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); + else + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); n_cached_elements = 0; row_in_cache = 0; @@ -360,10 +424,16 @@ namespace TrilinosWrappers for (unsigned int row=0; row - (new Epetra_FECrsGraph(Copy, row_map, - n_entries_per_row[input_row_map.MinMyGID()], - false)); + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); + else + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); Assert (graph->NumGlobalRows() == (int)sp.n_rows(), ExcDimensionMismatch (graph->NumGlobalRows(), @@ -417,10 +487,16 @@ namespace TrilinosWrappers for (unsigned int row=0; row - (new Epetra_FECrsGraph(Copy, row_map, - n_entries_per_row[input_row_map.MinMyGID()], - false)); + if (row_map.Comm().NumProc() > 1) + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); + else + graph = std::auto_ptr + (new Epetra_FECrsGraph(Copy, row_map, col_map, + n_entries_per_row[input_row_map.MinMyGID()], + false)); Assert (graph->NumGlobalRows() == (int)sp.n_rows(), ExcDimensionMismatch (graph->NumGlobalRows(), -- 2.39.5