1 /*
  2  *   Copyright (c) 2011, Michael Lehn
  3  *
  4  *   All rights reserved.
  5  *
  6  *   Redistribution and use in source and binary forms, with or without
  7  *   modification, are permitted provided that the following conditions
  8  *   are met:
  9  *
 10  *   1) Redistributions of source code must retain the above copyright
 11  *      notice, this list of conditions and the following disclaimer.
 12  *   2) Redistributions in binary form must reproduce the above copyright
 13  *      notice, this list of conditions and the following disclaimer in
 14  *      the documentation and/or other materials provided with the
 15  *      distribution.
 16  *   3) Neither the name of the FLENS development group nor the names of
 17  *      its contributors may be used to endorse or promote products derived
 18  *      from this software without specific prior written permission.
 19  *
 20  *   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 21  *   "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 22  *   LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 23  *   A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 24  *   OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 25  *   SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 26  *   LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 27  *   DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 28  *   THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 29  *   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 30  *   OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 31  */
 32 
 33 /* Based on
 34        INTEGER FUNCTION IPARMQ( spec, NAME, OPTS, N, ILO, IHI, LWORK )
 35  *
 36  *  -- LAPACK auxiliary routine (version 3.2) --
 37  *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 38  *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 39  *     November 2006
 40  */
 41 
 42 #ifndef FLENS_LAPACK_AUX_IPARMQ_TCC
 43 #define FLENS_LAPACK_AUX_IPARMQ_TCC 1
 44 
 45 #include <complex>
 46 #include <string>
 47 
 48 #include <flens/aux/issame.h>
 49 #include <flens/lapack/lapack.h>
 50 
 51 namespace flens { namespace lapack {
 52 
 53 //== publich interface / generic lapack implementation =========================
 54 
 55 template <typename T>
 56 int
 57 iparmq(int spec, const char *, const char *, intint iLo, int iHi, int)
 58 {
 59     using std::log;
 60     using std::max;
 61 
 62     const int INMIN = 12;
 63     const int INWIN = 13;
 64     const int INIBL = 14;
 65     const int ISHFTS = 15;
 66     const int IACC22 = 16;
 67 
 68     const float Two(2);
 69 
 70     const int nMin = 75;
 71     const int k22Min = 14;
 72     const int kacMin = 14;
 73     const int nibble = 14;
 74     const int knwSwp = 500;
 75 
 76     int nh = -1,
 77         ns = -1,
 78         result = -1;
 79 
 80     if ((spec==ISHFTS) || (spec==INWIN) || (spec==IACC22)) {
 81 //
 82 //      ==== Set the number simultaneous shifts ====
 83 //
 84         nh = iHi - iLo + 1;
 85         ns = 2;
 86         if (nh>=30) {
 87             ns = 4;
 88         }
 89         if (nh>=60) {
 90             ns = 10;
 91         }
 92         if (nh>=150) {
 93             ns = max(10int(nh / nint(log(float(nh)) / log(Two))));
 94         }
 95         if (nh>=590) {
 96             ns = 64;
 97         }
 98         if (nh>=3000) {
 99             ns = 128;
100         }
101         if (nh>=6000) {
102             ns = 256;
103         }
104         ns = max(2, ns - (ns % 2));
105     }
106 
107     if (spec==INMIN) {
108 //
109 //
110 //      ===== Matrices of order smaller than NMIN get sent
111 //      .     to xLAHQR, the classic double shift algorithm.
112 //      .     This must be at least 11. ====
113 //
114         result = nMin;
115 
116     } else if (spec==INIBL) {
117 //
118 //      ==== INIBL: skip a multi-shift qr iteration and
119 //      .    whenever aggressive early deflation finds
120 //      .    at least (NIBBLE*(window size)/100) deflations. ====
121 //
122         result = nibble;
123 
124     } else if (spec==ISHFTS) {
125 //
126 //      ==== NSHFTS: The number of simultaneous shifts =====
127 //
128         result = ns;
129 
130     } else if (spec==INWIN) {
131 //
132 //      ==== NW: deflation window size.  ====
133 //
134         if (nh<=knwSwp) {
135             result = ns;
136         } else {
137             result = 3*ns / 2;
138         }
139 
140     } else if (spec==IACC22) {
141 //
142 //      ==== IACC22: Whether to accumulate reflections
143 //      .     before updating the far-from-diagonal elements
144 //      .     and whether to use 2-by-2 block structure while
145 //      .     doing it.  A small amount of work could be saved
146 //      .     by making this choice dependent also upon the
147 //      .     NH=IHI-ILO+1.
148 //
149         result = 0;
150         if (ns>=kacMin) {
151             result = 1;
152         }
153         if (ns>=k22Min) {
154             result = 2;
155         }
156 
157     } else {
158 //
159 //      ===== invalid value of ispec =====
160 //
161         result = -1;
162 
163     }
164 
165     return result;
166 }
167 
168 } } // namespace lapack, flens
169 
170 #endif // FLENS_LAPACK_AUX_IPARMQ_TCC