Home
Browse Files
  1 //#define CXXBLAS_DEBUG_OUT(x)      std::cerr << x << std::endl;
  2 
  3 #define STR(x)      #x
  4 #define STRING(x)   STR(x)
  5 
  6 #define FLENS_DEFAULT_INDEXTYPE int
  7 
  8 #include <flens/lapack/interface/include/config.h>
  9 
 10 
 11 namespace flens { namespace lapack {
 12 
 13 extern "C" {
 14 
 15 //-- dgesvx --------------------------------------------------------------------
 16 void
 17 LAPACK_DECL(dgesvx)(const char       *FACT,
 18                     const char       *TRANS,
 19                     const INTEGER    *N,
 20                     const INTEGER    *NRHS,
 21                     DOUBLE           *A,
 22                     const INTEGER    *LDA,
 23                     DOUBLE           *AF,
 24                     const INTEGER    *LDAF,
 25                     INTEGER          *IPIV,
 26                     char             *EQUED,
 27                     DOUBLE           *R,
 28                     DOUBLE           *C,
 29                     DOUBLE           *B,
 30                     const INTEGER    *LDB,
 31                     DOUBLE           *X,
 32                     const INTEGER    *LDX,
 33                     DOUBLE           *RCOND,
 34                     DOUBLE           *FERR,
 35                     DOUBLE           *BERR,
 36                     DOUBLE           *WORK,
 37                     INTEGER          *IWORK,
 38                     INTEGER          *INFO)
 39 {
 40     DEBUG_FLENS_LAPACK("dgesvx");
 41 //
 42 //  Test the input parameters so that we pass LAPACK error checks
 43 //
 44     typedef DOUBLE   T;
 45     typedef INTEGER  IndexType;
 46 
 47     const T  Zero(0), One(1);
 48     const T bigNum = One / lamch<T>(SafeMin);
 49 
 50     bool rowEqu, colEqu;
 51 
 52     if (*FACT=='N' || *FACT=='E') {
 53         *EQUED = 'N';
 54         rowEqu = false;
 55         colEqu = false;
 56     } else {
 57         rowEqu = (*EQUED=='R' || *EQUED=='B');
 58         colEqu = (*EQUED=='C' || *EQUED=='B');
 59     }
 60 
 61     *INFO = 0;
 62     if (*FACT!='F' && *FACT!='N' && *FACT!='E') {
 63         *INFO = -1;
 64     } else if (*TRANS!='N' && *TRANS!='T' && *TRANS!='C') {
 65         *INFO = -2;
 66     } else if (*N<0) {
 67         *INFO = -3;
 68     } else if (*NRHS<0) {
 69         *INFO = -4;
 70     } else if (*LDA<std::max(INTEGER(1), *N)) {
 71         *INFO = -6;
 72     } else if (*LDAF<std::max(INTEGER(1), *N)) {
 73         *INFO = -8;
 74     } else if (*FACT=='F' && !(rowEqu || colEqu || *EQUED=='N')) {
 75         *INFO = -10;
 76     } else {
 77         if (rowEqu) {
 78             T rcMin = bigNum;
 79             for (IndexType j=0; j<*N; ++j) {
 80                 rcMin = min(rcMin, R[j]);
 81             }
 82             if (rcMin<=Zero) {
 83                 *INFO = -11;
 84             }
 85         }
 86         if (colEqu && *INFO==0) {
 87             T rcMin = bigNum;
 88             for (IndexType j=0; j<*N; ++j) {
 89                 rcMin = min(rcMin, C[j]);
 90             }
 91             if (rcMin<=Zero) {
 92                 *INFO = -12;
 93             }
 94         }
 95         if (*INFO==0) {
 96             if (*LDB<std::max(INTEGER(1), *N)) {
 97                 *INFO = -14;
 98             } else if (*LDX<std::max(INTEGER(1), *N)) {
 99                 *INFO = -16;
100             }
101         }
102     }
103     if (*INFO!=0) {
104         *INFO = -(*INFO);
105         LAPACK_ERROR("DGESVX", INFO);
106         *INFO = -(*INFO);
107         return;
108     }
109 //
110 //  Call FLENS implementation
111 //
112     SVX::Fact          _FACT  = SVX::Fact(*FACT);
113     Transpose          _TRANS = cxxblas::getCxxBlasEnum<Transpose>(*TRANS);
114     DGeMatrixView      _A     = DFSView(*N, *N, A, *LDA);
115     DGeMatrixView      _AF    = DFSView(*N, *N, AF, *LDAF);
116     IDenseVectorView   _IPIV  = IArrayView(*N, IPIV, 1);
117     SVX::Equilibration _EQUED = SVX::Equilibration(*EQUED);
118     DDenseVectorView   _R     = DArrayView(*N, R, 1);
119     DDenseVectorView   _C     = DArrayView(*N, C, 1);
120     DGeMatrixView      _B     = DFSView(*N, *NRHS, B, *LDB);
121     DGeMatrixView      _X     = DFSView(*N, *NRHS, X, *LDX);
122     DDenseVectorView   _FERR  = DArrayView(*NRHS, FERR, 1);
123     DDenseVectorView   _BERR  = DArrayView(*NRHS, BERR, 1);
124     DDenseVectorView   _WORK  = DArrayView(*N*4, WORK, 1);
125     IDenseVectorView   _IWORK = IArrayView(*N, IWORK, 1);
126 
127     *INFO = svx(_FACT, _TRANS, _A, _AF, _IPIV, _EQUED, _R, _C, _B, _X, *RCOND,
128                 _FERR, _BERR, _WORK, _IWORK);
129 
130     *EQUED = char(_EQUED);
131 }
132 
133 } // extern "C"
134 
135 } } // namespace lapack, flens
Home