1 // Copyright (C) 2010 Davis E. King (davis@dlib.net) 2 // License: Boost Software License See LICENSE.txt for the full license. 3 #ifndef DLIB_LAPACk_GETRF_Hh_ 4 #define DLIB_LAPACk_GETRF_Hh_ 5 6 #include "fortran_id.h" 7 #include "../matrix.h" 8 9 namespace dlib 10 { 11 namespace lapack 12 { 13 namespace binding 14 { 15 extern "C" 16 { 17 void DLIB_FORTRAN_ID(dgetrf) (const integer *m, const integer *n, double *a, 18 const integer *lda, integer *ipiv, integer *info); 19 20 void DLIB_FORTRAN_ID(sgetrf) (const integer *m, const integer *n, float *a, 21 const integer *lda, integer *ipiv, integer *info); 22 23 } 24 getrf(integer m,integer n,double * a,integer lda,integer * ipiv)25 inline int getrf (integer m, integer n, double *a, 26 integer lda, integer *ipiv) 27 { 28 integer info = 0; 29 DLIB_FORTRAN_ID(dgetrf)(&m, &n, a, &lda, ipiv, &info); 30 return info; 31 } 32 getrf(integer m,integer n,float * a,integer lda,integer * ipiv)33 inline int getrf (integer m, integer n, float *a, 34 integer lda, integer *ipiv) 35 { 36 integer info = 0; 37 DLIB_FORTRAN_ID(sgetrf)(&m, &n, a, &lda, ipiv, &info); 38 return info; 39 } 40 41 42 } 43 44 // ------------------------------------------------------------------------------------ 45 46 47 /* -- LAPACK routine (version 3.1) -- */ 48 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ 49 /* November 2006 */ 50 51 /* .. Scalar Arguments .. */ 52 /* .. */ 53 /* .. Array Arguments .. */ 54 /* .. */ 55 56 /* Purpose */ 57 /* ======= */ 58 59 /* DGETRF computes an LU factorization of a general M-by-N matrix A */ 60 /* using partial pivoting with row interchanges. */ 61 62 /* The factorization has the form */ 63 /* A = P * L * U */ 64 /* where P is a permutation matrix, L is lower triangular with unit */ 65 /* diagonal elements (lower trapezoidal if m > n), and U is upper */ 66 /* triangular (upper trapezoidal if m < n). */ 67 68 /* This is the right-looking Level 3 BLAS version of the algorithm. */ 69 70 /* Arguments */ 71 /* ========= */ 72 73 /* M (input) INTEGER */ 74 /* The number of rows of the matrix A. M >= 0. */ 75 76 /* N (input) INTEGER */ 77 /* The number of columns of the matrix A. N >= 0. */ 78 79 /* A (input/output) DOUBLE PRECISION array, dimension (LDA,N) */ 80 /* On entry, the M-by-N matrix to be factored. */ 81 /* On exit, the factors L and U from the factorization */ 82 /* A = P*L*U; the unit diagonal elements of L are not stored. */ 83 84 /* LDA (input) INTEGER */ 85 /* The leading dimension of the array A. LDA >= max(1,M). */ 86 87 /* IPIV (output) INTEGER array, dimension (min(M,N)) */ 88 /* The pivot indices; for 1 <= i <= min(M,N), row i of the */ 89 /* matrix was interchanged with row IPIV(i). */ 90 91 /* INFO (output) INTEGER */ 92 /* = 0: successful exit */ 93 /* < 0: if INFO = -i, the i-th argument had an illegal value */ 94 /* > 0: if INFO = i, U(i,i) is exactly zero. The factorization */ 95 /* has been completed, but the factor U is exactly */ 96 /* singular, and division by zero will occur if it is used */ 97 /* to solve a system of equations. */ 98 99 100 // ------------------------------------------------------------------------------------ 101 102 template < 103 typename T, 104 long NR1, long NR2, 105 long NC1, long NC2, 106 typename MM, 107 typename layout 108 > getrf(matrix<T,NR1,NC1,MM,column_major_layout> & a,matrix<integer,NR2,NC2,MM,layout> & ipiv)109 int getrf ( 110 matrix<T,NR1,NC1,MM,column_major_layout>& a, 111 matrix<integer,NR2,NC2,MM,layout>& ipiv 112 ) 113 { 114 const long m = a.nr(); 115 const long n = a.nc(); 116 117 ipiv.set_size(std::min(m,n), 1); 118 119 // compute the actual decomposition 120 return binding::getrf(m, n, &a(0,0), a.nr(), &ipiv(0,0)); 121 } 122 123 // ------------------------------------------------------------------------------------ 124 125 } 126 127 } 128 129 // ---------------------------------------------------------------------------------------- 130 131 #endif // DLIB_LAPACk_GETRF_Hh_ 132 133