From 7c308960d31a2e31e26e444c07584577426af814 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Fri, 26 Oct 2018 16:18:18 +0200 Subject: [PATCH] add DynamicSparsityPattern::column_index() --- .../deal.II/lac/dynamic_sparsity_pattern.h | 8 +++ source/lac/dynamic_sparsity_pattern.cc | 33 +++++++++++ tests/sparsity/dynamic_sparsity_pattern_17.cc | 57 +++++++++++++++++++ .../dynamic_sparsity_pattern_17.output | 2 + 4 files changed, 100 insertions(+) create mode 100644 tests/sparsity/dynamic_sparsity_pattern_17.cc create mode 100644 tests/sparsity/dynamic_sparsity_pattern_17.output diff --git a/include/deal.II/lac/dynamic_sparsity_pattern.h b/include/deal.II/lac/dynamic_sparsity_pattern.h index c33fe8edcb..4e21e2a9e4 100644 --- a/include/deal.II/lac/dynamic_sparsity_pattern.h +++ b/include/deal.II/lac/dynamic_sparsity_pattern.h @@ -505,6 +505,14 @@ public: size_type column_number(const size_type row, const size_type index) const; + /** + * Return index of column @p col in row @p row. If the column does not + * exist in this sparsity pattern, the returned value will be + * 'numbers::invalid_size_type'. + */ + size_type + column_index(const size_type row, const size_type col) const; + /** * @name Iterators */ diff --git a/source/lac/dynamic_sparsity_pattern.cc b/source/lac/dynamic_sparsity_pattern.cc index 62ca57ff73..0c5c17fd77 100644 --- a/source/lac/dynamic_sparsity_pattern.cc +++ b/source/lac/dynamic_sparsity_pattern.cc @@ -14,6 +14,7 @@ // --------------------------------------------------------------------- #include +#include #include #include @@ -577,6 +578,38 @@ DynamicSparsityPattern::memory_consumption() const } + +types::global_dof_index +DynamicSparsityPattern::column_index( + const DynamicSparsityPattern::size_type row, + const DynamicSparsityPattern::size_type col) const +{ + Assert(row < n_rows(), + ExcIndexRangeType(row, + 0, + n_rows())); + Assert(col < n_cols(), + ExcIndexRangeType(row, + 0, + n_cols())); + Assert(rowset.size() == 0 || rowset.is_element(row), ExcInternalError()); + + const DynamicSparsityPattern::size_type local_row = + rowset.size() ? rowset.index_within_set(row) : row; + + // now we need to do a binary search. Note that col indices are assumed to + // be sorted. + const auto &cols = lines[local_row].entries; + auto it = Utilities::lower_bound(cols.begin(), cols.end(), col); + + if ((it != cols.end()) && (*it == col)) + return (it - cols.begin()); + else + return numbers::invalid_size_type; +} + + + // explicit instantiations template void DynamicSparsityPattern::Line::add_entries(size_type *, size_type *, const bool); diff --git a/tests/sparsity/dynamic_sparsity_pattern_17.cc b/tests/sparsity/dynamic_sparsity_pattern_17.cc new file mode 100644 index 0000000000..0ddc547002 --- /dev/null +++ b/tests/sparsity/dynamic_sparsity_pattern_17.cc @@ -0,0 +1,57 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check DynamicSparsityPattern::column_index. +// this test is based on dynamic_sparsity_pattern_04.cc + +#include + +#include "../tests.h" + + +void +test() +{ + const unsigned int N = 100; + DynamicSparsityPattern csp(N, N); + for (unsigned int i = 0; i < N; ++i) + for (unsigned int j = 0; j < 10; ++j) + csp.add(i, (i + (i + 1) * (j * j + i)) % N); + + for (unsigned int i = 0; i < N; ++i) + for (unsigned int ind = 0; ind < csp.row_length(i); ++ind) + { + const auto j = csp.column_number(i, ind); + const auto val = csp.column_index(i, j); + AssertThrow(val == ind, + ExcMessage(std::to_string(val) + + "!=" + std::to_string(ind))); + } + + deallog << "Ok" << std::endl; +} + + + +int +main() +{ + initlog(); + + test(); + return 0; +} diff --git a/tests/sparsity/dynamic_sparsity_pattern_17.output b/tests/sparsity/dynamic_sparsity_pattern_17.output new file mode 100644 index 0000000000..6ea1f897a5 --- /dev/null +++ b/tests/sparsity/dynamic_sparsity_pattern_17.output @@ -0,0 +1,2 @@ + +DEAL::Ok -- 2.39.5