1 //=================================================================================================
2 /*!
3 // \file blaze/math/lapack/clapack/geqlf.h
4 // \brief Header file for the CLAPACK geqlf wrapper functions
5 //
6 // Copyright (C) 2012-2020 Klaus Iglberger - All Rights Reserved
7 //
8 // This file is part of the Blaze library. You can redistribute it and/or modify it under
9 // the terms of the New (Revised) BSD License. Redistribution and use in source and binary
10 // forms, with or without modification, are permitted provided that the following conditions
11 // are met:
12 //
13 // 1. Redistributions of source code must retain the above copyright notice, this list of
14 // conditions and the following disclaimer.
15 // 2. Redistributions in binary form must reproduce the above copyright notice, this list
16 // of conditions and the following disclaimer in the documentation and/or other materials
17 // provided with the distribution.
18 // 3. Neither the names of the Blaze development group nor the names of its contributors
19 // may be used to endorse or promote products derived from this software without specific
20 // prior written permission.
21 //
22 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
23 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
24 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
25 // SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
26 // INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
27 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
28 // BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
29 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30 // ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
31 // DAMAGE.
32 */
33 //=================================================================================================
34
35 #ifndef _BLAZE_MATH_LAPACK_CLAPACK_GEQLF_H_
36 #define _BLAZE_MATH_LAPACK_CLAPACK_GEQLF_H_
37
38
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42
43 #include <blaze/math/blas/Types.h>
44 #include <blaze/util/Complex.h>
45 #include <blaze/util/StaticAssert.h>
46
47
48 //=================================================================================================
49 //
50 // LAPACK FORWARD DECLARATIONS
51 //
52 //=================================================================================================
53
54 //*************************************************************************************************
55 /*! \cond BLAZE_INTERNAL */
56 #if !defined(INTEL_MKL_VERSION)
57 extern "C" {
58
59 void sgeqlf_( blaze::blas_int_t* m, blaze::blas_int_t* n, float* A, blaze::blas_int_t* lda,
60 float* tau, float* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info );
61 void dgeqlf_( blaze::blas_int_t* m, blaze::blas_int_t* n, double* A, blaze::blas_int_t* lda,
62 double* tau, double* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info );
63 void cgeqlf_( blaze::blas_int_t* m, blaze::blas_int_t* n, float* A, blaze::blas_int_t* lda,
64 float* tau, float* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info );
65 void zgeqlf_( blaze::blas_int_t* m, blaze::blas_int_t* n, double* A, blaze::blas_int_t* lda,
66 double* tau, double* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info );
67
68 }
69 #endif
70 /*! \endcond */
71 //*************************************************************************************************
72
73
74
75
76 namespace blaze {
77
78 //=================================================================================================
79 //
80 // LAPACK QL DECOMPOSITION FUNCTIONS (GEQLF)
81 //
82 //=================================================================================================
83
84 //*************************************************************************************************
85 /*!\name LAPACK QL decomposition functions (geqlf) */
86 //@{
87 void geqlf( blas_int_t m, blas_int_t n, float* A, blas_int_t lda,
88 float* tau, float* work, blas_int_t lwork, blas_int_t* info );
89
90 void geqlf( blas_int_t m, blas_int_t n, double* A, blas_int_t lda,
91 double* tau, double* work, blas_int_t lwork, blas_int_t* info );
92
93 void geqlf( blas_int_t m, blas_int_t n, complex<float>* A, blas_int_t lda,
94 complex<float>* tau, complex<float>* work, blas_int_t lwork, blas_int_t* info );
95
96 void geqlf( blas_int_t m, blas_int_t n, complex<double>* A, blas_int_t lda,
97 complex<double>* tau, complex<double>* work, blas_int_t lwork, blas_int_t* info );
98 //@}
99 //*************************************************************************************************
100
101
102 //*************************************************************************************************
103 /*!\brief LAPACK kernel for the QL decomposition of the given dense single precision column-major
104 // matrix.
105 // \ingroup lapack_decomposition
106 //
107 // \param m The number of rows of the given matrix \f$[0..\infty)\f$.
108 // \param n The number of columns of the given matrix \f$[0..\infty)\f$.
109 // \param A Pointer to the first element of the single precision column-major matrix.
110 // \param lda The total number of elements between two columns of the matrix \f$[0..\infty)\f$.
111 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
112 // \param work Auxiliary array; size >= max( 1, \a lwork ).
113 // \param lwork The dimension of the array \a work; size >= max( 1, \a n ).
114 // \param info Return code of the function call.
115 // \return void
116 //
117 // This function performs the dense matrix QL decomposition of a general \a m-by-\a n single
118 // precision column-major matrix based on the LAPACK sgeqlf() function. The resulting decomposition
119 // has the form
120
121 \f[ A = Q \cdot L, \f]
122
123 // where the \c Q is represented as a product of elementary reflectors
124
125 \f[ Q = H(k) ... H(2) H(1) \texttt{, with k = min(m,n).} \f]
126
127 // Each H(i) has the form
128
129 \f[ H(i) = I - tau \cdot v \cdot v^T, \f]
130
131 // where \c tau is a real scalar, and \c v is a real vector with <tt>v(m-k+i+1:m) = 0</tt> and
132 // <tt>v(m-k+i) = 1</tt>. <tt>v(1:m-k+i-1)</tt> is stored on exit in <tt>A(1:m-k+i-1,n-k+i)</tt>,
133 // and \c tau in \c tau(i). Thus in case \a m >= \a n, the lower triangle of the subarray
134 // A(m-n+1:m,1:n) contains the \a n-by-\a n lower triangular matrix \c L and in case \a m <= \a n,
135 // the elements on and below the (\a n-\a m)-th subdiagonal contain the \a m-by-\a n lower
136 // trapezoidal matrix \c L; the remaining elements in combination with the array \c tau represent
137 // the orthogonal matrix \c Q as a product of min(\a m,\a n) elementary reflectors.
138 //
139 // The \a info argument provides feedback on the success of the function call:
140 //
141 // - = 0: The decomposition finished successfully.
142 // - < 0: The i-th argument had an illegal value.
143 //
144 // For more information on the sgeqlf() function, see the LAPACK online documentation browser:
145 //
146 // http://www.netlib.org/lapack/explore-html/
147 //
148 // \note This function can only be used if a fitting LAPACK library, which supports this function,
149 // is available and linked to the executable. Otherwise a call to this function will result in a
150 // linker error.
151 */
geqlf(blas_int_t m,blas_int_t n,float * A,blas_int_t lda,float * tau,float * work,blas_int_t lwork,blas_int_t * info)152 inline void geqlf( blas_int_t m, blas_int_t n, float* A, blas_int_t lda,
153 float* tau, float* work, blas_int_t lwork, blas_int_t* info )
154 {
155 #if defined(INTEL_MKL_VERSION)
156 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
157 #endif
158
159 sgeqlf_( &m, &n, A, &lda, tau, work, &lwork, info );
160 }
161 //*************************************************************************************************
162
163
164 //*************************************************************************************************
165 /*!\brief LAPACK kernel for the QL decomposition of the given dense double precision column-major
166 // matrix.
167 // \ingroup lapack_decomposition
168 //
169 // \param m The number of rows of the given matrix \f$[0..\infty)\f$.
170 // \param n The number of columns of the given matrix \f$[0..\infty)\f$.
171 // \param A Pointer to the first element of the double precision column-major matrix.
172 // \param lda The total number of elements between two columns of the matrix \f$[0..\infty)\f$.
173 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
174 // \param work Auxiliary array; size >= max( 1, \a lwork ).
175 // \param lwork The dimension of the array \a work; size >= max( 1, \a n ).
176 // \param info Return code of the function call.
177 // \return void
178 //
179 // This function performs the dense matrix QL decomposition of a general \a m-by-\a n double
180 // precision column-major matrix based on the LAPACK sgeqlf() function. The resulting decomposition
181 // has the form
182
183 \f[ A = Q \cdot L, \f]
184
185 // where the \c Q is represented as a product of elementary reflectors
186
187 \f[ Q = H(k) ... H(2) H(1) \texttt{, with k = min(m,n).} \f]
188
189 // Each H(i) has the form
190
191 \f[ H(i) = I - tau \cdot v \cdot v^T, \f]
192
193 // where \c tau is a real scalar, and \c v is a real vector with <tt>v(m-k+i+1:m) = 0</tt> and
194 // <tt>v(m-k+i) = 1</tt>. <tt>v(1:m-k+i-1)</tt> is stored on exit in <tt>A(1:m-k+i-1,n-k+i)</tt>,
195 // and \c tau in \c tau(i). Thus in case \a m >= \a n, the lower triangle of the subarray
196 // A(m-n+1:m,1:n) contains the \a n-by-\a n lower triangular matrix \c L and in case \a m <= \a n,
197 // the elements on and below the (\a n-\a m)-th subdiagonal contain the \a m-by-\a n lower
198 // trapezoidal matrix \c L; the remaining elements in combination with the array \c tau represent
199 // the orthogonal matrix \c Q as a product of min(\a m,\a n) elementary reflectors.
200 //
201 // The \a info argument provides feedback on the success of the function call:
202 //
203 // - = 0: The decomposition finished successfully.
204 // - < 0: The i-th argument had an illegal value.
205 //
206 // For more information on the sgeqlf() function, see the LAPACK online documentation browser:
207 //
208 // http://www.netlib.org/lapack/explore-html/
209 //
210 // \note This function can only be used if a fitting LAPACK library, which supports this function,
211 // is available and linked to the executable. Otherwise a call to this function will result in a
212 // linker error.
213 */
geqlf(blas_int_t m,blas_int_t n,double * A,blas_int_t lda,double * tau,double * work,blas_int_t lwork,blas_int_t * info)214 inline void geqlf( blas_int_t m, blas_int_t n, double* A, blas_int_t lda,
215 double* tau, double* work, blas_int_t lwork, blas_int_t* info )
216 {
217 #if defined(INTEL_MKL_VERSION)
218 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
219 #endif
220
221 dgeqlf_( &m, &n, A, &lda, tau, work, &lwork, info );
222 }
223 //*************************************************************************************************
224
225
226 //*************************************************************************************************
227 /*!\brief LAPACK kernel for the QL decomposition of the given dense single precision complex
228 // column-major matrix.
229 // \ingroup lapack_decomposition
230 //
231 // \param m The number of rows of the given matrix \f$[0..\infty)\f$.
232 // \param n The number of columns of the given matrix \f$[0..\infty)\f$.
233 // \param A Pointer to the first element of the single precision complex column-major matrix.
234 // \param lda The total number of elements between two columns of the matrix \f$[0..\infty)\f$.
235 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
236 // \param work Auxiliary array; size >= max( 1, \a lwork ).
237 // \param lwork The dimension of the array \a work; size >= max( 1, \a n ).
238 // \param info Return code of the function call.
239 // \return void
240 //
241 // This function performs the dense matrix QL decomposition of a general \a m-by-\a n single
242 // precision complex column-major matrix based on the LAPACK sgeqlf() function. The resulting
243 // decomposition has the form
244
245 \f[ A = Q \cdot L, \f]
246
247 // where the \c Q is represented as a product of elementary reflectors
248
249 \f[ Q = H(k) ... H(2) H(1) \texttt{, with k = min(m,n).} \f]
250
251 // Each H(i) has the form
252
253 \f[ H(i) = I - tau \cdot v \cdot v^T, \f]
254
255 // where \c tau is a real scalar, and \c v is a real vector with <tt>v(m-k+i+1:m) = 0</tt> and
256 // <tt>v(m-k+i) = 1</tt>. <tt>v(1:m-k+i-1)</tt> is stored on exit in <tt>A(1:m-k+i-1,n-k+i)</tt>,
257 // and \c tau in \c tau(i). Thus in case \a m >= \a n, the lower triangle of the subarray
258 // A(m-n+1:m,1:n) contains the \a n-by-\a n lower triangular matrix \c L and in case \a m <= \a n,
259 // the elements on and below the (\a n-\a m)-th subdiagonal contain the \a m-by-\a n lower
260 // trapezoidal matrix \c L; the remaining elements in combination with the array \c tau represent
261 // the orthogonal matrix \c Q as a product of min(\a m,\a n) elementary reflectors.
262 //
263 // The \a info argument provides feedback on the success of the function call:
264 //
265 // - = 0: The decomposition finished successfully.
266 // - < 0: The i-th argument had an illegal value.
267 //
268 // For more information on the sgeqlf() function, see the LAPACK online documentation browser:
269 //
270 // http://www.netlib.org/lapack/explore-html/
271 //
272 // \note This function can only be used if a fitting LAPACK library, which supports this function,
273 // is available and linked to the executable. Otherwise a call to this function will result in a
274 // linker error.
275 */
geqlf(blas_int_t m,blas_int_t n,complex<float> * A,blas_int_t lda,complex<float> * tau,complex<float> * work,blas_int_t lwork,blas_int_t * info)276 inline void geqlf( blas_int_t m, blas_int_t n, complex<float>* A, blas_int_t lda,
277 complex<float>* tau, complex<float>* work, blas_int_t lwork, blas_int_t* info )
278 {
279 BLAZE_STATIC_ASSERT( sizeof( complex<float> ) == 2UL*sizeof( float ) );
280
281 #if defined(INTEL_MKL_VERSION)
282 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
283 BLAZE_STATIC_ASSERT( sizeof( MKL_Complex8 ) == sizeof( complex<float> ) );
284 using ET = MKL_Complex8;
285 #else
286 using ET = float;
287 #endif
288
289 cgeqlf_( &m, &n, reinterpret_cast<ET*>( A ), &lda, reinterpret_cast<ET*>( tau ),
290 reinterpret_cast<ET*>( work ), &lwork, info );
291 }
292 //*************************************************************************************************
293
294
295 //*************************************************************************************************
296 /*!\brief LAPACK kernel for the QL decomposition of the given dense double precision complex
297 // column-major matrix.
298 // \ingroup lapack_decomposition
299 //
300 // \param m The number of rows of the given matrix \f$[0..\infty)\f$.
301 // \param n The number of columns of the given matrix \f$[0..\infty)\f$.
302 // \param A Pointer to the first element of the double precision complex column-major matrix.
303 // \param lda The total number of elements between two columns of the matrix \f$[0..\infty)\f$.
304 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
305 // \param work Auxiliary array; size >= max( 1, \a lwork ).
306 // \param lwork The dimension of the array \a work; size >= max( 1, \a n ).
307 // \param info Return code of the function call.
308 // \return void
309 //
310 // This function performs the dense matrix QL decomposition of a general \a m-by-\a n double
311 // precision complex column-major matrix based on the LAPACK sgeqlf() function. The resulting
312 // decomposition has the form
313
314 \f[ A = Q \cdot L, \f]
315
316 // where the \c Q is represented as a product of elementary reflectors
317
318 \f[ Q = H(k) ... H(2) H(1) \texttt{, with k = min(m,n).} \f]
319
320 // Each H(i) has the form
321
322 \f[ H(i) = I - tau \cdot v \cdot v^T, \f]
323
324 // where \c tau is a real scalar, and \c v is a real vector with <tt>v(m-k+i+1:m) = 0</tt> and
325 // <tt>v(m-k+i) = 1</tt>. <tt>v(1:m-k+i-1)</tt> is stored on exit in <tt>A(1:m-k+i-1,n-k+i)</tt>,
326 // and \c tau in \c tau(i). Thus in case \a m >= \a n, the lower triangle of the subarray
327 // A(m-n+1:m,1:n) contains the \a n-by-\a n lower triangular matrix \c L and in case \a m <= \a n,
328 // the elements on and below the (\a n-\a m)-th subdiagonal contain the \a m-by-\a n lower
329 // trapezoidal matrix \c L; the remaining elements in combination with the array \c tau represent
330 // the orthogonal matrix \c Q as a product of min(\a m,\a n) elementary reflectors.
331 //
332 // The \a info argument provides feedback on the success of the function call:
333 //
334 // - = 0: The decomposition finished successfully.
335 // - < 0: The i-th argument had an illegal value.
336 //
337 // For more information on the sgeqlf() function, see the LAPACK online documentation browser:
338 //
339 // http://www.netlib.org/lapack/explore-html/
340 //
341 // \note This function can only be used if a fitting LAPACK library, which supports this function,
342 // is available and linked to the executable. Otherwise a call to this function will result in a
343 // linker error.
344 */
geqlf(blas_int_t m,blas_int_t n,complex<double> * A,blas_int_t lda,complex<double> * tau,complex<double> * work,blas_int_t lwork,blas_int_t * info)345 inline void geqlf( blas_int_t m, blas_int_t n, complex<double>* A, blas_int_t lda,
346 complex<double>* tau, complex<double>* work, blas_int_t lwork, blas_int_t* info )
347 {
348 BLAZE_STATIC_ASSERT( sizeof( complex<double> ) == 2UL*sizeof( double ) );
349
350 #if defined(INTEL_MKL_VERSION)
351 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
352 BLAZE_STATIC_ASSERT( sizeof( MKL_Complex16 ) == sizeof( complex<double> ) );
353 using ET = MKL_Complex16;
354 #else
355 using ET = double;
356 #endif
357
358 zgeqlf_( &m, &n, reinterpret_cast<ET*>( A ), &lda, reinterpret_cast<ET*>( tau ),
359 reinterpret_cast<ET*>( work ), &lwork, info );
360 }
361 //*************************************************************************************************
362
363 } // namespace blaze
364
365 #endif
366