From 20ff8d4e14f748a025ad013a7b34bdaadb05266f Mon Sep 17 00:00:00 2001 From: kronbichler Date: Tue, 9 Sep 2008 16:32:30 +0000 Subject: [PATCH] A few updates for parallel run. At least square matrices seem to work in parallel as well now. git-svn-id: https://svn.dealii.org/trunk@16771 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/trilinos_vector.h | 22 ++--- deal.II/lac/source/trilinos_sparse_matrix.cc | 84 +++++++------------- deal.II/lac/source/trilinos_vector.cc | 12 ++- 3 files changed, 53 insertions(+), 65 deletions(-) diff --git a/deal.II/lac/include/lac/trilinos_vector.h b/deal.II/lac/include/lac/trilinos_vector.h index de66fabbaa..94275380cf 100755 --- a/deal.II/lac/include/lac/trilinos_vector.h +++ b/deal.II/lac/include/lac/trilinos_vector.h @@ -341,10 +341,14 @@ namespace TrilinosWrappers typedef const internal::VectorReference const_reference; /** - * Default constructor. It - * doesn't do anything, derived - * classes will have to - * initialize the data. + * Default constructor that + * generates an empty (zero + * size) vector. The function + * reinit() will have + * to give the vector the + * correct size and + * distribution among processes + * in case of an MPI run. */ Vector (); @@ -353,9 +357,9 @@ namespace TrilinosWrappers * Epetra_Map that already * knows how to distribute the * individual components among - * the MPI processors, - * including the size of the - * vector. + * the MPI processors. It also + * includes information about + * the size of the vector. */ Vector (const Epetra_Map &InputMap); @@ -379,7 +383,7 @@ namespace TrilinosWrappers */ virtual ~Vector (); - /** + /** * Reinit functionality. This function * destroys the old vector content * and generates a new one based on @@ -387,7 +391,7 @@ namespace TrilinosWrappers */ void reinit (const Epetra_Map &input_map); - /** + /** * Reinit functionality. This function * copies the vector v to the current * one. diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc index 1f9af2c7dd..0d00c4fb09 100755 --- a/deal.II/lac/source/trilinos_sparse_matrix.cc +++ b/deal.II/lac/source/trilinos_sparse_matrix.cc @@ -186,25 +186,20 @@ namespace TrilinosWrappers // ExcDimensionMismatch (matrix->NumGlobalCols(), // sparsity_pattern.n_cols())); - std::vector n_entries_per_row(n_rows); - - for (unsigned int row=0; row values; std::vector row_indices; for (unsigned int row=0; rowInsertGlobalValues(row, row_length, - &values[0], &row_indices[0]); + matrix->InsertGlobalValues (row, row_length, + &values[0], &row_indices[0]); } } @@ -222,40 +217,7 @@ namespace TrilinosWrappers SparseMatrix::reinit (const Epetra_Map &input_map, const SparsityPattern &sparsity_pattern) { - // TODO: There seems to be problem - // in Epetra when a quadratic matrix - // is initialized with both row and - // column map. Maybe find something - // more out about this... reinit (input_map, input_map, sparsity_pattern); - - /* matrix.reset(); - - unsigned int n_rows = sparsity_pattern.n_rows(); - - if (row_map.Comm().MyPID() == 0) - { - Assert (input_map.NumGlobalElements() == (int)sparsity_pattern.n_rows(), - ExcDimensionMismatch (input_map.NumGlobalElements(), - sparsity_pattern.n_rows())); - Assert (input_map.NumGlobalElements() == (int)sparsity_pattern.n_cols(), - ExcDimensionMismatch (input_map.NumGlobalElements(), - sparsity_pattern.n_cols())); - } - - row_map = input_map; - col_map = row_map; - - std::vector n_entries_per_row(n_rows); - - for (unsigned int row=0; row - (new Epetra_FECrsMatrix(Copy, row_map, &n_entries_per_row[0], - false)); - - reinit (sparsity_pattern);*/ } @@ -287,9 +249,19 @@ namespace TrilinosWrappers for (unsigned int row=0; row - (new Epetra_FECrsMatrix(Copy, row_map, col_map, - &n_entries_per_row[0], false)); + // TODO: There seems to be problem + // in Epetra when a quadratic matrix + // is initialized with both row and + // column map. Maybe find something + // more out about this... + if (row_map.SameAs(col_map) == true) + matrix = std::auto_ptr + (new Epetra_FECrsMatrix(Copy, row_map, + &n_entries_per_row[0], false)); + else + matrix = std::auto_ptr + (new Epetra_FECrsMatrix(Copy, row_map, col_map, + &n_entries_per_row[0], false)); reinit (sparsity_pattern); } @@ -403,9 +375,9 @@ namespace TrilinosWrappers { // flush buffers int ierr; - //if (row_map.SameAs(col_map)) - //ierr = matrix->GlobalAssemble (); - //else + if (row_map.SameAs(col_map)) + ierr = matrix->GlobalAssemble (); + else ierr = matrix->GlobalAssemble (col_map, row_map); AssertThrow (ierr == 0, ExcTrilinosError(ierr)); @@ -569,11 +541,15 @@ namespace TrilinosWrappers // If the data is not on the // present processor, we can't - // continue. + // continue. Just print out + // zero. + + // TODO: Is this reasonable? if ((trilinos_i == -1 ) || (trilinos_j == -1)) { - Assert (false, ExcAccessToNonLocalElement(i, j, local_range().first, - local_range().second)); + return 0.; + //Assert (false, ExcAccessToNonLocalElement(i, j, local_range().first, + // local_range().second)); } else { @@ -639,7 +615,7 @@ namespace TrilinosWrappers unsigned int SparseMatrix::m () const { - int n_rows = matrix->NumGlobalRows(); + int n_rows = matrix -> NumGlobalRows(); return n_rows; } @@ -669,8 +645,8 @@ namespace TrilinosWrappers SparseMatrix::local_range () const { int begin, end; - begin = matrix->RowMap().MinMyGID(); - end = matrix->RowMap().MaxMyGID()+1; + begin = matrix -> RowMap().MinMyGID(); + end = matrix -> RowMap().MaxMyGID()+1; return std::make_pair (begin, end); } diff --git a/deal.II/lac/source/trilinos_vector.cc b/deal.II/lac/source/trilinos_vector.cc index c10a0a7f8e..623788c09a 100755 --- a/deal.II/lac/source/trilinos_vector.cc +++ b/deal.II/lac/source/trilinos_vector.cc @@ -220,10 +220,18 @@ namespace TrilinosWrappers // the vector. int trilinos_i = vector->Map().LID(index); TrilinosScalar value = 0.; + + // If the element is not + // present on the current + // processor, we can't + // continue. Just print out 0. + + // TODO: Is this reasonable? if (trilinos_i == -1 ) { - Assert (false, ExcAccessToNonlocalElement(index, local_range().first, - local_range().second)); + return 0.; + //Assert (false, ExcAccessToNonlocalElement(index, local_range().first, + // local_range().second)); } else value = (*vector)[0][trilinos_i]; -- 2.39.5