1 //=================================================================================================
2 /*!
3 // \file blaze/math/lapack/clapack/ormlq.h
4 // \brief Header file for the CLAPACK ormlq 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_ORMLQ_H_
36 #define _BLAZE_MATH_LAPACK_CLAPACK_ORMLQ_H_
37
38
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42
43 #include <blaze/math/blas/Types.h>
44 #include <blaze/util/StaticAssert.h>
45 #include <blaze/util/Types.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 sormlq_( char* side, char* trans, blaze::blas_int_t* m, blaze::blas_int_t* n,
60 blaze::blas_int_t* k, float* A, blaze::blas_int_t* lda, float* tau, float* C,
61 blaze::blas_int_t* ldc, float* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info,
62 blaze::fortran_charlen_t nside, blaze::fortran_charlen_t ntrans );
63 void dormlq_( char* side, char* trans, blaze::blas_int_t* m, blaze::blas_int_t* n,
64 blaze::blas_int_t* k, double* A, blaze::blas_int_t* lda, double* tau, double* C,
65 blaze::blas_int_t* ldc, double* work, blaze::blas_int_t* lwork, blaze::blas_int_t* info,
66 blaze::fortran_charlen_t nside, blaze::fortran_charlen_t ntrans );
67
68 }
69 #endif
70 /*! \endcond */
71 //*************************************************************************************************
72
73
74
75
76 namespace blaze {
77
78 //=================================================================================================
79 //
80 // LAPACK FUNCTIONS TO MULTIPLY Q FROM A LQ DECOMPOSITION WITH A MATRIX (ORMLQ)
81 //
82 //=================================================================================================
83
84 //*************************************************************************************************
85 /*!\name LAPACK functions to multiply Q from a LQ decomposition with a matrix (ormlq) */
86 //@{
87 void ormlq( char side, char trans, blas_int_t m, blas_int_t n, blas_int_t k,
88 const float* A, blas_int_t lda, const float* tau, float* C, blas_int_t ldc,
89 float* work, blas_int_t lwork, blas_int_t* info );
90
91 void ormlq( char side, char trans, blas_int_t m, blas_int_t n, blas_int_t k,
92 const double* A, blas_int_t lda, const double* tau, double* C, blas_int_t ldc,
93 double* work, blas_int_t lwork, blas_int_t* info );
94 //@}
95 //*************************************************************************************************
96
97
98 //*************************************************************************************************
99 /*!\brief LAPACK kernel for the multiplication of the single precision Q from a LQ decomposition
100 // with another matrix.
101 // \ingroup lapack_decomposition
102 //
103 // \param side \c 'L' to apply \f$ Q \f$ or \f$ Q^T \f$ from the left, \c 'R' to apply from the right.
104 // \param trans \c 'N' for \f$ Q \f$, \c 'T' for \f$ Q^T \f$.
105 // \param m The number of rows of the matrix \a C \f$[0..\infty)\f$.
106 // \param n The number of columns of the matrix \a C \f$[0..\infty)\f$.
107 // \param k The number of elementary reflectors, whose product defines the matrix \a Q.
108 // \param A Pointer to the first element of the LQ decomposed single precision column-major matrix.
109 // \param lda The total number of elements between two columns of the matrix \a A \f$[0..\infty)\f$.
110 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
111 // \param C Pointer to the first element of the single precision column-major matrix multiplicator.
112 // \param ldc The total number of elements between two columns of the matrix \a C \f$[0..\infty)\f$.
113 // \param work Auxiliary array; size >= max( 1, \a lwork ).
114 // \param lwork The dimension of the array \a work.
115 // \param info Return code of the function call.
116 // \return void
117 //
118 // This function multiplies the \a Q matrix resulting from the LQ decomposition of the sgelqf()
119 // function with the given general single precision \a m-by-\a n matrix \a C. Depending on the
120 // settings of \a side and \a trans it overwrites \a C with
121
122 \code
123 | side = 'L' | side = 'R'
124 -------------|--------------|--------------
125 trans = 'N': | Q * C | C * Q
126 trans = 'T': | trans(Q) * C | C * trans(Q)
127 \endcode
128
129 // In case \a side is specified as \c 'L', the \a Q matrix is expected to be of size \a m-by-\a m,
130 // in case \a side is set to \c 'R', \a Q is expected to be of size \a n-by-\a n.
131 //
132 // The \a info argument provides feedback on the success of the function call:
133 //
134 // - = 0: The decomposition finished successfully.
135 // - < 0: The i-th argument had an illegal value.
136 //
137 // For more information on the sormlq() function, see the LAPACK online documentation browser:
138 //
139 // http://www.netlib.org/lapack/explore-html/
140 //
141 // \note This function can only be used if a fitting LAPACK library, which supports this function,
142 // is available and linked to the executable. Otherwise a call to this function will result in a
143 // linker error.
144 */
ormlq(char side,char trans,blas_int_t m,blas_int_t n,blas_int_t k,const float * A,blas_int_t lda,const float * tau,float * C,blas_int_t ldc,float * work,blas_int_t lwork,blas_int_t * info)145 inline void ormlq( char side, char trans, blas_int_t m, blas_int_t n, blas_int_t k,
146 const float* A, blas_int_t lda, const float* tau, float* C, blas_int_t ldc,
147 float* work, blas_int_t lwork, blas_int_t* info )
148 {
149 #if defined(INTEL_MKL_VERSION)
150 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
151 #endif
152
153 sormlq_( &side, &trans, &m, &n, &k, const_cast<float*>( A ), &lda,
154 const_cast<float*>( tau ), C, &ldc, work, &lwork, info
155 #if !defined(INTEL_MKL_VERSION)
156 , blaze::fortran_charlen_t(1), blaze::fortran_charlen_t(1)
157 #endif
158 );
159 }
160 //*************************************************************************************************
161
162
163 //*************************************************************************************************
164 /*!\brief LAPACK kernel for the multiplication of the double precision Q from a LQ decomposition
165 // with another matrix.
166 // \ingroup lapack_decomposition
167 //
168 // \param side \c 'L' to apply \f$ Q \f$ or \f$ Q^T \f$ from the left, \c 'R' to apply from the right.
169 // \param trans \c 'N' for \f$ Q \f$, \c 'T' for \f$ Q^T \f$.
170 // \param m The number of rows of the matrix \a C \f$[0..\infty)\f$.
171 // \param n The number of columns of the matrix \a C \f$[0..\infty)\f$.
172 // \param k The number of elementary reflectors, whose product defines the matrix \a Q.
173 // \param A Pointer to the first element of the LQ decomposed double precision column-major matrix.
174 // \param lda The total number of elements between two columns of the matrix \a A \f$[0..\infty)\f$.
175 // \param tau Array for the scalar factors of the elementary reflectors; size >= min( \a m, \a n ).
176 // \param C Pointer to the first element of the double precision column-major matrix multiplicator.
177 // \param ldc The total number of elements between two columns of the matrix \a C \f$[0..\infty)\f$.
178 // \param work Auxiliary array; size >= max( 1, \a lwork ).
179 // \param lwork The dimension of the array \a work.
180 // \param info Return code of the function call.
181 // \return void
182 //
183 // This function multiplies the \a Q matrix resulting from the LQ decomposition of the sgelqf()
184 // function with the given general double precision \a m-by-\a n matrix \a C. Depending on the
185 // settings of \a side and \a trans it overwrites \a C with
186
187 \code
188 | side = 'L' | side = 'R'
189 -------------|--------------|--------------
190 trans = 'N': | Q * C | C * Q
191 trans = 'T': | trans(Q) * C | C * trans(Q)
192 \endcode
193
194 // In case \a side is specified as \c 'L', the \a Q matrix is expected to be of size \a m-by-\a m,
195 // in case \a side is set to \c 'R', \a Q is expected to be of size \a n-by-\a n.
196 //
197 // The \a info argument provides feedback on the success of the function call:
198 //
199 // - = 0: The decomposition finished successfully.
200 // - < 0: The i-th argument had an illegal value.
201 //
202 // For more information on the dormlq() function, see the LAPACK online documentation browser:
203 //
204 // http://www.netlib.org/lapack/explore-html/
205 //
206 // \note This function can only be used if a fitting LAPACK library, which supports this function,
207 // is available and linked to the executable. Otherwise a call to this function will result in a
208 // linker error.
209 */
ormlq(char side,char trans,blas_int_t m,blas_int_t n,blas_int_t k,const double * A,blas_int_t lda,const double * tau,double * C,blas_int_t ldc,double * work,blas_int_t lwork,blas_int_t * info)210 inline void ormlq( char side, char trans, blas_int_t m, blas_int_t n, blas_int_t k,
211 const double* A, blas_int_t lda, const double* tau, double* C, blas_int_t ldc,
212 double* work, blas_int_t lwork, blas_int_t* info )
213 {
214 #if defined(INTEL_MKL_VERSION)
215 BLAZE_STATIC_ASSERT( sizeof( MKL_INT ) == sizeof( blas_int_t ) );
216 #endif
217
218 dormlq_( &side, &trans, &m, &n, &k, const_cast<double*>( A ), &lda,
219 const_cast<double*>( tau ), C, &ldc, work, &lwork, info
220 #if !defined(INTEL_MKL_VERSION)
221 , blaze::fortran_charlen_t(1), blaze::fortran_charlen_t(1)
222 #endif
223 );
224 }
225 //*************************************************************************************************
226
227 } // namespace blaze
228
229 #endif
230