From d51c66dfe87af104515deffdf57a310d82bc24c8 Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 15 Mar 2004 22:51:28 +0000 Subject: [PATCH] Purge the use of short ints. Noone knows how many entries we are going to store in a row. Fix various other bugs. git-svn-id: https://svn.dealii.org/trunk@8753 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/petsc_matrix_base.h | 56 ++++++++++++++------- deal.II/lac/source/petsc_matrix_base.cc | 40 ++++++++++++--- 2 files changed, 72 insertions(+), 24 deletions(-) diff --git a/deal.II/lac/include/lac/petsc_matrix_base.h b/deal.II/lac/include/lac/petsc_matrix_base.h index 4878b4eb70..6519013071 100644 --- a/deal.II/lac/include/lac/petsc_matrix_base.h +++ b/deal.II/lac/include/lac/petsc_matrix_base.h @@ -63,7 +63,7 @@ namespace PETScWrappers */ Accessor (const MatrixBase *matrix, const unsigned int row, - const unsigned short index); + const unsigned int index); /** * Row number of the element @@ -77,7 +77,7 @@ namespace PETScWrappers * represented by this * object. */ - unsigned short index() const; + unsigned int index() const; /** * Column number of the @@ -100,7 +100,7 @@ namespace PETScWrappers /** * The matrix accessed. */ - const MatrixBase *matrix; + mutable MatrixBase *matrix; /** * Current row number. @@ -167,9 +167,9 @@ namespace PETScWrappers * into the matrix @p matrix for the * given row and the index within it. */ - const_iterator (const MatrixBase *matrix, - const unsigned int row, - const unsigned short index); + const_iterator (const MatrixBase *matrix, + const unsigned int row, + const unsigned int index); /** * Prefix increment. @@ -358,6 +358,23 @@ namespace PETScWrappers PetscScalar el (const unsigned int i, const unsigned int j) const; + /** + * Return the main diagonal + * element in the ith + * row. This function throws an + * error if the matrix is not + * quadratic. + * + * Since we do not have direct access + * to the underlying data structure, + * this function is no faster than the + * elementwise access using the el() + * function. However, we provide this + * function for compatibility with the + * SparseMatrix class. + */ + PetscScalar diag_element (const unsigned int i) const; + /** * Return the number of rows in this * matrix. @@ -592,6 +609,10 @@ namespace PETScWrappers * Exception */ DeclException0 (ExcSourceEqualsDestination); + /** + * Exception + */ + DeclException0 (ExcMatrixNotSquare); protected: /** @@ -645,16 +666,17 @@ namespace PETScWrappers inline const_iterator::Accessor:: - Accessor (const MatrixBase *matrix, - const unsigned int row, - const unsigned short index) + Accessor (const MatrixBase *matrix, + const unsigned int row, + const unsigned int index) : - matrix(matrix), + matrix(const_cast(matrix)), a_row(row), a_index(index) { visit_present_row (); - Assert (index < colnum_cache->size(), + Assert ((row == matrix->m()) || + (index < colnum_cache->size()), ExcInvalidIndexWithinRow (index, row)); } @@ -678,7 +700,7 @@ namespace PETScWrappers inline - unsigned short + unsigned int const_iterator::Accessor::index() const { Assert (a_row < matrix->m(), ExcBeyondEndOfMatrix()); @@ -697,9 +719,9 @@ namespace PETScWrappers inline const_iterator:: - const_iterator(const MatrixBase *matrix, - const unsigned int row, - const unsigned short index) + const_iterator(const MatrixBase *matrix, + const unsigned int row, + const unsigned int index) : accessor(matrix, row, index) {} @@ -763,8 +785,8 @@ namespace PETScWrappers const_iterator:: operator == (const const_iterator& other) const { - return (accessor.row() == other.accessor.row() && - accessor.index() == other.accessor.index()); + return (accessor.a_row == other.accessor.a_row && + accessor.a_index == other.accessor.a_index); } diff --git a/deal.II/lac/source/petsc_matrix_base.cc b/deal.II/lac/source/petsc_matrix_base.cc index 4f240b3a54..8251719562 100644 --- a/deal.II/lac/source/petsc_matrix_base.cc +++ b/deal.II/lac/source/petsc_matrix_base.cc @@ -26,24 +26,37 @@ namespace PETScWrappers MatrixBase::const_iterator::Accessor:: visit_present_row () { - const int a_row_int = this->a_row; // need signed int for PETSc? + // if we are asked to visit the + // past-the-end line, then simply + // release all our caches and go on + // with life + if (this->a_row == matrix->m()) + { + colnum_cache.reset (); + value_cache.reset (); + + return; + } + + // otherwise first flush PETSc caches + matrix->compress (); + + const signed int a_row_int = this->a_row; // get a representation of the present // row int ncols; int *colnums; PetscScalar *values; + int ierr; ierr = MatGetRow(*matrix, a_row_int, &ncols, &colnums, &values); AssertThrow (ierr == 0, MatrixBase::ExcPETScError(ierr)); // copy it into our caches - colnum_cache = - boost::shared_ptr > - ( new std::vector (colnums, colnums+ncols)); - value_cache = - boost::shared_ptr > - ( new std::vector (values, values+ncols)); + colnum_cache.reset (new std::vector (colnums, + colnums+ncols)); + value_cache.reset (new std::vector (values, values+ncols)); // and finally restore the matrix ierr = MatRestoreRow(*matrix, a_row_int, &ncols, &colnums, &values); @@ -149,6 +162,19 @@ namespace PETScWrappers } + + PetscScalar + MatrixBase::diag_element (const unsigned int i) const + { + Assert (m() == n(), ExcMatrixNotSquare()); + + // this doesn't seem to work any + // different than any other element + return el(i,i); + } + + + void MatrixBase::compress () { -- 2.39.5