From be365be1bfc3e2907930edfca4f7afec42474d84 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Fri, 2 Feb 2018 11:15:04 +0100 Subject: [PATCH] Change size_type for FullMatrix --- include/deal.II/lac/full_matrix.h | 37 ++++---- include/deal.II/lac/full_matrix.templates.h | 100 +++++++++----------- source/lac/full_matrix.inst.in | 12 +-- 3 files changed, 70 insertions(+), 79 deletions(-) diff --git a/include/deal.II/lac/full_matrix.h b/include/deal.II/lac/full_matrix.h index 9b3e315a59..f76bdde858 100644 --- a/include/deal.II/lac/full_matrix.h +++ b/include/deal.II/lac/full_matrix.h @@ -65,12 +65,9 @@ class FullMatrix : public Table<2,number> { public: /** - * A type of used to index into this container. Because we can not expect to - * store matrices bigger than what can be indexed by a regular unsigned - * integer, unsigned int is completely sufficient as an index - * type. + * A type of used to index into this container. */ - typedef unsigned int size_type; + typedef std::size_t size_type; /** * Type of matrix entries. This typedef is analogous to value_type @@ -317,10 +314,10 @@ public: template void copy_from (const Tensor<2,dim> &T, - const size_type src_r_i=0, - const size_type src_r_j=dim-1, - const size_type src_c_i=0, - const size_type src_c_j=dim-1, + const unsigned int src_r_i=0, + const unsigned int src_r_j=dim-1, + const unsigned int src_c_i=0, + const unsigned int src_c_j=dim-1, const size_type dst_r=0, const size_type dst_c=0); @@ -338,8 +335,8 @@ public: const size_type src_r_j=dim-1, const size_type src_c_i=0, const size_type src_c_j=dim-1, - const size_type dst_r=0, - const size_type dst_c=0) const; + const unsigned int dst_r=0, + const unsigned int dst_c=0) const; /** * Copy a subset of the rows and columns of another matrix into the current @@ -731,7 +728,7 @@ public: */ template void add (const size_type row, - const unsigned int n_cols, + const size_type n_cols, const index_type *col_indices, const number2 *values, const bool elide_zero_values = true, @@ -1475,7 +1472,7 @@ inline typename FullMatrix::const_iterator FullMatrix::begin (const size_type r) const { - Assert (r::const_iterator FullMatrix::end (const size_type r) const { - Assert (r template inline void -FullMatrix::add (const size_type row, - const unsigned int n_cols, - const index_type *col_indices, - const number2 *values, +FullMatrix::add (const size_type row, + const size_type n_cols, + const index_type *col_indices, + const number2 *values, const bool, const bool) { @@ -1536,8 +1533,8 @@ FullMatrix::print (StreamType &s, Assert (!this->empty(), ExcEmptyMatrix()); // save the state of out stream - const unsigned int old_precision = s.precision (p); - const unsigned int old_width = s.width (w); + const std::streamsize old_precision = s.precision (p); + const std::streamsize old_width = s.width (w); for (size_type i=0; im(); ++i) { diff --git a/include/deal.II/lac/full_matrix.templates.h b/include/deal.II/lac/full_matrix.templates.h index 06375b3af4..8e54331910 100644 --- a/include/deal.II/lac/full_matrix.templates.h +++ b/include/deal.II/lac/full_matrix.templates.h @@ -301,7 +301,7 @@ void FullMatrix::backward (Vector &dst, size_type j; size_type nu = (m()=0; --i) + for (std::make_signed::type i=nu-1; i>=0; --i) { typename ProductType::type s = src(i); for (j=i+1; j::fill (const FullMatrix &src, const size_type src_offset_i, const size_type src_offset_j) { - Assert (dst_offset_i < m(), - ExcIndexRange (dst_offset_i, 0, m())); - Assert (dst_offset_j < n(), - ExcIndexRange (dst_offset_j, 0, n())); - Assert (src_offset_i < src.m(), - ExcIndexRange (src_offset_i, 0, src.m())); - Assert (src_offset_j < src.n(), - ExcIndexRange (src_offset_j, 0, src.n())); + AssertIndexRange (dst_offset_i, m()); + AssertIndexRange (dst_offset_j, n()); + AssertIndexRange (src_offset_i, src.m()); + AssertIndexRange (src_offset_j, src.n()); // Compute maximal size of copied block const size_type rows = std::min (m() - dst_offset_i, @@ -1085,14 +1081,10 @@ void FullMatrix::add (const FullMatrix &src, const size_type src_offset_i, const size_type src_offset_j) { - Assert (dst_offset_i < m(), - ExcIndexRange (dst_offset_i, 0, m())); - Assert (dst_offset_j < n(), - ExcIndexRange (dst_offset_j, 0, n())); - Assert (src_offset_i < src.m(), - ExcIndexRange (src_offset_i, 0, src.m())); - Assert (src_offset_j < src.n(), - ExcIndexRange (src_offset_j, 0, src.n())); + AssertIndexRange (dst_offset_i, m()); + AssertIndexRange (dst_offset_j, n()); + AssertIndexRange (src_offset_i, src.m()); + AssertIndexRange (src_offset_j, src.n()); // Compute maximal size of copied block const size_type rows = std::min (m() - dst_offset_i, src.m() - src_offset_i); @@ -1115,14 +1107,10 @@ void FullMatrix::Tadd (const FullMatrix &src, const size_type src_offset_i, const size_type src_offset_j) { - Assert (dst_offset_i < m(), - ExcIndexRange (dst_offset_i, 0, m())); - Assert (dst_offset_j < n(), - ExcIndexRange (dst_offset_j, 0, n())); - Assert (src_offset_i < src.n(), - ExcIndexRange (src_offset_i, 0, src.n())); - Assert (src_offset_j < src.m(), - ExcIndexRange (src_offset_j, 0, src.m())); + AssertIndexRange (dst_offset_i, m()); + AssertIndexRange (dst_offset_j, n()); + AssertIndexRange (src_offset_i, src.n()); + AssertIndexRange (src_offset_j, src.m()); // Compute maximal size of copied block const size_type rows = std::min (m() - dst_offset_i, src.n() - src_offset_j); @@ -1193,7 +1181,11 @@ namespace internal #ifdef DEAL_II_WITH_LAPACK static number value (const FullMatrix &A) { - LAPACKFullMatrix lp_A (A.m(), A.n()); + using s_type = typename LAPACKFullMatrix::size_type; + AssertIndexRange (A.m(), std::numeric_limits::max()+1); + AssertIndexRange (A.n(), std::numeric_limits::max()+1); + LAPACKFullMatrix lp_A (static_cast(A.m()), + static_cast(A.n())); lp_A = A; lp_A.compute_lu_factorization(); return lp_A.determinant(); @@ -1618,27 +1610,28 @@ template template void FullMatrix::copy_from (const Tensor<2,dim> &T, - const size_type src_r_i, - const size_type src_r_j, - const size_type src_c_i, - const size_type src_c_j, + const unsigned int src_r_i, + const unsigned int src_r_j, + const unsigned int src_c_i, + const unsigned int src_c_j, const size_type dst_r, const size_type dst_c) { - Assert (!this->empty(), ExcEmptyMatrix()); - Assert(this->m()-dst_r>src_r_j-src_r_i, - ExcIndexRange(this->m()-dst_r,0,src_r_j-src_r_i)); - Assert(this->n()-dst_c>src_c_j-src_c_i, - ExcIndexRange(this->n()-dst_c,0,src_c_j-src_c_i)); - Assert(dim>src_r_j, ExcIndexRange(dim,0,src_r_j)); - Assert(dim>src_c_j, ExcIndexRange(dim,0,src_r_j)); - Assert(src_r_j>=src_r_i, ExcIndexRange(src_r_j,0,src_r_i)); - Assert(src_c_j>=src_c_i, ExcIndexRange(src_r_j,0,src_r_i)); + AssertIndexRange(src_r_j-src_r_i, this->m()-dst_r); + AssertIndexRange(src_c_j-src_c_i, this->n()-dst_c); + AssertIndexRange(src_r_j,dim); + AssertIndexRange(src_c_j,dim); + AssertIndexRange(src_r_i,src_r_j+1); + AssertIndexRange(src_c_i,src_c_j+1); for (size_type i=0; i(i+src_r_i); + const unsigned int src_c_index = static_cast(j+src_c_i); + (*this)(i+dst_r,j+dst_c) = number(T[src_r_index][src_c_index]); + } } @@ -1651,23 +1644,24 @@ FullMatrix::copy_to (Tensor<2,dim> &T, const size_type src_r_j, const size_type src_c_i, const size_type src_c_j, - const size_type dst_r, - const size_type dst_c) const + const unsigned int dst_r, + const unsigned int dst_c) const { Assert (!this->empty(), ExcEmptyMatrix()); - Assert(dim-dst_r>src_r_j-src_r_i, - ExcIndexRange(dim-dst_r,0,src_r_j-src_r_i)); - Assert(dim-dst_c>src_c_j-src_c_i, - ExcIndexRange(dim-dst_c,0,src_c_j-src_c_i)); - Assert(this->m()>src_r_j, ExcIndexRange(dim,0,src_r_j)); - Assert(this->n()>src_c_j, ExcIndexRange(dim,0,src_r_j)); - Assert(src_r_j>=src_r_i, ExcIndexRange(src_r_j,0,src_r_i)); - Assert(src_c_j>=src_c_i, ExcIndexRange(src_r_j,0,src_r_i)); - + AssertIndexRange(src_r_j-src_r_i,dim-dst_r); + AssertIndexRange(src_c_j-src_c_i,dim-dst_c); + AssertIndexRange(src_r_j,this->m()); + AssertIndexRange(src_r_j,this->n()); + AssertIndexRange(src_r_i,src_r_j+1); + AssertIndexRange(src_c_j,src_c_j+1); for (size_type i=0; i(i+dst_r); + const unsigned int dst_c_index = static_cast(j+dst_c); + T[dst_r_index][dst_c_index] = double ((*this)(i+src_r_i,j+src_c_i)); + } } @@ -1712,7 +1706,7 @@ FullMatrix::print_formatted ( // set output format, but store old // state std::ios::fmtflags old_flags = out.flags(); - unsigned int old_precision = out.precision (precision); + std::streamsize old_precision = out.precision (precision); if (scientific) { diff --git a/source/lac/full_matrix.inst.in b/source/lac/full_matrix.inst.in index 389b72e8fe..0bb473eb43 100644 --- a/source/lac/full_matrix.inst.in +++ b/source/lac/full_matrix.inst.in @@ -25,26 +25,26 @@ for (S : REAL_AND_COMPLEX_SCALARS) std::ostream&, const unsigned int, const unsigned int) const; template void FullMatrix::copy_from<1>( - const Tensor<2,1>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type); + const Tensor<2,1>&, const unsigned int, const unsigned int, unsigned int, const unsigned int, const size_type, const size_type); template void FullMatrix::copy_from<2>( - const Tensor<2,2>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type); + const Tensor<2,2>&, const unsigned int, const unsigned int, unsigned int, const unsigned int, const size_type, const size_type); template void FullMatrix::copy_from<3>( - const Tensor<2,3>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type); + const Tensor<2,3>&, const unsigned int, const unsigned int, unsigned int, const unsigned int, const size_type, const size_type); } for (S : REAL_SCALARS) { template void FullMatrix::copy_to<1>( - Tensor<2,1>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type) const; + Tensor<2,1>&, const size_type, const size_type, const size_type, const size_type, const unsigned int, const unsigned int) const; template void FullMatrix::copy_to<2>( - Tensor<2,2>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type) const; + Tensor<2,2>&, const size_type, const size_type, const size_type, const size_type, const unsigned int, const unsigned int) const; template void FullMatrix::copy_to<3>( - Tensor<2,3>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type) const; + Tensor<2,3>&, const size_type, const size_type, const size_type, const size_type, const unsigned int, const unsigned int) const; } -- 2.39.5