1! Signatures for f2py-wrappers of FORTRAN LAPACK General Matrix functions. 2! 3 4subroutine <prefix>gebal(scale,permute,n,a,m,lo,hi,pivscale,info) 5 ! 6 ! ba,lo,hi,pivscale,info = gebal(a,scale=0,permute=0,overwrite_a=0) 7 ! Balance general matrix a. 8 ! hi,lo are such that ba[i][j]==0 if i>j and j=0...lo-1 or i=hi+1..n-1 9 ! pivscale([0:lo], [lo:hi+1], [hi:n+1]) = (p1,d,p2) where (p1,p2)[j] is 10 ! the index of the row and column interchanged with row and column j. 11 ! d[j] is the scaling factor applied to row and column j. 12 ! The order in which the interchanges are made is n-1 to hi+1, then 0 to lo-1. 13 ! 14 ! P * A * P = [[T1,X,Y],[0,B,Z],[0,0,T2]] 15 ! BA = [[T1,X*D,Y],[0,inv(D)*B*D,ind(D)*Z],[0,0,T2]] 16 ! where D = diag(d), T1,T2 are upper triangular matrices. 17 ! lo,hi mark the starting and ending columns of submatrix B. 18 19 callstatement { (*f2py_func)((permute?(scale?"B":"P"):(scale?"S":"N")),&n,a,&m,&lo,&hi,pivscale,&info); hi--; lo--; } 20 callprotoargument char*,F_INT*,<ctype>*,F_INT*,F_INT*,F_INT*,<ctypereal>*,F_INT* 21 integer intent(in),optional :: permute = 0 22 integer intent(in),optional :: scale = 0 23 integer intent(hide),depend(a,n) :: m = shape(a,0) 24 integer intent(hide),depend(a) :: n = shape(a,1) 25 check(m>=n) m 26 integer intent(out) :: hi,lo 27 <ftypereal> dimension(n),intent(out),depend(n) :: pivscale 28 <ftype> dimension(m,n),intent(in,out,copy,out=ba) :: a 29 integer intent(out) :: info 30 31end subroutine <prefix>gebal 32 33subroutine <prefix>gehrd(n,lo,hi,a,tau,work,lwork,info) 34 ! 35 ! hq,tau,info = gehrd(a,lo=0,hi=n-1,lwork=n,overwrite_a=0) 36 ! Reduce general matrix A to upper Hessenberg form H by unitary similarity 37 ! transform Q^H * A * Q = H 38 ! 39 ! Q = H(lo) * H(lo+1) * ... * H(hi-1) 40 ! H(i) = I - tau * v * v^H 41 ! v[0:i+1] = 0, v[i+1]=1, v[hi+1:n] = 0 42 ! v[i+2:hi+1] is stored in hq[i+2:hi+i,i] 43 ! tau is tau[i] 44 ! 45 ! hq for n=7,lo=1,hi=5: 46 ! [a a h h h h a 47 ! a h h h h a 48 ! h h h h h h 49 ! v2h h h h h 50 ! v2v3h h h h 51 ! v2v3v4h h h 52 ! a] 53 ! 54 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,a,&n,tau,work,&lwork,&info); } 55 callprotoargument F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 56 integer intent(hide),depend(a) :: n = shape(a,0) 57 <ftype> dimension(n,n),intent(in,out,copy,out=ht,aligned8),check(shape(a,0)==shape(a,1)) :: a 58 integer intent(in),optional :: lo = 0 59 integer intent(in),optional,depend(n) :: hi = n-1 60 <ftype> dimension(n-1),intent(out),depend(n) :: tau 61 <ftype> dimension(lwork),intent(cache,hide),depend(lwork) :: work 62 integer intent(in),optional,depend(n),check(lwork>=MAX(n,1)) :: lwork = MAX(n,1) 63 integer intent(out) :: info 64 65end subroutine <prefix>gehrd 66 67subroutine <prefix>gehrd_lwork(n,lo,hi,a,tau,work,lwork,info) 68 ! LWORK computation for GEHRD 69 fortranname <prefix>gehrd 70 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,&a,&n,&tau,&work,&lwork,&info); } 71 callprotoargument F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 72 integer intent(in) :: n 73 <ftype> intent(hide) :: a 74 integer intent(in),optional :: lo = 0 75 integer intent(in),optional,depend(n) :: hi = n-1 76 <ftype> intent(hide) :: tau 77 <ftype> intent(out) :: work 78 integer intent(hide) :: lwork = -1 79 integer intent(out) :: info 80end subroutine <prefix>gehrd_lwork 81 82subroutine <prefix>gesv(n,nrhs,a,piv,b,info) 83 ! lu,piv,x,info = gesv(a,b,overwrite_a=0,overwrite_b=0) 84 ! Solve A * X = B. 85 ! A = P * L * U 86 ! U is upper diagonal triangular, L is unit lower triangular, 87 ! piv pivots columns. 88 89 callstatement {F_INT i;(*f2py_func)(&n,&nrhs,a,&n,piv,b,&n,&info);for(i=0;i\<n;--piv[i++]);} 90 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* 91 92 integer depend(a),intent(hide):: n = shape(a,0) 93 integer depend(b),intent(hide):: nrhs = shape(b,1) 94 <ftype> dimension(n,n),check(shape(a,0)==shape(a,1)) :: a 95 integer dimension(n),depend(n),intent(out) :: piv 96 <ftype> dimension(n,nrhs),check(shape(a,0)==shape(b,0)),depend(n) :: b 97 integer intent(out)::info 98 intent(in,out,copy,out=x) b 99 intent(in,out,copy,out=lu) a 100end subroutine <prefix>gesv 101 102subroutine <prefix2>gesvx(fact,trans,n,nrhs,a,lda,af,ldaf,ipiv,equed,r,c,b,ldb,x,ldx,rcond,ferr,berr,work,iwork,info) 103 ! Solve A * X = B using LU decomposition 104 ! The expert driver of ?GESV with condition number, backward/forward error estimates, and iterative refinement 105 ! This part takes care of the data types, single and double reals (sgesvx and dgesvx) 106 threadsafe 107 callstatement {F_INT i;(*f2py_func)(fact,trans,&n,&nrhs,a,&lda,af,&ldaf,ipiv,equed,r,c,b,&ldb,x,&ldx,&rcond,ferr,berr,work,iwork,&info);for(i=0;i\<n;--ipiv[i++]);} 108 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,char*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 109 110 character optional,intent(in):: trans = "N" 111 character optional,intent(in):: fact = "E" 112 integer depend(a),intent(hide):: n = shape(a,0) 113 integer depend(b),intent(hide):: nrhs = shape(b,1) 114 <ftype2> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=as):: a 115 integer depend(a),intent(hide):: lda = shape(a,0) 116 <ftype2> optional,dimension(n,n),intent(in,out,out=lu):: af 117 integer optional,depend(af),intent(hide):: ldaf = shape(af,0) 118 integer optional,dimension(n),depend(n),intent(in,out):: ipiv 119 character optional,intent(in,out):: equed = "B" 120 <ftype2> optional,dimension(n),depend(n),intent(in,out,out=rs):: r 121 <ftype2> optional,dimension(n),depend(n),intent(in,out,out=cs):: c 122 <ftype2> depend(n),dimension(n,nrhs),intent(in,out,copy,out=bs):: b 123 integer depend(b),intent(hide):: ldb = shape(b,0) 124 <ftype2> dimension(n,nrhs),depend(n,nrhs),intent(out):: x 125 integer depend(n),intent(hide):: ldx = n 126 <ftype2> intent(out):: rcond 127 <ftype2> intent(out),dimension(nrhs),depend(nrhs):: ferr 128 <ftype2> intent(out),dimension(nrhs),depend(nrhs):: berr 129 <ftype2> dimension(4*n),depend(n),intent(hide,cache):: work 130 integer intent(hide,cache),dimension(n),depend(n) :: iwork 131 integer intent(out):: info 132 133end subroutine <prefix2>gesvx 134 135subroutine <prefix2c>gesvx(fact,trans,n,nrhs,a,lda,af,ldaf,ipiv,equed,r,c,b,ldb,x,ldx,rcond,ferr,berr,work,rwork,info) 136 ! Solve A * X = B using LU decomposition 137 ! The expert driver of ?GESV with condition number, backward/forward error estimates, and iterative refinement 138 ! This part takes care of the data types, complex and double complex (cgesvx and zgesvx) 139 threadsafe 140 callstatement {F_INT i;(*f2py_func)(fact,trans,&n,&nrhs,a,&lda,af,&ldaf,ipiv,equed,r,c,b,&ldb,x,&ldx,&rcond,ferr,berr,work,rwork,&info);for(i=0;i\<n;--ipiv[i++]);} 141 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,char*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,<ctype2>*,F_INT* 142 143 character optional,intent(in):: trans = "N" 144 character optional,intent(in):: fact = "E" 145 integer depend(a),intent(hide):: n = shape(a,0) 146 integer depend(b),intent(hide):: nrhs = shape(b,1) 147 <ftype2c> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=as):: a 148 integer depend(a),intent(hide):: lda = shape(a,0) 149 <ftype2c> optional,dimension(n,n),depend(n),intent(in,out,out=lu):: af 150 integer optional,depend(af),intent(hide):: ldaf = shape(af,0) 151 integer optional,dimension(n),depend(n),intent(in,out):: ipiv 152 character optional,intent(in,out):: equed = "B" 153 <ftype2> optional,dimension(n),depend(n),intent(in,out,out=rs):: r 154 <ftype2> optional,dimension(n),depend(n),intent(in,out,out=cs):: c 155 <ftype2c> depend(n),dimension(n,nrhs),intent(in,out,copy,out=bs):: b 156 integer depend(b),intent(hide):: ldb = shape(b,0) 157 <ftype2c> dimension(n,nrhs),depend(n,nrhs),intent(out):: x 158 integer depend(n),intent(hide):: ldx = n 159 <ftype2> intent(out):: rcond 160 <ftype2> intent(out),dimension(nrhs),depend(nrhs):: ferr 161 <ftype2> intent(out),dimension(nrhs),depend(nrhs):: berr 162 <ftype2c> dimension(2*n),depend(n),intent(hide,cache):: work 163 <ftype2> intent(hide,cache),dimension(2*n),depend(n) :: rwork 164 integer intent(out):: info 165 166end subroutine <prefix2c>gesvx 167 168subroutine <prefix>gecon(norm,n,a,lda,anorm,rcond,work,irwork,info) 169 ! Computes the 1- or inf- norm reciprocal condition number estimate. 170 threadsafe 171 callstatement (*f2py_func)(norm,&n,a,&lda,&anorm,&rcond,work,irwork,&info) 172 callprotoargument char*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctype>*,<F_INT,F_INT,float,double>*,F_INT* 173 174 character optional,intent(in):: norm = '1' 175 integer depend(a),intent(hide):: n = shape(a,0) 176 <ftype> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in):: a 177 integer depend(a),intent(hide):: lda = shape(a,0) 178 <ftypereal> intent(in):: anorm 179 <ftypereal> intent(out):: rcond 180 <ftype> depend(n),dimension(<4*n,4*n,2*n,2*n>),intent(hide,cache):: work 181 <integer,integer,real, double precision> depend(n),dimension(<n,n,2*n,2*n>),intent(hide,cache):: irwork 182 integer intent(out):: info 183 184end subroutine <prefix>gecon 185 186subroutine <prefix>getrf(m,n,a,piv,info) 187 ! lu,piv,info = getrf(a,overwrite_a=0) 188 ! Compute an LU factorization of a general M-by-N matrix A. 189 ! A = P * L * U 190 threadsafe 191 callstatement {F_INT i;(*f2py_func)(&m,&n,a,&m,piv,&info);for(i=0,n=MIN(m,n);i\<n;--piv[i++]);} 192 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,F_INT*,F_INT* 193 194 integer depend(a),intent(hide):: m = shape(a,0) 195 integer depend(a),intent(hide):: n = shape(a,1) 196 <ftype> dimension(m,n),intent(in,out,copy,out=lu) :: a 197 integer dimension(MIN(m,n)),depend(m,n),intent(out) :: piv 198 integer intent(out):: info 199 200end subroutine <prefix>getrf 201 202subroutine <prefix>getrs(n,nrhs,lu,piv,b,info,trans) 203 ! x,info = getrs(lu,piv,b,trans=0,overwrite_b=0) 204 ! Solve A * X = B if trans=0 205 ! Solve A^T * X = B if trans=1 206 ! Solve A^H * X = B if trans=2 207 ! A = P * L * U 208 threadsafe 209 callstatement {F_INT i;for(i=0;i\<n;++piv[i++]);(*f2py_func)((trans?(trans==2?"C":"T"):"N"),&n,&nrhs,lu,&n,piv,b,&n,&info);for(i=0;i\<n;--piv[i++]);} 210 callprotoargument char*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* 211 212 integer optional,intent(in),check(trans>=0 && trans <=2) :: trans = 0 213 214 integer depend(lu),intent(hide):: n = shape(lu,0) 215 integer depend(b),intent(hide):: nrhs = shape(b,1) 216 <ftype> dimension(n,n),intent(in) :: lu 217 check(shape(lu,0)==shape(lu,1)) :: lu 218 integer dimension(n),intent(in),depend(n) :: piv 219 <ftype> dimension(n,nrhs),intent(in,out,copy,out=x),depend(n),check(shape(lu,0)==shape(b,0)) :: b 220 integer intent(out):: info 221 222end subroutine <prefix>getrs 223 224subroutine <prefix>getc2(n,a,lda,ipiv,jpiv,info) 225 ! lu,ipiv,jpiv,info = getc2(a,overwrite_a=0) 226 ! Compute an LU factorization of with complete pivoting of a general n-by-n matrix. 227 ! A = P * L * U * Q 228 threadsafe 229 callstatement {F_INT i;(*f2py_func)(&n,a,&lda,ipiv,jpiv,&info);for(i=0;i\<n;--ipiv[i],--jpiv[i++]);} 230 callprotoargument F_INT*,<ctype>*,F_INT*,F_INT*,F_INT*,F_INT* 231 232 integer depend(a),intent(hide):: n = shape(a,0) 233 <ftype> dimension(n,n),intent(in,out,copy,out=lu) :: a 234 check(shape(a,0)==shape(a,1)) :: a 235 integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) 236 integer dimension(n),depend(n),intent(out) :: ipiv 237 integer dimension(n),depend(n),intent(out) :: jpiv 238 integer intent(out):: info 239 240end subroutine <prefix>getc2 241 242subroutine <prefix>gesc2(n,lu,lda,rhs,ipiv,jpiv,scale) 243 ! x,scale = gesc2(lu,rhs,ipiv,jpiv,overwrite_rhs=0) 244 ! Solve A * X = scale * RHS 245 ! A = P * L * U * Q 246 threadsafe 247 callstatement {F_INT i;for(i=0;i\<n;++ipiv[i],++jpiv[i++]);(*f2py_func)(&n,lu,&lda,rhs,ipiv,jpiv,&scale);for(i=0;i\<n;--ipiv[i],--jpiv[i++]);} 248 callprotoargument F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctypereal>* 249 250 integer depend(lu),intent(hide):: n = shape(lu,0) 251 <ftype> dimension(n,n),intent(in) :: lu 252 check(shape(lu,0)==shape(lu,1)) :: lu 253 integer intent(hide),depend(lu):: lda = MAX(1,shape(lu,0)) 254 <ftype> dimension(n),intent(in,out,copy,out=x),depend(n),check(shape(lu,0)==len(rhs)) :: rhs 255 integer dimension(n),intent(in),depend(n) :: ipiv 256 integer dimension(n),intent(in),depend(n) :: jpiv 257 <ftypereal> intent(out):: scale 258 259end subroutine <prefix>gesc2 260 261subroutine <prefix>getri(n,lu,piv,work,lwork,info) 262 ! inv_a,info = getri(lu,piv,lwork=3*n,overwrite_lu=0) 263 ! Find A inverse A^-1. 264 ! A = P * L * U 265 266 callstatement {F_INT i;for(i=0;i\<n;++piv[i++]);(*f2py_func)(&n,lu,&n,piv,work,&lwork,&info);for(i=0;i\<n;--piv[i++]);} 267 callprotoargument F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* 268 269 integer depend(lu),intent(hide):: n = shape(lu,0) 270 <ftype> dimension(n,n),intent(in,out,copy,out=inv_a) :: lu 271 check(shape(lu,0)==shape(lu,1)) :: lu 272 integer dimension(n),intent(in),depend(n) :: piv 273 integer intent(out):: info 274 integer optional,intent(in),depend(n),check(lwork>=n) :: lwork=max(3*n,1) 275 <ftype> dimension(lwork),intent(hide,cache),depend(lwork) :: work 276 277end subroutine <prefix>getri 278 279subroutine <prefix>getri_lwork(n,lu,piv,work,lwork,info) 280 ! *GETRI LWORK query 281 fortranname <prefix>getri 282 callstatement (*f2py_func)(&n,&lu,&n,&piv,&work,&lwork,&info) 283 callprotoargument F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* 284 285 integer intent(in):: n 286 <ftype> intent(hide) :: lu 287 integer intent(hide) :: piv 288 integer intent(out):: info 289 integer intent(hide) :: lwork=-1 290 <ftype> intent(out) :: work 291end subroutine <prefix>getri_lwork 292 293subroutine <prefix2>gesdd(m,n,minmn,u0,u1,vt0,vt1,a,compute_uv,full_matrices,u,s,vt,work,lwork,iwork,info) 294 ! u,s,vt,info = gesdd(a,compute_uv=1,lwork=..,overwrite_a=0) 295 ! Compute the singular value decomposition (SVD) using divide and conquer: 296 ! A = U * SIGMA * transpose(V) 297 ! A - M x N matrix 298 ! U - M x M matrix or min(M,N) x N if full_matrices=False 299 ! SIGMA - M x N zero matrix with a main diagonal filled with min(M,N) 300 ! singular values 301 ! transpose(V) - N x N matrix or N x min(M,N) if full_matrices=False 302 303 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,a,&m,s,u,&u0,vt,&vt0,work,&lwork,iwork,&info) 304 callprotoargument char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 305 306 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 307 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 308 integer intent(hide),depend(a):: m = shape(a,0) 309 integer intent(hide),depend(a):: n = shape(a,1) 310 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 311 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 312 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: u1 = (compute_uv?(full_matrices?m:minmn):1) 313 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 314 integer intent(hide),depend(compute_uv,minmn) :: vt1 = (compute_uv?n:1) 315 <ftype2> dimension(m,n),intent(in,copy,aligned8) :: a 316 <ftype2> dimension(minmn),intent(out),depend(minmn) :: s 317 <ftype2> dimension(u0,u1),intent(out),depend(u0, u1) :: u 318 <ftype2> dimension(vt0,vt1),intent(out),depend(vt0, vt1) :: vt 319 <ftype2> dimension(lwork),intent(hide,cache),depend(lwork) :: work 320 integer optional,intent(in),depend(minmn,compute_uv) & 321 :: lwork = max((compute_uv?4*minmn*minmn+MAX(m,n)+9*minmn:MAX(14*minmn+4,10*minmn+2+25*(25+8))+MAX(m,n)),1) 322 integer intent(hide,cache),dimension(8*minmn),depend(minmn) :: iwork 323 integer intent(out)::info 324 325end subroutine <prefix2>gesdd 326 327subroutine <prefix2>gesdd_lwork(m,n,minmn,u0,vt0,a,compute_uv,full_matrices,u,s,vt,work,lwork,iwork,info) 328 ! LWORK computation for (S/D)GESDD 329 330 fortranname <prefix2>gesdd 331 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,&a,&m,&s,&u,&u0,&vt,&vt0,&work,&lwork,&iwork,&info) 332 callprotoargument char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 333 334 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 335 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 336 integer intent(in) :: m 337 integer intent(in) :: n 338 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 339 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 340 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 341 <ftype2> intent(hide) :: a 342 <ftype2> intent(hide) :: s 343 <ftype2> intent(hide) :: u 344 <ftype2> intent(hide) :: vt 345 <ftype2> intent(out) :: work 346 integer intent(hide) :: lwork = -1 347 integer intent(hide) :: iwork 348 integer intent(out) :: info 349 350end subroutine <prefix2>gesdd_lwork 351 352subroutine <prefix2c>gesdd(m,n,minmn,u0,u1,vt0,vt1,a,compute_uv,full_matrices,u,s,vt,work,rwork,lwork,iwork,info) 353 ! u,s,vt,info = gesdd(a,compute_uv=1,lwork=..,overwrite_a=0) 354 ! Compute the singular value decomposition (SVD) using divide and conquer: 355 ! A = U * SIGMA * conjugate-transpose(V) 356 ! A - M x N matrix 357 ! U - M x M matrix or min(M,N) x N if full_matrices=False 358 ! SIGMA - M x N zero matrix with a main diagonal filled with min(M,N) 359 ! singular values 360 ! transpose(V) - N x N matrix or N x min(M,N) if full_matrices=False 361 362 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,a,&m,s,u,&u0,vt,&vt0,work,&lwork,rwork,iwork,&info) 363 callprotoargument char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT* 364 365 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 366 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 367 integer intent(hide),depend(a):: m = shape(a,0) 368 integer intent(hide),depend(a):: n = shape(a,1) 369 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 370 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 371 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: u1 = (compute_uv?(full_matrices?m:minmn):1) 372 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 373 integer intent(hide),depend(compute_uv,minmn) :: vt1 = (compute_uv?n:1) 374 <ftype2c> dimension(m,n),intent(in,copy) :: a 375 <ftype2> dimension(minmn),intent(out),depend(minmn) :: s 376 <ftype2c> dimension(u0,u1),intent(out),depend(u0,u1) :: u 377 <ftype2c> dimension(vt0,vt1),intent(out),depend(vt0,vt1) :: vt 378 <ftype2c> dimension(lwork),intent(hide,cache),depend(lwork) :: work 379 <ftype2> dimension((compute_uv?minmn*MAX(5*minmn+7, 2*MAX(m,n)+2*minmn+1):7*minmn)),intent(hide,cache),depend(minmn,compute_uv) :: rwork 380 integer optional,intent(in),depend(minmn,compute_uv) & 381 :: lwork = max((compute_uv?2*minmn*minmn+MAX(m,n)+2*minmn:2*minmn+MAX(m,n)),1) 382 integer intent(hide,cache),dimension(8*minmn),depend(minmn) :: iwork 383 integer intent(out)::info 384 385end subroutine <prefix2c>gesdd 386 387subroutine <prefix2c>gesdd_lwork(m,n,minmn,u0,vt0,a,compute_uv,full_matrices,u,s,vt,work,rwork,lwork,iwork,info) 388 ! (C/Z)GESDD call with LWORK=-1 -- copypaste of above gesdd with dummy arrays 389 390 fortranname <prefix2c>gesdd 391 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,&a,&m,&s,&u,&u0,&vt,&vt0,&work,&lwork,&rwork,&iwork,&info) 392 callprotoargument char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT* 393 394 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 395 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 396 integer intent(in),depend(a):: m 397 integer intent(in),depend(a):: n 398 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 399 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 400 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 401 <ftype2c> intent(hide) :: a 402 <ftype2> intent(hide) :: s 403 <ftype2c> intent(hide) :: u 404 <ftype2c> intent(hide) :: vt 405 <ftype2c> intent(out) :: work 406 <ftype2> intent(hide) :: rwork 407 integer intent(hide) :: lwork = -1 408 integer intent(hide) :: iwork 409 integer intent(out) :: info 410 411end subroutine <prefix2c>gesdd_lwork 412 413subroutine <prefix2>gesvd(m,n,minmn,u0,u1,vt0,vt1,a,compute_uv,full_matrices,u,s,vt,work,lwork,info) 414 ! u,s,vt,info = gesvd(a,compute_uv=1,lwork=..,overwrite_a=0) 415 ! Compute the singular value decomposition (SVD): 416 ! A = U * SIGMA * transpose(V) 417 ! A - M x N matrix 418 ! U - M x M matrix or min(M,N) x N if full_matrices=False 419 ! SIGMA - M x N zero matrix with a main diagonal filled with min(M,N) 420 ! singular values 421 ! transpose(V) - N x N matrix or N x min(M,N) if full_matrices=False 422 423 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),(compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,a,&m,s,u,&u0,vt,&vt0,work,&lwork,&info) 424 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 425 426 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 427 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 428 integer intent(hide),depend(a):: m = shape(a,0) 429 integer intent(hide),depend(a):: n = shape(a,1) 430 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 431 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 432 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: u1 = (compute_uv?(full_matrices?m:minmn):1) 433 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 434 integer intent(hide),depend(compute_uv,minmn) :: vt1 = (compute_uv?n:1) 435 <ftype2> dimension(m,n),intent(in,copy,aligned8) :: a 436 <ftype2> dimension(minmn),intent(out),depend(minmn) :: s 437 <ftype2> dimension(u0,u1),intent(out),depend(u0, u1) :: u 438 <ftype2> dimension(vt0,vt1),intent(out),depend(vt0, vt1) :: vt 439 <ftype2> dimension(lwork),intent(hide,cache),depend(lwork) :: work 440 integer optional,intent(in),depend(minmn) :: lwork = max(MAX(3*minmn+MAX(m,n),5*minmn),1) 441 integer intent(out) :: info 442 443end subroutine <prefix2>gesvd 444 445subroutine <prefix2>gesvd_lwork(m,n,minmn,u0,vt0,a,compute_uv,full_matrices,u,s,vt,work,lwork,info) 446 ! LWORK computation for (S/D)GESVD 447 448 fortranname <prefix2>gesvd 449 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),(compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,&a,&m,&s,&u,&u0,&vt,&vt0,&work,&lwork,&info) 450 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 451 452 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 453 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 454 integer intent(in) :: m 455 integer intent(in) :: n 456 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 457 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 458 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 459 <ftype2> intent(hide) :: a 460 <ftype2> intent(hide) :: s 461 <ftype2> intent(hide) :: u 462 <ftype2> intent(hide) :: vt 463 integer intent(hide) :: lwork = -1 464 <ftype2> intent(out) :: work 465 integer intent(out) :: info 466 467end subroutine <prefix2>gesvd_lwork 468 469subroutine <prefix2c>gesvd(m,n,minmn,u0,u1,vt0,vt1,a,compute_uv,full_matrices,u,s,vt,work,rwork,lwork,info) 470 ! u,s,vt,info = gesvd(a,compute_uv=1,lwork=..,overwrite_a=0) 471 ! Compute the singular value decomposition (SVD): 472 ! A = U * SIGMA * conjugate-transpose(V) 473 ! A - M x N matrix 474 ! U - M x M matrix or min(M,N) x N if full_matrices=False 475 ! SIGMA - M x N zero matrix with a main diagonal filled with min(M,N) 476 ! singular values 477 ! transpose(V) - N x N matrix or N x min(M,N) if full_matrices=False 478 479 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),(compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,a,&m,s,u,&u0,vt,&vt0,work,&lwork,rwork,&info) 480 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 481 482 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 483 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 484 integer intent(hide),depend(a):: m = shape(a,0) 485 integer intent(hide),depend(a):: n = shape(a,1) 486 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 487 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 488 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: u1 = (compute_uv?(full_matrices?m:minmn):1) 489 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 490 integer intent(hide),depend(compute_uv,minmn) :: vt1 = (compute_uv?n:1) 491 <ftype2c> dimension(m,n),intent(in,copy) :: a 492 <ftype2> dimension(minmn),intent(out),depend(minmn) :: s 493 <ftype2c> dimension(u0,u1),intent(out),depend(u0,u1) :: u 494 <ftype2c> dimension(vt0,vt1),intent(out),depend(vt0,vt1) :: vt 495 <ftype2c> dimension(lwork),intent(hide,cache),depend(lwork) :: work 496 <ftype2> dimension((MAX(1,5*minmn))),intent(hide,cache),depend(minmn) :: rwork 497 integer optional,intent(in),depend(minmn) :: lwork = MAX(2*minmn+MAX(m,n),1) 498 integer intent(out) :: info 499 500end subroutine <prefix2c>gesvd 501 502subroutine <prefix2c>gesvd_lwork(m,n,minmn,u0,vt0,a,compute_uv,full_matrices,u,s,vt,work,rwork,lwork,info) 503 ! (C/Z)GESVD call with LWORK=-1 -- copypaste of above gesvd with dummy arrays 504 505 fortranname <prefix2c>gesvd 506 callstatement (*f2py_func)((compute_uv?(full_matrices?"A":"S"):"N"),(compute_uv?(full_matrices?"A":"S"):"N"),&m,&n,&a,&m,&s,&u,&u0,&vt,&vt0,&work,&lwork,&rwork,&info) 507 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 508 509 integer intent(in),optional,check(compute_uv==0||compute_uv==1):: compute_uv = 1 510 integer intent(in),optional,check(full_matrices==0||full_matrices==1):: full_matrices = 1 511 integer intent(in),depend(a):: m 512 integer intent(in),depend(a):: n 513 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 514 integer intent(hide),depend(compute_uv,minmn) :: u0 = (compute_uv?m:1) 515 integer intent(hide),depend(compute_uv,minmn, full_matrices) :: vt0 = (compute_uv?(full_matrices?n:minmn):1) 516 integer intent(hide) :: lwork = -1 517 <ftype2c> intent(hide) :: a 518 <ftype2> intent(hide) :: s 519 <ftype2c> intent(hide) :: u 520 <ftype2c> intent(hide) :: vt 521 <ftype2c> intent(out) :: work 522 <ftype2> intent(hide) :: rwork 523 integer intent(out) :: info 524 525end subroutine <prefix2c>gesvd_lwork 526 527subroutine <prefix>gels(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info) 528 ! lqr,x,info = gels(a,b,lwork=..,overwrite_a=False,overwrite_b=False) 529 ! Solve Minimize 2-norm(A * X - B). 530 531 callstatement (*f2py_func)(trans,&m,&n,&nrhs,a,&lda,b,&ldb,work,&lwork,&info) 532 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 533 534 character optional,intent(in),check(*trans=='N'||*trans==<'T',\0,'C',\2>):: trans = 'N' 535 integer intent(hide),depend(a):: m = shape(a,0) 536 integer intent(hide),depend(a):: n = shape(a,1) 537 integer intent(hide),depend(b):: nrhs = shape(b,1) 538 <ftype> dimension(m,n),intent(in,out,copy,out=lqr):: a 539 integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) 540 <ftype> intent(in,out,copy,out=x),depend(trans,m,n),dimension(MAX(m,n),nrhs),check(shape(b,0)==MAX(m,n)) :: b 541 integer depend(b),intent(hide):: ldb = MAX(1,MAX(m,n)) 542 integer optional,intent(in),depend(nrhs,m,n),check(lwork>=1||lwork==-1)::lwork=MAX(MIN(m,n)+MAX(MIN(m,n),nrhs),1) 543 <ftype> depend(lwork),dimension(MAX(1,lwork)),intent(hide,cache):: work 544 integer intent(out)::info 545 546end subroutine <prefix>gels 547 548subroutine <prefix>gels_lwork(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info) 549 ! ?GELS LWORK Query for optimal block size 550 551 fortranname <prefix>gels 552 callstatement (*f2py_func)(trans,&m,&n,&nrhs,&a,&lda,&b,&ldb,&work,&lwork,&info) 553 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 554 555 character optional,intent(in),check(*trans=='N'||*trans==<'T',\0,'C',\2>):: trans = 'N' 556 integer intent(in),check(m>=0) :: m 557 integer intent(in),check(n>=0) :: n 558 integer intent(in),check(nrhs>=0) :: nrhs 559 560 <ftype> intent(hide):: a 561 integer intent(hide):: lda = MAX(1,m) 562 <ftype> intent(hide):: b 563 integer intent(hide):: ldb = MAX(1,MAX(m,n)) 564 integer intent(hide):: lwork=-1 565 566 <ftype> intent(out):: work 567 integer intent(out)::info 568 569end subroutine <prefix>gels_lwork 570 571subroutine <prefix2>gelss(m,n,minmn,maxmn,nrhs,a,b,s,cond,r,work,lwork,info) 572 ! v,x,s,rank,work,info = gelss(a,b,cond=-1.0,overwrite_a=0,overwrite_b=0) 573 ! Solve Minimize 2-norm(A * X - B). 574 575 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,s,&cond,&r,work,&lwork,&info) 576 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 577 578 integer intent(hide),depend(a):: m = shape(a,0) 579 integer intent(hide),depend(a):: n = shape(a,1) 580 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 581 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 582 <ftype2> dimension(m,n),intent(in,out,copy,out=v) :: a 583 584 integer depend(b),intent(hide):: nrhs = shape(b,1) 585 <ftype2> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 586 intent(in,out,copy,out=x) b 587 588 <ftype2> intent(in),optional :: cond = -1.0 589 integer intent(out,out=rank) :: r 590 <ftype2> intent(out),dimension(minmn),depend(minmn) :: s 591 592 integer optional,intent(in),depend(nrhs,minmn,maxmn),& 593 check(lwork>=1||lwork==-1) & 594 :: lwork=max(3*minmn+MAX(2*minmn,MAX(maxmn,nrhs)),1) 595 !check(lwork>=3*minmn+MAX(2*minmn,MAX(maxmn,nrhs))) 596 <ftype2> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 597 integer intent(out)::info 598 599end subroutine <prefix2>gelss 600 601subroutine <prefix2>gelss_lwork(m,n,maxmn,nrhs,a,b,s,cond,r,work,lwork,info) 602 ! work,info = gelss_lwork(m,n,nrhs,cond=-1.0) 603 ! Query for optimal lwork size 604 fortranname <prefix2>gelss 605 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&s,&cond,&r,&work,&lwork,&info) 606 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 607 608 integer intent(in):: m 609 integer intent(in):: n 610 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 611 <ftype2> intent(hide) :: a 612 613 integer intent(in):: nrhs 614 <ftype2> intent(hide) :: b 615 616 <ftype2> intent(in),optional :: cond = -1.0 617 integer intent(hide) :: r 618 <ftype2> intent(hide) :: s 619 620 integer optional,intent(in) :: lwork=-1 621 <ftype2> intent(out) :: work 622 integer intent(out)::info 623 624end subroutine <prefix2>gelss_lwork 625 626subroutine <prefix2c>gelss(m,n,minmn,maxmn,nrhs,a,b,s,cond,r,work,rwork,lwork,info) 627 ! v,x,s,rank,work,info = gelss(a,b,cond=-1.0,overwrite_a=0,overwrite_b=0) 628 ! Solve Minimize 2-norm(A * X - B). 629 630 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,s,&cond,&r,work,&lwork,rwork,&info) 631 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 632 633 integer intent(hide),depend(a):: m = shape(a,0) 634 integer intent(hide),depend(a):: n = shape(a,1) 635 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 636 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 637 <ftype2c> dimension(m,n),intent(in,out,copy,out=v) :: a 638 639 integer depend(b),intent(hide):: nrhs = shape(b,1) 640 <ftype2c> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 641 intent(in,out,copy,out=x) b 642 643 <ftype2> intent(in),optional :: cond = -1.0 644 integer intent(out,out=rank) :: r 645 <ftype2> intent(out),dimension(minmn),depend(minmn) :: s 646 647 integer optional,intent(in),depend(nrhs,minmn,maxmn),& 648 check(lwork>=1||lwork==-1) & 649 :: lwork=max(2*minmn+MAX(maxmn,nrhs),1) 650 ! check(lwork>=2*minmn+MAX(maxmn,nrhs)) 651 <ftype2c> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 652 <ftype2> dimension(5*minmn),intent(hide),depend(lwork) :: rwork 653 integer intent(out)::info 654 655end subroutine <prefix2c>gelss 656 657subroutine <prefix2c>gelss_lwork(m,n,maxmn,nrhs,a,b,s,cond,r,work,rwork,lwork,info) 658 ! work,info = gelss_lwork(m,n,nrhs,cond=-1.0) 659 ! Query for optimal lwork size 660 661 fortranname <prefix2c>gelss 662 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&s,&cond,&r,&work,&lwork,&rwork,&info) 663 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 664 665 integer intent(in):: m 666 integer intent(in):: n 667 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 668 <ftype2c> intent(hide) :: a 669 670 integer intent(in):: nrhs 671 <ftype2c> intent(hide) :: b 672 673 <ftype2> intent(in),optional :: cond = -1.0 674 integer intent(hide) :: r 675 <ftype2> intent(hide) :: s 676 677 integer optional,intent(in) :: lwork=-1 678 <ftype2c> intent(out) :: work 679 <ftype2> intent(hide) :: rwork 680 integer intent(out)::info 681 682end subroutine <prefix2>gelss_lwork 683 684subroutine <prefix2>gelsy(m,n,maxmn,minmn,nrhs,a,b,jptv,cond,r,work,lwork,info) 685 ! v,x,j,rank,info = dgelsy(a,b,jptv,cond,lwork,overwrite_a=True,overwrite_b=True) 686 ! Solve Minimize 2-norm(A * X - B). 687 688 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,jptv,&cond,&r,work,&lwork,&info) 689 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 690 691 integer intent(hide),depend(a):: m = shape(a,0) 692 integer intent(hide),depend(a):: n = shape(a,1) 693 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 694 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 695 <ftype2> dimension(m,n),intent(in,out,copy,out=v) :: a 696 697 integer depend(b),intent(hide):: nrhs = shape(b,1) 698 <ftype2> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 699 intent(in,out,copy,out=x) b 700 701 <ftype2> intent(in) :: cond 702 integer intent(out,out=rank) :: r 703 integer intent(in,out,out=j),dimension(n),depend(n) :: jptv 704 705 ! LWORK is obtained by the query call 706 integer intent(in),depend(nrhs,m,n,minmn) :: lwork 707 check(lwork>=MAX(minmn+3*n+1, 2*minmn+nrhs)) :: lwork 708 <ftype2> dimension(lwork),intent(cache,hide),depend(lwork) :: work 709 integer intent(out)::info 710 711end subroutine <prefix2>gelsy 712 713subroutine <prefix2>gelsy_lwork(m,n,maxmn,nrhs,a,b,jptv,cond,r,work,lwork,info) 714 ! work,info = dgelsy_lwork(m,n,nrhs,cond) 715 ! Query for optimal lwork size 716 717 fortranname <prefix2>gelsy 718 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&jptv,&cond,&r,&work,&lwork,&info) 719 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 720 721 integer intent(in) :: m 722 integer intent(in) :: n 723 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 724 <ftype2> intent(hide) :: a 725 726 integer intent(in):: nrhs 727 <ftype2> intent(hide):: b 728 729 <ftype2> intent(in) :: cond 730 integer intent(hide) :: r 731 integer intent(hide):: jptv 732 733 integer intent(in),optional :: lwork = -1 734 <ftype2> intent(out) :: work 735 integer intent(out)::info 736 737end subroutine <prefix2>gelsy_lwork 738 739subroutine <prefix2c>gelsy(m,n,maxmn,minmn,nrhs,a,b,jptv,cond,r,work,lwork,rwork,info) 740 ! v,x,j,rank,info = zgelsy(a,b,jptv,cond,lwork,overwrite_a=True,overwrite_b=True) 741 ! Solve Minimize 2-norm(A * X - B). 742 743 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,jptv,&cond,&r,work,&lwork,rwork,&info) 744 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 745 746 integer intent(hide),depend(a):: m = shape(a,0) 747 integer intent(hide),depend(a):: n = shape(a,1) 748 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 749 integer intent(hide), depend(m,n) :: minmn = MIN(m,n) 750 <ftype2c> dimension(m,n),intent(in,out,copy,out=v) :: a 751 752 integer depend(b),intent(hide):: nrhs = shape(b,1) 753 <ftype2c> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 754 intent(in,out,copy,out=x) b 755 756 <ftype2> intent(in) :: cond 757 integer intent(out,out=rank) :: r 758 integer intent(in,out,out=j),dimension(n),depend(n) :: jptv 759 760 ! LWORK is obtained by the query call 761 integer intent(in),depend(nrhs,m,n,minmn) :: lwork 762 check(lwork>=minmn+MAX(MAX(2*minmn, n+1), minmn+nrhs)) :: lwork 763 <ftype2c> dimension(lwork),intent(hide,cache),depend(lwork) :: work 764 <ftype2> dimension(2*n),intent(hide,cache),depend(n) :: rwork 765 integer intent(out)::info 766 767end subroutine <prefix2c>gelsy 768 769subroutine <prefix2c>gelsy_lwork(m,n,maxmn,nrhs,a,b,jptv,cond,r,work,lwork,rwork,info) 770 ! work,info = zgelsy_lwork(m,n,nrhs,cond) 771 ! Query for optimal lwork size 772 773 fortranname <prefix2c>gelsy 774 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&jptv,&cond,&r,&work,&lwork,&rwork,&info) 775 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 776 777 integer intent(in) :: m 778 integer intent(in) :: n 779 integer intent(hide) :: maxmn = MAX(m,n) 780 <ftype2c> intent(hide) :: a 781 782 integer intent(in):: nrhs 783 <ftype2c> intent(hide) :: b 784 785 <ftype2> intent(in) :: cond 786 integer intent(hide) :: r 787 integer intent(hide) :: jptv 788 789 integer intent(in),optional :: lwork = -1 790 <ftype2c> intent(out) :: work 791 <ftype2> intent(hide) :: rwork 792 integer intent(out)::info 793 794end subroutine <prefix2c>gelsy_lwork 795 796subroutine <prefix2>gelsd(m,n,minmn,maxmn,nrhs,a,b,s,cond,r,work,lwork,size_iwork,iwork,info) 797 ! x,s,rank,info = dgelsd(a,b,lwork,size_iwork,cond=-1.0,overwrite_a=True,overwrite_b=True) 798 ! Solve Minimize 2-norm(A * X - B). 799 800 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,s,&cond,&r,work,&lwork,iwork,&info) 801 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 802 803 integer intent(hide),depend(a):: m = shape(a,0) 804 integer intent(hide),depend(a):: n = shape(a,1) 805 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 806 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 807 <ftype2> dimension(m,n),intent(in,copy) :: a 808 809 integer depend(b),intent(hide):: nrhs = shape(b,1) 810 <ftype2> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 811 intent(in,out,copy,out=x) b 812 813 <ftype2> intent(in),optional :: cond=-1.0 814 integer intent(out,out=rank) :: r 815 <ftype2> intent(out),dimension(minmn),depend(minmn) :: s 816 817 integer intent(in),check(lwork>=1) :: lwork 818 ! Impossible to calculate lwork explicitly, need to obtain it from query call first 819 ! Same for size_iwork 820 <ftype2> dimension(lwork),intent(cache,hide),depend(lwork) :: work 821 822 integer intent(in) :: size_iwork 823 integer intent(cache,hide),dimension(MAX(1,size_iwork)),depend(size_iwork) :: iwork 824 integer intent(out)::info 825 826end subroutine <prefix2>gelsd 827 828subroutine <prefix2>gelsd_lwork(m,n,maxmn,nrhs,a,b,s,cond,r,work,lwork,iwork,info) 829 ! work,iwork,info = dgelsd_lwork(m,n,nrhs,cond=-1.0) 830 ! Query for optimal lwork size 831 832 fortranname <prefix2>gelsd 833 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&s,&cond,&r,&work,&lwork,&iwork,&info) 834 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 835 836 integer intent(in) :: m 837 integer intent(in) :: n 838 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 839 <ftype2> intent(hide) :: a 840 841 integer intent(in):: nrhs 842 <ftype2> intent(hide) :: b 843 844 <ftype2> intent(in),optional :: cond=-1.0 845 integer intent(hide) :: r 846 <ftype2> intent(hide) :: s 847 848 integer intent(in),optional :: lwork = -1 849 <ftype2> intent(out) :: work 850 851 integer intent(out) :: iwork 852 integer intent(out)::info 853 854end subroutine <prefix2>gelsd_lwork 855 856subroutine <prefix2c>gelsd(m,n,minmn,maxmn,nrhs,a,b,s,cond,r,work,lwork,size_rwork,rwork, size_iwork,iwork,info) 857 ! x,s,rank,info = zgelsd(a,b,lwork,size_rwork,size_iwork,cond=-1.0,overwrite_a=True,overwrite_b=True) 858 ! Solve Minimize 2-norm(A * X - B). 859 860 callstatement (*f2py_func)(&m,&n,&nrhs,a,&m,b,&maxmn,s,&cond,&r,work,&lwork, rwork, iwork,&info) 861 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*, <ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT* 862 863 integer intent(hide),depend(a):: m = shape(a,0) 864 integer intent(hide),depend(a):: n = shape(a,1) 865 integer intent(hide),depend(m,n):: minmn = MIN(m,n) 866 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 867 <ftype2c> dimension(m,n),intent(in,copy) :: a 868 869 integer depend(b),intent(hide):: nrhs = shape(b,1) 870 <ftype2c> dimension(maxmn,nrhs),check(maxmn==shape(b,0)),depend(maxmn) :: b 871 intent(in,out,copy,out=x) b 872 873 <ftype2> intent(in),optional :: cond=-1.0 874 integer intent(out,out=rank) :: r 875 <ftype2> intent(out),dimension(minmn),depend(minmn) :: s 876 877 integer intent(in),check(lwork>=1||lwork==-1) :: lwork 878 ! Impossible to calculate lwork explicitly, need to obtain it from query call first 879 ! Same for size_rwork, size_iwork 880 <ftype2c> dimension(MAX(lwork,1)),intent(cache,hide),depend(lwork) :: work 881 882 integer intent(in) :: size_rwork 883 <ftype2> intent(cache,hide),dimension(MAX(1,size_rwork)),depend(size_rwork) :: rwork 884 885 integer intent(in) :: size_iwork 886 integer intent(cache,hide),dimension(MAX(1,size_iwork)),depend(size_iwork) :: iwork 887 integer intent(out)::info 888 889end subroutine <prefix2c>gelsd 890 891subroutine <prefix2c>gelsd_lwork(m,n,maxmn,nrhs,a,b,s,cond,r,work,lwork,rwork,iwork,info) 892 ! work,rwork,iwork,info = zgelsd_lwork(m,n,nrhs,lwork=-1.0,cond=-1.0) 893 ! Query for optimal lwork size 894 895 fortranname <prefix2c>gelsd 896 callstatement (*f2py_func)(&m,&n,&nrhs,&a,&m,&b,&maxmn,&s,&cond,&r,&work,&lwork, &rwork, &iwork,&info) 897 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*, <ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT* 898 899 integer intent(in) :: m 900 integer intent(in) :: n 901 integer intent(hide),depend(m,n):: maxmn = MAX(m,n) 902 <ftype2c> intent(hide) :: a 903 904 integer intent(in):: nrhs 905 <ftype2c> intent(hide) :: b 906 907 <ftype2> intent(in),optional :: cond=-1.0 908 integer intent(hide) :: r 909 <ftype2> intent(hide) :: s 910 911 integer intent(in),optional :: lwork = -1 912 <ftype2c> intent(out) :: work 913 <ftype2> intent(out) :: rwork 914 915 integer intent(out) :: iwork 916 integer intent(out)::info 917 918end subroutine <prefix2c>gelsd_lwork 919 920subroutine <prefix2>geqp3(m,n,a,jpvt,tau,work,lwork,info) 921 ! qr_a,jpvt,tau,work,info = geqp3(a,lwork=3*(n+1),overwrite_a=0) 922 ! Compute a QR factorization of a real M-by-N matrix A with column pivoting: 923 ! A * P = Q * R. 924 925 threadsafe 926 callstatement (*f2py_func)(&m,&n,a,&m,jpvt,tau,work,&lwork,&info) 927 callprotoargument F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 928 integer intent(hide),depend(a):: m = shape(a,0) 929 integer intent(hide),depend(a):: n = shape(a,1) 930 <ftype2> dimension(m,n),intent(in,out,copy,out=qr,aligned8) :: a 931 integer dimension(n),intent(out) :: jpvt 932 <ftype2> dimension(MIN(m,n)),intent(out) :: tau 933 934 integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) :: lwork=max(3*(n+1),1) 935 <ftype2> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 936 integer intent(out) :: info 937end subroutine <prefix2>geqp3 938 939subroutine <prefix2c>geqp3(m,n,a,jpvt,tau,work,lwork,rwork,info) 940 ! qr_a,jpvt,tau,work,info = geqp3(a,lwork,overwrite_a=0) 941 ! Compute a QR factorization of a complex M-by-N matrix A with column pivoting: 942 ! A * P = Q * R. 943 944 threadsafe 945 callstatement (*f2py_func)(&m,&n,a,&m,jpvt,tau,work,&lwork,rwork,&info) 946 callprotoargument F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 947 948 integer intent(hide),depend(a):: m = shape(a,0) 949 integer intent(hide),depend(a):: n = shape(a,1) 950 <ftype2c> dimension(m,n),intent(in,out,copy,out=qr,aligned8) :: a 951 integer dimension(n),intent(out) :: jpvt 952 <ftype2c> dimension(MIN(m,n)),intent(out) :: tau 953 954 integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) :: lwork=max(3*(n+1),1) 955 <ftype2c> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 956 <ftype2> dimension(2*n),intent(hide),depend(n) :: rwork 957 integer intent(out) :: info 958 959end subroutine <prefix2c>geqp3 960 961subroutine <prefix>geqrf(m,n,a,tau,work,lwork,info) 962 ! qr_a,tau,work,info = geqrf(a,lwork=3*n,overwrite_a=0) 963 ! Compute a QR factorization of a real M-by-N matrix A: 964 ! A = Q * R. 965 966 threadsafe 967 callstatement (*f2py_func)(&m,&n,a,&m,tau,work,&lwork,&info) 968 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 969 970 integer intent(hide),depend(a):: m = shape(a,0) 971 integer intent(hide),depend(a):: n = shape(a,1) 972 <ftype> dimension(m,n),intent(in,out,copy,out=qr,aligned8) :: a 973 <ftype> dimension(MIN(m,n)),intent(out) :: tau 974 975 integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) :: lwork=max(3*n,1) 976 <ftype> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 977 integer intent(out) :: info 978 979end subroutine <prefix>geqrf 980 981subroutine <prefix>geqrf_lwork(m,n,a,tau,work,lwork,info) 982 ! work, info = geqrf_lwork(m, n) 983 ! Calculate the optimal size of the ?geqrf work array. 984 985 fortranname <prefix>geqrf 986 callstatement (*f2py_func)(&m,&n,&a,&m,&tau,&work,&lwork,&info) 987 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 988 989 integer intent(in), check(m > 0) :: m 990 integer intent(in), check(n > 0) :: n 991 <ftype> intent(hide) :: a 992 <ftype> intent(hide) :: tau 993 <ftype> intent(out) :: work 994 integer intent(hide) :: lwork = -1 995 integer intent(out) :: info 996 997end subroutine <prefix>geqrf_lwork 998 999subroutine <prefix>geqrfp(m,n,a,lda,tau,work,lwork,info) 1000 ! qr,tau,work,info = geqrfp(a,lwork=3*n,overwrite_a=0) 1001 ! DGEQR2P computes a QR factorization of a real M-by-N matrix A: 1002 ! 1003 ! A = Q * ( R ), 1004 ! ( 0 ) 1005 ! 1006 ! where: 1007 ! 1008 ! Q is a M-by-M orthogonal matrix; 1009 ! R is an upper-triangular N-by-N matrix with nonnegative diagonal 1010 ! entries; 1011 ! 0 is a (M-N)-by-N zero matrix, if M > N. 1012 1013 callstatement (*f2py_func)(&m,&n,a,&lda,tau,work,&lwork,&info) 1014 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1015 1016 integer intent(hide), check(m > 0), depend(a) :: m = shape(a, 0) 1017 integer intent(hide), check(n > 0), depend(a) :: n = shape(a, 1) 1018 <ftype> intent(in,out,copy,out=qr), dimension(m,n) :: a 1019 integer intent(hide), depend(a) :: lda = max(1, shape(a, 0)) 1020 <ftype> intent(out), depend(m,n), dimension(MIN(m,n)) :: tau 1021 <ftype> intent(hide), depend(lwork), dimension(lwork) :: work 1022 integer intent(in), check(lwork>=n||lwork==-1) :: lwork = MAX(1, n) 1023 integer intent(out) :: info 1024 1025end subroutine <prefix>geqrfp 1026 1027subroutine <prefix>geqrfp_lwork(m,n,a,lda,tau,work,lwork,info) 1028 ! work, info = geqrfp_lwork(m, n) 1029 1030 fortranname <prefix>geqrfp 1031 callstatement (*f2py_func)(&m,&n,&a,&lda,&tau,&work,&lwork,&info) 1032 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1033 1034 integer intent(in), check(m > 0) :: m 1035 integer intent(in), check(n > 0) :: n 1036 <ftype> intent(hide) :: a 1037 integer intent(hide), depend(m) :: lda = m 1038 <ftype> intent(hide) :: tau 1039 <ftype> intent(out) :: work 1040 integer intent(hide) :: lwork = -1 1041 integer intent(out) :: info 1042 1043end subroutine <prefix>geqrfp_lwork 1044 1045subroutine <prefix>gerqf(m,n,a,tau,work,lwork,info) 1046 ! rq_a,tau,work,info = gerqf(a,lwork=3*n,overwrite_a=0) 1047 ! Compute an RQ factorization of a real M-by-N matrix A: 1048 ! A = R * Q. 1049 1050 threadsafe 1051 callstatement (*f2py_func)(&m,&n,a,&m,tau,work,&lwork,&info) 1052 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1053 1054 integer intent(hide),depend(a):: m = shape(a,0) 1055 integer intent(hide),depend(a):: n = shape(a,1) 1056 <ftype> dimension(m,n),intent(in,out,copy,out=qr,aligned8) :: a 1057 <ftype> dimension(MIN(m,n)),intent(out) :: tau 1058 1059 integer optional,intent(in),depend(m),check(lwork>=m||lwork==-1) :: lwork=max(3*m,1) 1060 <ftype> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 1061 integer intent(out) :: info 1062 1063end subroutine <prefix>gerqf 1064 1065subroutine <prefix2>geev(compute_vl,compute_vr,n,a,wr,wi,vl,ldvl,vr,ldvr,work,lwork,info) 1066 ! wr,wi,vl,vr,info = geev(a,compute_vl=1,compute_vr=1,lwork=4*n,overwrite_a=0) 1067 1068 callstatement {(*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,wr,wi,vl,&ldvl,vr,&ldvr,work,&lwork,&info);} 1069 callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 1070 1071 integer optional,intent(in):: compute_vl = 1 1072 check(compute_vl==1||compute_vl==0) compute_vl 1073 integer optional,intent(in):: compute_vr = 1 1074 check(compute_vr==1||compute_vr==0) compute_vr 1075 1076 integer intent(hide),depend(a) :: n = shape(a,0) 1077 <ftype2> dimension(n,n),intent(in,copy,aligned8) :: a 1078 check(shape(a,0)==shape(a,1)) :: a 1079 1080 <ftype2> dimension(n),intent(out),depend(n) :: wr 1081 <ftype2> dimension(n),intent(out),depend(n) :: wi 1082 1083 <ftype2> dimension(ldvl,n),intent(out) :: vl 1084 integer intent(hide),depend(n,compute_vl) :: ldvl=(compute_vl?n:1) 1085 1086 <ftype2> dimension(ldvr,n),intent(out) :: vr 1087 integer intent(hide),depend(n,compute_vr) :: ldvr=(compute_vr?n:1) 1088 1089 integer optional,intent(in),depend(n,compute_vl,compute_vr) :: lwork=max(4*n,1) 1090 check(lwork>=((compute_vl||compute_vr)?4*n:3*n)) :: lwork 1091 <ftype2> dimension(lwork),intent(hide,cache),depend(lwork) :: work 1092 1093 integer intent(out):: info 1094 1095end subroutine <prefix2>geev 1096 1097subroutine <prefix2>geev_lwork(compute_vl,compute_vr,n,a,wr,wi,vl,ldvl,vr,ldvr,work,lwork,info) 1098 ! LWORK=-1 call for (S/D)GEEV --- keep in sync with above (S/D)GEEV definition 1099 1100 fortranname <prefix2>geev 1101 callstatement {(*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,&a,&n,&wr,&wi,&vl,&ldvl,&vr,&ldvr,&work,&lwork,&info);} 1102 callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 1103 1104 integer optional,intent(in):: compute_vl = 1 1105 check(compute_vl==1||compute_vl==0) compute_vl 1106 integer optional,intent(in):: compute_vr = 1 1107 check(compute_vr==1||compute_vr==0) compute_vr 1108 1109 integer intent(in) :: n 1110 <ftype2> intent(hide) :: a 1111 1112 <ftype2> intent(hide) :: wr 1113 <ftype2> intent(hide) :: wi 1114 1115 <ftype2> intent(hide) :: vl 1116 integer intent(hide),depend(n,compute_vl) :: ldvl=(compute_vl?n:1) 1117 1118 <ftype2> intent(hide) :: vr 1119 integer intent(hide),depend(n,compute_vr) :: ldvr=(compute_vr?n:1) 1120 1121 integer intent(hide) :: lwork = -1 1122 <ftype2> intent(out) :: work 1123 1124 integer intent(out):: info 1125 1126end subroutine <prefix2>geev_lwork 1127 1128subroutine <prefix2c>geev(compute_vl,compute_vr,n,a,w,vl,ldvl,vr,ldvr,work,lwork,rwork,info) 1129 ! w,vl,vr,info = geev(a,compute_vl=1,compute_vr=1,lwork=2*n,overwrite_a=0) 1130 1131 callstatement (*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,w,vl,&ldvl,vr,&ldvr,work,&lwork,rwork,&info) 1132 callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 1133 1134 integer optional,intent(in):: compute_vl = 1 1135 check(compute_vl==1||compute_vl==0) compute_vl 1136 integer optional,intent(in):: compute_vr = 1 1137 check(compute_vr==1||compute_vr==0) compute_vr 1138 1139 integer intent(hide),depend(a) :: n = shape(a,0) 1140 <ftype2c> dimension(n,n),intent(in,copy) :: a 1141 check(shape(a,0)==shape(a,1)) :: a 1142 1143 <ftype2c> dimension(n),intent(out),depend(n) :: w 1144 1145 <ftype2c> dimension(ldvl,n),depend(ldvl),intent(out) :: vl 1146 integer intent(hide),depend(compute_vl,n) :: ldvl=(compute_vl?n:1) 1147 1148 <ftype2c> dimension(ldvr,n),depend(ldvr),intent(out) :: vr 1149 integer intent(hide),depend(compute_vr,n) :: ldvr=(compute_vr?n:1) 1150 1151 integer optional,intent(in),depend(n) :: lwork=max(2*n,1) 1152 check(lwork>=2*n) :: lwork 1153 <ftype2c> dimension(lwork),intent(hide),depend(lwork) :: work 1154 <ftype2> dimension(2*n),intent(hide,cache),depend(n) :: rwork 1155 1156 integer intent(out):: info 1157 1158end subroutine <prefix2c>geev 1159 1160subroutine <prefix2c>geev_lwork(compute_vl,compute_vr,n,a,w,vl,ldvl,vr,ldvr,work,lwork,rwork,info) 1161 ! LWORK=-1 call for (C/Z)GEEV --- keep in sync with above (C/Z)GEEV definition 1162 1163 fortranname <prefix2c>geev 1164 callstatement (*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,&a,&n,&w,&vl,&ldvl,&vr,&ldvr,&work,&lwork,&rwork,&info) 1165 callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 1166 1167 integer optional,intent(in):: compute_vl = 1 1168 check(compute_vl==1||compute_vl==0) compute_vl 1169 integer optional,intent(in):: compute_vr = 1 1170 check(compute_vr==1||compute_vr==0) compute_vr 1171 1172 integer intent(in) :: n 1173 <ftype2c> intent(hide) :: a 1174 1175 <ftype2c> intent(hide) :: w 1176 1177 <ftype2c> intent(hide) :: vl 1178 integer intent(hide),depend(compute_vl,n) :: ldvl=(compute_vl?n:1) 1179 1180 <ftype2c> intent(hide) :: vr 1181 integer intent(hide),depend(compute_vr,n) :: ldvr=(compute_vr?n:1) 1182 1183 integer intent(hide) :: lwork = -1 1184 <ftype2c> intent(out) :: work 1185 <ftype2> intent(hide) :: rwork 1186 1187 integer intent(out):: info 1188 1189end subroutine <prefix2c>geev_lwork 1190 1191subroutine <prefix2>gegv(compute_vl,compute_vr,n,a,b,alphar,alphai,beta,vl,ldvl,vr,ldvr,work,lwork,info) 1192 ! Compute the generalized eigenvalues (alphar +/- alphai*i, beta) 1193 ! of the real nonsymmetric matrices A and B: det(A-w*B)=0 where w=alpha/beta. 1194 ! Optionally, compute the left and/or right generalized eigenvectors: 1195 ! (A - w B) r = 0, l^H * (A - w B) = 0 1196 ! 1197 ! alphar,alphai,beta,vl,vr,info = gegv(a,b,compute_vl=1,compute_vr=1,lwork=8*n,overwrite_a=0,overwrite_b=0) 1198 1199 callstatement (*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,b,&n,alphar,alphai,beta,vl,&ldvl,vr,&ldvr,work,&lwork,&info) 1200 callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 1201 1202 integer optional,intent(in):: compute_vl = 1 1203 check(compute_vl==1||compute_vl==0) compute_vl 1204 integer optional,intent(in):: compute_vr = 1 1205 check(compute_vr==1||compute_vr==0) compute_vr 1206 1207 integer intent(hide),depend(a) :: n = shape(a,0) 1208 <ftype2> dimension(n,n),intent(in,copy) :: a 1209 check(shape(a,0)==shape(a,1)) :: a 1210 1211 <ftype2> dimension(n,n),depend(n),intent(in,copy) :: b 1212 check(shape(b,0)==shape(b,1) && shape(b,0)==n) :: b 1213 1214 <ftype2> dimension(n),depend(n),intent(out) :: alphar 1215 <ftype2> dimension(n),depend(n),intent(out) :: alphai 1216 <ftype2> dimension(n),depend(n),intent(out) :: beta 1217 1218 <ftype2> dimension(ldvl,n),intent(out),depend(ldvl) :: vl 1219 integer intent(hide),depend(compute_vl,n) :: ldvl=(compute_vl?n:1) 1220 1221 <ftype2> dimension(ldvr,n),intent(out),depend(ldvr) :: vr 1222 integer intent(hide),depend(compute_vr,n) :: ldvr=(compute_vr?n:1) 1223 1224 integer optional,intent(in),depend(n) :: lwork=max(8*n,1) 1225 check(lwork>=8*n) :: lwork 1226 <ftype2> dimension(lwork),intent(hide),depend(lwork) :: work 1227 1228 integer intent(out):: info 1229 1230end subroutine <prefix2>gegv 1231 1232subroutine <prefix2c>gegv(compute_vl,compute_vr,n,a,b,alpha,beta,vl,ldvl,vr,ldvr,work,lwork,rwork,info) 1233 ! Compute the generalized eigenvalues (alpha, beta) 1234 ! of the comples nonsymmetric matrices A and B: det(A-w*B)=0 where w=alpha/beta. 1235 ! Optionally, compute the left and/or right generalized eigenvectors: 1236 ! (A - w B) r = 0, l^H * (A - w B) = 0 1237 ! 1238 ! alpha,beta,vl,vr,info = gegv(a,b,compute_vl=1,compute_vr=1,lwork=2*n,overwrite_a=0,overwrite_b=0) 1239 1240 callstatement (*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,b,&n,alpha,beta,vl,&ldvl,vr,&ldvr,work,&lwork,rwork,&info) 1241 callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 1242 1243 integer optional,intent(in):: compute_vl = 1 1244 check(compute_vl==1||compute_vl==0) compute_vl 1245 integer optional,intent(in):: compute_vr = 1 1246 check(compute_vr==1||compute_vr==0) compute_vr 1247 1248 integer intent(hide),depend(a) :: n = shape(a,0) 1249 <ftype2c> dimension(n,n),intent(in,copy) :: a 1250 check(shape(a,0)==shape(a,1)) :: a 1251 1252 <ftype2c> dimension(n,n),depend(n),intent(in,copy) :: b 1253 check(shape(b,0)==shape(b,1) && shape(b,0)==n) :: b 1254 1255 <ftype2c> dimension(n),depend(n),intent(out) :: alpha 1256 <ftype2c> dimension(n),depend(n),intent(out) :: beta 1257 1258 <ftype2c> dimension(ldvl,n),intent(out),depend(ldvl) :: vl 1259 integer intent(hide),depend(compute_vl,n) :: ldvl=(compute_vl?n:1) 1260 1261 <ftype2c> dimension(ldvr,n),intent(out),depend(ldvr) :: vr 1262 integer intent(hide),depend(compute_vr,n) :: ldvr=(compute_vr?n:1) 1263 1264 integer optional,intent(in),depend(n) :: lwork=max(2*n,1) 1265 check(lwork>=2*n) :: lwork 1266 <ftype2c> dimension(lwork),intent(hide),depend(lwork) :: work 1267 <ftype2> dimension(8*n),intent(hide),depend(n) :: rwork 1268 1269 integer intent(out):: info 1270 1271end subroutine <prefix2c>gegv 1272 1273subroutine <prefix2c>gees(compute_v,sort_t,<prefix2c>select,n,a,nrows,sdim,w,vs,ldvs,work,lwork,rwork,bwork,info) 1274 ! t,sdim,w,vs,work,info=gees(compute_v=1,sort_t=0,select,a,lwork=3*n) 1275 ! For an NxN matrix compute the eigenvalues, the schur form T, and optionally 1276 ! the matrix of Schur vectors Z. This gives the Schur factorization 1277 ! A = Z * T * Z^H -- a complex matrix is in Schur form if it is upper 1278 ! triangular 1279 1280 callstatement (*f2py_func)((compute_v?"V":"N"),(sort_t?"S":"N"),cb_<prefix2c>select_in_gees__user__routines,&n,a,&nrows,&sdim,w,vs,&ldvs,work,&lwork,rwork,bwork,&info,1,1) 1281 callprotoargument char*,char*,F_INT(*)(<ctype2c>*),F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT,F_INT 1282 1283 use gees__user__routines 1284 1285 integer optional,intent(in),check(compute_v==0||compute_v==1) :: compute_v = 1 1286 integer optional,intent(in),check(sort_t==0||sort_t==1) :: sort_t = 0 1287 external <prefix2c>select 1288 integer intent(hide),depend(a) :: n = shape(a,1) 1289 <ftype2c> intent(in,out,copy,out=t),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a 1290 integer intent(hide),depend(a) :: nrows=shape(a,0) 1291 integer intent(out) :: sdim=0 1292 <ftype2c> intent(out),dimension(n) :: w 1293 <ftype2c> intent(out),depend(ldvs,n),dimension(ldvs,n) :: vs 1294 integer intent(hide),depend(compute_v,n) :: ldvs=((compute_v==1)?n:1) 1295 <ftype2c> intent(out),depend(lwork),dimension(MAX(lwork,1)) :: work 1296 integer optional,intent(in),check((lwork==-1)||(lwork >= MAX(1,2*n))),depend(n) :: lwork = max(3*n,1) 1297 <ftype2> optional,intent(hide),depend(n),dimension(n) :: rwork 1298 logical optional,intent(hide),depend(n),dimension(n) :: bwork 1299 integer intent(out) :: info 1300 1301end subroutine <prefix2c>gees 1302 1303subroutine <prefix2>gees(compute_v,sort_t,<prefix2>select,n,a,nrows,sdim,wr,wi,vs,ldvs,work,lwork,bwork,info) 1304 ! t,sdim,w,vs,work,info=gees(compute_v=1,sort_t=0,select,a,lwork=3*n) 1305 ! For an NxN matrix compute the eigenvalues, the schur form T, and optionally 1306 ! the matrix of Schur vectors Z. This gives the Schur factorization 1307 ! A = Z * T * Z^H -- a real matrix is in Schur form if it is upper quasi- 1308 ! triangular with 1x1 and 2x2 blocks. 1309 1310 callstatement (*f2py_func)((compute_v?"V":"N"),(sort_t?"S":"N"),cb_<prefix2>select_in_gees__user__routines,&n,a,&nrows,&sdim,wr,wi,vs,&ldvs,work,&lwork,bwork,&info,1,1) 1311 callprotoargument char*,char*,F_INT(*)(<ctype2>*,<ctype2>*),F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT,F_INT 1312 1313 use gees__user__routines 1314 1315 integer optional,intent(in),check(compute_v==0||compute_v==1) :: compute_v = 1 1316 integer optional,intent(in),check(sort_t==0||sort_t==1) :: sort_t = 0 1317 external <prefix2>select 1318 integer intent(hide),depend(a) :: n = shape(a,1) 1319 <ftype2> intent(in,out,copy,out=t,aligned8),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a 1320 integer intent(hide),depend(a) :: nrows=shape(a,0) 1321 integer intent(out) :: sdim=0 1322 <ftype2> intent(out),dimension(n) :: wr 1323 <ftype2> intent(out),dimension(n) :: wi 1324 <ftype2> intent(out),depend(ldvs,n),dimension(ldvs,n) :: vs 1325 integer intent(hide),depend(compute_v,n) :: ldvs=((compute_v==1)?n:1) 1326 <ftype2> intent(out),depend(lwork),dimension(MAX(lwork,1)) :: work 1327 integer optional,intent(in),check((lwork==-1)||(lwork >= MAX(1,2*n))),depend(n) :: lwork = max(3*n,1) 1328 <ftype2> optional,intent(hide),depend(n),dimension(n) :: rwork 1329 logical optional,intent(hide),depend(n),dimension(n) :: bwork 1330 integer intent(out) :: info 1331 1332end subroutine <prefix2>gees 1333 1334subroutine <prefix2>gges(jobvsl,jobvsr,sort_t,<prefix2>select,n,a,lda,b,ldb,sdim,alphar,alphai,beta,vsl,ldvsl,vsr,ldvsr,work,lwork,bwork,info) 1335 ! For a pair of N-by-N real nonsymmetric matrices (A,B) computes 1336 ! the generalized eigenvalues, the generalized real Schur form (S,T), 1337 ! optionally, the left and/or right matrices of Schur vectors (VSL 1338 ! and VSR). 1339 ! (A,B) = ( (VSL)*S*(VSR)**T, (VSL)*T*(VSR)**T ) 1340 1341 callstatement (*f2py_func)((jobvsl?"V":"N"),(jobvsr?"V":"N"),(sort_t?"S":"N"),cb_<prefix2>select_in_gges__user__routines,&n,a,&lda,b,&ldb,&sdim,alphar,alphai,beta,vsl,&ldvsl,vsr,&ldvsr,work,&lwork,bwork,&info) 1342 callprotoargument char*,char*,char*,F_INT(*)(<ctype2>*,<ctype2>*,<ctype2>*),F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 1343 1344 use gges__user__routines 1345 1346 integer optional,intent(in),check(jobvsl==0||jobvsl==1) :: jobvsl=1 1347 integer optional,intent(in),check(jobvsr==0||jobvsr==1) :: jobvsr=1 1348 integer optional,intent(in),check(sort_t==0||sort_t==1) :: sort_t=0 1349 external <prefix2>select 1350 integer intent(hide),depend(a) :: n=shape(a,1) 1351 <ftype2> intent(in,out,copy),dimension(lda,n) :: a 1352 integer intent(hide),depend(a) :: lda=shape(a,0) 1353 <ftype2> intent(in,out,copy),dimension(ldb,n) :: b 1354 integer intent(hide),depend(b) :: ldb=shape(b,0) 1355 integer intent(out) :: sdim=0 1356 <ftype2> intent(out),dimension(n) :: alphar 1357 <ftype2> intent(out),dimension(n) :: alphai 1358 <ftype2> intent(out),dimension(n) :: beta 1359 <ftype2> intent(out),depend(ldvsl,n),dimension(ldvsl,n) :: vsl 1360 integer optional,intent(in),depend(jobvsl,n) :: ldvsl=((jobvsl==1)?n:1) 1361 <ftype2> intent(out),depend(ldvsr,n),dimension(ldvsr,n) :: vsr 1362 integer optional,intent(in),depend(jobvsr,n) :: ldvsr=((jobvsr==1)?n:1) 1363 <ftype2> intent(out),depend(lwork),dimension(MAX(lwork,1)) :: work 1364 integer optional,intent(in),depend(n),check(lwork>=MAX(1,MAX(8*n, 6*n+16))||lwork==-1):: lwork=max(8*n+16,1) 1365 logical optional,intent(hide),depend(n),dimension(n) :: bwork 1366 integer intent(out):: info 1367 1368end subroutine <prefix2>gges 1369 1370subroutine <prefix2c>gges(jobvsl,jobvsr,sort_t,<prefix2c>select,n,a,lda,b,ldb,sdim,alpha,beta,vsl,ldvsl,vsr,ldvsr,work,lwork,rwork,bwork,info) 1371 ! For a pair of N-by-N complex nonsymmetric matrices (A,B) computes 1372 ! the generalized eigenvalues, the generalized real Schur form (S,T), 1373 ! optionally, the left and/or right matrices of Schur vectors (VSL 1374 ! and VSR). 1375 ! (A,B) = ( (VSL)*S*(VSR)**T, (VSL)*T*(VSR)**H ) 1376 1377 callstatement (*f2py_func)((jobvsl?"V":"N"),(jobvsr?"V":"N"),(sort_t?"S":"N"),cb_<prefix2c>select_in_gges__user__routines,&n,a,&lda,b,&ldb,&sdim,alpha,beta,vsl,&ldvsl,vsr,&ldvsr,work,&lwork,rwork,bwork,&info) 1378 callprotoargument char*,char*,char*,F_INT(*)(<ctype2c>*,<ctype2c>*),F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT* 1379 1380 use gges__user__routines 1381 1382 integer optional,intent(in),check(jobvsl==0||jobvsl==1) :: jobvsl=1 1383 integer optional,intent(in),check(jobvsr==0||jobvsr==1) :: jobvsr=1 1384 integer optional,intent(in),check(sort_t==0||sort_t==1) :: sort_t=0 1385 external <prefix2c>select 1386 integer intent(hide),depend(a) :: n=shape(a,1) 1387 <ftype2c> intent(in,out,copy),dimension(lda,n) :: a 1388 integer intent(hide),depend(a) :: lda=shape(a,0) 1389 <ftype2c> intent(in,out,copy),dimension(ldb,n) :: b 1390 integer intent(hide),depend(b) :: ldb=shape(b,0) 1391 integer intent(out) :: sdim=0 1392 <ftype2c> intent(out),dimension(n) :: alpha 1393 <ftype2c> intent(out),dimension(n) :: beta 1394 <ftype2c> intent(out),depend(ldvsl,n),dimension(ldvsl,n) :: vsl 1395 integer optional,intent(in),depend(jobvsl,n) :: ldvsl=((jobvsl==1)?n:1) 1396 <ftype2c> intent(out),depend(ldvsr,n),dimension(ldvsr,n) :: vsr 1397 integer optional,intent(in),depend(jobvsr,n) :: ldvsr=((jobvsr==1)?n:1) 1398 <ftype2c> intent(out),depend(lwork),dimension(MAX(lwork,1)) :: work 1399 integer optional,intent(in),depend(n),check(lwork>=MAX(1,2*n)||lwork==-1):: lwork=max(2*n,1) 1400 <ftype2> intent(hide),dimension(8*n) :: rwork 1401 logical optional,intent(hide),depend(n),dimension(n) :: bwork 1402 integer intent(out):: info 1403 1404end subroutine <prefix2c>gges 1405 1406subroutine <prefix2>ggev(compute_vl,compute_vr,n,a,b,alphar,alphai,beta,vl,ldvl,vr,ldvr,work,lwork,info) 1407 1408 callstatement {(*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,b,&n,alphar,alphai,beta,vl,&ldvl,vr,&ldvr,work,&lwork,&info);} 1409 callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 1410 1411 integer optional,intent(in):: compute_vl = 1 1412 check(compute_vl==1||compute_vl==0) compute_vl 1413 integer optional,intent(in):: compute_vr = 1 1414 check(compute_vr==1||compute_vr==0) compute_vr 1415 1416 integer intent(hide),depend(a) :: n = shape(a,0) 1417 <ftype2> dimension(n,n),intent(in,copy) :: a 1418 check(shape(a,0)==shape(a,1)) :: a 1419 1420 <ftype2> intent(in,copy), dimension(n,n) :: b 1421 check(shape(b,0)==shape(b,1)) :: b 1422 1423 <ftype2> intent(out), dimension(n), depend(n) :: alphar 1424 <ftype2> intent(out), dimension(n), depend(n) :: alphai 1425 <ftype2> intent(out), dimension(n), depend(n) :: beta 1426 1427 <ftype2> depend(ldvl,n), dimension(ldvl,n),intent(out) :: vl 1428 integer intent(hide),depend(n,compute_vl) :: ldvl=(compute_vl?n:1) 1429 1430 <ftype2> depend(ldvr,n), dimension(ldvr,n),intent(out) :: vr 1431 integer intent(hide),depend(n,compute_vr) :: ldvr=(compute_vr?n:1) 1432 1433 integer optional,intent(in),depend(n,compute_vl,compute_vr) :: lwork=max(8*n,1) 1434 check((lwork==-1) || (lwork>=MAX(1,8*n))) :: lwork 1435 <ftype2> intent(out), dimension(MAX(lwork,1)), depend(lwork) :: work 1436 1437 integer intent(out):: info 1438 1439end subroutine <prefix2>ggev 1440 1441subroutine <prefix2c>ggev(compute_vl,compute_vr,n,a,b,alpha,beta,vl,ldvl,vr,ldvr,work,lwork,rwork,info) 1442 1443 callstatement {(*f2py_func)((compute_vl?"V":"N"),(compute_vr?"V":"N"),&n,a,&n,b,&n,alpha,beta,vl,&ldvl,vr,&ldvr,work,&lwork,rwork,&info);} 1444 callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* 1445 1446 integer optional,intent(in):: compute_vl = 1 1447 check(compute_vl==1||compute_vl==0) compute_vl 1448 integer optional,intent(in):: compute_vr = 1 1449 check(compute_vr==1||compute_vr==0) compute_vr 1450 1451 integer intent(hide),depend(a) :: n = shape(a,0) 1452 <ftype2c> dimension(n,n),intent(in,copy) :: a 1453 check(shape(a,0)==shape(a,1)) :: a 1454 1455 <ftype2c> intent(in,copy), dimension(n,n) :: b 1456 check(shape(b,0)==shape(b,1)) :: b 1457 1458 <ftype2c> intent(out), dimension(n), depend(n) :: alpha 1459 <ftype2c> intent(out), dimension(n), depend(n) :: beta 1460 1461 <ftype2c> depend(ldvl,n), dimension(ldvl,n),intent(out) :: vl 1462 integer intent(hide),depend(n,compute_vl) :: ldvl=(compute_vl?n:1) 1463 1464 <ftype2c> depend(ldvr,n), dimension(ldvr,n),intent(out) :: vr 1465 integer intent(hide),depend(n,compute_vr) :: ldvr=(compute_vr?n:1) 1466 1467 integer optional,intent(in),depend(n,compute_vl,compute_vr) :: lwork=max(2*n,1) 1468 check((lwork==-1) || (lwork>=MAX(1,2*n))) :: lwork 1469 <ftype2c> intent(out), dimension(MAX(lwork,1)), depend(lwork) :: work 1470 <ftype2> intent(hide), dimension(8*n), depend(n) :: rwork 1471 1472 integer intent(out):: info 1473 1474end subroutine <prefix2>ggev 1475 1476subroutine <prefix>geequ(m,n,a,lda,r,c,rowcnd,colcnd,amax,info) 1477 1478 callstatement (*f2py_func)(&m,&n,a,&lda,r,c,&rowcnd,&colcnd,&amax,&info) 1479 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctypereal>*,F_INT* 1480 1481 integer intent(hide),depend(a) :: m = shape(a,0) 1482 integer intent(hide),depend(a) :: n = shape(a,1) 1483 integer intent(hide),depend(a) :: lda = MAX(1,shape(a,0)) 1484 1485 <ftype> intent(in), dimension(m,n) :: a 1486 <ftypereal> intent(out),dimension(m),depend(m) :: r 1487 <ftypereal> intent(out),dimension(n),depend(n) :: c 1488 <ftypereal> intent(out) :: rowcnd 1489 <ftypereal> intent(out) :: colcnd 1490 <ftypereal> intent(out) :: amax 1491 integer intent(out) :: info 1492 1493end subroutine <prefix>geequ 1494 1495subroutine <prefix>geequb(m,n,a,lda,r,c,rowcnd,colcnd,amax,info) 1496 1497 callstatement (*f2py_func)(&m,&n,a,&lda,r,c,&rowcnd,&colcnd,&amax,&info) 1498 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctypereal>*,F_INT* 1499 1500 integer intent(hide),depend(a) :: m = shape(a,0) 1501 integer intent(hide),depend(a) :: n = shape(a,1) 1502 integer intent(hide),depend(a) :: lda = MAX(1,shape(a,0)) 1503 1504 <ftype> intent(in),dimension(m,n) :: a 1505 <ftypereal> intent(out),dimension(m),depend(m) :: r 1506 <ftypereal> intent(out),dimension(n),depend(n) :: c 1507 <ftypereal> intent(out) :: rowcnd 1508 <ftypereal> intent(out) :: colcnd 1509 <ftypereal> intent(out) :: amax 1510 integer intent(out) :: info 1511 1512end subroutine <prefix>geequb 1513