vendor/gonum.org/v1/gonum/lapack/gonum/iparmq.go

   1 // Copyright ©2016 The Gonum Authors. All rights reserved.
   2 // Use of this source code is governed by a BSD-style
   3 // license that can be found in the LICENSE file.
   4
   5 package gonum
   6
   7 import "math"
   8
   9 // Iparmq returns problem and machine dependent parameters useful for Dhseqr and
  10 // related subroutines for eigenvalue problems.
  11 //
  12 // ispec specifies the parameter to return:
  13 //  12: Crossover point between Dlahqr and Dlaqr0. Will be at least 11.
  14 //  13: Deflation window size.
  15 //  14: Nibble crossover point. Determines when to skip a multi-shift QR sweep.
  16 //  15: Number of simultaneous shifts in a multishift QR iteration.
  17 //  16: Select structured matrix multiply.
  18 // For other values of ispec Iparmq will panic.
  19 //
  20 // name is the name of the calling function. name must be in uppercase but this
  21 // is not checked.
  22 //
  23 // opts is not used and exists for future use.
  24 //
  25 // n is the order of the Hessenberg matrix H.
  26 //
  27 // ilo and ihi specify the block [ilo:ihi+1,ilo:ihi+1] that is being processed.
  28 //
  29 // lwork is the amount of workspace available.
  30 //
  31 // Except for ispec input parameters are not checked.
  32 //
  33 // Iparmq is an internal routine. It is exported for testing purposes.
  34 func (Implementation) Iparmq(ispec int, name, opts string, n, ilo, ihi, lwork int) int {
  35         nh := ihi - ilo + 1
  36         ns := 2
  37         switch {
  38         case nh >= 30:
  39                 ns = 4
  40         case nh >= 60:
  41                 ns = 10
  42         case nh >= 150:
  43                 ns = max(10, nh/int(math.Log(float64(nh))/math.Ln2))
  44         case nh >= 590:
  45                 ns = 64
  46         case nh >= 3000:
  47                 ns = 128
  48         case nh >= 6000:
  49                 ns = 256
  50         }
  51         ns = max(2, ns-(ns%2))
  52
  53         switch ispec {
  54         default:
  55                 panic("lapack: bad ispec")
  56
  57         case 12:
  58                 // Matrices of order smaller than nmin get sent to Dlahqr, the
  59                 // classic double shift algorithm. This must be at least 11.
  60                 const nmin = 75
  61                 return nmin
  62
  63         case 13:
  64                 const knwswp = 500
  65                 if nh <= knwswp {
  66                         return ns
  67                 }
  68                 return 3 * ns / 2
  69
  70         case 14:
  71                 // Skip a computationally expensive multi-shift QR sweep with
  72                 // Dlaqr5 whenever aggressive early deflation finds at least
  73                 // nibble*(window size)/100 deflations. The default, small,
  74                 // value reflects the expectation that the cost of looking
  75                 // through the deflation window with Dlaqr3 will be
  76                 // substantially smaller.
  77                 const nibble = 14
  78                 return nibble
  79
  80         case 15:
  81                 return ns
  82
  83         case 16:
  84                 if len(name) != 6 {
  85                         panic("lapack: bad name")
  86                 }
  87                 const (
  88                         k22min = 14
  89                         kacmin = 14
  90                 )
  91                 var acc22 int
  92                 switch {
  93                 case name[1:] == "GGHRD" || name[1:] == "GGHD3":
  94                         acc22 = 1
  95                         if nh >= k22min {
  96                                 acc22 = 2
  97                         }
  98                 case name[3:] == "EXC":
  99                         if nh >= kacmin {
 100                                 acc22 = 1
 101                         }
 102                         if nh >= k22min {
 103                                 acc22 = 2
 104                         }
 105                 case name[1:] == "HSEQR" || name[1:5] == "LAQR":
 106                         if ns >= kacmin {
 107                                 acc22 = 1
 108                         }
 109                         if ns >= k22min {
 110                                 acc22 = 2
 111                         }
 112                 }
 113                 return acc22
 114         }
 115 }