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 //-- dgetrs --------------------------------------------------------------------
 16 void
 17 LAPACK_DECL(dgejsv)(const char       *JOBA,
 18                     const char       *JOBU,
 19                     const char       *JOBV,
 20                     const char       *JOBR,
 21                     const char       *JOBT,
 22                     const char       *JOBP,
 23                     const INTEGER    *M,
 24                     const INTEGER    *N,
 25                     DOUBLE           *A,
 26                     const INTEGER    *LDA,
 27                     DOUBLE           *SVA,
 28                     DOUBLE           *U,
 29                     const INTEGER    *LDU,
 30                     DOUBLE           *V,
 31                     const INTEGER    *LDV,
 32                     DOUBLE           *WORK,
 33                     const INTEGER    *LWORK,
 34                     INTEGER          *IWORK,
 35                     INTEGER          *INFO)
 36 {
 37     DEBUG_FLENS_LAPACK("dgejsv");
 38 //
 39 //  Test the input parameters so that we pass LAPACK error checks
 40 //
 41     const bool lsvec  = (*JOBU=='U') || (*JOBU=='F');
 42     const bool jracc  = (*JOBV=='J');
 43     const bool rsvec  = (*JOBV=='V') || jracc;
 44     const bool rowpiv = (*JOBA=='F') || (*JOBA=='G');
 45     const bool l2rank = (*JOBA=='R');
 46     const bool l2aber = (*JOBA=='A');
 47     const bool errest = (*JOBA=='E') || (*JOBA=='G');
 48     const bool l2tran = (*JOBT=='T');
 49     const bool l2kill = (*JOBR=='R');
 50     const bool defr   = (*JOBR=='N');
 51     const bool l2pert = (*JOBP=='P');
 52 
 53     const INTEGER m = *M;
 54     const INTEGER n = *N;
 55 
 56     *INFO = 0;
 57     if (!(rowpiv || l2rank || l2aber || errest || *JOBA=='C')) {
 58         *INFO = -1;
 59     } else if (!(lsvec || *JOBU=='N' || *JOBU=='W')) {
 60         *INFO = -2;
 61     } else if (!(rsvec || *JOBV=='N' || *JOBV=='W') || (jracc && !lsvec)) {
 62         *INFO = -3;
 63     } else if (!(l2kill || defr)) {
 64         *INFO = -4;
 65     } else if (!(l2tran || *JOBT=='N')) {
 66         *INFO = -5;
 67     } else if (!(l2pert || *JOBP=='N')) {
 68         *INFO = -6;
 69     } else if (m<0) {
 70         *INFO = -7;
 71     } else if ((n<0) || (n>m)) {
 72         *INFO = -8;
 73     } else if (*LDA<m) {
 74         *INFO = -10;
 75     } else if (lsvec && *LDU<m) {
 76         *INFO = -13;
 77     } else if (rsvec && *LDV<n) {
 78         *INFO = -14;
 79     } else if ((!(lsvec || rsvec || errest) || (*LWORK < max(7,4*n+1,2*m+n)))
 80         || (!(lsvec || rsvec) || errest || (*LWORK < max(7,4*n+n*n,2*m+n)))
 81         || (lsvec || (!rsvec) || (*LWORK < max(7,2*m+n,4*n+1)))
 82         || (rsvec || (!lsvec) || (*LWORK < max(7,2*m+n,4*n+1)))
 83         || (lsvec || rsvec || (!jracc) || (*LWORK<max(2*m+n,6*n+2*n*n)))
 84         || (lsvec || rsvec || jracc || *LWORK<max(2*m+n,4*n+n*n,2*n+n*n+6)))
 85     {
 86         *INFO = -17;
 87     }
 88 
 89     if (*INFO!=0) {
 90         *INFO = -(*INFO);
 91         LAPACK_ERROR("DGEJSV", INFO);
 92         *INFO = -(*INFO);
 93         return;
 94     }
 95 
 96 //
 97 //  Call FLENS implementation
 98 //
 99     JSV::Accuracy  accuracy    = JSV::Accuracy(*JOBA);
100     JSV::JobU  jobU            = JSV::JobU(*JOBU);
101     JSV::JobV  jobV            = JSV::JobV(*JOBV);
102     const bool restrictedRange = (*JOBR=='R');
103     const bool considerTransA  = (*JOBT=='T');
104     const bool perturb         = (*JOBR=='P');
105 
106     DGeMatrixView     _A       = DFSView(m, n, A, *LDA);
107     DDenseVectorView  _sva     = DArrayView(n, SVA, INTEGER(1));
108     DGeMatrixView     _U       = DFSView(m, n, U, *LDU);
109     DGeMatrixView     _V       = DFSView(n, n, V, *LDV);
110     DDenseVectorView  _work    = DArrayView(*LWORK, WORK, INTEGER(1));
111     IDenseVectorView  _iwork   = IArrayView(m+3*n, IWORK, INTEGER(1));
112 
113     *INFO = jsv(accuracy, jobU, jobV,
114                 restrictedRange, considerTransA, perturb,
115                 _A, _sva, _U, _V, _work, _iwork);
116 
117     if (*INFO<0) {
118         *INFO = -(*INFO);
119         LAPACK_ERROR("DGEJSV", INFO);
120         *INFO = -(*INFO);
121         return;
122     }
123 }
124 
125 } // extern "C"
126 
127 } } // namespace lapack, flens
Home