From b8c06ec67cc4d52dc8b2b5dcaf3a1231f8751bcf Mon Sep 17 00:00:00 2001 From: kronbichler Date: Wed, 6 Feb 2013 11:27:10 +0000 Subject: [PATCH] Replace deprecated functions. git-svn-id: https://svn.dealii.org/trunk@28249 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/sparse_matrix.h | 2 - .../include/deal.II/lac/sparsity_pattern.h | 18 ++++ deal.II/source/lac/sparse_direct.cc | 100 +++++++----------- deal.II/source/lac/sparsity_tools.cc | 98 +++++++---------- 4 files changed, 96 insertions(+), 122 deletions(-) diff --git a/deal.II/include/deal.II/lac/sparse_matrix.h b/deal.II/include/deal.II/lac/sparse_matrix.h index e061abc4f7..af8c02dc53 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix.h +++ b/deal.II/include/deal.II/lac/sparse_matrix.h @@ -2320,8 +2320,6 @@ namespace SparseMatrixIterators Assert (&accessor.get_matrix() == &other.accessor.get_matrix(), ExcInternalError()); - const SparsityPattern &sparsity = accessor.get_matrix().get_sparsity_pattern(); - return (*this)->a_index - other->a_index; } diff --git a/deal.II/include/deal.II/lac/sparsity_pattern.h b/deal.II/include/deal.II/lac/sparsity_pattern.h index 4334b3adf0..252feb11fb 100644 --- a/deal.II/include/deal.II/lac/sparsity_pattern.h +++ b/deal.II/include/deal.II/lac/sparsity_pattern.h @@ -264,6 +264,14 @@ namespace SparsityPatternIterators */ bool operator < (const Iterator &) const; + /** + * Return the distance between the current iterator and the argument. + * The distance is given by how many times one has to apply operator++ + * to the current iterator to get the argument (for a positive return + * value), or operator-- (for a negative return value). + */ + int operator - (const Iterator &p) const; + private: /** * Store an object of the accessor class. @@ -1588,6 +1596,16 @@ namespace SparsityPatternIterators return accessor < other.accessor; } + + inline + int + Iterator::operator - (const Iterator &other) const + { + Assert (accessor.sparsity_pattern == other.accessor.sparsity_pattern, + ExcInternalError()); + + return (*this)->a_index - other->a_index; + } } diff --git a/deal.II/source/lac/sparse_direct.cc b/deal.II/source/lac/sparse_direct.cc index 9a7cf69508..4da14da2ad 100644 --- a/deal.II/source/lac/sparse_direct.cc +++ b/deal.II/source/lac/sparse_direct.cc @@ -540,47 +540,36 @@ SparseDirectMA27::initialize (const SparsityPattern &sp) const unsigned int n_rows = sparsity_pattern->n_rows(); - const std::size_t *const - rowstart_indices = sparsity_pattern->get_rowstart_indices(); - const unsigned int *const - col_nums = sparsity_pattern->get_column_numbers(); - // first count number of nonzero - // elements in the upper right - // part. the matrix is symmetric, - // so this suffices + // first count number of nonzero elements in the upper right part. the + // matrix is symmetric, so this suffices n_nonzero_elements = 0; for (unsigned int row=0; rowbegin(row); + col < sparsity_pattern->end(row); ++col) + if (row <= col->column()) ++n_nonzero_elements; - // fill the row numbers and column - // numbers arrays from the sparsity - // pattern. note that we have - // Fortran convention, i.e. indices - // need to be 1-base, as opposed to - // C's 0-based convention! + // fill the row numbers and column numbers arrays from the sparsity + // pattern. note that we have Fortran convention, i.e. indices need to be + // 1-base, as opposed to C's 0-based convention! row_numbers.resize (n_nonzero_elements); column_numbers.resize (n_nonzero_elements); unsigned int global_index = 0; for (unsigned int row=0; rowbegin(row); + col < sparsity_pattern->end(row); ++col) // note that the matrix must be // symmetric, so only treat the // upper right part - if (row <= *col) + if (row <= col->column()) { Assert (global_index < n_nonzero_elements, ExcInternalError()); row_numbers[global_index] = row+1; - column_numbers[global_index] = *col+1; + column_numbers[global_index] = col->column()+1; ++global_index; }; Assert (global_index == n_nonzero_elements, ExcInternalError()); @@ -870,40 +859,37 @@ SparseDirectMA27::fill_A (const SparseMatrix &matrix) const SparsityPattern &sparsity_pattern = matrix.get_sparsity_pattern (); const unsigned int n_rows = sparsity_pattern.n_rows(); - const std::size_t *rowstart_indices = sparsity_pattern.get_rowstart_indices(); - const unsigned int *col_nums = sparsity_pattern.get_column_numbers(); unsigned int global_index = 0; for (unsigned int row=0; row::const_iterator col=matrix.begin(row); + col < matrix.end(row); ++col) // note that the matrix must be // symmetric, so only treat the // upper right part - if (row <= *col) + if (row <= col->column()) { Assert (global_index < n_nonzero_elements, ExcInternalError()); - A[global_index] = matrix(row,*col); + A[global_index] = col->value(); ++global_index; // make sure that the symmetric // entry exists and has the same // value, unless this one is zero - Assert ((matrix(row,*col) == 0) + Assert ((col->value() == 0) || - (std::fabs(matrix(row,*col) - matrix(*col,row)) - <= 1e-15 * std::fabs (matrix(row,*col))), + (std::fabs(col->value() - matrix(col->column(),row)) + <= 1e-15 * std::fabs (col->value())), ExcMatrixNotSymmetric()); } else // lower left part. just check // symmetry - Assert ((matrix(row,*col) == 0) + Assert ((col->value() == 0) || - (std::fabs(matrix(row,*col) - matrix(*col,row)) - <= 1e-15 * std::fabs (matrix(row,*col))), + (std::fabs(col->value() - matrix(col->column(),row)) + <= 1e-15 * std::fabs (col->value())), ExcMatrixNotSymmetric()); Assert (global_index == n_nonzero_elements, ExcInternalError()); @@ -1117,10 +1103,6 @@ SparseDirectMA47::initialize (const SparseMatrix &m) const unsigned int n_rows = sparsity_pattern.n_rows(); - const std::size_t *const - rowstart_indices = sparsity_pattern.get_rowstart_indices(); - const unsigned int *const - col_nums = sparsity_pattern.get_column_numbers(); // first count number of nonzero // elements in the upper right @@ -1128,12 +1110,10 @@ SparseDirectMA47::initialize (const SparseMatrix &m) // so this suffices n_nonzero_elements = 0; for (unsigned int row=0; row::const_iterator col = m.begin(row); + col < m.end(row); ++col) + // skip zero elements, as required by the docs of MA47 + if (row <= col->column() && col->value() != 0) ++n_nonzero_elements; @@ -1148,18 +1128,17 @@ SparseDirectMA47::initialize (const SparseMatrix &m) unsigned int global_index = 0; for (unsigned int row=0; row::const_iterator col = m.begin(row); + col < m.end(row); ++col) // note that the matrix must be // symmetric, so only treat the // upper right part - if ((row <= *col) && (m(row,*col) != 0)) + if ((row <= col->column()) && (col->value() != 0)) { Assert (global_index < n_nonzero_elements, ExcInternalError()); row_numbers[global_index] = row+1; - column_numbers[global_index] = *col+1; + column_numbers[global_index] = col->column()+1; ++global_index; }; Assert (global_index == n_nonzero_elements, ExcInternalError()); @@ -1389,38 +1368,35 @@ SparseDirectMA47::fill_A (const SparseMatrix &matrix) const SparsityPattern &sparsity_pattern = matrix.get_sparsity_pattern (); const unsigned int n_rows = sparsity_pattern.n_rows(); - const std::size_t *rowstart_indices = sparsity_pattern.get_rowstart_indices(); - const unsigned int *col_nums = sparsity_pattern.get_column_numbers(); unsigned int global_index = 0; for (unsigned int row=0; row::const_iterator col=matrix.begin(row); + col < matrix.end(row); ++col) // note that the matrix must be // symmetric, so only treat the // upper right part - if ((row <= *col) && (matrix(row,*col) != 0)) + if ((row <= col->column()) && (col->value() != 0)) { Assert (global_index < n_nonzero_elements, ExcInternalError()); - A[global_index] = matrix(row,*col); + A[global_index] = col->value(); ++global_index; // make sure that the symmetric // entry exists and has the same // value, unless this one is zero - Assert ((matrix(row,*col) == 0) + Assert ((col->value() == 0) || - (matrix(row,*col) == matrix(*col,row)), + (col->value() == matrix(col->column(),row)), ExcMatrixNotSymmetric()); } else // lower left part. just check // symmetry - Assert ((matrix(row,*col) == 0) + Assert ((col->value() == 0) || - (matrix(row,*col) == matrix(*col,row)), + (col->value() == matrix(col->column(),row)), ExcMatrixNotSymmetric()); Assert (global_index == n_nonzero_elements, ExcInternalError()); diff --git a/deal.II/source/lac/sparsity_tools.cc b/deal.II/source/lac/sparsity_tools.cc index f0752975d5..c86cb0c60c 100644 --- a/deal.II/source/lac/sparsity_tools.cc +++ b/deal.II/source/lac/sparsity_tools.cc @@ -87,12 +87,17 @@ namespace SparsityTools // one more nuisance: we have to copy our // own data to arrays that store signed // integers :-( - std::vector int_rowstart (sparsity_pattern.get_rowstart_indices(), - sparsity_pattern.get_rowstart_indices() + - sparsity_pattern.n_rows()+1); - std::vector int_colnums (sparsity_pattern.get_column_numbers(), - sparsity_pattern.get_column_numbers()+ - int_rowstart[sparsity_pattern.n_rows()]); + std::vector int_rowstart(1); + int_rowstart.reserve(sparsity_pattern.n_rows()+1); + std::vector int_colnums; + int_colnums.reserve(sparsity_pattern.n_nonzero_elements()); + for (unsigned int row=0; rowcolumn()); + int_rowstart.push_back(int_colnums.size()); + } std::vector int_partition_indices (sparsity_pattern.n_rows()); @@ -135,13 +140,9 @@ namespace SparsityTools namespace internal { /** - * Given a connectivity graph and - * a list of indices (where - * invalid_unsigned_int indicates - * that a node has not been - * numbered yet), pick a valid - * starting index among the - * as-yet unnumbered one. + * Given a connectivity graph and a list of indices (where + * invalid_unsigned_int indicates that a node has not been numbered yet), + * pick a valid starting index among the as-yet unnumbered one. */ unsigned int find_unnumbered_starting_index (const SparsityPattern &sparsity, @@ -151,46 +152,31 @@ namespace SparsityTools unsigned int starting_point = numbers::invalid_unsigned_int; unsigned int min_coordination = sparsity.n_rows(); for (unsigned int row=0; rowis_valid_entry() == false) break; - // post-condition after loop: - // coordination, i.e. the number - // of entries in this row is now - // j-rowstart[row] - if (j-sparsity.get_rowstart_indices()[row] < min_coordination) + // post-condition after loop: coordination, i.e. the number of + // entries in this row is now j-rowstart[row] + if (j-sparsity.begin(row) < min_coordination) { - min_coordination = j-sparsity.get_rowstart_indices()[row]; + min_coordination = j-sparsity.begin(row); starting_point = row; } } - // now we still have to care - // for the case that no - // unnumbered dof has a - // coordination number less - // than - // sparsity.n_rows(). this - // rather exotic case only - // happens if we only have - // one cell, as far as I can - // see, but there may be - // others as well. + // now we still have to care for the case that no unnumbered dof has a + // coordination number less than sparsity.n_rows(). this rather exotic + // case only happens if we only have one cell, as far as I can see, + // but there may be others as well. // - // if that should be the - // case, we can chose an - // arbitrary dof as starting - // point, e.g. the first - // unnumbered one + // if that should be the case, we can chose an arbitrary dof as + // starting point, e.g. the first unnumbered one if (starting_point == numbers::invalid_unsigned_int) { for (unsigned int i=0; i last_round_dofs (starting_indices); - // initialize the new_indices array with - // invalid values + // initialize the new_indices array with invalid values std::fill (new_indices.begin(), new_indices.end(), numbers::invalid_unsigned_int); @@ -245,9 +229,7 @@ namespace SparsityTools std::bind2nd(std::equal_to(), numbers::invalid_unsigned_int)); - // now if no valid points remain: - // find dof with lowest coordination - // number + // now if no valid points remain: find dof with lowest coordination number if (last_round_dofs.empty()) last_round_dofs .push_back (internal::find_unnumbered_starting_index (sparsity, @@ -272,12 +254,12 @@ namespace SparsityTools // dofs numbered in the last // round for (unsigned int i=0; iis_valid_entry() == false) break; else - next_round_dofs.push_back (sparsity.get_column_numbers()[j]); + next_round_dofs.push_back (j->column()); // sort dof numbers std::sort (next_round_dofs.begin(), next_round_dofs.end()); @@ -354,9 +336,9 @@ namespace SparsityTools s!=next_round_dofs.end(); ++s) { unsigned int coordination = 0; - for (unsigned int j=sparsity.get_rowstart_indices()[*s]; - jis_valid_entry() == false) break; else ++coordination; -- 2.39.5