From ec704621922b5118a0244b54b308260bddb5b9c2 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Wed, 6 Feb 2013 10:45:52 +0000 Subject: [PATCH] Fix a few deprecation warnings. git-svn-id: https://svn.dealii.org/trunk@28248 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/lac/constraint_matrix.templates.h | 95 ++++++++----------- .../lac/sparse_decomposition.templates.h | 59 +++++------- .../deal.II/lac/sparse_ilu.templates.h | 69 +++++--------- .../deal.II/lac/sparse_mic.templates.h | 16 +--- deal.II/source/lac/constraint_matrix.cc | 77 ++++++--------- 5 files changed, 123 insertions(+), 193 deletions(-) diff --git a/deal.II/include/deal.II/lac/constraint_matrix.templates.h b/deal.II/include/deal.II/lac/constraint_matrix.templates.h index c95e72b0a0..56154db949 100644 --- a/deal.II/include/deal.II/lac/constraint_matrix.templates.h +++ b/deal.II/include/deal.II/lac/constraint_matrix.templates.h @@ -1729,28 +1729,18 @@ namespace internals { template static inline - void add_value (const double value, - const unsigned int row, - const unsigned int column, - const unsigned int *col_ptr, - const bool are_on_diagonal, - unsigned int &counter, - SparseMatrixIterator val_ptr) + void add_value (const double value, + const unsigned int row, + const unsigned int column, + SparseMatrixIterator &matrix_values) { if (value != 0.) { - Assert (col_ptr != 0, - typename SparseMatrix::ExcInvalidIndex (row, column)); - if (are_on_diagonal) - { - val_ptr->value() += value; - return; - } - while (col_ptr[counter] < column) - ++counter; - Assert (col_ptr[counter] == column, + while (matrix_values->column() < column) + ++matrix_values; + Assert (matrix_values->column() == column, typename SparseMatrix::ExcInvalidIndex(row, column)); - (val_ptr+counter)->value() += value; + matrix_values->value() += value; } } } @@ -1780,20 +1770,12 @@ namespace internals if (sparsity.n_nonzero_elements() == 0) return; - const std::size_t *row_start = sparsity.get_rowstart_indices(); - const unsigned int *sparsity_struct = sparsity.get_column_numbers(); - const unsigned int row = global_rows.global_row(i); const unsigned int loc_row = global_rows.local_row(i); - const unsigned int *col_ptr = sparsity.row_length(row) == 0 ? 0 : - &sparsity_struct[row_start[row]]; typename SparseMatrix::iterator - val_ptr = (sparsity.row_length(row) == 0 ? - sparse_matrix->end() : - sparse_matrix->begin(row)); - const bool optimize_diagonal = sparsity.optimize_diagonal(); - unsigned int counter = optimize_diagonal; + matrix_values = sparse_matrix->begin(row); + const bool optimize_diagonal = sparsity.n_rows() == sparsity.n_cols(); // distinguish three cases about what can // happen for checking whether the diagonal is @@ -1814,8 +1796,7 @@ namespace internals const double col_val = matrix_ptr[loc_col]; dealiiSparseMatrix::add_value (col_val, row, global_rows.global_row(j), - col_ptr, false, counter, - val_ptr); + matrix_values); } } else @@ -1825,55 +1806,57 @@ namespace internals double col_val = resolve_matrix_entry (global_rows, global_rows, i, j, loc_row, local_matrix); dealiiSparseMatrix::add_value (col_val, row, - global_rows.global_row(j), col_ptr, - false, counter, val_ptr); + global_rows.global_row(j), + matrix_values); } } } else if (i>=column_start && ibegin(row)->value() += matrix_ptr[loc_row]; for (unsigned int j=column_start; jvalue() += matrix_ptr[loc_row]; for (unsigned int j=i+1; jbegin(row)->value() += + resolve_matrix_entry (global_rows, global_rows, i, i, + loc_row, local_matrix); for (unsigned int j=column_start; jvalue() += resolve_matrix_entry (global_rows, global_rows, i, i, loc_row, - local_matrix); for (unsigned int j=i+1; jbegin(row)->value() += col_val; + else + dealiiSparseMatrix::add_value(col_val, row, + global_rows.global_row(j), + matrix_values); } } else { + ++matrix_values; // jump over diagonal element for (unsigned int j=column_start; jbegin(row)->value() += col_val; + else + dealiiSparseMatrix::add_value (col_val, row, + global_rows.global_row(j), + matrix_values); } } } diff --git a/deal.II/include/deal.II/lac/sparse_decomposition.templates.h b/deal.II/include/deal.II/lac/sparse_decomposition.templates.h index 0f4ca7017f..33a47c6e22 100644 --- a/deal.II/include/deal.II/lac/sparse_decomposition.templates.h +++ b/deal.II/include/deal.II/lac/sparse_decomposition.templates.h @@ -125,7 +125,7 @@ void SparseLUDecomposition::initialize ( } // now use this sparsity pattern - Assert (sparsity_pattern_to_use->optimize_diagonal(), + Assert (sparsity_pattern_to_use->n_rows()==sparsity_pattern_to_use->n_cols(), typename SparsityPattern::ExcDiagonalNotOptimized()); decomposed = false; { @@ -156,7 +156,7 @@ decompose (const SparseMatrix &matrix, template void SparseLUDecomposition::reinit (const SparsityPattern &sparsity) { - Assert (sparsity.optimize_diagonal(), + Assert (sparsity.n_rows() == sparsity.n_cols(), typename SparsityPattern::ExcDiagonalNotOptimized()); decomposed = false; { @@ -173,19 +173,19 @@ void SparseLUDecomposition::prebuild_lower_bound() { const unsigned int *const - column_numbers = this->get_sparsity_pattern().get_column_numbers(); + column_numbers = this->get_sparsity_pattern().colnums; const std::size_t *const - rowstart_indices = this->get_sparsity_pattern().get_rowstart_indices(); + rowstart_indices = this->get_sparsity_pattern().rowstart; const unsigned int N = this->m(); - - prebuilt_lower_bound.resize (N); - - for (unsigned int row=0; row::copy_from (const SparseMatrix &matrix // awkward way so that we find the corresponding function in the base class. SparseMatrix::operator= (number(0)); - // note: pointers to the sparsity - // pattern of the old matrix! - const std::size_t *const in_rowstart_indices - = matrix.get_sparsity_pattern().get_rowstart_indices(); - const unsigned int *const in_cols - = matrix.get_sparsity_pattern().get_column_numbers(); - const unsigned int *cols = this->get_sparsity_pattern().get_column_numbers(); - const std::size_t *rowstart_indices = - this->get_sparsity_pattern().get_rowstart_indices(); - - // both allow more and less entries - // in the new matrix + // both allow more and less entries in the new matrix std::size_t in_index, index; for (unsigned int row=0; rowm(); ++row) { - index = rowstart_indices[row]; - in_index = in_rowstart_indices[row]; - this->val[index++] = matrix.val[in_index++]; - while (in_index < in_rowstart_indices[row+1] && - index < rowstart_indices[row+1]) + typename SparseMatrix::iterator index = this->begin(row); + typename SparseMatrix::const_iterator + in_index = matrix.begin(row); + index->value() = in_index->value(); + ++index, ++in_index; + while (index < this->end(row) && in_index < matrix.end(row)) { - while (cols[index] < in_cols[in_index] && index < rowstart_indices[row+1]) + while (index->column() < in_index->column() && index < this->end(row)) ++index; - while (in_cols[in_index] < cols[index] && in_index < in_rowstart_indices[row+1]) + while (in_index->column() < index->column() && in_index < matrix.end(row)) ++in_index; - this->val[index++] = matrix.val[in_index++]; + index->value() = in_index->value(); + ++index, ++in_index; } } } @@ -254,7 +245,7 @@ SparseLUDecomposition::strengthen_diagonal_impl () { // get the global index of the first // non-diagonal element in this row - Assert (this->cols->optimize_diagonal(), ExcNotImplemented()); + Assert (this->m() == this->n(), ExcNotImplemented()); typename SparseMatrix::iterator diagonal_element = this->begin(row); diff --git a/deal.II/include/deal.II/lac/sparse_ilu.templates.h b/deal.II/include/deal.II/lac/sparse_ilu.templates.h index 6bb71fa990..f80470bc3c 100644 --- a/deal.II/include/deal.II/lac/sparse_ilu.templates.h +++ b/deal.II/include/deal.II/lac/sparse_ilu.templates.h @@ -66,15 +66,12 @@ void SparseILU::decompose (const SparseMatrix &matrix, if (strengthen_diagonal>0) this->strengthen_diagonal_impl(); - // in the following, we implement - // algorithm 10.4 in the book by - // Saad by translating in essence - // the algorithm given at the end - // of section 10.3.2, using the - // names of variables used there + // in the following, we implement algorithm 10.4 in the book by Saad by + // translating in essence the algorithm given at the end of section 10.3.2, + // using the names of variables used there const SparsityPattern &sparsity = this->get_sparsity_pattern(); - const std::size_t *const ia = sparsity.get_rowstart_indices(); - const unsigned int *const ja = sparsity.get_column_numbers(); + const std::size_t *const ia = sparsity.rowstart; + const unsigned int *const ja = sparsity.colnums; number *luval = this->SparseMatrix::val; @@ -91,22 +88,14 @@ void SparseILU::decompose (const SparseMatrix &matrix, for (unsigned int j=j1; j<=j2; ++j) iw[ja[j]] = j; - // the algorithm in the book - // works on the elements of row - // k left of the - // diagonal. however, since we - // store the diagonal element - // at the first position, start - // at the element after the - // diagonal and run as long as - // we don't walk into the right - // half + // the algorithm in the book works on the elements of row k left of the + // diagonal. however, since we store the diagonal element at the first + // position, start at the element after the diagonal and run as long as + // we don't walk into the right half unsigned int j = j1+1; - // pathological case: the current row - // of the matrix has only the - // diagonal entry. then we have - // nothing to do. + // pathological case: the current row of the matrix has only the + // diagonal entry. then we have nothing to do. if (j > j2) goto label_200; @@ -121,9 +110,7 @@ label_150: number t1 = luval[j] * luval[ia[jrow]]; luval[j] = t1; - // jj runs from just right of - // the diagonal to the end of - // the row + // jj runs from just right of the diagonal to the end of the row unsigned int jj = ia[jrow]+1; while (ja[jj] < jrow) ++jj; @@ -141,25 +128,15 @@ label_150: label_200: - // in the book there is an - // assertion that we have hit - // the diagonal element, - // i.e. that jrow==k. however, - // we store the diagonal - // element at the front, so - // jrow must actually be larger - // than k or j is already in + // in the book there is an assertion that we have hit the diagonal + // element, i.e. that jrow==k. however, we store the diagonal element at + // the front, so jrow must actually be larger than k or j is already in // the next row Assert ((jrow > k) || (j==ia[k+1]), ExcInternalError()); - // now we have to deal with the - // diagonal element. in the - // book it is located at - // position 'j', but here we - // use the convention of - // storing the diagonal element - // first, so instead of j we - // use uptr[k]=ia[k] + // now we have to deal with the diagonal element. in the book it is + // located at position 'j', but here we use the convention of storing + // the diagonal element first, so instead of j we use uptr[k]=ia[k] Assert (luval[ia[k]] != 0, ExcInternalError()); luval[ia[k]] = 1./luval[ia[k]]; @@ -181,10 +158,10 @@ void SparseILU::vmult (Vector &dst, Assert (dst.size() == this->m(), ExcDimensionMismatch(dst.size(), this->m())); const unsigned int N=dst.size(); - const std::size_t *const rowstart_indices - = this->get_sparsity_pattern().get_rowstart_indices(); + const std::size_t *const rowstart_indices + = this->get_sparsity_pattern().rowstart; const unsigned int *const column_numbers - = this->get_sparsity_pattern().get_column_numbers(); + = this->get_sparsity_pattern().colnums; // solve LUx=b in two steps: // first Ly = b, then @@ -256,9 +233,9 @@ void SparseILU::Tvmult (Vector &dst, const unsigned int N=dst.size(); const std::size_t *const rowstart_indices - = this->get_sparsity_pattern().get_rowstart_indices(); + = this->get_sparsity_pattern().rowstart; const unsigned int *const column_numbers - = this->get_sparsity_pattern().get_column_numbers(); + = this->get_sparsity_pattern().colnums; // solve (LU)'x=b in two steps: // first U'y = b, then diff --git a/deal.II/include/deal.II/lac/sparse_mic.templates.h b/deal.II/include/deal.II/lac/sparse_mic.templates.h index fc322ef102..f760b38490 100644 --- a/deal.II/include/deal.II/lac/sparse_mic.templates.h +++ b/deal.II/include/deal.II/lac/sparse_mic.templates.h @@ -132,12 +132,9 @@ void SparseMIC::decompose (const SparseMatrix &matrix, for (unsigned int row=0; rowm(); row++) inner_sums[row] = get_rowsum(row); - const unsigned int *const col_nums = this->get_sparsity_pattern().get_column_numbers(); - const std::size_t *const rowstarts = this->get_sparsity_pattern().get_rowstart_indices(); - for (unsigned int row=0; rowm(); row++) { - const number temp = this->diag_element(row); + const number temp = this->begin(row)->value(); number temp1 = 0; // work on the lower left part of the matrix. we know @@ -187,15 +184,10 @@ SparseMIC::vmult (Vector &dst, Assert (dst.size() == this->m(), ExcDimensionMismatch(dst.size(), this->m())); const unsigned int N=dst.size(); - const std::size_t *const rowstart_indices = this->get_sparsity_pattern().get_rowstart_indices(); - const unsigned int *const column_numbers = this->get_sparsity_pattern().get_column_numbers(); - // We assume the underlying matrix A is: - // A = X - L - U, where -L and -U are - // strictly lower- and upper- diagonal - // parts of the system. + // We assume the underlying matrix A is: A = X - L - U, where -L and -U are + // strictly lower- and upper- diagonal parts of the system. // - // Solve (X-L)X{-1}(X-U) x = b - // in 3 steps: + // Solve (X-L)X{-1}(X-U) x = b in 3 steps: dst = src; for (unsigned int row=0; rowcolumn()]) + condensed.add (new_line[row], new_line[j->column()]); else { - // let c point to the constraint - // of this column + // let c point to the constraint of this column std::vector::const_iterator c = lines.begin(); - while (c->line != uncondensed.get_column_numbers()[j]) + while (c->line != j->column()) ++c; for (unsigned int q=0; q!=c->entries.size(); ++q) @@ -841,14 +834,14 @@ void ConstraintMatrix::condense (const SparsityPattern &uncondensed, else // line must be distributed { - for (unsigned int j=uncondensed.get_rowstart_indices()[row]; - jcolumn()] != -1) // column is not constrained for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) condensed.add (new_line[next_constraint->entries[q].first], - new_line[uncondensed.get_column_numbers()[j]]); + new_line[j->column()]); else // not only this line but @@ -857,7 +850,7 @@ void ConstraintMatrix::condense (const SparsityPattern &uncondensed, // let c point to the constraint // of this column std::vector::const_iterator c = lines.begin(); - while (c->line != uncondensed.get_column_numbers()[j]) ++c; + while (c->line != j->column()) ++c; for (unsigned int p=0; p!=c->entries.size(); ++p) for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) @@ -877,8 +870,7 @@ void ConstraintMatrix::condense (SparsityPattern &sparsity) const { Assert (sorted == true, ExcMatrixNotClosed()); Assert (sparsity.is_compressed() == false, ExcMatrixIsClosed()); - Assert (sparsity.n_rows() == sparsity.n_cols(), - ExcNotQuadratic()); + Assert (sparsity.n_rows() == sparsity.n_cols(), ExcNotQuadratic()); // store for each index whether it must be // distributed or not. If entry is @@ -898,28 +890,17 @@ void ConstraintMatrix::condense (SparsityPattern &sparsity) const { if (distribute[row] == numbers::invalid_unsigned_int) { - // regular line. loop over cols all - // valid cols. note that this - // changes the line we are - // presently working on: we add - // additional entries. these are - // put to the end of the - // row. however, as constrained - // nodes cannot be constrained to - // other constrained nodes, nothing - // will happen if we run into these - // added nodes, as they can't be - // distributed further. we might - // store the position of the last - // old entry and stop work there, - // but since operating on the newly - // added ones only takes two - // comparisons (column index valid, - // distribute[column] necessarily - // ==numbers::invalid_unsigned_int), - // it is cheaper to not do so and - // run right until the end of the - // line + // regular line. loop over cols all valid cols. note that this + // changes the line we are presently working on: we add additional + // entries. these are put to the end of the row. however, as + // constrained nodes cannot be constrained to other constrained + // nodes, nothing will happen if we run into these added nodes, as + // they can't be distributed further. we might store the position of + // the last old entry and stop work there, but since operating on + // the newly added ones only takes two comparisons (column index + // valid, distribute[column] necessarily + // ==numbers::invalid_unsigned_int), it is cheaper to not do so and + // run right until the end of the line for (SparsityPattern::iterator entry = sparsity.begin(row); ((entry != sparsity.end(row)) && entry->is_valid_entry()); -- 2.39.5