gonum.org/v1/gonum@v0.14.0/lapack/gonum/iparmq.go (about)

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