]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change LAPACKFullMatrix::size_type to unsigned int. 5652/head
authorDavid Wells <wellsd2@rpi.edu>
Wed, 20 Dec 2017 15:54:51 +0000 (10:54 -0500)
committerDavid Wells <wellsd2@rpi.edu>
Thu, 21 Dec 2017 01:56:41 +0000 (20:56 -0500)
This fixes compilation with 64 bit indices and also makes things more consistent
with base class functions (e.g., TransposeTable::n_rows() returns an unsigned
int).

doc/news/changes/incompatibilities/20171220DavidWells [new file with mode: 0644]
include/deal.II/lac/lapack_full_matrix.h
source/lac/lapack_full_matrix.cc

diff --git a/doc/news/changes/incompatibilities/20171220DavidWells b/doc/news/changes/incompatibilities/20171220DavidWells
new file mode 100644 (file)
index 0000000..57124b8
--- /dev/null
@@ -0,0 +1,3 @@
+Changed: The type <code>LAPACKFullMatrix::size_type</code> is now <code>unsigned int</code>, which matches <code>FullMatrix::size_type</code>.
+<br>
+(David Wells, 2017/12/20)
index 518fede4b3cf1d89a76cd11975e26baa49a3700e..f0aa974f29f61744163339217bd5381204b463f9 100644 (file)
@@ -57,7 +57,7 @@ public:
   /**
    * Declare type for container size.
    */
-  typedef types::global_dof_index size_type;
+  typedef unsigned int size_type;
 
   /**
    * Constructor. Initialize the matrix as a square matrix with dimension
@@ -217,14 +217,14 @@ public:
    *
    * @note The matrix is of dimension $m \times n$.
    */
-  unsigned int m () const;
+  size_type m () const;
 
   /**
    * Return the dimension of the domain space.
    *
    * @note The matrix is of dimension $m \times n$.
    */
-  unsigned int n () const;
+  size_type n () const;
 
   /**
    * Fill rectangular block.
@@ -859,7 +859,7 @@ private:
 
 template <typename number>
 inline
-unsigned int
+typename LAPACKFullMatrix<number>::size_type
 LAPACKFullMatrix<number>::m () const
 {
   return this->n_rows ();
@@ -867,7 +867,7 @@ LAPACKFullMatrix<number>::m () const
 
 template <typename number>
 inline
-unsigned int
+typename LAPACKFullMatrix<number>::size_type
 LAPACKFullMatrix<number>::n () const
 {
   return this->n_cols ();
index d4155d510fbce89740152aa3fef57807064a7c64..c217d22381bfc5bb8fc263c900c78d5f2b983d4b 100644 (file)
@@ -88,8 +88,8 @@ LAPACKFullMatrix<number>::grow_or_shrink (const size_type n)
   TransposeTable<number> copy(*this);
   const size_type s = std::min(std::min(this->m(), n), this->n());
   this->TransposeTable<number>::reinit (n, n);
-  for (unsigned int i = 0; i < s; ++i)
-    for (unsigned int j = 0; j < s; ++j)
+  for (size_type i = 0; i < s; ++i)
+    for (size_type j = 0; j < s; ++j)
       (*this)(i,j) = copy(i,j);
 }
 
@@ -162,8 +162,8 @@ LAPACKFullMatrix<number>::operator*= (const number factor)
          state == LAPACKSupport::inverse_matrix,
          ExcState(state));
 
-  for (unsigned int column = 0; column<this->n(); ++column)
-    for (unsigned int row = 0; row<this->m(); ++row)
+  for (size_type column = 0; column<this->n(); ++column)
+    for (size_type row = 0; row<this->m(); ++row)
       (*this)(row,column) *= factor;
 
   return *this;
@@ -181,8 +181,8 @@ LAPACKFullMatrix<number>::operator/= (const number factor)
   AssertIsFinite(factor);
   Assert (factor != number(0.), ExcZero() );
 
-  for (unsigned int column = 0; column<this->n(); ++column)
-    for (unsigned int row = 0; row<this->m(); ++row)
+  for (size_type column = 0; column<this->n(); ++column)
+    for (size_type row = 0; row<this->m(); ++row)
       (*this)(row,column) /= factor;
 
   return *this;
@@ -218,7 +218,7 @@ namespace
                       const number a,
                       const Vector<number> &v)
   {
-    const unsigned int N = A.n();
+    const typename LAPACKFullMatrix<number>::size_type N = A.n();
     Vector<number> z(v);
     // Cholesky update / downdate, see
     // 6.5.4 Cholesky Updating and Downdating, Golub 2013 Matrix computations
@@ -238,11 +238,11 @@ namespace
         // rotations to make the matrix lower-triangular
         // Also see LINPACK's dchud http://www.netlib.org/linpack/dchud.f
         z *= std::sqrt(a);
-        for (unsigned int k = 0; k < N; ++k)
+        for (typename LAPACKFullMatrix<number>::size_type k = 0; k < N; ++k)
           {
             const std::array<number,3> csr = Utilities::LinearAlgebra::givens_rotation(A(k,k),z(k));
             A(k,k) = csr[2];
-            for (unsigned int i = k+1; i < N; ++i)
+            for (typename LAPACKFullMatrix<number>::size_type i = k+1; i < N; ++i)
               {
                 const number t = A(i,k);
                 A(i,k) =  csr[0] * A(i,k) + csr[1] * z(i);
@@ -278,11 +278,11 @@ namespace
         // https://infoscience.epfl.ch/record/161468/files/cholupdate.pdf and
         // "Analysis of a recursive Least Squares Hyperbolic Rotation Algorithm for Signal Processing", Alexander, Pan, Plemmons, 1988.
         z *= std::sqrt(-a);
-        for (unsigned int k = 0; k < N; ++k)
+        for (typename LAPACKFullMatrix<number>::size_type k = 0; k < N; ++k)
           {
             const std::array<number,3> csr = Utilities::LinearAlgebra::hyperbolic_rotation(A(k,k),z(k));
             A(k,k) = csr[2];
-            for (unsigned int i = k+1; i < N; ++i)
+            for (typename LAPACKFullMatrix<number>::size_type i = k+1; i < N; ++i)
               {
                 const number t = A(i,k);
                 A(i,k) =  csr[0] * A(i,k) - csr[1] * z(i);
@@ -319,18 +319,21 @@ LAPACKFullMatrix<number>::rank1_update(const number a,
 
   if (state == LAPACKSupport::matrix)
     {
-      const int N = this->n_rows();
-      const char uplo = LAPACKSupport::U;
-      const int lda = N;
-      const int incx=1;
+      {
+        const int N = this->n_rows();
+        const char uplo = LAPACKSupport::U;
+        const int lda = N;
+        const int incx=1;
 
-      syr(&uplo, &N, &a, v.begin(), &incx, this->values.begin(), &lda);
+        syr(&uplo, &N, &a, v.begin(), &incx, this->values.begin(), &lda);
+      }
 
+      const size_type N = this->n_rows();
       // FIXME: we should really only update upper or lower triangular parts
       // of symmetric matrices and make sure the interface is consistent,
       // for example operator(i,j) gives correct results regardless of storage.
-      for (unsigned int i = 0; i < N; ++i)
-        for (unsigned int j = 0; j < i; ++j)
+      for (size_type i = 0; i < N; ++i)
+        for (size_type j = 0; j < i; ++j)
           (*this)(i,j) = (*this)(j,i);
     }
   else if (state == LAPACKSupport::cholesky)

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.