From 5898a384aecaeeeb0e08f26b3995dc384d25de89 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Mon, 23 Nov 2009 02:35:34 +0000 Subject: [PATCH] Fix some compiling problems. Fix sparse matrix iterator. git-svn-id: https://svn.dealii.org/trunk@20147 0785d39b-7218-0410-832d-ea1e28bc413d --- .../lac/include/lac/trilinos_sparse_matrix.h | 6 +-- deal.II/lac/include/lac/trilinos_vector.h | 47 ++++++++++++------- deal.II/lac/source/trilinos_sparse_matrix.cc | 28 +++++++---- 3 files changed, 52 insertions(+), 29 deletions(-) diff --git a/deal.II/lac/include/lac/trilinos_sparse_matrix.h b/deal.II/lac/include/lac/trilinos_sparse_matrix.h index e40e8caace..2b64416e3a 100644 --- a/deal.II/lac/include/lac/trilinos_sparse_matrix.h +++ b/deal.II/lac/include/lac/trilinos_sparse_matrix.h @@ -167,13 +167,13 @@ namespace TrilinosWrappers * than one accessor can access * this data if necessary. */ - std_cxx1x::shared_ptr > colnum_cache; + std_cxx1x::shared_ptr > colnum_cache; /** * Similar cache for the values * of this row. */ - std_cxx1x::shared_ptr > value_cache; + std_cxx1x::shared_ptr > value_cache; /** * Discard the old row caches @@ -3235,7 +3235,7 @@ namespace TrilinosWrappers temp_vector.reinit(dst, true); - vmult (temp_vector, src); + Tvmult (temp_vector, src); dst += temp_vector; } diff --git a/deal.II/lac/include/lac/trilinos_vector.h b/deal.II/lac/include/lac/trilinos_vector.h index 5bff7053eb..ec0a4e82e7 100644 --- a/deal.II/lac/include/lac/trilinos_vector.h +++ b/deal.II/lac/include/lac/trilinos_vector.h @@ -367,6 +367,15 @@ namespace TrilinosWrappers explicit Vector (const Epetra_Map ¶llel_partitioning, const VectorBase &v); + /** + * Reinitialize from a deal.II + * vector. The Epetra_Map specifies the + * %parallel partitioning. + */ + template + void reinit (const Epetra_Map ¶llel_partitioner, + const dealii::Vector &v); + /** * Reinit functionality. This * function destroys the old @@ -469,19 +478,7 @@ namespace TrilinosWrappers Vector::Vector (const Epetra_Map &input_map, const dealii::Vector &v) { - vector = std::auto_ptr (new Epetra_FEVector(input_map)); - - const int min_my_id = input_map.MinMyGID(); - const int size = input_map.NumMyElements(); - - Assert (input_map.MaxLID() == size-1, - ExcDimensionMismatch(input_map.MaxLID(), size-1)); - - // Need to copy out values, since the - // deal.II might not use doubles, so - // that a direct access is not possible. - for (int i=0; i + void Vector::reinit (const Epetra_Map ¶llel_partitioner, + const dealii::Vector &v) + { + if (&*vector != 0 && vector->Map().SameAs(parallel_partitioner)) + vector = std::auto_ptr + (new Epetra_FEVector(parallel_partitioner)); + + const int size = parallel_partitioner.NumMyElements(); + + // Need to copy out values, since the + // deal.II might not use doubles, so + // that a direct access is not possible. + for (int i=0; i + vector = std::auto_ptr (new Epetra_FEVector(Epetra_Map (v.size(), 0, #ifdef DEAL_II_COMPILER_SUPPORTS_MPI Epetra_MpiComm(MPI_COMM_SELF) @@ -763,8 +779,7 @@ namespace TrilinosWrappers vector = std::auto_ptr (new Epetra_FEVector(map)); } - Epetra_Map & map = vector_partitioner(); - const int min_my_id = map.MinMyGID(); + const Epetra_Map & map = vector_partitioner(); const int size = map.NumMyElements(); Assert (map.MaxLID() == size-1, diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc index 94de3a0f70..363c8a1788 100644 --- a/deal.II/lac/source/trilinos_sparse_matrix.cc +++ b/deal.II/lac/source/trilinos_sparse_matrix.cc @@ -57,12 +57,24 @@ namespace TrilinosWrappers // row int ncols; int colnums = matrix->n(); - TrilinosScalar *values = new TrilinosScalar(colnums); + if (value_cache.get() == 0) + { + value_cache.reset (new std::vector (matrix->n())); + colnum_cache.reset (new std::vector (matrix->n())); + } + else + { + value_cache->resize (matrix->n()); + colnum_cache->resize (matrix->n()); + } - int ierr; - ierr = matrix->trilinos_matrix().ExtractGlobalRowCopy((int)this->a_row, - colnums, - ncols, &(values[0])); + int ierr = matrix->trilinos_matrix(). + ExtractGlobalRowCopy((int)this->a_row, + colnums, + ncols, &((*value_cache)[0]), + reinterpret_cast(&((*colnum_cache)[0]))); + value_cache->resize (ncols); + colnum_cache->resize (ncols); AssertThrow (ierr == 0, ExcTrilinosError(ierr)); // copy it into our caches if the @@ -71,10 +83,6 @@ namespace TrilinosWrappers // we shouldn't have initialized an // iterator for an empty line (what // would it point to?) - Assert (ncols != 0, ExcInternalError()); - colnum_cache.reset (new std::vector (colnums, - colnums+ncols)); - value_cache.reset (new std::vector (values, values+ncols)); } } @@ -1397,7 +1405,7 @@ namespace TrilinosWrappers { matrix->ExtractMyRowView (i, num_entries, values, indices); for (int j=0; jGRID(j)] << ") " + out << "(" << matrix->GRID(i) << "," << indices[matrix->GCID(j)] << ") " << values[j] << std::endl; } } -- 2.39.5