]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Corrections 6062/head
authorBenjamin Brands <benjamin.brands@fau.de>
Mon, 19 Mar 2018 15:57:42 +0000 (16:57 +0100)
committerbenjamin <benjamin.brands@fau.de>
Fri, 6 Apr 2018 20:51:14 +0000 (22:51 +0200)
include/deal.II/lac/lapack_templates.h
include/deal.II/lac/scalapack.h
source/lac/scalapack.cc
tests/scalapack/scalapack_15.cc
tests/scalapack/scalapack_15_b.cc

index cab270705828f4a7327376940b55a795c0b24ecc..e3c11a310e06508c8ebf8b681653bb28ad5cd473 100644 (file)
@@ -1403,36 +1403,7 @@ syevr(const char *jobz,
 {
   dsyevr_(jobz,range,uplo,n,A,lda,vl,vu,il,iu,abstol,m,w,z,ldz,isuppz,work,lwork,iwork,liwork,info);
 }
-#else
-inline void
-syevr(const char * /*jobz*/,
-      const char * /*range*/,
-      const char * /*uplo*/,
-      const types::blas_int * /*n*/,
-      double * /*A*/,
-      const types::blas_int * /*lda*/,
-      const double * /*vl*/,
-      const double * /*vu*/,
-      const types::blas_int * /*il*/,
-      const types::blas_int * /*iu*/,
-      const double * /*abstol*/,
-      types::blas_int * /*m*/,
-      double * /*w*/,
-      double * /*z*/,
-      const types::blas_int * /*ldz*/,
-      types::blas_int * /*isuppz*/,
-      double * /*work*/,
-      types::blas_int * /*lwork*/,
-      types::blas_int * /*iwork*/,
-      types::blas_int * /*liwork*/,
-      types::blas_int * /*info*/)
-{
-  Assert (false, LAPACKSupport::ExcMissing("dsyevr"));
-}
-#endif
 
-
-#ifdef DEAL_II_WITH_LAPACK
 inline void
 syevr(const char *jobz,
       const char *range,
@@ -1458,32 +1429,6 @@ syevr(const char *jobz,
 {
   ssyevr_(jobz,range,uplo,n,A,lda,vl,vu,il,iu,abstol,m,w,z,ldz,isuppz,work,lwork,iwork,liwork,info);
 }
-#else
-inline void
-syevr(const char * /*jobz*/,
-      const char * /*range*/,
-      const char * /*uplo*/,
-      const types::blas_int * /*n*/,
-      float * /*A*/,
-      const types::blas_int * /*lda*/,
-      const float * /*vl*/,
-      const float * /*vu*/,
-      const types::blas_int * /*il*/,
-      const types::blas_int * /*iu*/,
-      const float * /*abstol*/,
-      types::blas_int * /*m*/,
-      float * /*w*/,
-      float * /*z*/,
-      const types::blas_int * /*ldz*/,
-      types::blas_int * /*isuppz*/,
-      float * /*work*/,
-      types::blas_int * /*lwork*/,
-      types::blas_int * /*iwork*/,
-      types::blas_int * /*liwork*/,
-      types::blas_int * /*info*/)
-{
-  Assert (false, LAPACKSupport::ExcMissing("ssyevr"));
-}
 #endif
 
 
@@ -1795,26 +1740,12 @@ lamch(const char *chmach, double &precision)
 {
   precision = dlamch_(chmach);
 }
-#else
-inline void
-lamch(const char * /*chmach*/, double & /*precision*/)
-{
-  Assert(false, LAPACKSupport::ExcMissing("dlamch"));
-}
-#endif
 
-#ifdef DEAL_II_WITH_LAPACK
 inline void
 lamch(const char *chmach, float &precision)
 {
   precision = slamch_(chmach);
 }
-#else
-inline void
-lamch(const char * /*chmach*/, float & /*precision*/)
-{
-  Assert(false, LAPACKSupport::ExcMissing("slamch"));
-}
 #endif
 
 DEAL_II_NAMESPACE_CLOSE
index 5a4aa99ac545fc51b1ea9abcec965eb612240fb3..0f0a1bf3e73764a936b6d1814d60d677acfc3027 100644 (file)
@@ -423,7 +423,7 @@ public:
 
   /**
    * Computing selected eigenvalues and, optionally, the eigenvectors of the real symmetric
-   * matrix $A \in \mathbb{R}^{M \times M}$ using the MRRR algorithm.
+   * matrix $\mathbf{A} \in \mathbb{R}^{M \times M}$ using the MRRR algorithm.
    *
    * The eigenvalues/eigenvectors are selected by prescribing a range of indices @p index_limits.
    *
@@ -439,7 +439,8 @@ public:
                                                              const bool compute_eigenvectors);
 
   /**
-   * Computing selected eigenvalues and, optionally, the eigenvectors using the MRRR algorithm.
+   * Computing selected eigenvalues and, optionally, the eigenvectors of the real symmetric
+   * matrix $\mathbf{A} \in \mathbb{R}^{M \times M}$ using the MRRR algorithm.
    * The eigenvalues/eigenvectors are selected by prescribing a range of values @p value_limits for the eigenvalues.
    *
    * If successful, the computed eigenvalues are arranged in ascending order.
@@ -635,14 +636,19 @@ private:
                                                  std::make_pair(std::numeric_limits<NumberType>::quiet_NaN(),std::numeric_limits<NumberType>::quiet_NaN()));
 
   /**
-   * Computing selected eigenvalues and, optionally, the eigenvectors using the MRRR algorithm.
+   * Computing selected eigenvalues and, optionally, the eigenvectors of the real symmetric
+   * matrix $\mathbf{A} \in \mathbb{R}^{M \times M}$ using the MRRR algorithm.
    * The eigenvalues/eigenvectors are selected by either prescribing a range of indices @p index_limits
    * or a range of values @p value_limits for the eigenvalues. The function will throw an exception
    * if both ranges are prescribed (meaning that both ranges differ from the default value)
    * as this ambiguity is prohibited.
+   *
+   * By calling this function the original content of the matrix will be overwritten.
+   * If requested, the eigenvectors are stored in the columns of the matrix.
+   * Also in the case that just the eigenvalues are required,
+   * the content of the matrix will be overwritten.
+   *
    * If successful, the computed eigenvalues are arranged in ascending order.
-   * The eigenvectors are stored in the columns of the matrix, thereby
-   * overwriting the original content of the matrix.
    */
   std::vector<NumberType> eigenpairs_symmetric_MRRR(const bool compute_eigenvectors,
                                                     const std::pair<unsigned int,unsigned int> &index_limits=
index 6ab1b4439ad294f2e5ddb37f681b5045a567f2d1..97ad95f20a09a1cf27fc20890a25c497b58ff115 100644 (file)
@@ -915,14 +915,14 @@ std::vector<NumberType> ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_by_ind
     const bool compute_eigenvectors)
 {
   // Check validity of index limits.
-  Assert (index_limits.first < (unsigned int)n_rows,ExcIndexRange(index_limits.first,0,n_rows));
-  Assert (index_limits.second < (unsigned int)n_rows,ExcIndexRange(index_limits.second,0,n_rows));
+  AssertIndexRange(index_limits.first,static_cast<unsigned int>(n_rows));
+  AssertIndexRange(index_limits.second,static_cast<unsigned int>(n_rows));
 
-  std::pair<unsigned int,unsigned int> idx = std::make_pair(std::min(index_limits.first,index_limits.second),
-                                                            std::max(index_limits.first,index_limits.second));
+  const std::pair<unsigned int,unsigned int> idx = std::make_pair(std::min(index_limits.first,index_limits.second),
+                                                   std::max(index_limits.first,index_limits.second));
 
   // Compute all eigenvalues/eigenvectors.
-  if (idx.first==0 && idx.second==(unsigned int)n_rows-1)
+  if (idx.first==0 && idx.second==static_cast<unsigned int>(n_rows-1))
     return eigenpairs_symmetric_MRRR(compute_eigenvectors);
   else
     return eigenpairs_symmetric_MRRR(compute_eigenvectors,idx);
@@ -934,10 +934,10 @@ template <typename NumberType>
 std::vector<NumberType> ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_by_value_MRRR(const std::pair<NumberType,NumberType> &value_limits,
     const bool compute_eigenvectors)
 {
-  Assert (!std::isnan(value_limits.first),ExcMessage("value_limits.first is NaN"));
-  Assert (!std::isnan(value_limits.second),ExcMessage("value_limits.second is NaN"));
+  AssertIsFinite(value_limits.first);
+  AssertIsFinite(value_limits.second);
 
-  std::pair<unsigned int,unsigned int> indices = std::make_pair(numbers::invalid_unsigned_int,numbers::invalid_unsigned_int);
+  const std::pair<unsigned int,unsigned int> indices = std::make_pair(numbers::invalid_unsigned_int,numbers::invalid_unsigned_int);
 
   return eigenpairs_symmetric_MRRR(compute_eigenvectors,indices,value_limits);
 }
@@ -960,7 +960,8 @@ ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_MRRR(const bool compute_eigenv
   const bool use_values = (std::isnan(eigenvalue_limits.first) || std::isnan(eigenvalue_limits.second)) ? false : true;
   const bool use_indices = ((eigenvalue_idx.first==numbers::invalid_unsigned_int) || (eigenvalue_idx.second==numbers::invalid_unsigned_int)) ? false : true;
 
-  Assert(!(use_values && use_indices),ExcMessage("Prescribing both the index and value range for the eigenvalues is ambiguous"));
+  Assert(!(use_values && use_indices),
+         ExcMessage("Prescribing both the index and value range for the eigenvalues is ambiguous"));
 
   // If computation of eigenvectors is not required, use a sufficiently small distributed matrix.
   std::unique_ptr<ScaLAPACKMatrix<NumberType>> eigenvectors = compute_eigenvectors ?
@@ -983,9 +984,8 @@ ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_MRRR(const bool compute_eigenv
        * For jobz==N only eigenvalues are computed, for jobz='V' also the eigenvectors of the matrix are computed.
        */
       char jobz = compute_eigenvectors ? 'V' : 'N';
-      char range='A';
       // Default value is to compute all eigenvalues and optionally eigenvectors.
-      bool all_eigenpairs=true;
+      char range='A';
       NumberType vl=NumberType(),vu=NumberType();
       int il=1,iu=1;
 
@@ -996,12 +996,10 @@ ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_MRRR(const bool compute_eigenv
           if (!use_values)
             {
               range = 'A';
-              all_eigenpairs = true;
             }
           else
             {
               range = 'V';
-              all_eigenpairs = false;
               vl = std::min(eigenvalue_limits.first,eigenvalue_limits.second);
               vu = std::max(eigenvalue_limits.first,eigenvalue_limits.second);
             }
@@ -1009,7 +1007,6 @@ ScaLAPACKMatrix<NumberType>::eigenpairs_symmetric_MRRR(const bool compute_eigenv
       else
         {
           range = 'I';
-          all_eigenpairs = false;
           // As Fortran starts counting/indexing from 1 unlike C/C++, where it starts from 0.
           il = std::min(eigenvalue_idx.first,eigenvalue_idx.second) + 1;
           iu = std::max(eigenvalue_idx.first,eigenvalue_idx.second) + 1;
index 3b1ef5ceb57421d3ed6f6c3720f896d40e9af41c..af60c275d765992241182345d4a081f42de30e59 100644 (file)
@@ -132,7 +132,7 @@ void test(const unsigned int size, const unsigned int block_size, const NumberTy
   for (unsigned int i=0; i<max_n_eigenvalues; ++i)
     {
       const NumberType product = p_eigenvectors_[i] * s_eigenvectors_[i];
-      //the requirement for alignment of the eigenvectors has to be released (primarily for floats)
+      //the requirement for alignment of the eigenvectors has to be released
       AssertThrow (std::abs(std::abs(product)-1) < tol*10,
                    ExcInternalError());
     }
index c8b675db821df36a340ce76cd59f541864d1a1a7..20d51187944b879c307cd3c3953be5554bf805e2 100644 (file)
@@ -132,7 +132,7 @@ void test(const unsigned int size, const unsigned int block_size, const NumberTy
   for (unsigned int i=0; i<max_n_eigenvalues; ++i)
     {
       const NumberType product = p_eigenvectors_[i] * s_eigenvectors_[i];
-      //the requirement for alignment of the eigenvectors has to be released (primarily for floats)
+      //the requirement for alignment of the eigenvectors has to be released
       AssertThrow (std::abs(std::abs(product)-1) < tol*10,
                    ExcInternalError());
     }

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.