1! Signatures for f2py-wrappers of FORTRAN LAPACK Other Auxillary and Computational functions. 2! 3 4subroutine <prefix2>gejsv(joba,jobu,jobv,jobr,jobt,jobp,m,n,a,lda,sva,u,ldu,v,ldv,work,workout,lwork,iwork,iworkout,info) 5 ! ?GEJSV computes the singular value decomposition (SVD) of a complex 6 ! M-by-N matrix [A], where M >= N. The SVD of [A] is written as 7 ! 8 ! [A] = [U] * [SIGMA] * [V]^*, 9 ! 10 ! where [SIGMA] is an N-by-N (M-by-N) matrix which is zero except for its N 11 ! diagonal elements, [U] is an M-by-N (or M-by-M) unitary matrix, and 12 ! [V] is an N-by-N unitary matrix. The diagonal elements of [SIGMA] are 13 ! the singular values of [A]. The columns of [U] and [V] are the left and 14 ! the right singular vectors of [A], respectively. The matrices [U] and [V] 15 ! are computed and stored in the arrays U and V, respectively. The diagonal 16 ! of [SIGMA] is computed and stored in the array SVA. 17 callstatement {F_INT i;(*f2py_func)(&"CEFGAR"[joba],&"UFWN"[jobu],&"VJWN"[jobv],(jobr?"R":"N"),(jobt?"T":"N"),(jobp?"P":"N"),&m,&n,a,&lda,sva,u,&ldu,v,&ldv,work,&lwork,iwork,&info);for(i=0;i\<7;i++){workout[i] = work[i];}for(i=0;i\<3;i++){iworkout[i] = iwork[i];}} 18 callprotoargument char*,char*,char*,char*,char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 19 20 integer intent(in, optional), check((0 <= joba) && (joba < 6)) :: joba = 4 21 integer intent(in, optional), check((0 <= jobu) && (jobu < 4)) :: jobu = 0 22 integer intent(in, optional), check((0 <= jobv) && (jobv < 4) && ((jobv < 1) || (1 < jobv) || ((jobv == 1) && (jobu < 2)))), depend(jobu) :: jobv = 0 23 integer intent(in, optional), check((jobr == 0) || (jobr == 1)) :: jobr = 1 24 integer intent(in, optional), check((jobt == 0) || (jobt == 1)) :: jobt = 0 25 integer intent(in, optional), check((jobp == 0) || (jobp == 1)) :: jobp = 1 26 27 integer intent(hide), depend(a), check(m>=n), depend(n) :: m = shape(a, 0) 28 integer intent(hide), depend(a) :: n = shape(a, 1) 29 <ftype2> intent(in,copy), dimension(lda, n) :: a 30 integer intent(hide), depend(a) :: lda = max(1, shape(a, 0)) 31 <ftype2> intent(out), depend(n), dimension(n) :: sva 32 <ftype2> intent(out), depend(ldu, n, jobt, jobu, m), dimension(((jobt == 0)&&(jobu == 3)?0:m), ((jobt == 0)&&(jobu == 3)?0:(jobu == 1?m:n))) :: u 33 integer intent(hide), depend(m, jobu) :: ldu = max(1, (jobu < 3?m:1)) 34 <ftype2> intent(out), depend(ldv, n, jobt, jobv), dimension(((jobt == 0)&&(jobv == 3)?0:ldv),((jobt == 0)&&(jobv == 3)?0:n)) :: v 35 integer intent(hide), depend(n, jobv) :: ldv = max(1, (jobv < 3?n:1)) 36 <ftype2> intent(hide), depend(lwork), dimension(lwork) :: work 37 <ftype2> intent(out), dimension(7) :: workout 38 integer intent(in, optional), depend(m, n), check(lwork>=7) :: lwork = max(6*n+2*n*n, max(2*m+n, max(4*n+n*n, max(2*n+n*n+6, 7)))) 39 integer intent(hide), depend(m ,n), dimension(MAX(3, m+3*n)) :: iwork 40 integer intent(out), dimension(3) :: iworkout 41 integer intent(out) :: info 42 43end subroutine <prefix2>gejsv 44 45subroutine <prefix2>tgexc(wantq,wantz,n,a,lda,b,ldb,q,ldq,z,ldz,ifst,ilst,work,lwork,info) 46 ! Reorder the generalized Schur decomposition of a real matrix 47 ! pair using an orthogonal or unitary equivalence transformation. 48 49 callstatement { ifst++; ilst++; (*f2py_func)(&wantq,&wantz,&n,a,&lda,b,&ldb,q,&ldq,z,&ldz,&ifst,&ilst,work,&lwork,&info); } 50 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT* 51 52 integer intent(hide),check(wantq==0||wantq==1) :: wantq=1 53 integer intent(hide),check(wantz==0||wantz==1) :: wantz=1 54 integer intent(hide),depend(a) :: n=shape(a,1) 55 <ftype2> intent(in,out,copy),dimension(lda,n) :: a 56 integer intent(hide),depend(a) :: lda=shape(a,0) 57 <ftype2> intent(in,out,copy),dimension(ldb,n) :: b 58 integer intent(hide),depend(b) :: ldb=shape(b,0) 59 <ftype2> intent(in,out,copy),dimension(ldq,n) :: q 60 integer intent(hide),depend(q) :: ldq=shape(q,0) 61 <ftype2> intent(in,out,copy),dimension(ldz,n) :: z 62 integer intent(hide),depend(z) :: ldz=shape(z,0) 63 integer intent(in) :: ifst 64 integer intent(in) :: ilst 65 <ftype2> intent(out),dimension(MAX(lwork,1)) :: work 66 integer optional,intent(in),depend(n),check(lwork == -1 || lwork >= 4*n+16) :: lwork=max(4*n+16,1) 67 integer intent(out) :: info 68 69end subroutine <prefix2>tgexc 70 71subroutine <prefix2c>tgexc(wantq,wantz,n,a,lda,b,ldb,q,ldq,z,ldz,ifst,ilst,info) 72 ! Reorder the generalized Schur decomposition of a complex matrix 73 ! pair using an orthogonal or unitary equivalence transformation. 74 75 callstatement { ifst++; ilst++; (*f2py_func)(&wantq,&wantz,&n,a,&lda,b,&ldb,q,&ldq,z,&ldz,&ifst,&ilst,&info); } 76 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,F_INT*,F_INT* 77 78 integer intent(hide),check(wantq==0||wantq==1) :: wantq=1 79 integer intent(hide),check(wantz==0||wantz==1) :: wantz=1 80 integer intent(hide),depend(a) :: n=shape(a,1) 81 <ftype2c> intent(in,out,copy),dimension(lda,n) :: a 82 integer intent(hide),depend(a) :: lda=shape(a,0) 83 <ftype2c> intent(in,out,copy),dimension(ldb,n) :: b 84 integer intent(hide),depend(b) :: ldb=shape(b,0) 85 <ftype2c> intent(in,out,copy),dimension(ldq,n) :: q 86 integer intent(hide),depend(q) :: ldq=shape(q,0) 87 <ftype2c> intent(in,out,copy),dimension(ldz,n) :: z 88 integer intent(hide),depend(z) :: ldz=shape(z,0) 89 integer intent(in) :: ifst 90 integer intent(in) :: ilst 91 integer intent(out) :: info 92 93end subroutine <prefix2c>tgexc 94 95subroutine <prefix2>tgsen(ijob,wantq,wantz,select,n,a,lda,b,ldb,alphar,alphai,beta,q,ldq,z,ldz,m,pl,pr,dif,work,lwork,iwork,liwork,info) 96 97 callstatement (*f2py_func)(&ijob,&wantq,&wantz,select,&n,a,&lda,b,&ldb,alphar,alphai,beta,q,&ldq,z,&ldz,&m,&pl,&pr,dif,work,&lwork,iwork,&liwork,&info) 98 callprotoargument F_INT*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 99 100 integer optional, intent(in),check(ijob>=0&&ijob<=5) :: ijob=4 101 logical optional, intent(in),check(wantq==0||wantq==1) :: wantq=1 102 logical optional, intent(in),check(wantz==0||wantz==1) :: wantz=1 103 logical intent(in),dimension(n),depend(n) :: select 104 integer intent(hide),depend(a) :: n=shape(a,0) 105 <ftype2> intent(in,out,copy,out=as),dimension(n,n) :: a 106 integer intent(hide),depend(a) :: lda=MAX(1,shape(a,0)) 107 <ftype2> intent(in,out,copy,out=bs),dimension(n,n) :: b 108 integer intent(hide),depend(b) :: ldb = MAX(1, shape(b,0)) 109 <ftype2> intent(out),dimension(n),depend(n) :: alphar 110 <ftype2> intent(out),dimension(n),depend(n) :: alphai 111 <ftype2> intent(out),dimension(n),depend(n) :: beta 112 <ftype2> intent(in,out,copy,out=qs),dimension(n,n),depend(n) :: q 113 integer intent(hide),depend(q) :: ldq=MAX(1,shape(q,0)) 114 <ftype2> intent(in,out,copy,out=zs),dimension(n,n),depend(n) :: z 115 integer intent(hide),depend(z) :: ldz=MAX(1,shape(z,0)) 116 integer intent(out) :: m 117 <ftype2> intent(out) :: pl 118 <ftype2> intent(out) :: pr 119 <ftype2> intent(out),dimension(2) :: dif 120 <ftype2> intent(hide),dimension(MAX(lwork,1)) :: work 121 ! these lwork and liwork values are bare minimum estiamates only for cases ijob == 1,2,4 122 ! a separate lwork query is a prerequisite due to m dependence. 123 ! Depending on the m value lwork can go upto n**2 / 4 for n = 2m hence it is not 124 ! possible to give a minimal value without a potential excessive memory waste 125 integer optional,intent(in),depend(n),check(lwork == -1 || lwork >= 1) :: lwork=4*n+16 126 integer intent(hide),dimension(MAX(1,liwork)) :: iwork 127 integer optional,intent(in),depend(n),check(liwork == -1 || liwork >= 1) :: liwork=n+6 128 integer intent(out) :: info 129 130end subroutine <prefix2>tgsen 131 132subroutine <prefix2>tgsen_lwork(ijob,wantq,wantz,select,n,a,lda,b,ldb,alphar,alphai,beta,q,ldq,z,ldz,m,pl,pr,dif,work,lwork,iwork,liwork,info) 133 134 fortranname <prefix2>tgsen 135 callstatement (*f2py_func)(&ijob,&wantq,&wantz,select,&n,a,&lda,&b,&ldb,&alphar,&alphai,&beta,&q,&ldq,&z,&ldz,&m,&pl,&pr,&dif,&work,&lwork,&iwork,&liwork,&info) 136 callprotoargument F_INT*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 137 138 integer optional, intent(in),check(ijob>=0&&ijob<=5) :: ijob = 4 139 logical intent(in),dimension(n),depend(n) :: select 140 <ftype2> intent(in),dimension(n,n) :: a 141 142 logical intent(hide) :: wantq = 0 143 logical intent(hide) :: wantz = 0 144 integer intent(hide),depend(a) :: n = shape(a,0) 145 integer intent(hide),depend(n) :: lda = MAX(1, n) 146 <ftype2> intent(hide) :: b 147 integer intent(hide),depend(n) :: ldb = MAX(1, n) 148 <ftype2> intent(hide) :: alphar 149 <ftype2> intent(hide) :: alphai 150 <ftype2> intent(hide) :: beta 151 <ftype2> intent(hide) :: q 152 integer intent(hide),depend(n) :: ldq = MAX(1, n) 153 <ftype2> intent(hide) :: z 154 integer intent(hide),depend(n) :: ldz = MAX(1, n) 155 integer intent(hide) :: m 156 <ftype2> intent(hide) :: pl 157 <ftype2> intent(hide) :: pr 158 <ftype2> intent(hide) :: dif 159 integer intent(hide):: lwork=-1 160 integer intent(hide) :: liwork=-1 161 162 <ftype2> intent(out) :: work 163 integer intent(out) :: iwork 164 integer intent(out) :: info 165 166end subroutine <prefix2>tgsen_lwork 167 168subroutine <prefix2c>tgsen(ijob,wantq,wantz,select,n,a,lda,b,ldb,alpha,beta,q,ldq,z,ldz,m,pl,pr,dif,work,lwork,iwork,liwork,info) 169 callstatement (*f2py_func)(&ijob,&wantq,&wantz,select,&n,a,&lda,b,&ldb,alpha,beta,q,&ldq,z,&ldz,&m,&pl,&pr,dif,work,&lwork,iwork,&liwork,&info) 170 callprotoargument F_INT*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT*,F_INT*,F_INT*,F_INT* 171 172 integer optional, intent(in),check(ijob>=0&&ijob<=5) :: ijob=4 173 logical optional, intent(in),check(wantq==0||wantq==1) :: wantq=1 174 logical optional, intent(in),check(wantz==0||wantz==1) :: wantz=1 175 logical intent(in),dimension(n),depend(n) :: select 176 integer intent(hide),depend(a) :: n=shape(a,0) 177 <ftype2c> intent(in,out,copy,out=as),dimension(n,n) :: a 178 integer intent(hide),depend(a) :: lda=MAX(1,shape(a,0)) 179 <ftype2c> intent(in,out,copy,out=bs),dimension(n,n) :: b 180 integer intent(hide),depend(b) :: ldb = MAX(1, shape(b,0)) 181 <ftype2c> intent(out),dimension(n),depend(n) :: alpha 182 <ftype2c> intent(out),dimension(n),depend(n) :: beta 183 <ftype2c> intent(in,out,copy,out=qs),dimension(n,n),depend(n) :: q 184 integer intent(hide),depend(q) :: ldq=MAX(1,shape(q,0)) 185 <ftype2c> intent(in,out,copy,out=zs),dimension(n,n),depend(n) :: z 186 integer intent(hide),depend(z) :: ldz=MAX(1,shape(z,0)) 187 integer intent(out) :: m 188 <ftype2> intent(out) :: pl 189 <ftype2> intent(out) :: pr 190 <ftype2> intent(out),dimension(2) :: dif 191 <ftype2c> intent(hide),dimension(MAX(lwork,1)) :: work 192 ! these lwork and liwork values are bare minimum estiamates only for cases ijob ==0,1,2,4 193 ! a separate lwork query is a prerequisite due to m dependence. 194 ! Depending on the m value lwork can go upto n**2 / 4 for n = 2m hence it is not 195 ! possible to give a minimal value without a potential excessive memory waste 196 integer optional,intent(in),depend(n,ijob),check(lwork == -1 || lwork >= 1) :: lwork=(ijob==0?1:n+2) 197 integer intent(hide),dimension(liwork):: iwork 198 integer optional,intent(in),depend(n,ijob),check(liwork == -1 || liwork>=1) :: liwork=(ijob==0?1:n+2) 199 integer intent(out) :: info 200 201end subroutine <prefix2c>tgsen 202 203subroutine <prefix2c>tgsen_lwork(ijob,wantq,wantz,select,n,a,lda,b,ldb,alpha,beta,q,ldq,z,ldz,m,pl,pr,dif,work,lwork,iwork,liwork,info) 204 205 fortranname <prefix2c>tgsen 206 callstatement (*f2py_func)(&ijob,&wantq,&wantz,select,&n,a,&lda,b,&ldb,alpha,beta,&q,&ldq,&z,&ldz,&m,&pl,&pr,&dif,&work,&lwork,&iwork,&liwork,&info) 207 callprotoargument F_INT*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT*,F_INT*,F_INT*,F_INT* 208 209 integer optional, intent(in),check(ijob>=0&&ijob<=5) :: ijob = 4 210 logical intent(in),dimension(n),depend(n) :: select 211 <ftype2c> intent(in),dimension(n,n) :: a 212 <ftype2c> intent(in),dimension(n,n) :: b 213 214 logical intent(hide) :: wantq = 0 215 logical intent(hide) :: wantz = 0 216 integer intent(hide),depend(a) :: n = shape(a,0) 217 integer intent(hide),depend(n) :: lda = MAX(1, n) 218 integer intent(hide),depend(n) :: ldb = MAX(1, n) 219 <ftype2c> intent(hide),dimension(n) :: alpha 220 <ftype2c> intent(hide),dimension(n) :: beta 221 <ftype2c> intent(hide) :: q 222 integer intent(hide),depend(n) :: ldq = MAX(1, n) 223 <ftype2c> intent(hide) :: z 224 integer intent(hide),depend(n) :: ldz = MAX(1, n) 225 integer intent(hide) :: m 226 <ftype2> intent(hide) :: pl 227 <ftype2> intent(hide) :: pr 228 <ftype2> intent(hide) :: dif 229 integer intent(hide):: lwork=-1 230 integer intent(hide) :: liwork=-1 231 232 <ftype2c> intent(out) :: work 233 integer intent(out) :: iwork 234 integer intent(out) :: info 235 236end subroutine <prefix2c>tgsen_lwork 237 238subroutine <prefix2>pbtrf(lower,n,kd,ab,ldab,info) 239 ! Compute Cholesky decomposition of banded symmetric positive definite 240 ! matrix: 241 ! A = U^T * U, C = U if lower = 0 242 ! A = L * L^T, C = L if lower = 1 243 ! C is triangular matrix of the corresponding Cholesky decomposition. 244 245 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,ab,&ldab,&info); 246 callprotoargument char*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT* 247 248 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 249 integer intent(hide),depend(ab) :: n=shape(ab,1) 250 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 251 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 252 253 <ftype2> dimension(ldab,n),intent(in,out,copy,out=c) :: ab 254 integer intent(out) :: info 255 256end subroutine <prefix2>pbtrf 257 258subroutine <prefix2c>pbtrf(lower,n,kd,ab,ldab,info) 259 ! Compute Cholesky decomposition of banded symmetric positive definite 260 ! matrix: 261 ! A = U^H * U, C = U if lower = 0 262 ! A = L * L^H, C = L if lower = 1 263 ! C is triangular matrix of the corresponding Cholesky decomposition. 264 265 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,ab,&ldab,&info); 266 callprotoargument char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT* 267 268 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 269 integer intent(hide),depend(ab) :: n=shape(ab,1) 270 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 271 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 272 273 <ftype2c> dimension(ldab,n),intent(in,out,copy,out=c) :: ab 274 integer intent(out) :: info 275 276end subroutine <prefix2c>pbtrf 277 278subroutine <prefix2>pbtrs(lower, n, kd, nrhs, ab, ldab, b, ldb, info) 279 280 ! Solve a system of linear equations A*X = B with a symmetric 281 ! positive definite band matrix A using the Cholesky factorization. 282 ! AB is the triangular factur U or L from the Cholesky factorization 283 ! previously computed with *PBTRF. 284 ! A = U^T * U, AB = U if lower = 0 285 ! A = L * L^T, AB = L if lower = 1 286 287 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,&nrhs,ab,&ldab,b,&ldb,&info); 288 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 289 290 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 291 integer intent(hide),depend(ab) :: n=shape(ab,1) 292 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 293 integer intent(hide),depend(b) :: ldb=shape(b,0) 294 integer intent(hide),depend(b) :: nrhs=shape(b,1) 295 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 296 297 <ftype2> dimension(ldb, nrhs),intent(in,out,copy,out=x) :: b 298 <ftype2> dimension(ldab,n),intent(in) :: ab 299 integer intent(out) :: info 300 301end subroutine <tchar=s,d>pbtrs 302 303subroutine <prefix2c>pbtrs(lower, n, kd, nrhs, ab, ldab, b, ldb, info) 304 ! Solve a system of linear equations A*X = B with a symmetric 305 ! positive definite band matrix A using the Cholesky factorization. 306 ! AB is the triangular factur U or L from the Cholesky factorization 307 ! previously computed with *PBTRF. 308 ! A = U^T * U, AB = U if lower = 0 309 ! A = L * L^T, AB = L if lower = 1 310 311 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,&nrhs,ab,&ldab,b,&ldb,&info); 312 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* 313 314 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 315 integer intent(hide),depend(ab) :: n=shape(ab,1) 316 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 317 integer intent(hide),depend(b) :: ldb=shape(b,0) 318 integer intent(hide),depend(b) :: nrhs=shape(b,1) 319 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 320 321 <ftype2c> dimension(ldb, nrhs),intent(in,out,copy,out=x) :: b 322 <ftype2c> dimension(ldab,n),intent(in) :: ab 323 integer intent(out) :: info 324 325end subroutine <prefix2c>pbtrs 326 327subroutine <prefix>trtrs(lower, trans, unitdiag, n, nrhs, a, lda, b, ldb, info) 328 329 ! Solve a system of linear equations A*X = B with a triangular 330 ! matrix A. 331 332 callstatement (*f2py_func)((lower?"L":"U"),(trans?(trans==2?"C":"T"):"N"),(unitdiag?"U":"N"),&n,&nrhs,a,&lda,b,&ldb,&info); 333 callprotoargument char*,char*,char*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 334 335 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 336 integer optional,intent(in),check(trans>=0 && trans <=2) :: trans = 0 337 integer optional,intent(in),check(unitdiag==0||unitdiag==1) :: unitdiag = 0 338 integer optional,check(shape(a,0)==lda),depend(a) :: lda=shape(a,0) 339 integer intent(hide),depend(a) :: n=shape(a,1) 340 integer intent(hide),depend(b) :: ldb=shape(b,0) 341 integer intent(hide),depend(b) :: nrhs=shape(b,1) 342 343 <ftype> dimension(ldb, nrhs),intent(in,out,copy,out=x) :: b 344 <ftype> dimension(lda,n),intent(in) :: a 345 integer intent(out) :: info 346 347end subroutine <prefix>trtrs 348 349 350subroutine <prefix>tbtrs(uplo,trans,diag,n,kd,nrhs,ab,ldab,b,ldb,info) 351 ! x, info = tbtrs(ab, b, uplo="U", trans="N", diag="N", overwrite_b=0) 352 ! ?TBTRS solves a triangular system of the form 353 ! 354 ! A * X = B or A**T * X = B, 355 ! 356 ! where: 357 ! * A is a triangular band matrix of order N, 358 ! * B is an N-by NRHS matrix. 359 ! A check is made to verify that A is nonsingular. 360 361 callstatement (*f2py_func)(uplo,trans,diag,&n,&kd,&nrhs,ab,&ldab,b,&ldb,&info) 362 callprotoargument char*,char*,char*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 363 364 <ftype> intent(in), dimension(ldab, n) :: ab 365 <ftype> intent(in,out,copy,out=x), dimension(ldb, nrhs) :: b 366 integer intent(out) :: info 367 368 character optional intent(in), check(*uplo=='U'||*uplo=='L') :: uplo = 'U' 369 character optional intent(in), check(*trans=='N'||*trans=='T'||*trans=='C') :: trans = 'N' 370 character optional intent(in), check(*diag=='N'||*diag=='U') :: diag = 'N' 371 372 integer intent(hide), depend(ab) :: ldab = MAX(1, shape(ab, 0)) 373 integer intent(hide), depend(ab) :: n = MAX(1, shape(ab, 1)) 374 integer intent(hide), depend(ldab) :: kd = ldab - 1 375 integer intent(hide), depend(b, n), check(ldb >= n) :: ldb = MAX(1, shape(b, 0)) 376 integer intent(hide), depend(b) :: nrhs = MAX(1, shape(b, 1)) 377 378end subroutine <prefix>tbtrs 379 380subroutine <prefix2>pbsv(lower,n,kd,nrhs,ab,ldab,b,ldb,info) 381 ! 382 ! Computes the solution to a real system of linear equations 383 ! A * X = B, 384 ! where A is an N-by-N symmetric positive definite band matrix and X 385 ! and B are N-by-NRHS matrices. 386 ! 387 ! The Cholesky decomposition is used to factor A as 388 ! A = U**T * U, if lower=1, or 389 ! A = L * L**T, if lower=0 390 ! where U is an upper triangular band matrix, and L is a lower 391 ! triangular band matrix, with the same number of superdiagonals or 392 ! subdiagonals as A. The factored form of A is then used to solve the 393 ! system of equations A * X = B. 394 395 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,&nrhs,ab,&ldab,b,&ldb,&info); 396 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 397 398 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 399 integer intent(hide),depend(ab) :: n=shape(ab,1) 400 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 401 integer intent(hide),depend(b) :: ldb=shape(b,0) 402 integer intent(hide),depend(b) :: nrhs=shape(b,1) 403 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 404 405 <ftype2> dimension(ldb, nrhs),intent(in,out,copy,out=x) :: b 406 <ftype2> dimension(ldab,n),intent(in,out,copy,out=c) :: ab 407 integer intent(out) :: info 408 409end subroutine <prefix2>pbsv 410 411subroutine <prefix2c>pbsv(lower,n,kd,nrhs,ab,ldab,b,ldb,info) 412 ! Computes the solution to a real system of linear equations 413 ! A * X = B, 414 ! where A is an N-by-N Hermitian positive definite band matrix and X 415 ! and B are N-by-NRHS matrices. 416 ! 417 ! The Cholesky decomposition is used to factor A as 418 ! A = U**H * U, if lower=1, or 419 ! A = L * L**H, if lower=0 420 ! where U is an upper triangular band matrix, and L is a lower 421 ! triangular band matrix, with the same number of superdiagonals or 422 ! subdiagonals as A. The factored form of A is then used to solve the 423 ! system of equations A * X = B. 424 425 callstatement (*f2py_func)((lower?"L":"U"),&n,&kd,&nrhs,ab,&ldab,b,&ldb,&info); 426 callprotoargument char*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* 427 428 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 429 integer intent(hide),depend(ab) :: n=shape(ab,1) 430 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 431 integer intent(hide),depend(b) :: ldb=shape(b,0) 432 integer intent(hide),depend(b) :: nrhs=shape(b,1) 433 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 434 435 <ftype2c> dimension(ldb, nrhs),intent(in,out,copy,out=x) :: b 436 <ftype2c> dimension(ldab,n),intent(in,out,copy,out=c) :: ab 437 integer intent(out) :: info 438 439end subroutine <prefix2c>pbsv 440 441subroutine <prefix2>orcsd(compute_u1,compute_u2,compute_v1t,compute_v2t,trans,signs,m,p,q,x11,ldx11,x12,ldx12,x21,ldx21,x22,ldx22,theta,u1,ldu1,u2,ldu2,v1t,ldv1t,v2t,ldv2t,work,lwork,iwork,info,mmp,mmq) 442 ! DORCSD computes the CS decomposition of an M-by-M partitioned 443 ! unitary matrix X: 444 ! 445 ! [ I11 0 0 | 0 0 0 ] 446 ! [ 0 C 0 | 0 -S 0 ] 447 ! [ X11 | X12 ] [ U1 | ] [ 0 0 0 | 0 0 -I12 ] [ V1 | ]**T 448 ! [-----------] = [---------] [---------------------------] [---------] . 449 ! [ X21 | X22 ] [ | U2 ] [ 0 0 0 | I22 0 0 ] [ | V2 ] 450 ! [ 0 S 0 | 0 C 0 ] 451 ! [ 0 0 I21 | 0 0 0 ] 452 ! 453 ! U1, U2, V1, V2 are square orthogonal matrices of 454 ! dimensions (p,p), (m-p,m-p), (q,q), (m-q,m-q), respectively, 455 ! and C and S are (r, r) nonnegative diagonal matrices satisfying 456 ! C^2 + S^2 = I where r = min(p, m-p, q, m-q). 457 ! 458 ! Moreover, the rank of the identity matrices are min(p, q) - r, 459 ! min(p, m - q) - r, min(m - p, q) - r, and min(m - p, m - q) - r 460 ! respectively. 461 462 callstatement (*f2py_func)((compute_u1?"Y":"N"),(compute_u2?"Y":"N"),(compute_v1t?"Y":"N"),(compute_v2t?"Y":"N"),(trans?"T":"N"),(signs?"O":"D"),&m,&p,&q,x11,&ldx11,x12,&ldx12,x21,&ldx21,x22,&ldx22,theta,u1,&ldu1,u2,&ldu2,v1t,&ldv1t,v2t,&ldv2t,work,&lwork,iwork,&info) 463 callprotoargument char*,char*,char*,char*,char*,char*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 464 465 integer optional,intent(in),check(compute_u1==0||compute_u1==1) :: compute_u1 = 1 466 integer optional,intent(in),check(compute_u2==0||compute_u2==1) :: compute_u2 = 1 467 integer optional,intent(in),check(compute_v1t==0||compute_v1t==1) :: compute_v1t = 1 468 integer optional,intent(in),check(compute_v2t==0||compute_v2t==1) :: compute_v2t = 1 469 integer optional,intent(in),check(trans==0||trans==1) :: trans = 0 470 integer optional,intent(in),check(signs==0||signs==1) :: signs = 0 471 472 integer intent(hide),check(p+mmp==q+mmq),depend(p,q,mmp,mmq) :: m = p + mmp 473 474 <ftype2> intent(in,out,copy,out=cs11),dimension(p,q) :: x11 475 <ftype2> intent(in,out,copy,out=cs22),dimension(mmp,mmq) :: x22 476 <ftype2> intent(in,out,copy,out=cs12),dimension(p,mmq),check(mmq==shape(x12,1)||p==shape(x12,0)),depend(p,mmq) :: x12 477 <ftype2> intent(in,out,copy,out=cs21),dimension(mmp,q),check(mmp==shape(x21,0)||q==shape(x21,1)),depend(mmp,q) :: x21 478 479 integer intent(hide),depend(x11) :: p = shape(x11,0) 480 integer intent(hide),depend(x11) :: q = shape(x11,1) 481 integer intent(hide),depend(x22) :: mmp = shape(x22,0) 482 integer intent(hide),depend(x22) :: mmq = shape(x22,1) 483 484 integer intent(hide),depend(x11) :: ldx11 = MAX(1,shape(x11,0)) 485 integer intent(hide),depend(x12) :: ldx12 = MAX(1,shape(x12,0)) 486 integer intent(hide),depend(x21) :: ldx21 = MAX(1,shape(x21,0)) 487 integer intent(hide),depend(x22) :: ldx22 = MAX(1,shape(x22,0)) 488 489 <ftype2> intent(out),dimension(min(min(p,mmp),min(q,mmq))),depend(p,q,mmp,mmq) :: theta 490 <ftype2> intent(out),dimension((compute_u1?p:0),(compute_u1?p:0)),depend(p) :: u1 491 <ftype2> intent(out),dimension((compute_u2?mmp:0),(compute_u2?mmp:0)),depend(mmp) :: u2 492 <ftype2> intent(out),dimension((compute_v1t?q:0),(compute_v1t?q:0)),depend(q) :: v1t 493 <ftype2> intent(out),dimension((compute_v2t?mmq:0),(compute_v2t?mmq:0)),depend(mmq) :: v2t 494 495 integer intent(hide),depend(p) :: ldu1 = MAX(1,p) 496 integer intent(hide),depend(mmp) :: ldu2 = MAX(1,mmp) 497 integer intent(hide),depend(q) :: ldv1t = MAX(1,q) 498 integer intent(hide),depend(mmq) :: ldv2t = MAX(1,mmq) 499 500 <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work 501 integer intent(hide),dimension(p+mmp-MIN(MIN(p,mmp),MIN(q,mmq))),depend(p,q,mmp,mmq) :: iwork 502 503 integer optional,intent(in),check(lwork==-1||lwork>0),depend(m,mmp,mmq) :: lwork = 2 + 2*m + 5*MAX(1,q-1) + 4*MAX(1,q) + 8*q 504 integer intent(out) :: info 505 506end subroutine <prefix2>orcsd 507 508subroutine <prefix2>orcsd_lwork(m,p,q,x11,ldx11,x12,ldx12,x21,ldx21,x22,ldx22,theta,u1,ldu1,u2,ldu2,v1t,ldv1t,v2t,ldv2t,work,lwork,iwork,info,mmp,mmq) 509 ! LWORK computation for (S/D)ORCSD 510 511 fortranname <prefix2>orcsd 512 callstatement (*f2py_func)("Y","Y","Y","Y","N","D",&m,&p,&q,&x11,&ldx11,&x12,&ldx12,&x21,&ldx21,&x22,&ldx22,&theta,&u1,&ldu1,&u2,&ldu2,&v1t,&ldv1t,&v2t,&ldv2t,&work,&lwork,&iwork,&info) 513 callprotoargument char*,char*,char*,char*,char*,char*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 514 515 integer intent(in) :: m 516 integer intent(in) :: p 517 integer intent(in) :: q 518 519 <ftype2> intent(hide) :: x11 520 <ftype2> intent(hide) :: x22 521 <ftype2> intent(hide) :: x12 522 <ftype2> intent(hide) :: x21 523 integer intent(hide),depend(m,p) :: mmp = m - p 524 integer intent(hide),depend(m,q) :: mmq = m - q 525 integer intent(hide),depend(p) :: ldx11 = MAX(1,p) 526 integer intent(hide),depend(p) :: ldx12 = MAX(1,p) 527 integer intent(hide),depend(mmp) :: ldx21 = MAX(1,mmp) 528 integer intent(hide),depend(mmp) :: ldx22 = MAX(1,mmp) 529 <ftype2> intent(hide) :: theta 530 <ftype2> intent(hide) :: u1 531 <ftype2> intent(hide) :: u2 532 <ftype2> intent(hide) :: v1t 533 <ftype2> intent(hide) :: v2t 534 integer intent(hide),depend(p) :: ldu1 = MAX(1,p) 535 integer intent(hide),depend(mmp) :: ldu2 = MAX(1,mmp) 536 integer intent(hide),depend(q) :: ldv1t = MAX(1,q) 537 integer intent(hide),depend(mmq) :: ldv2t = MAX(1,mmq) 538 integer intent(hide) :: iwork 539 integer intent(hide) :: lwork = -1 540 541 <ftype2> intent(out) :: work 542 integer intent(out) :: info 543 544end subroutine <prefix2>orcsd_lwork 545 546subroutine <prefix2c>uncsd(compute_u1,compute_u2,compute_v1t,compute_v2t,trans,signs,m,p,q,x11,ldx11,x12,ldx12,x21,ldx21,x22,ldx22,theta,u1,ldu1,u2,ldu2,v1t,ldv1t,v2t,ldv2t,work,lwork,rwork,lrwork,iwork,info,mmp,mmq) 547 ! ZUNCSD computes the CS decomposition of an M-by-M partitioned 548 ! unitary matrix X: 549 ! 550 ! [ I11 0 0 | 0 0 0 ] 551 ! [ 0 C 0 | 0 -S 0 ] 552 ! [ X11 | X12 ] [ U1 | ] [ 0 0 0 | 0 0 -I12 ] [ V1 | ]* 553 ! [-----------] = [---------] [---------------------------] [---------] . 554 ! [ X21 | X22 ] [ | U2 ] [ 0 0 0 | I22 0 0 ] [ | V2 ] 555 ! [ 0 S 0 | 0 C 0 ] 556 ! [ 0 0 I21 | 0 0 0 ] 557 ! 558 ! U1, U2, V1, V2 are square orthogonal matrices of 559 ! dimensions (p,p), (m-p,m-p), (q,q), (m-q,m-q), respectively, 560 ! and C and S are (r, r) nonnegative diagonal matrices satisfying 561 ! C^2 + S^2 = I where r = min(p, m-p, q, m-q). 562 ! 563 ! Moreover, the rank of the identity matrices are min(p, q) - r, 564 ! min(p, m - q) - r, min(m - p, q) - r, and min(m - p, m - q) - r 565 ! respectively. 566 567 568 callstatement (*f2py_func)((compute_u1?"Y":"N"),(compute_u2?"Y":"N"),(compute_v1t?"Y":"N"),(compute_v2t?"Y":"N"),(trans?"T":"N"),(signs?"O":"D"),&m,&p,&q,x11,&ldx11,x12,&ldx12,x21,&ldx21,x22,&ldx22,theta,u1,&ldu1,u2,&ldu2,v1t,&ldv1t,v2t,&ldv2t,work,&lwork,rwork,&lrwork,iwork,&info) 569 callprotoargument char*,char*,char*,char*,char*,char*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 570 571 integer optional,intent(in),check(compute_u1==0||compute_u1==1) :: compute_u1 = 1 572 integer optional,intent(in),check(compute_u2==0||compute_u2==1) :: compute_u2 = 1 573 integer optional,intent(in),check(compute_v1t==0||compute_v1t==1) :: compute_v1t = 1 574 integer optional,intent(in),check(compute_v2t==0||compute_v2t==1) :: compute_v2t = 1 575 integer optional,intent(in),check(trans==0||trans==1) :: trans = 0 576 integer optional,intent(in),check(signs==0||signs==1) :: signs = 0 577 578 integer intent(hide),check(p+mmp==q+mmq),depend(p,q,mmp,mmq) :: m = p + mmp 579 580 <ftype2c> intent(in,out,copy,out=cs11),dimension(p,q) :: x11 581 <ftype2c> intent(in,out,copy,out=cs22),dimension(mmp,mmq) :: x22 582 <ftype2c> intent(in,out,copy,out=cs12),dimension(p,mmq),check(mmq==shape(x12,1)||p==shape(x12,0)),depend(p,mmq) :: x12 583 <ftype2c> intent(in,out,copy,out=cs21),dimension(mmp,q),check(mmp==shape(x21,0)||q==shape(x21,1)),depend(mmp,q) :: x21 584 585 integer intent(hide),depend(x11) :: p = shape(x11,0) 586 integer intent(hide),depend(x11) :: q = shape(x11,1) 587 integer intent(hide),depend(x22) :: mmp = shape(x22,0) 588 integer intent(hide),depend(x22) :: mmq = shape(x22,1) 589 590 integer intent(hide),depend(x11) :: ldx11 = MAX(1,shape(x11,0)) 591 integer intent(hide),depend(x12) :: ldx12 = MAX(1,shape(x12,0)) 592 integer intent(hide),depend(x21) :: ldx21 = MAX(1,shape(x21,0)) 593 integer intent(hide),depend(x22) :: ldx22 = MAX(1,shape(x22,0)) 594 595 <ftype2> intent(out),dimension(min(min(p,mmp),min(q,mmq))),depend(p,q,mmp,mmq) :: theta 596 <ftype2c> intent(out),dimension((compute_u1?p:0),(compute_u1?p:0)),depend(p) :: u1 597 <ftype2c> intent(out),dimension((compute_u2?mmp:0),(compute_u2?mmp:0)),depend(mmp) :: u2 598 <ftype2c> intent(out),dimension((compute_v1t?q:0),(compute_v1t?q:0)),depend(q) :: v1t 599 <ftype2c> intent(out),dimension((compute_v2t?mmq:0),(compute_v2t?mmq:0)),depend(mmq) :: v2t 600 601 integer intent(hide),depend(p) :: ldu1 = MAX(1,p) 602 integer intent(hide),depend(mmp) :: ldu2 = MAX(1,mmp) 603 integer intent(hide),depend(q) :: ldv1t = MAX(1,q) 604 integer intent(hide),depend(mmq) :: ldv2t = MAX(1,mmq) 605 606 <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work 607 <ftype2> intent(hide),dimension(lrwork),depend(lrwork) :: rwork 608 integer intent(hide),dimension(p+mmp-MIN(MIN(p,mmp),MIN(q,mmq))),depend(p,q,mmp,mmq) :: iwork 609 610 integer optional,intent(in),check(lwork==-1||lwork>0),depend(m,mmp,mmq) :: lwork = 2*m + MAX(1,MAX(mmp,mmq)) + 1 611 integer optional,intent(in),check(lrwork==-1||lrwork>0),depend(q) :: lrwork = 5*MAX(1,q-1) + 4*MAX(1,q) + 8*q + 1 612 integer intent(out) :: info 613 614end subroutine <prefix2c>uncsd 615 616subroutine <prefix2c>uncsd_lwork(m,p,q,x11,ldx11,x12,ldx12,x21,ldx21,x22,ldx22,theta,u1,ldu1,u2,ldu2,v1t,ldv1t,v2t,ldv2t,work,lwork,rwork,lrwork,iwork,info,mmp,mmq) 617 ! LWORK computation for (C/Z)UNCSD 618 619 fortranname <prefix2c>uncsd 620 callstatement (*f2py_func)("Y","Y","Y","Y","N","D",&m,&p,&q,&x11,&ldx11,&x12,&ldx12,&x21,&ldx21,&x22,&ldx22,&theta,&u1,&ldu1,&u2,&ldu2,&v1t,&ldv1t,&v2t,&ldv2t,&work,&lwork,&rwork,&lrwork,&iwork,&info) 621 callprotoargument char*,char*,char*,char*,char*,char*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 622 623 integer intent(in) :: m 624 integer intent(in) :: p 625 integer intent(in) :: q 626 627 <ftype2c> intent(hide) :: x11 628 <ftype2c> intent(hide) :: x22 629 <ftype2c> intent(hide) :: x12 630 <ftype2c> intent(hide) :: x21 631 integer intent(hide),depend(m,p) :: mmp = m - p 632 integer intent(hide),depend(m,q) :: mmq = m - q 633 integer intent(hide),depend(p) :: ldx11 = MAX(1,p) 634 integer intent(hide),depend(p) :: ldx12 = MAX(1,p) 635 integer intent(hide),depend(mmp) :: ldx21 = MAX(1,mmp) 636 integer intent(hide),depend(mmp) :: ldx22 = MAX(1,mmp) 637 <ftype2> intent(hide) :: theta 638 <ftype2c> intent(hide) :: u1 639 <ftype2c> intent(hide) :: u2 640 <ftype2c> intent(hide) :: v1t 641 <ftype2c> intent(hide) :: v2t 642 integer intent(hide),depend(p) :: ldu1 = MAX(1,p) 643 integer intent(hide),depend(mmp) :: ldu2 = MAX(1,mmp) 644 integer intent(hide),depend(q) :: ldv1t = MAX(1,q) 645 integer intent(hide),depend(mmq) :: ldv2t = MAX(1,mmq) 646 integer intent(hide) :: iwork 647 integer intent(hide) :: lwork = -1 648 integer intent(hide) :: lrwork = -1 649 650 <ftype2c> intent(out) :: work 651 <ftype2> intent(out) :: rwork 652 integer intent(out) :: info 653 654end subroutine <prefix2c>uncsd_lwork 655 656subroutine <prefix2>orghr(n,lo,hi,a,tau,work,lwork,info) 657 ! 658 ! q,info = orghr(a,tau,lo=0,hi=n-1,lwork=n,overwrite_a=0) 659 ! Compute orthogonal matrix Q for Hessenberg reduction from the matrix 660 ! that was computed by gehrd 661 ! 662 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,a,&n,tau,work,&lwork,&info); } 663 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 664 integer intent(hide),depend(a) :: n = shape(a,0) 665 <ftype2> dimension(n,n),intent(in,out,copy,out=ht,aligned8),check(shape(a,0)==shape(a,1)) :: a 666 integer intent(in),optional :: lo = 0 667 integer intent(in),optional,depend(n) :: hi = n-1 668 <ftype2> dimension(n-1),intent(in),depend(n) :: tau 669 <ftype2> dimension(lwork),intent(cache,hide),depend(lwork) :: work 670 integer intent(in),optional,depend(n),check(lwork>=hi-lo) :: lwork = max(hi-lo,1) 671 integer intent(out) :: info 672 673end subroutine <prefix2>orghr 674 675subroutine <prefix2>orghr_lwork(n,lo,hi,a,tau,work,lwork,info) 676 ! LWORK computation ofr orghr 677 fortranname <prefix2>orghr 678 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,&a,&n,&tau,&work,&lwork,&info); } 679 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 680 integer intent(in) :: n 681 <ftype2> intent(hide) :: a 682 integer intent(in), optional :: lo = 0 683 integer intent(in), optional, depend(n) :: hi = n-1 684 <ftype2> intent(hide) :: tau 685 <ftype2> intent(out) :: work 686 integer intent(hide) :: lwork = -1 687 integer intent(out) :: info 688 689end subroutine <prefix2>orghr_lwork 690 691subroutine <prefix2c>unghr(n,lo,hi,a,tau,work,lwork,info) 692 ! q,info = orghr(a,tau,lo=0,hi=n-1,lwork=n,overwrite_a=0) 693 ! Compute orthogonal matrix Q for Hessenberg reduction from the matrix 694 ! that was computed by gehrd 695 ! 696 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,a,&n,tau,work,&lwork,&info); } 697 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* 698 integer intent(hide),depend(a) :: n = shape(a,0) 699 <ftype2c> dimension(n,n),intent(in,out,copy,out=ht,aligned8),check(shape(a,0)==shape(a,1)) :: a 700 integer intent(in),optional :: lo = 0 701 integer intent(in),optional,depend(n) :: hi = n-1 702 <ftype2c> dimension(n-1),intent(in),depend(n) :: tau 703 <ftype2c> dimension(lwork),intent(cache,hide),depend(lwork) :: work 704 integer intent(in),optional,depend(n),check(lwork>=hi-lo) :: lwork = max(hi-lo,1) 705 integer intent(out) :: info 706 707end subroutine <prefix2c>unghr 708 709subroutine <prefix2c>unghr_lwork(n,lo,hi,a,tau,work,lwork,info) 710 ! LWORK computation for orghr 711 fortranname <prefix2c>unghr 712 callstatement { hi++; lo++; (*f2py_func)(&n,&lo,&hi,&a,&n,&tau,&work,&lwork,&info); } 713 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* 714 integer intent(in) :: n 715 <ftype2c> intent(hide) :: a 716 integer intent(in), optional :: lo = 0 717 integer intent(in), optional, depend(n) :: hi = n-1 718 <ftype2c> intent(hide) :: tau 719 <ftype2c> intent(out) :: work 720 integer intent(hide) :: lwork = -1 721 integer intent(out) :: info 722 723end subroutine <prefix2c>unghr_lwork 724 725subroutine <prefix2>orgqr(m,n,k,a,tau,work,lwork,info) 726 ! q,work,info = orgqr(a,lwork=3*n,overwrite_a=0) 727 ! Generates an M-by-N real matrix Q with orthonormal columns, 728 ! which is defined as the first N columns of a product of K elementary 729 ! reflectors of order M (e.g. output of geqrf) 730 731 threadsafe 732 callstatement (*f2py_func)(&m,&n,&k,a,&m,tau,work,&lwork,&info) 733 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 734 735 integer intent(hide),depend(a):: m = shape(a,0) 736 integer intent(hide),depend(a):: n = shape(a,1) 737 integer intent(hide),depend(tau):: k = shape(tau,0) 738 <ftype2> dimension(m,n),intent(in,out,copy,out=q) :: a 739 <ftype2> dimension(k),intent(in) :: tau 740 741 integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) :: lwork=max(3*n,1) 742 <ftype2> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 743 integer intent(out) :: info 744 745end subroutine <prefix2>orgqr 746 747subroutine <prefix2c>ungqr(m,n,k,a,tau,work,lwork,info) 748 ! q,work,info = ungqr(a,lwork=3*n,overwrite_a=0) 749 ! Generates an M-by-N complex matrix Q with unitary columns, 750 ! which is defined as the first N columns of a product of K elementary 751 ! reflectors of order M (e.g. output of geqrf) 752 753 threadsafe 754 callstatement (*f2py_func)(&m,&n,&k,a,&m,tau,work,&lwork,&info) 755 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* 756 757 integer intent(hide),depend(a):: m = shape(a,0) 758 integer intent(hide),depend(a):: n = shape(a,1) 759 integer intent(hide),depend(tau):: k = shape(tau,0) 760 <ftype2c> dimension(m,n),intent(in,out,copy,out=q) :: a 761 <ftype2c> dimension(k),intent(in) :: tau 762 763 integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) :: lwork=max(3*n,1) 764 <ftype2c> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 765 integer intent(out) :: info 766 767end subroutine <prefix2c>ungqr 768 769subroutine <prefix2>ormqr(side,trans,m,n,k,a,lda,tau,c,ldc,work,lwork,info) 770 ! cq,work,info = ormqr(side,trans,a,tau,c,lwork) 771 ! multiplies the real matrix C with the real orthogonal matrix Q, 772 ! which is defined as the first N columns of a product of K elementary 773 ! reflectors of order M (e.g. output of geqrf) 774 775 threadsafe 776 callstatement (*f2py_func)(side,trans,&m,&n,&k,a,&lda,tau,c,&ldc,work,&lwork,&info) 777 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* 778 779 character intent(in),check(*side=='L'||*side=='R'):: side 780 character intent(in),check(*trans=='N'||*trans=='T'):: trans 781 integer intent(hide),depend(c):: m = shape(c,0) 782 integer intent(hide),depend(c):: n = shape(c,1) 783 integer intent(hide),depend(a):: k = shape(a,1) 784 <ftype2> dimension(lda,k),intent(in):: a 785 integer intent(hide),depend(a):: lda = shape(a, 0) 786 <ftype2> dimension(k),intent(in):: tau 787 <ftype2> dimension(ldc,n),intent(in,out,copy,out=cq):: c 788 integer intent(hide),depend(c):: ldc = shape(c, 0) 789 <ftype2> dimension(MAX(lwork,1)),intent(out):: work 790 integer intent(in):: lwork 791 integer intent(out) :: info 792 793end subroutine <prefix2>ormqr 794 795subroutine <prefix2c>unmqr(side,trans,m,n,k,a,lda,tau,c,ldc,work,lwork,info) 796 ! cq,work,info = unmqr(side,trans,a,tau,c,lwork) 797 ! multiplies the complex matrix C with the complex unitary matrix Q, 798 ! which is defined as the first N columns of a product of K elementary 799 ! reflectors of order M (e.g. output of geqrf) 800 801 threadsafe 802 callstatement (*f2py_func)(side,trans,&m,&n,&k,a,&lda,tau,c,&ldc,work,&lwork,&info) 803 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* 804 805 character intent(in),check(*side=='L'||*side=='R'):: side 806 character intent(in),check(*trans=='N'||*trans=='C'):: trans 807 integer intent(hide),depend(c):: m = shape(c,0) 808 integer intent(hide),depend(c):: n = shape(c,1) 809 integer intent(hide),depend(a):: k = shape(a,1) 810 <ftype2c> dimension(lda,k),intent(in):: a 811 integer intent(hide),depend(a):: lda = shape(a, 0) 812 <ftype2c> dimension(k),intent(in):: tau 813 <ftype2c> dimension(ldc,n),intent(in,out,copy,out=cq):: c 814 integer intent(hide),depend(c):: ldc = shape(c, 0) 815 <ftype2c> dimension(MAX(lwork,1)),intent(out):: work 816 integer intent(in):: lwork 817 integer intent(out) :: info 818 819end subroutine <prefix2c>unmqr 820 821subroutine <prefix>geqrt(m,n,nb,a,lda,t,ldt,work,info) 822 ! a,t,info = geqrt(nb,a,[overwrite_a=0]) 823 ! 824 ! Computes a QR factorization with block size nb of a general matrix A, 825 ! using the compact WY representation for Q. T stores the upper triangular 826 ! block reflectors. 827 828 callstatement (*f2py_func)(&m,&n,&nb,a,&lda,t,&ldt,work,&info) 829 callprotoargument F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 830 831 integer intent(hide),depend(a):: m = shape(a,0) 832 integer intent(hide),depend(a):: n = shape(a,1) 833 integer intent(in),depend(m,n),check(MIN(m,n)>=nb&&nb>=1):: nb 834 <ftype> dimension(m,n),intent(in,out,copy):: a 835 integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) 836 <ftype> dimension(nb,MIN(m,n)),intent(out):: t 837 integer intent(hide),depend(t):: ldt = MAX(1,shape(t,0)) 838 <ftype> dimension(nb,n),intent(hide,cache):: work 839 integer intent(out):: info 840 841end subroutine <prefix>geqrt 842 843subroutine <prefix>gemqrt(side,trans,m,n,k,nb,v,ldv,t,ldt,c,ldc,work,info) 844 ! c,info = gemqrt(side,trans,v,t,c,[overwrite_c=0]) 845 ! 846 ! Multiplies a general matrix C by the orthogonal matrix Q defined by the 847 ! elementary reflectors stored in matrix V and the upper triangular block 848 ! reflectors in matrix T. C may be multiplied by Q, its transpose (for real 849 ! matrices), or its adjoint (for complex matrices) from the left or right. 850 851 callstatement (*f2py_func)(side,trans,&m,&n,&k,&nb,v,&ldv,t,&ldt,c,&ldc,work,&info) 852 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 853 854 character optional,intent(in),check(*side=='L'||*side=='R'):: side = 'L' 855 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 856 integer intent(hide),depend(c):: m = shape(c,0) 857 integer intent(hide),depend(c):: n = shape(c,1) 858 integer intent(hide),depend(m,n,v),check((*side=='L'?m:n)>=k&&k>=0):: k = shape(v,1) 859 integer intent(hide),depend(k,t),check(k>=nb&&nb>=1):: nb = shape(t,0) 860 <ftype> dimension((side[0]=='L'?m:n),k),intent(in):: v 861 integer intent(hide),depend(v):: ldv = MAX(1,shape(v,0)) 862 <ftype> dimension(nb,k),intent(in):: t 863 integer intent(hide),depend(t):: ldt = MAX(1,shape(t,0)) 864 <ftype> dimension(m,n),intent(in,out,copy):: c 865 integer intent(hide),depend(c):: ldc = MAX(1,shape(c,0)) 866 <ftype> dimension((side[0]=='L'?n:m)*nb),intent(hide,cache):: work 867 integer intent(out):: info 868 869end subroutine <prefix>gemqrt 870 871subroutine <prefix>tpqrt(m,n,l,nb,a,lda,b,ldb,t,ldt,work,info) 872 ! a,b,t,info = tpqrt(l,nb,a,b,[overwrite_a=0,overwrite_b=0]) 873 ! 874 ! Computes a QR factorization with block size nb of a triangular-pentagonal 875 ! matrix consisting of square upper triangular matrix A and pentagonal 876 ! matrix B, in the compact WY representation. L is the order of the 877 ! trapezoidal part of matrix B. T stores the upper triangular block 878 ! reflectors. 879 880 callstatement (*f2py_func)(&m,&n,&l,&nb,a,&lda,b,&ldb,t,&ldt,work,&info) 881 callprotoargument F_INT*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 882 883 integer intent(hide),depend(b):: m = shape(b,0) 884 integer intent(hide),depend(b):: n = shape(b,1) 885 integer intent(in),depend(m,n),check(MIN(m,n)>=l&&l>=0):: l 886 integer intent(in),depend(n),check(n>=nb&&nb>=1):: nb 887 <ftype> dimension(n,n),intent(in,out,copy):: a 888 integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) 889 <ftype> dimension(m,n),intent(in,out,copy):: b 890 integer intent(hide),depend(b):: ldb = MAX(1,shape(b,0)) 891 <ftype> dimension(nb,n),intent(out):: t 892 integer intent(hide),depend(t):: ldt = MAX(1,shape(t,0)) 893 <ftype> dimension(nb,n),intent(hide,cache):: work 894 integer intent(out):: info 895 896end subroutine <prefix>tpqrt 897 898subroutine <prefix>tpmqrt(side,trans,m,n,k,l,nb,v,ldv,t,ldt,a,lda,b,ldb,work,info) 899 ! a,b,info = tpmqrt(side,trans,l,v,t,a,b,[overwrite_a=0,overwrite_b=0]) 900 ! 901 ! Multiplies a general matrix C by the orthogonal matrix Q defined by the 902 ! elementary reflectors stored in the pentagonal matrix V and the upper 903 ! triangular block reflectors in matrix T. L is the order of the trapezoidal 904 ! part of matrix V. Matrix C consists of blocks A and B, and may be 905 ! multiplied by Q, its transpose (for real matrices), or its adjoint (for 906 ! complex matrices) from the left or right. 907 908 callstatement (*f2py_func)(side,trans,&m,&n,&k,&l,&nb,v,&ldv,t,&ldt,a,&lda,b,&ldb,work,&info) 909 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 910 911 character optional,intent(in),check(*side=='L'||*side=='R'):: side = 'L' 912 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 913 integer intent(hide),depend(b):: m = shape(b,0) 914 integer intent(hide),depend(b):: n = shape(b,1) 915 integer intent(hide),depend(t):: k = shape(t,1) 916 integer intent(in),depend(k),check(k>=l&&l>=0):: l 917 integer intent(hide),depend(k,t),check(k>=nb&&nb>=1):: nb = shape(t,0) 918 <ftype> dimension((side[0]=='L'?m:n),k),intent(in):: v 919 integer intent(hide),depend(v):: ldv = MAX(1,shape(v,0)) 920 <ftype> dimension(nb,k),intent(in):: t 921 integer intent(hide),depend(t):: ldt = MAX(1,shape(t,0)) 922 <ftype> dimension((side[0]=='L'?k:m),(side[0]=='L'?n:k)),intent(in,out,copy):: a 923 integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) 924 <ftype> dimension(m,n),intent(in,out,copy):: b 925 integer intent(hide),depend(b):: ldb = MAX(1,shape(b,0)) 926 <ftype> dimension((side[0]=='L'?n:m)*nb),intent(hide,cache):: work 927 integer intent(out):: info 928 929end subroutine <prefix>tpmqrt 930 931subroutine <prefix><or,or,un,un>mrz(side,trans,m,n,k,l,a,lda,nt,tau,c,ldc,work,lwork,info) 932 ! 933 ! ?OR/UNMRZ overwrites the general complex M-by-N matrix C with 934 ! 935 ! SIDE = 'L' SIDE = 'R' 936 ! TRANS = 'N': Q * C C * Q 937 ! TRANS = 'C': Q**H * C C * Q**H 938 ! 939 ! where Q is a complex unitary matrix defined as the product of k 940 ! elementary reflectors 941 ! 942 ! Q = H(1) H(2) . . . H(k) 943 ! 944 ! as returned by ?TZRZF. Q is of order M if SIDE = 'L' and of order N 945 ! if SIDE = 'R'. 946 ! 947 callstatement (*f2py_func)(side,trans,&m,&n,&k,&l,a,&lda,tau,c,&ldc,work,&lwork,&info) 948 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 949 950 character optional,intent(in),check(*side=='L'||*side=='R'):: side = 'L' 951 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 952 integer intent(hide),depend(c),check(m>=0):: m=shape(c,0) 953 integer intent(hide),depend(c),check(n>=0):: n=shape(c,1) 954 integer intent(hide),depend(a):: k=shape(a,0) 955 integer intent(hide),depend(a):: l=shape(a,1)-shape(a,0) 956 <ftype> dimension(k,nt),intent(in),check(shape(a,1)>=shape(a,0)):: a 957 integer intent(hide),depend(a,side,n,m),check((*side=='L'?m:n)==nt):: nt=shape(a,1) 958 integer intent(hide),depend(a):: lda=MAX(shape(a,0),1) 959 <ftype> dimension(k),depend(k),check(rank(tau)==1),intent(in):: tau 960 <ftype> dimension(m,n),intent(in,out,copy,out=cq):: c 961 integer intent(hide),depend(c):: ldc = MAX(shape(c,0),1) 962 <ftype> dimension(lwork),depend(lwork),intent(hide,cache):: work 963 integer optional,intent(in),depend(side,m,n),check(lwork>=(*side=='L'?n:m)||lwork==-1):: lwork=MAX((side[0]=='L'?n:m),1) 964 integer intent(out):: info 965 966end subroutine <prefix><or,or,un,un>mrz 967 968subroutine <prefix><or,or,un,un>mrz_lwork(side,trans,m,n,k,l,a,lda,tau,c,ldc,work,lwork,info) 969 ! 970 ! ?OR/UNMRZ overwrites the general complex M-by-N matrix C with 971 ! 972 ! SIDE = 'L' SIDE = 'R' 973 ! TRANS = 'N': Q * C C * Q 974 ! TRANS = 'C': Q**H * C C * Q**H 975 ! 976 ! where Q is a complex unitary matrix defined as the product of k 977 ! elementary reflectors 978 ! 979 ! Q = H(1) H(2) . . . H(k) 980 ! 981 ! as returned by ?TZRZF. Q is of order M if SIDE = 'L' and of order N 982 ! if SIDE = 'R'. 983 ! 984 fortranname <prefix><or,or,un,un>mrz 985 callstatement (*f2py_func)(side,trans,&m,&n,&k,&l,&a,&lda,&tau,&c,&ldc,&work,&lwork,&info) 986 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* 987 988 character optional,intent(in),check(*side=='L'||*side=='R'):: side = 'L' 989 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 990 integer intent(in):: m 991 integer intent(in):: n 992 integer intent(hide),depend(side,m,n):: k=(*side=='L'?m:n) 993 integer intent(hide):: l 994 <ftype> intent(hide):: a 995 integer intent(hide),depend(k):: lda=k 996 <ftype> intent(hide):: c 997 integer intent(hide),depend(m):: ldc=m 998 <ftype> intent(hide):: tau 999 <ftype> intent(out):: work 1000 integer intent(hide):: lwork=-1 1001 integer intent(out):: info 1002 1003end subroutine <prefix><or,or,un,un>mrz_lwork 1004 1005subroutine <prefix2>orgrq(m,n,k,a,tau,work,lwork,info) 1006 ! q,work,info = orgrq(a,lwork=3*n,overwrite_a=0) 1007 ! Generates an M-by-N real matrix Q with orthonormal columns, 1008 ! which is defined as the first N columns of a product of K elementary 1009 ! reflectors of order M (e.g. output of gerqf) 1010 1011 threadsafe 1012 callstatement (*f2py_func)(&m,&n,&k,a,&m,tau,work,&lwork,&info) 1013 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* 1014 1015 integer intent(hide),depend(a):: m = shape(a,0) 1016 integer intent(hide),depend(a):: n = shape(a,1) 1017 integer intent(hide),depend(tau):: k = shape(tau,0) 1018 <ftype2> dimension(m,n),intent(in,out,copy,out=q) :: a 1019 <ftype2> dimension(k),intent(in) :: tau 1020 1021 integer optional,intent(in),depend(m),check(lwork>=m||lwork==-1) :: lwork=max(3*m,1) 1022 <ftype2> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 1023 integer intent(out) :: info 1024 1025end subroutine <prefix2>orgrq 1026 1027subroutine <prefix2c>ungrq(m,n,k,a,tau,work,lwork,info) 1028 ! q,work,info = ungrq(a,lwork=3*n,overwrite_a=0) 1029 ! Generates an M-by-N complex matrix Q with unitary columns, 1030 ! which is defined as the first N columns of a product of K elementary 1031 ! reflectors of order M (e.g. output of gerqf) 1032 1033 threadsafe 1034 callstatement (*f2py_func)(&m,&n,&k,a,&m,tau,work,&lwork,&info) 1035 callprotoargument F_INT*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* 1036 1037 integer intent(hide),depend(a):: m = shape(a,0) 1038 integer intent(hide),depend(a):: n = shape(a,1) 1039 integer intent(hide),depend(tau):: k = shape(tau,0) 1040 <ftype2c> dimension(m,n),intent(in,out,copy,out=q) :: a 1041 <ftype2c> dimension(k),intent(in) :: tau 1042 1043 integer optional,intent(in),depend(m),check(lwork>=m||lwork==-1) :: lwork=max(3*m,1) 1044 <ftype2c> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work 1045 integer intent(out) :: info 1046 1047end subroutine <prefix2c>ungrq 1048 1049subroutine <prefix>trtri(n,c,info,lower,unitdiag) 1050 1051 ! inv_c,info = trtri(c,lower=0,unitdiag=1,overwrite_c=0) 1052 ! Compute C inverse C^-1 where 1053 ! C = U if lower = 0 1054 ! C = L if lower = 1 1055 ! C is non-unit triangular matrix if unitdiag = 0 1056 ! C is unit triangular matrix if unitdiag = 1 1057 1058 callstatement (*f2py_func)((lower?"L":"U"),(unitdiag?"U":"N"),&n,c,&n,&info) 1059 callprotoargument char*,char*,F_INT*,<ctype>*,F_INT*,F_INT* 1060 1061 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1062 integer optional,intent(in),check(unitdiag==0||unitdiag==1) :: unitdiag = 0 1063 1064 integer depend(c),intent(hide):: n = shape(c,0) 1065 <ftype> dimension(n,n),intent(in,out,copy,out=inv_c) :: c 1066 check(shape(c,0)==shape(c,1)) :: c 1067 integer intent(out) :: info 1068 1069end subroutine <prefix>trtri 1070 1071subroutine <prefix>trsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info) 1072 ! x,scale,info = trsyl(trana='N', tranb='N', isgn, a, b, c) 1073 ! 1074 ! Solves the real Sylvester matrix equation: 1075 ! 1076 ! op(A)*X + X*op(B) = scale*C or op(A)*X - X*op(B) = scale*C 1077 ! 1078 ! where A and B are both quasi-triangular matrices. A and B must be in 1079 ! Schur canonical form. op(A) and op(B) are specified via trana and tranb 1080 ! respectively, and may take the forms 'N' (no transpose), 'T' (transpose), 1081 ! or 'C' (conjugate transpose, where applicable) to indicate the operation 1082 ! to be performed. The value of isgn (1 or -1) specifies the sign of the 1083 ! X*op(B) term in the equation. 1084 ! 1085 ! Upon exit, x contains the solution, scale represnets scale factor, set 1086 ! <= 1 to avoid overflow in the solution, and info contains the exit 1087 ! status: 1088 ! 1089 ! 0: success 1090 ! < 0: if info = -i, the i-th argument had an illegal value 1091 ! 1: A and B have common or very close eigenvalues; perturbed values 1092 ! were used to solve the equation 1093 1094 callstatement (*f2py_func)(trana,tranb,&isgn,&m,&n,a,&lda,b,&ldb,c,&ldc,&scale,&info) 1095 callprotoargument char*,char*,F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,F_INT* 1096 1097 character optional,intent(in),check(*trana=='N'||*trana=='T'||*trana=='C'):: trana='N' 1098 character optional,intent(in),check(*tranb=='N'||*tranb=='T'||*tranb=='C'):: tranb='N' 1099 1100 integer optional,intent(in),check(isgn==1||isgn==-1)::isgn=1 1101 1102 integer depend(a),intent(hide):: m = shape(a,0) 1103 integer depend(b),intent(hide):: n = shape(b,0) 1104 1105 <ftype> dimension(m,m),intent(in) :: a 1106 check(shape(a,0)==shape(a,1)) :: a 1107 integer depend(a),intent(hide):: lda = shape(a,0) 1108 1109 <ftype> dimension(n,n),intent(in) :: b 1110 check(shape(b,0)==shape(b,1)) :: b 1111 integer depend(b),intent(hide):: ldb = shape(b,0) 1112 1113 <ftype> dimension(m,n),intent(in,out,copy,out=x) :: c 1114 integer depend(c),intent(hide):: ldc = shape(c,0) 1115 1116 <ftypereal> intent(out) :: scale 1117 1118 integer intent(out) :: info 1119 1120end subroutine <prefix>trsyl 1121 1122subroutine <prefix2c>hbevd(ab,compute_v,lower,n,ldab,kd,w,z,ldz,work,lwork,rwork,lrwork,iwork,liwork,info) 1123 ! in :Band:zubevd.f 1124 1125 callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,&kd,ab,&ldab,w,z,&ldz,work,&lwork,rwork,&lrwork,iwork,&liwork,&info) 1126 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 1127 1128 ! Remark: if ab is fortran contigous on input 1129 ! and overwrite_ab=1 ab will be overwritten. 1130 <ftype2c> dimension(ldab,n), intent(in, overwrite) :: ab 1131 1132 integer optional,intent(in):: compute_v = 1 1133 check( compute_v==1||compute_v==0) compute_v 1134 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1135 1136 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 1137 integer intent(hide),depend(ab) :: n=shape(ab,1) 1138 ! case n=0 is omitted in calculaton of lwork, lrwork, liwork 1139 ! so we forbid it 1140 check( n>0 ) n 1141 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 1142 1143 <ftype2> dimension(n),intent(out),depend(n) :: w 1144 1145 ! For compute_v=1 z is used and contains the eigenvectors 1146 integer intent(hide),depend(n) :: ldz=(compute_v?n:1) 1147 <ftype2c> dimension(ldz,ldz),intent(out),depend(ldz) :: z 1148 1149 integer intent(hide),depend(n) :: lwork=max((compute_v?2*n*n:n),1) 1150 <ftype2c> dimension(lwork),intent(hide),depend(lwork) :: work 1151 integer intent(out)::info 1152 1153 integer optional, check(lrwork>=(compute_v?1+5*n+2*n*n:n)),depend(n) :: lrwork=(compute_v?1+5*n+2*n*n:n) 1154 1155 <ftype2> intent(hide),dimension(lrwork),depend(lrwork) :: rwork 1156 1157 ! documentation says liwork >=2+5*n, but that crashes, +1 helps 1158 integer optional, check(liwork>=(compute_v?3+5*n:1)),depend(n) :: liwork=(compute_v?3+5*n:1) 1159 integer intent(hide),dimension(liwork),depend(liwork) :: iwork 1160 1161end subroutine <prefix2c>hbevd 1162 1163subroutine <prefix2c>hbevx(ab,ldab,compute_v,range,lower,n,kd,q,ldq,vl,vu,il,iu,abstol,w,z,m,mmax,ldz,work,rwork,iwork,ifail,info) ! in :Band:dsbevx.f 1164 1165 callstatement (*f2py_func)((compute_v?"V":"N"),(range>0?(range==1?"V":"I"):"A"),(lower?"L":"U"),&n,&kd,ab,&ldab,q,&ldq,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,rwork,iwork,ifail,&info) 1166 callprotoargument char*,char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,<ctype2>*,F_INT*,F_INT*,F_INT* 1167 1168 integer optional,intent(in):: compute_v = 1 1169 check(compute_v==1||compute_v==0) compute_v 1170 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1171 1172 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 1173 integer intent(hide),depend(ab) :: n=shape(ab,1) 1174 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 1175 1176 integer optional,intent(in):: range = 0 1177 check(range==2||range==1||range==0) range 1178 1179 ! Remark: if ab is fortran contigous on input 1180 ! and overwrite_ab=1 ab will be overwritten. 1181 <ftype2c> dimension(ldab,n),intent(in, overwrite) :: ab 1182 1183 ! FIXME: do we need to make q available for outside usage ??? 1184 ! If so: how to make this optional 1185 !* Q (output) DOUBLE PRECISION array, dimension (LDQ, N) 1186 !* If JOBZ = 'V', the N-by-N orthogonal matrix used in the 1187 !* reduction to tridiagonal form. 1188 !* If JOBZ = 'N', the array Q is not referenced. 1189 integer intent(hide),depend(n,compute_v) :: ldq=(compute_v?n:1) 1190 <ftype2c> dimension(ldq,ldq),intent(hide),depend(ldq) :: q 1191 1192 <ftype2> :: vl 1193 <ftype2> :: vu 1194 integer,check((il\>=1 && il\<=n)),depend(n) :: il 1195 integer,check((iu\>=1 && iu\<=n && iu\>=il)),depend(n,il) :: iu 1196 1197 ! Remark, we don't use python indexing here, because 1198 ! if someone uses ?sbevx directly, 1199 ! he should expect Fortran style indexing. 1200 !integer,check((il\>=0 && il\<n)),depend(n) :: il+1 1201 !integer,check((iu\>=0 && iu\<n && iu\>=il)),depend(n,il) :: iu+1 1202 1203 ! Remark: 1204 ! Eigenvalues will be computed most accurately when ABSTOL is 1205 ! set to twice the underflow threshold 2*DLAMCH('S'), not zero. 1206 ! 1207 ! The easiest is to wrap DLAMCH (done below) 1208 ! and let the user provide the value. 1209 <ftype2> optional,intent(in):: abstol=0.0 1210 1211 <ftype2> dimension(n),intent(out),depend(n) :: w 1212 1213 <ftype2c> dimension(ldz,mmax),depend(ldz,mmax),intent(out) :: z 1214 integer intent(hide),depend(n,compute_v) :: ldz=(compute_v?n:1) 1215 1216 ! We use the mmax parameter to fix the size of z 1217 ! (only if eigenvalues are requested) 1218 ! Otherwise we would allocate a (possibly) huge 1219 ! region of memory for the eigenvectors, even 1220 ! in cases where only a few are requested. 1221 ! If RANGE = 'V' (range=1) we a priori don't know the 1222 ! number of eigenvalues in the interval in advance. 1223 ! As default we use the maximum value 1224 ! but the user should use an appropriate mmax. 1225 integer intent(in),depend(n,iu,il,compute_v,range) :: mmax=(compute_v?(range==2?(iu-il+1):n):1) 1226 integer intent(out) :: m 1227 1228 <ftype2c> dimension(n),depend(n),intent(hide) :: work 1229 <ftype2> dimension(7*n),depend(n),intent(hide) :: rwork 1230 integer dimension(5*n),depend(n),intent(hide) :: iwork 1231 integer dimension((compute_v?n:1)),depend(n,compute_v),intent(out) :: ifail 1232 integer intent(out):: info 1233 1234end subroutine <prefix2c>hbevx 1235 1236subroutine <prefix>gglse(m,n,p,a,lda,b,ldb,c,d,x,work,lwork,info) 1237 ! Solves the linear equality-constrained least squares (LSE) 1238 ! problem: 1239 ! 1240 ! minimize || c - A*x ||_2 subject to B*x = d 1241 ! 1242 ! where A is an M-by-N matrix, B is a P-by-N matrix, c is a given 1243 ! M-vector, and d is a given P-vector. It is assumed that 1244 ! P <= N <= M+P, and 1245 ! 1246 ! rank(B) = P and rank( (A) ) = N. 1247 ! ( (B) ) 1248 ! 1249 ! These conditions ensure that the LSE problem has a unique solution, 1250 ! which is obtained using a generalized RQ factorization of the 1251 ! matrices (B, A) given by 1252 ! 1253 ! B = (0 R)*Q, A = Z*T*Q. 1254 1255 callstatement (*f2py_func)(&m,&n,&p,a,&lda,b,&ldb,c,d,x,work,&lwork,&info) 1256 callprotoargument F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,<ctype>*,<ctype>*,F_INT*,F_INT* 1257 1258 integer intent(hide),depend(a),check(m>=0) :: m = shape(a,0) 1259 integer intent(hide),depend(a),check(n>=0) :: n = shape(a,1) 1260 integer intent(out) :: info 1261 1262 <ftype> dimension(m,n),intent(in,out,copy,out=t) :: a 1263 integer intent(hide),depend(a) :: lda = MAX(shape(a,0),1) 1264 <ftype> dimension(p,n),depend(n),intent(in,out,copy,out=r) :: b 1265 integer intent(hide),depend(b) :: ldb = MAX(shape(b,0),1) 1266 integer intent(hide),depend(b,m,n),check((p>=n-m)&&(p>=0)) :: p = shape(b,0) 1267 1268 <ftype> dimension(m),depend(m),intent(in,out,copy,out=res) :: c 1269 <ftype> dimension(p),depend(p),intent(in,copy) :: d 1270 <ftype> dimension(n),depend(n),intent(out) :: x 1271 1272 integer optional,intent(in),depend(p,m,n),check((lwork==-1)||(lwork>=1)) :: lwork = max(m+n+p,1) 1273 <ftype> intent(hide),dimension(lwork),depend(lwork) :: work 1274 1275end subroutine <prefix>gglse 1276 1277 1278subroutine <prefix>gglse_lwork(m,n,p,a,lda,b,ldb,c,d,x,work,lwork,info) 1279 ! 1280 ! lwork routine for ?gglse 1281 ! 1282 fortranname <prefix>gglse 1283 callstatement (*f2py_func)(&m,&n,&p,&a,&lda,&b,&ldb,&c,&d,&x,&work,&lwork,&info) 1284 callprotoargument F_INT*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,<ctype>*,<ctype>*,F_INT*,F_INT* 1285 1286 integer intent(in),check(m>=0) :: m 1287 integer intent(in),check(n>=0) :: n 1288 integer intent(in),depend(m,n),check((p>=n-m)&&(p>=0)&&p<=n) :: p 1289 integer intent(out) :: info 1290 <ftype> intent(out) :: work 1291 1292 <ftype> intent(hide) :: a 1293 integer intent(hide),depend(m) :: lda = max(1,m) 1294 <ftype> intent(hide) :: b 1295 integer intent(hide),depend(p) :: ldb = max(1,p) 1296 <ftype> intent(hide) :: c 1297 <ftype> intent(hide) :: d 1298 <ftype> intent(hide) :: x 1299 integer intent(hide) :: lwork = -1 1300 1301end subroutine <prefix>gglse_lwork 1302 1303subroutine <prefix>ppcon(lower,n,ap,anorm,rcond,work,irwork,info,L) 1304 ! ?PPCON estimates the reciprocal of the condition number (in the 1305 ! 1-norm) of a symmetric/Hermitian positive definite packed matrix using 1306 ! the Cholesky factorization A = U**T*U or A = L*L**T computed by 1307 ! DPPTRF. 1308 ! 1309 ! An estimate is obtained for norm(inv(A)), and the reciprocal of the 1310 ! condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))). 1311 threadsafe 1312 callstatement (*f2py_func)((lower?"L":"U"),&n,ap,&anorm,&rcond,work,irwork,&info) 1313 callprotoargument char*,F_INT*,<ctype>*,<ctypereal>*,<ctypereal>*,<ctype>*,<F_INT,F_INT,float,double>*,F_INT* 1314 1315 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1316 integer intent(in),check(n>=0) :: n 1317 <ftype> dimension(L),intent(in) :: ap 1318 integer intent(hide),depend(ap,n),check(L>=(n*(n+1)/2)) :: L = len(ap) 1319 <ftypereal> intent(in):: anorm 1320 <ftypereal> intent(out):: rcond 1321 <ftype> depend(n),dimension(<3*n,3*n,2*n,2*n>),intent(hide,cache):: work 1322 <integer,integer,real, double precision> dimension(n), intent(hide,cache),depend(n) :: irwork 1323 integer intent(out):: info 1324 1325end subroutine <prefix>ppcon 1326 1327subroutine <prefix>ppsv(lower,n,nrhs,ap,b,ldb,info,L) 1328 ! DPPSV computes the solution to a real system of linear equations 1329 ! A * X = B, 1330 ! where A is an N-by-N symmetric positive definite matrix stored in 1331 ! packed format and X and B are N-by-NRHS matrices. 1332 ! 1333 ! The Cholesky decomposition is used to factor A as 1334 ! A = U**T* U, if UPLO = 'U', or 1335 ! A = L * L**T, if UPLO = 'L', 1336 ! where U is an upper triangular matrix and L is a lower triangular 1337 ! matrix. The factored form of A is then used to solve the system of 1338 ! equations A * X = B. 1339 threadsafe 1340 callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,ap,b,&ldb,&info) 1341 callprotoargument char*,F_INT*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1342 1343 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1344 integer intent(in),check(n>=0) :: n 1345 integer intent(hide),depend(b) :: ldb = max(1, shape(b,0)) 1346 integer intent(hide),depend(b) :: nrhs = shape(b,1) 1347 <ftype> dimension(L),intent(in) :: ap 1348 integer intent(hide),depend(ap,n),check(L>=(n*(n+1)/2)) :: L = len(ap) 1349 <ftype> dimension(ldb,nrhs),intent(in,out,copy,out=x) :: b 1350 integer intent(out) :: info 1351 1352end subroutine <prefix>ppsv 1353 1354subroutine <prefix>pptrf(lower,n,ap,info,L) 1355 ! ?PPTRF computes the Cholesky factorization of a symmetric/hermitian 1356 ! positive definite matrix A stored in packed format. 1357 ! 1358 ! The factorization has the form 1359 ! A = U**T * U, if UPLO = 'U', or 1360 ! A = L * L**T, if UPLO = 'L', 1361 ! where U is an upper triangular matrix and L is lower triangular. 1362 threadsafe 1363 callstatement (*f2py_func)((lower?"L":"U"),&n,ap,&info) 1364 callprotoargument char*,F_INT*,<ctype>*,F_INT* 1365 1366 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1367 integer intent(in),check(n>=0) :: n 1368 <ftype> dimension(L),intent(in,out,copy,out=ul) :: ap 1369 integer intent(hide),depend(ap,n),check(L>=(n*(n+1)/2)) :: L = len(ap) 1370 integer intent(out) :: info 1371 1372end subroutine<prefix>pptrf 1373 1374subroutine <prefix>pptri(lower,n,ap,info,L) 1375 ! ?PPTRI computes the inverse of a symmetric/Hermitian positive definite 1376 ! matrix A using the Cholesky factorization A = U**T*U or A = L*L**T 1377 ! computed by ?PPTRF. 1378 threadsafe 1379 callstatement (*f2py_func)((lower?"L":"U"),&n,ap,&info) 1380 callprotoargument char*,F_INT*,<ctype>*,F_INT* 1381 1382 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1383 integer intent(in),check(n>=0) :: n 1384 <ftype> dimension(L),intent(in,out,copy,out=uli) :: ap 1385 integer intent(hide),depend(ap,n),check(L>=(n*(n+1)/2)) :: L = len(ap) 1386 integer intent(out) :: info 1387 1388end subroutine<prefix>pptri 1389 1390subroutine <prefix>pptrs(lower,n,nrhs,ap,b,ldb,info,L) 1391 ! DPPTRS solves a system of linear equations A*X = B with a symmetric 1392 ! positive definite matrix A in packed storage using the Cholesky 1393 ! factorization A = U**T*U or A = L*L**T computed by DPPTRF. 1394 threadsafe 1395 callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,ap,b,&ldb,&info) 1396 callprotoargument char*,F_INT*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1397 1398 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1399 integer intent(in),check(n>=0) :: n 1400 integer intent(hide),depend(b) :: ldb = max(1, shape(b,0)) 1401 integer intent(hide),depend(b) :: nrhs = shape(b,1) 1402 <ftype> dimension(L),intent(in) :: ap 1403 integer intent(hide),depend(ap,n),check(L>=(n*(n+1)/2)) :: L = len(ap) 1404 <ftype> dimension(ldb,nrhs),intent(in,out,copy,out=x) :: b 1405 integer intent(out) :: info 1406 1407end subroutine <prefix>pptrs 1408 1409subroutine <prefix2>sbev(ab,compute_v,lower,n,ldab,kd,w,z,ldz,work,info) 1410 ! in :Band:dsbev.f 1411 ! principally <s,d>sbevd does the same, and are recommended for use. 1412 ! (see man dsbevd) 1413 1414 callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,&kd,ab,&ldab,w,z,&ldz,work,&info) 1415 1416 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT* 1417 1418 ! Remark: if ab is fortran contigous on input 1419 ! and overwrite_ab=1 ab will be overwritten. 1420 <ftype2> dimension(ldab,n), intent(in,overwrite) :: ab 1421 1422 integer optional,intent(in):: compute_v = 1 1423 check(compute_v==1||compute_v==0) compute_v 1424 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1425 1426 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 1427 integer intent(hide),depend(ab) :: n=shape(ab,1) 1428 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 1429 1430 <ftype2> dimension(n),intent(out),depend(n) :: w 1431 1432 ! For compute_v=1 z is used and contains the eigenvectors 1433 integer intent(hide),depend(n) :: ldz=(compute_v?n:1) 1434 <ftype2> dimension(ldz,ldz),intent(out),depend(ldz) :: z 1435 1436 <ftype2> dimension(MAX(1,3*n-1)),intent(hide),depend(n) :: work 1437 integer intent(out)::info 1438 1439end subroutine <prefix2>sbev 1440 1441subroutine <prefix2>sbevd(ab,compute_v,lower,n,ldab,kd,w,z,ldz,work,lwork,iwork,liwork,info) 1442 ! in :Band:dsbevd.f 1443 callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,&kd,ab,&ldab,w,z,&ldz,work,&lwork,iwork,&liwork,&info) 1444 callprotoargument char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 1445 1446 ! Remark: if ab is fortran contigous on input 1447 ! and overwrite_ab=1 ab will be overwritten. 1448 <ftype2> dimension(ldab,n), intent(in, overwrite) :: ab 1449 1450 integer optional,intent(in):: compute_v = 1 1451 check( compute_v==1||compute_v==0) compute_v 1452 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1453 1454 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 1455 integer intent(hide),depend(ab) :: n=shape(ab,1) 1456 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 1457 1458 <ftype2> dimension(n),intent(out),depend(n) :: w 1459 <ftype2> dimension(ldz,ldz),intent(out),depend(ldz) :: z 1460 1461 ! For compute_v=1 z is used and contains the eigenvectors 1462 integer intent(hide),depend(n) :: ldz=(compute_v?n:1) 1463 <ftype2> dimension(ldz,ldz),depend(ldz) :: z 1464 1465 integer intent(hide),depend(n) :: lwork=max((compute_v?1+5*n+2*n*n:2*n),1) 1466 <ftype2> dimension(lwork),intent(hide),depend(lwork) :: work 1467 integer intent(out)::info 1468 1469 integer optional,check(liwork>=(compute_v?3+5*n:1)),depend(n) :: liwork=(compute_v?3+5*n:1) 1470 integer intent(hide),dimension(liwork),depend(liwork) :: iwork 1471 1472end subroutine <prefix2>sbevd 1473 1474subroutine <prefix2>sbevx(ab,ldab,compute_v,range,lower,n,kd,q,ldq,vl,vu,il,iu,abstol,w,z,m,mmax,ldz,work,iwork,ifail,info) ! in :Band:dsbevx.f 1475 1476 callstatement (*f2py_func)((compute_v?"V":"N"),(range>0?(range==1?"V":"I"):"A"),(lower?"L":"U"),&n,&kd,ab,&ldab,q,&ldq,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,iwork,ifail,&info) 1477 callprotoargument char*,char*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*, F_INT*,F_INT*,F_INT* 1478 1479 integer optional,intent(in):: compute_v = 1 1480 check(compute_v==1||compute_v==0) compute_v 1481 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 1482 1483 integer optional,check(shape(ab,0)==ldab),depend(ab) :: ldab=shape(ab,0) 1484 integer intent(hide),depend(ab) :: n=shape(ab,1) 1485 integer intent(hide),depend(ab) :: kd=shape(ab,0)-1 1486 1487 integer optional,intent(in):: range = 0 1488 check(range==2||range==1||range==0) range 1489 1490 1491 ! Remark: if ab is fortran contigous on input 1492 ! and overwrite_ab=1 ab will be overwritten. 1493 <ftype2> dimension(ldab,n),intent(in, overwrite) :: ab 1494 1495 1496 ! FIXME: do we need to make q available for outside usage ??? 1497 ! If so: how to make this optional 1498 !* Q (output) DOUBLE PRECISION array, dimension (LDQ, N) 1499 !* If JOBZ = 'V', the N-by-N orthogonal matrix used in the 1500 !* reduction to tridiagonal form. 1501 !* If JOBZ = 'N', the array Q is not referenced. 1502 integer intent(hide),depend(n,compute_v) :: ldq=(compute_v?n:1) 1503 <ftype2> dimension(ldq,ldq),intent(hide),depend(ldq) :: q 1504 1505 1506 <ftype2> :: vl 1507 <ftype2> :: vu 1508 integer,check((il\>=1 && il\<=n)),depend(n) :: il 1509 integer,check((iu\>=1 && iu\<=n && iu\>=il)),depend(n,il) :: iu 1510 1511 ! Remark, we don't use python indexing here, because 1512 ! if someone uses ?sbevx directly, 1513 ! he should expect Fortran style indexing. 1514 !integer,check((il\>=0 && il\<n)),depend(n) :: il+1 1515 !integer,check((iu\>=0 && iu\<n && iu\>=il)),depend(n,il) :: iu+1 1516 1517 ! Remark: 1518 ! Eigenvalues will be computed most accurately when ABSTOL is 1519 ! set to twice the underflow threshold 2*DLAMCH('S'), not zero. 1520 ! 1521 ! The easiest is to wrap DLAMCH (done below) 1522 ! and let the user provide the value. 1523 <ftype2> optional,intent(in):: abstol=0.0 1524 1525 <ftype2> dimension(n),intent(out),depend(n) :: w 1526 1527 <ftype2> dimension(ldz,mmax),depend(ldz,mmax),intent(out) :: z 1528 integer intent(hide),depend(n,compute_v) :: ldz=(compute_v?n:1) 1529 1530 ! We use the mmax parameter to fix the size of z 1531 ! (only if eigenvalues are requested) 1532 ! Otherwise we would allocate a (possibly) huge 1533 ! region of memory for the eigenvectors, even 1534 ! in cases where only a few are requested. 1535 ! If RANGE = 'V' (range=1) we a priori don't know the 1536 ! number of eigenvalues in the interval in advance. 1537 ! As default we use the maximum value 1538 ! but the user should use an appropriate mmax. 1539 integer intent(in),depend(n,iu,il,compute_v,range) :: mmax=(compute_v?(range==2?(iu-il+1):n):1) 1540 integer intent(out) :: m 1541 1542 <ftype2> dimension(7*n),depend(n),intent(hide) :: work 1543 integer dimension(5*n),depend(n),intent(hide) :: iwork 1544 integer dimension((compute_v?n:1)),depend(n,compute_v),intent(out) :: ifail 1545 integer intent(out):: info 1546 1547end subroutine <prefix2>sbevx 1548 1549subroutine <prefix2>stebz(d,e,range,vl,vu,il,iu,tol,order,n,work,iwork,m,nsplit,w,iblock,isplit,info) 1550 ! computes all or selected eigenvalues of a real, symmetric tridiagonal 1551 ! matrix. 1552 1553 callstatement (*f2py_func)((range>0?(range==1?"V":"I"):"A"),order,&n,&vl,&vu,&il,&iu,&tol,d,e,&m,&nsplit,w,iblock,isplit,work,iwork,&info) 1554 callprotoargument char*,char*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT* 1555 1556 <ftype2> dimension(n),intent(in) :: d 1557 <ftype2> dimension(n-1),depend(n),intent(in) :: e 1558 <ftype2> intent(in) :: vl 1559 <ftype2> intent(in) :: vu 1560 character intent(in) :: order 1561 integer intent(in) :: il 1562 integer intent(in) :: iu 1563 <ftype2> intent(in) :: tol 1564 integer intent(in) :: range 1565 integer depend(d),intent(hide),check(n>0) :: n = shape(d,0) 1566 <ftype2> dimension(4*n),depend(n),intent(hide) :: work 1567 integer dimension(3*n),depend(n),intent(hide) :: iwork 1568 integer intent(out) :: m 1569 integer intent(hide) :: nsplit 1570 <ftype2> dimension(n),depend(n),intent(out) :: w 1571 integer dimension(n),depend(n),intent(out) :: iblock 1572 integer dimension(n),depend(n),intent(out) :: isplit 1573 integer intent(out) :: info 1574 1575end subroutine <prefix2>stebz 1576 1577subroutine <prefix2>sterf(d,e,n,info) 1578 ! computes all eigenvalues of a real, symmetric tridiagonal matrix. 1579 1580 callstatement (*f2py_func)(&n,d,e,&info) 1581 callprotoargument F_INT*,<ctype2>*,<ctype2>*,F_INT* 1582 1583 <ftype2> dimension(n),intent(in,out,copy,out=vals) :: d 1584 <ftype2> dimension(n-1),depend(n),intent(in,copy) :: e 1585 integer depend(d),intent(hide) :: n = shape(d,0) 1586 integer intent(out) :: info 1587 1588end subroutine <prefix2>sterf 1589 1590subroutine <prefix2>stein(d,e,w,iblock,isplit,m,n,z,ldz,work,iwork,ifail,info) 1591 ! computes eigenvectors corresponding to eigenvalues of a real, symmetric 1592 ! tridiagonal matrix. 1593 1594 callstatement (*f2py_func)(&n,d,e,&m,w,iblock,isplit,z,&ldz,work,iwork,ifail,&info) 1595 callprotoargument F_INT*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* 1596 1597 <ftype2> dimension(n),intent(in) :: d 1598 <ftype2> dimension(n-1),depend(n),intent(in) :: e 1599 <ftype2> dimension(m),intent(in) :: w 1600 integer depend(w),intent(hide) :: m = shape(w,0) 1601 integer depend(d),intent(hide),check(n>0) :: n = shape(d,0) 1602 integer dimension(n),depend(n),intent(in) :: iblock 1603 integer dimension(n),depend(n),intent(in) :: isplit 1604 <ftype2> dimension(ldz,m),intent(out) :: z 1605 integer depend(n),intent(hide) :: ldz = n 1606 <ftype2> dimension(5*n),intent(hide) :: work 1607 integer dimension(n),intent(hide) :: iwork 1608 integer dimension(m),depend(m),intent(hide) :: ifail 1609 integer intent(out) :: info 1610 1611end subroutine <prefix2>stein 1612 1613subroutine <prefix2>stemr(d,e,range,vl,vu,il,iu,compute_v,n,m,w,z,ldz,nzc,isuppz,tryrac,work,lwork,iwork,liwork,info) 1614 ! computes all eigenvalues of a real, symmetric tridiagonal matrix. 1615 1616 callstatement (*f2py_func)((compute_v?"V":"N"),(range>0?(range==1?"V":"I"):"A"),&n,d,e,&vl,&vu,&il,&iu,&m,w,z,&ldz,&nzc,isuppz,&tryrac,work,&lwork,iwork,&liwork,&info) 1617 callprotoargument char*,char*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 1618 1619 <ftype2> dimension(n),intent(in,copy) :: d 1620 <ftype2> dimension(n),intent(in) :: e 1621 integer intent(in) :: range 1622 <ftype2> intent(in) :: vl 1623 <ftype2> intent(in) :: vu 1624 integer intent(in) :: il 1625 integer intent(in) :: iu 1626 integer optional,intent(in) :: compute_v = 1 1627 integer intent(hide),depend(d),check(n>0) :: n = shape(d,0) 1628 integer intent(out) :: m 1629 <ftype2> dimension(n),depend(n),intent(out) :: w 1630 <ftype2> dimension(n,n),depend(n),intent(out) :: z 1631 integer depend(n),intent(hide) :: ldz = (compute_v?n:1) ! could be made more efficient for index queries 1632 integer depend(n),intent(hide) :: nzc = n ! can also be passed as -1 to do a query 1633 integer dimension((compute_v?2*n:1)),depend(n),intent(hide) :: isuppz 1634 integer intent(hide) :: tryrac = 1 1635 integer depend(n),optional,intent(in),check(lwork>=(compute_v?18*n:12*n)) :: lwork = max((compute_v?18*n:12*n),1) 1636 <ftype2> dimension(lwork),depend(lwork),intent(hide) :: work 1637 integer depend(n),optional,intent(in),check(liwork>=(compute_v?10*n:8*n)) :: liwork = (compute_v?10*n:8*n) 1638 integer dimension(liwork),depend(liwork),intent(hide) :: iwork 1639 integer intent(out) :: info 1640 1641end subroutine <prefix2>stemr 1642 1643subroutine <prefix2>stemr_lwork(d,e,range,vl,vu,il,iu,compute_v,n,m,w,z,ldz,nzc,isuppz,tryrac,work,lwork,iwork,liwork,info) 1644 ! LWORK=-1, LIWORK=-1 call for STEMR 1645 1646 fortranname <prefix2c>stemr 1647 callstatement (*f2py_func)((compute_v?"V":"N"),(range>0?(range==1?"V":"I"):"A"),&n,d,e,&vl,&vu,&il,&iu,&m,w,z,&ldz,&nzc,isuppz,&tryrac,&work,&lwork,&iwork,&liwork,&info) 1648 callprotoargument char*,char*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* 1649 1650 <ftype2> dimension(n),intent(in,copy) :: d 1651 <ftype2> dimension(n),intent(in,copy) :: e 1652 integer intent(in) :: range 1653 <ftype2> intent(in) :: vl 1654 <ftype2> intent(in) :: vu 1655 integer intent(in) :: il 1656 integer intent(in) :: iu 1657 integer optional,intent(in) :: compute_v = 1 1658 integer intent(hide),depend(d),check(n>0) :: n = shape(d,0) 1659 integer intent(hide) :: m 1660 <ftype2> dimension(n),depend(n),intent(hide) :: w 1661 <ftype2> dimension(n,n),depend(n),intent(hide) :: z 1662 integer depend(n),intent(hide) :: ldz = (compute_v?n:1) ! could be made more efficient for index queries 1663 integer depend(n),intent(hide) :: nzc = n ! can also be passed as -1 to do a query 1664 integer dimension((compute_v?2*n:1)),depend(n),intent(hide) :: isuppz 1665 integer intent(hide) :: tryrac = 1 1666 integer depend(n),intent(hide) :: lwork = -1 1667 <ftype2> intent(out) :: work 1668 integer intent(hide) :: liwork = -1 1669 integer intent(out) :: iwork 1670 integer intent(out) :: info 1671end subroutine <prefix2>stemr_lwork 1672 1673subroutine <prefix2>stev(d,e,compute_v,n,z,ldz,work,info) 1674 ! computes all eigenvalues, and, optionally eigvectors of a real, 1675 ! symmetric tridiagonal matrix. 1676 1677 callstatement (*f2py_func)((compute_v?"V":"N"),&n,d,e,z,&ldz,work,&info) 1678 callprotoargument char*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,<ctype2>*,F_INT* 1679 1680 integer optional,intent(in):: compute_v = 1 1681 <ftype2> dimension(n),intent(in,out,copy,out=vals) :: d 1682 integer depend(d),intent(hide),check(n>0) :: n = shape(d,0) 1683 <ftype2> depend(n),dimension(MAX(n-1,1)),intent(in,copy) :: e 1684 integer intent(hide),depend(n) :: ldz=(compute_v?n:1) 1685 <ftype2> dimension(ldz,(compute_v?n:1)),intent(out),depend(n,ldz) :: z 1686 <ftype2> dimension((compute_v?MAX(1,2*n-2):1)),depend(n),intent(hide) :: work 1687 integer intent(out) :: info 1688 1689end subroutine <prefix2>stev 1690 1691subroutine <prefix><s,s,h,h>frk(transr,uplo,trans,n,k,nt,ka,alpha,a,lda,beta,c) 1692 ! 1693 ! (S/D)SFRK - (C/Z)HFRK performs one of the Sym/Hermitian rank--k operations 1694 ! 1695 ! C := alpha*A*A**H + beta*C, or C := alpha*A**H*A + beta*C 1696 ! 1697 ! where alpha and beta are real scalars, C is an n--by--n Sym/Hermitian 1698 ! matrix and A is an n--by--k matrix in the first case and a k--by--n 1699 ! matrix in the second case. 1700 1701 callstatement (*f2py_func)(transr,uplo,trans,&n,&k,&alpha,a,&lda,&beta,c) 1702 callprotoargument char*,char*,char*,F_INT*,F_INT*,<ctypereal>*,<ctype>*,F_INT*,<ctypereal>*,<ctype>* 1703 1704 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1705 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1706 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 1707 integer intent(in), check(n>=0):: n 1708 integer intent(in), check(k>=0):: k 1709 <ftypereal> intent(in):: alpha 1710 <ftypereal> intent(in):: beta 1711 <ftype> dimension(lda,ka),intent(in):: a 1712 integer intent(hide),depend(a,trans,n,k),check(ka==(*trans=='N'?k:n)):: ka=shape(a,1) 1713 integer intent(hide),depend(trans,a,n,k):: lda = MAX((*trans=='N'?n:k),1) 1714 <ftype> dimension(nt),intent(in,out,copy,out=cout):: c 1715 integer intent(hide),depend(c),check(nt==(n*(n+1)/2)):: nt=shape(c,0) 1716 1717end subroutine <prefix><s,s,h,h>frk 1718 1719subroutine <prefix>tpttf(transr,uplo,n,nt,ap,arf,info) 1720 ! 1721 ! copies a triangular matrix from the standard packed format (TP) to the 1722 ! rectangular full packed format (TF). 1723 ! 1724 1725 callstatement (*f2py_func)(transr,uplo,&n,ap,arf,&info) 1726 callprotoargument char*,char*,F_INT*,<ctype>*,<ctype>*,F_INT* 1727 1728 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1729 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1730 integer intent(in), check(n>=0):: n 1731 <ftype> dimension(nt),intent(in):: ap 1732 integer intent(hide),depend(ap,n),check(nt==(n*(n+1)/2)):: nt=shape(ap,0) 1733 <ftype> dimension(nt),intent(out),depend(nt):: arf 1734 integer intent(out):: info 1735 1736end subroutine <prefix>tpttf 1737 1738subroutine <prefix>tpttr(uplo,n,nt,ap,a,lda,info) 1739 ! 1740 ! TPTTR copies a triangular matrix from the standard packed format (TP) 1741 ! to the standard full format (TR). 1742 ! 1743 1744 callstatement (*f2py_func)(uplo,&n,ap,a,&lda,&info) 1745 callprotoargument char*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1746 1747 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1748 integer intent(in), check(n>=0):: n 1749 <ftype> dimension(nt),intent(in):: ap 1750 integer intent(hide),depend(ap,n),check(nt==(n*(n+1)/2)):: nt=shape(ap,0) 1751 <ftype> dimension(n,n),intent(out),depend(n):: a 1752 integer intent(hide),depend(n):: lda=MAX(n,1) 1753 integer intent(out):: info 1754 1755end subroutine <prefix>tpttr 1756 1757subroutine <prefix>tfttp(transr,uplo,n,nt,ap,arf,info) 1758 ! 1759 ! copies a triangular matrix from the standard packed format (TP) to the 1760 ! rectangular full packed format (TF). 1761 ! 1762 1763 callstatement (*f2py_func)(transr,uplo,&n,arf,ap,&info) 1764 callprotoargument char*,char*,F_INT*,<ctype>*,<ctype>*,F_INT* 1765 1766 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1767 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1768 integer intent(in), check(n>=0):: n 1769 <ftype> dimension(nt),intent(in):: arf 1770 integer intent(hide),depend(arf,n),check(nt==(n*(n+1)/2)):: nt=shape(arf,0) 1771 <ftype> dimension(nt),intent(out),depend(nt):: ap 1772 integer intent(out):: info 1773 1774end subroutine <prefix>tfttp 1775 1776subroutine <prefix>tfttr(transr,uplo,n,nt,arf,a,lda,info) 1777 ! 1778 ! TFTTR copies a triangular matrix from the rectangular full packed 1779 ! format (TF) to the standard full format (TR). 1780 ! 1781 callstatement (*f2py_func)(transr,uplo,&n,arf,a,&lda,&info) 1782 callprotoargument char*,char*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1783 1784 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1785 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1786 integer intent(in), check(n>=0):: n 1787 <ftype> dimension(nt),intent(in):: arf 1788 integer intent(hide),depend(arf,n),check(nt==(n*(n+1)/2)):: nt=shape(arf,0) 1789 <ftype> dimension(lda,n),depend(lda,n),intent(out):: a 1790 integer intent(hide),depend(n):: lda = MAX(n,1) 1791 integer intent(out):: info 1792 1793end subroutine <prefix>tfttr 1794 1795subroutine <prefix>trttf(transr,uplo,n,a,lda,arf,info) 1796 ! 1797 ! TRTTF copies a triangular matrix A from standard full format (TR) 1798 ! to rectangular full packed format (TF). 1799 ! 1800 callstatement (*f2py_func)(transr,uplo,&n,a,&lda,arf,&info) 1801 callprotoargument char*,char*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 1802 1803 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1804 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1805 integer intent(hide),depend(a):: n = shape(a,1) 1806 integer intent(hide),depend(a):: lda = MAX(shape(a,0),1) 1807 <ftype> dimension(lda,n),intent(in),check(shape(a,0)==shape(a,1)):: a 1808 <ftype> dimension(n*(n+1)/2),intent(out),depend(n):: arf 1809 integer intent(out):: info 1810 1811end subroutine <prefix>trttf 1812 1813subroutine <prefix>trttp(uplo,n,a,lda,ap,info) 1814 ! 1815 ! TRTTP copies a triangular matrix from the standard full format (TR) to 1816 ! the standard packed format (TP). 1817 ! 1818 callstatement (*f2py_func)(uplo,&n,a,&lda,ap,&info) 1819 callprotoargument char*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT* 1820 1821 character optional, intent(F_INT),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1822 integer intent(hide),depend(a):: n = shape(a,1) 1823 integer intent(hide),depend(a):: lda = MAX(shape(a,0),1) 1824 <ftype> dimension(lda,n),intent(in),check(shape(a,0)==shape(a,1)):: a 1825 <ftype> dimension(n*(n+1)/2),intent(out),depend(n):: ap 1826 integer intent(out):: info 1827 1828end subroutine <prefix>trttp 1829 1830subroutine <prefix>tfsm(transr,side,uplo,trans,diag,m,n,nt,alpha,a,b,ldb) 1831 ! 1832 ! Level 3 BLAS like routine for A in RFP Format. 1833 ! 1834 ! TFSM solves the matrix equation 1835 ! 1836 ! op( A )*X = alpha*B or X*op( A ) = alpha*B 1837 ! 1838 ! where alpha is a scalar, X and B are m by n matrices, A is a unit, or 1839 ! non-unit, upper or lower triangular matrix and op( A ) is one of 1840 ! 1841 ! op( A ) = A or op( A ) = A**H. 1842 ! 1843 ! A is in Rectangular Full Packed (RFP) Format. The matrix X is overwritten on B. 1844 ! 1845 1846 callstatement (*f2py_func)(transr,side,uplo,trans,diag,&m,&n,&alpha,a,b,&ldb) 1847 callprotoargument char*,char*,char*,char*,char*,F_INT*,F_INT*,<ctype>*,<ctype>*,<ctype>*,F_INT* 1848 1849 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1850 character optional,intent(in),check(*side=='L'||*side=='R'):: side = 'L' 1851 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1852 character optional,intent(in),check(*trans=='N'||*trans==<'T','T','C','C'>):: trans = 'N' 1853 character optional,intent(in),check(*diag=='U'||*diag=='N'):: diag = 'N' 1854 integer intent(hide),depend(b):: m=shape(b,0) 1855 integer intent(hide),depend(b):: n=shape(b,1) 1856 <ftype> intent(in):: alpha 1857 <ftype> dimension(nt),intent(in):: a 1858 integer intent(hide),depend(a,m,n,side),check(*side=='L'?nt==(m*(m+1)/2):nt==(n*(n+1)/2)):: nt=shape(a,0) 1859 <ftype> dimension(m,n),intent(in,out,copy,out=x):: b 1860 integer intent(hide),depend(b):: ldb=MAX(shape(b,0),1) 1861 1862end subroutine <prefix>tfsm 1863 1864subroutine <prefix>pftrf(transr,uplo,n,nt,a,info) 1865 ! 1866 ! Computes the Cholesky factorization of a complex Hermitian positive definite matrix A. 1867 ! The factorization has the form 1868 ! 1869 ! A = U**H * U, if UPLO = 'U', or A = L * L**H, if UPLO = 'L', 1870 ! 1871 ! where U is an upper triangular matrix and L is lower triangular. 1872 ! This is the block version of the algorithm, calling Level 3 BLAS. 1873 ! 1874 callstatement (*f2py_func)(transr,uplo,&n,a,&info) 1875 callprotoargument char*,char*,F_INT*,<ctype>*,F_INT* 1876 1877 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1878 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1879 integer intent(in), check(n>=0):: n 1880 <ftype> dimension(nt),intent(in,out,copy,out=achol):: a 1881 integer intent(hide),depend(a,n),check(nt==(n*(n+1)/2)):: nt=shape(a,0) 1882 integer intent(out):: info 1883 1884end subroutine <prefix>pftrf 1885 1886subroutine <prefix>pftri(transr,uplo,n,nt,a,info) 1887 ! 1888 ! Computes the inverse of a real/complex Sym/Hermitian positive definite 1889 ! matrix A using the Cholesky factorization A = U**H*U or A = L*L**H 1890 ! computed by ?PFTRF. 1891 ! 1892 1893 callstatement (*f2py_func)(transr,uplo,&n,a,&info) 1894 callprotoargument char*,char*,F_INT*,<ctype>*,F_INT* 1895 1896 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1897 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1898 integer intent(in), check(n>=0):: n 1899 <ftype> dimension(nt),intent(in,out,copy,out=ainv):: a 1900 integer intent(hide),depend(a,n),check(nt==(n*(n+1)/2)):: nt=shape(a,0) 1901 integer intent(out):: info 1902 1903end subroutine <prefix>pftri 1904 1905subroutine <prefix>pftrs(transr,uplo,n,nhrs,nt,a,b,ldb,info) 1906 ! 1907 ! ?PFTRS solves a system of linear equations A*X = B with a Sym/Hermitian 1908 ! positive definite matrix A using the Cholesky factorization 1909 ! A = U**H*U or A = L*L**H computed by ?PFTRF. 1910 ! 1911 callstatement (*f2py_func)(transr,uplo,&n,&nhrs,a,b,&ldb,&info) 1912 callprotoargument char*,char*,F_INT*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1913 1914 character optional,intent(in),check(*transr=='N'||*transr==<'T','T','C','C'>):: transr = 'N' 1915 character optional,intent(in),check(*uplo=='U'||*uplo=='L'):: uplo = 'U' 1916 integer intent(in), check(n>=0):: n 1917 <ftype> dimension(nt),intent(in):: a 1918 integer intent(hide),depend(a,n),check(nt==(n*(n+1)/2)):: nt=shape(a,0) 1919 <ftype> dimension(ldb, nhrs),intent(in,out,copy,out=x),depend(n),check(shape(b,0)>=n):: b 1920 integer intent(hide),depend(b):: nhrs = shape(b,1) 1921 integer intent(hide),depend(b):: ldb = MAX(shape(b,0),1) 1922 integer intent(out):: info 1923 1924end subroutine<prefix>pftrs 1925 1926subroutine <prefix>tzrzf(m,n,a,lda,tau,work,lwork,info) 1927 ! TZRZF reduces the M-by-N ( M<=N ) complex upper trapezoidal matrix A 1928 ! to upper triangular form by means of unitary transformations. 1929 ! 1930 ! The upper trapezoidal matrix A is factored as 1931 ! 1932 ! A = ( R 0 ) * Z, 1933 ! 1934 ! where Z is an N-by-N unitary matrix and R is an M-by-M upper 1935 ! triangular matrix. 1936 ! 1937 callstatement (*f2py_func)(&m,&n,a,&lda,tau,work,&lwork,&info) 1938 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1939 1940 integer intent(hide),depend(a):: m=shape(a,0) 1941 integer intent(hide),depend(a):: n=shape(a,1) 1942 integer intent(hide),depend(a):: lda=MAX(shape(a,0),1) 1943 <ftype> dimension(m,n),intent(in,out,copy,out=rz),check(shape(a,1)>=shape(a,0)):: a 1944 <ftype> dimension(m),intent(out),depend(m):: tau 1945 <ftype> dimension(MAX(lwork,1)),depend(lwork),intent(hide,cache):: work 1946 integer optional, intent(in),depend(m),check(lwork>=m):: lwork = MAX(m,1) 1947 integer intent(out):: info 1948 1949end subroutine<prefix>tzrzf 1950 1951subroutine <prefix>tzrzf_lwork(m,n,a,lda,tau,work,lwork,info) 1952 ! lwork computation for tzrzf 1953 fortranname <prefix>tzrzf 1954 callstatement (*f2py_func)(&m,&n,&a,&lda,&tau,&work,&lwork,&info) 1955 callprotoargument F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,<ctype>*,F_INT*,F_INT* 1956 1957 integer intent(in):: m 1958 integer intent(in):: n 1959 integer intent(hide),depend(m):: lda=MAX(m,1) 1960 <ftype> intent(hide) :: a 1961 <ftype> intent(hide) :: tau 1962 integer intent(hide) :: lwork = -1 1963 <ftype> intent(out) :: work 1964 integer intent(out) :: info 1965 1966end subroutine<prefix>tzrzf_lwork 1967 1968subroutine <prefix2>lasd4( n, i, d, z, delta, rho, sigma, work, info ) 1969 ! sigma, delta, work, info = lasd4(d,z,i,rho=1.0) 1970 ! Computes i-th square root of eigenvalue of rank one augmented diagonal matrix. Needed by SVD update procedure 1971 1972 callstatement { i++; (*f2py_func)( &n, &i, d, z, delta, &rho, &sigma, work, &info); } 1973 callprotoargument F_INT*, F_INT*, <ctype2>*, <ctype2>*, <ctype2>*, <ctype2>*, <ctype2>*, <ctype2>*, F_INT* 1974 1975 integer intent(hide),depend(d):: n = shape(d,0) 1976 integer intent(in),depend(d),check(i>=0 && i<=(shape(d,0)-1)):: i 1977 1978 <ftype2> dimension(n),intent(in) :: d 1979 <ftype2> dimension(n),intent(in),depend(n) :: z 1980 1981 <ftype2> intent(out) :: sigma 1982 <ftype2> dimension(n),intent(out),depend(n) :: delta 1983 <ftype2> intent(in),optional:: rho = 1.0 1984 1985 <ftype2> dimension(n),intent(out),depend(n) :: work 1986 integer intent(out) :: info 1987 1988end subroutine <prefix2>lasd4 1989 1990subroutine <prefix>lauum(n,c,info,lower) 1991 ! a,info = lauum(c,lower=0,overwrite_c=0) 1992 ! Compute product 1993 ! U^T * U, C = U if lower = 0 1994 ! L * L^T, C = L if lower = 1 1995 ! C is triangular matrix of the corresponding Cholesky decomposition. 1996 1997 callstatement (*f2py_func)((lower?"L":"U"),&n,c,&n,&info) 1998 callprotoargument char*,F_INT*,<ctype>*,F_INT*,F_INT* 1999 2000 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 2001 2002 integer depend(c),intent(hide):: n = shape(c,0) 2003 <ftype> dimension(n,n),intent(in,out,copy,out=a) :: c 2004 check(shape(c,0)==shape(c,1)) :: c 2005 integer intent(out) :: info 2006 2007end subroutine <prefix>lauum 2008 2009subroutine <prefix>laswp(n,a,nrows,k1,k2,piv,off,inc,m,npiv) 2010 ! a = laswp(a,piv,k1=0,k2=len(piv)-1,off=0,inc=1,overwrite_a=0) 2011 ! Perform row interchanges on the matrix A for each of row k1 through k2 2012 ! 2013 ! piv pivots rows. 2014 2015 callstatement {F_INT i;m=len(piv);for(i=0;i<m;++piv[i++]);++k1;++k2; (*f2py_func)(&n,a,&nrows,&k1,&k2,piv+off,&inc); for(i=0;i<m;--piv[i++]);} 2016 callprotoargument F_INT*,<ctype>*,F_INT*,F_INT*,F_INT*,F_INT*,F_INT* 2017 2018 integer depend(a),intent(hide):: nrows = shape(a,0) 2019 integer depend(a),intent(hide):: n = shape(a,1) 2020 <ftype> dimension(nrows,n),intent(in,out,copy) :: a 2021 integer dimension(npiv),intent(in) :: piv 2022 integer intent(hide),depend(piv,nrows),check(npiv<=nrows) :: npiv = len(piv) 2023 !XXX: how to check that all elements in piv are < n? 2024 2025 integer optional,intent(in),check(0<=k1) :: k1 = 0 2026 integer optional,intent(in),depend(k1,npiv,off),check(k1<=k2 && k2<npiv-off) :: k2 = npiv-1 2027 2028 integer optional, intent(in),check(inc>0||inc<0) :: inc = 1 2029 integer optional,intent(in),depend(npiv),check(off>=0 && off<len(piv)) :: off=0 2030 integer intent(hide),depend(npiv,inc,off),check(npiv-off>(m-1)*abs(inc)) :: m = (len(piv)-off)/abs(inc) 2031 2032end subroutine <prefix>laswp 2033 2034! dlamch = dlamch(cmach) 2035! 2036! determine double precision machine parameters 2037! CMACH (input) CHARACTER*1 2038! Specifies the value to be returned by DLAMCH: 2039! = 'E' or 'e', DLAMCH := eps 2040! = 'S' or 's , DLAMCH := sfmin 2041! = 'B' or 'b', DLAMCH := base 2042! = 'P' or 'p', DLAMCH := eps*base 2043! = 'N' or 'n', DLAMCH := t 2044! = 'R' or 'r', DLAMCH := rnd 2045! = 'M' or 'm', DLAMCH := emin 2046! = 'U' or 'u', DLAMCH := rmin 2047! = 'L' or 'l', DLAMCH := emax 2048! = 'O' or 'o', DLAMCH := rmax 2049! 2050! where 2051! 2052! eps = relative machine precision 2053! sfmin = safe minimum, such that 1/sfmin does not overflow 2054! base = base of the machine 2055! prec = eps*base 2056! t = number of (base) digits in the mantissa 2057! rnd = 1.0 when rounding occurs in addition, 0.0 otherwise 2058! emin = minimum exponent before (gradual) underflow 2059! rmin = underflow threshold - base**(emin-1) 2060! emax = largest exponent before overflow 2061! rmax = overflow threshold - (base**emax)*(1-eps) 2062function dlamch(cmach) 2063 character :: cmach 2064 double precision intent(out):: dlamch 2065end function dlamch 2066 2067function slamch(cmach) 2068 character :: cmach 2069 real intent(out):: slamch 2070end function slamch 2071 2072function <prefix2>lange(norm,m,n,a,lda,work) result(n2) 2073 ! the one norm, or the Frobenius norm, or the infinity norm, or the 2074 ! element of largest absolute value of a real matrix A. 2075 <ftype2> <prefix2>lange, n2 2076 callstatement (*f2py_func)(&<prefix2>lange,norm,&m,&n,a,&lda,work) 2077 callprotoargument <ctype2>*,char*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>* 2078 2079 character intent(in),check(*norm=='M'||*norm=='m'||*norm=='1'||*norm=='O'||*norm=='o'||*norm=='I'||*norm=='i'||*norm=='F'||*norm=='f'||*norm=='E'||*norm=='e'):: norm 2080 integer intent(hide),depend(a,n) :: m = shape(a,0) 2081 integer intent(hide),depend(m) :: lda = max(1,m) 2082 integer intent(hide),depend(a) :: n = shape(a,1) 2083 <ftype2> dimension(m,n),intent(in) :: a 2084 <ftype2> dimension(m+1),intent(cache,hide) :: work 2085end function <prefix2>lange 2086 2087function <prefix2c>lange(norm,m,n,a,lda,work) result(n2) 2088 ! the one norm, or the Frobenius norm, or the infinity norm, or the 2089 ! element of largest absolute value of a complex matrix A. 2090 <ftype2> <prefix2c>lange, n2 2091 callstatement (*f2py_func)(&<prefix2c>lange,norm,&m,&n,a,&lda,work) 2092 callprotoargument <ctype2>*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>* 2093 2094 character intent(in),check(*norm=='M'||*norm=='m'||*norm=='1'||*norm=='O'||*norm=='o'||*norm=='I'||*norm=='i'||*norm=='F'||*norm=='f'||*norm=='E'||*norm=='e'):: norm 2095 integer intent(hide),depend(a,n) :: m = shape(a,0) 2096 integer intent(hide),depend(m) :: lda = max(1,m) 2097 integer intent(hide),depend(a) :: n = shape(a,1) 2098 <ftype2c> dimension(m,n),intent(in) :: a 2099 <ftype2> dimension(m+1),intent(cache,hide) :: work 2100end function <prefix2c>lange 2101 2102subroutine <prefix>larfg(n, alpha, x, incx, tau, lx) 2103 integer intent(in), check(n>=1) :: n 2104 <ftype> intent(in,out) :: alpha 2105 <ftype> intent(in,copy,out), dimension(lx) :: x 2106 integer intent(in), check(incx>0||incx<0) :: incx = 1 2107 <ftype> intent(out) :: tau 2108 integer intent(hide),depend(x,n,incx),check(lx > (n-2)*incx) :: lx = len(x) 2109end subroutine <prefix>larfg 2110 2111subroutine <prefix>larf(side,m,n,v,incv,tau,c,ldc,work,lwork) 2112 character intent(in), check(side[0]=='L'||side[0]=='R') :: side = 'L' 2113 integer intent(in,hide), depend(c) :: m = shape(c,0) 2114 integer intent(in,hide), depend(c) :: n = shape(c,1) 2115 <ftype> intent(in),dimension((side[0]=='L'?(1 + (m-1)*abs(incv)):(1 + (n-1)*abs(incv)))),depend(n,m,side,incv) :: v 2116 integer intent(in), check(incv>0||incv<0) :: incv = 1 2117 <ftype> intent(in) :: tau 2118 <ftype> dimension(m,n), intent(in,copy,out) :: c 2119 integer intent(in,hide) :: ldc = max(1,shape(c,0)) 2120 ! FIXME: work should not have been an input argument but kept here for backwards compatibility! 2121 <ftype> intent(in),dimension(lwork),depend(side,m,n) :: work 2122 integer intent(hide),depend(work),check(lwork >= (side[0]=='L'?n:m)) :: lwork = len(work) 2123end subroutine <prefix>larf 2124 2125subroutine <prefix>lartg(f,g,cs,sn,r) 2126 <ftype> intent(in) :: f 2127 <ftype> intent(in) :: g 2128 <real,double precision,\0,\1> intent(out) :: cs 2129 <ftype> intent(out) :: sn 2130 <ftype> intent(out) :: r 2131end subroutine <prefix>lartg 2132 2133subroutine <prefix2c>rot(n,x,offx,incx,y,offy,incy,c,s,lx,ly) 2134 callstatement (*f2py_func)(&n,x+offx,&incx,y+offy,&incy,&c,&s) 2135 callprotoargument F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<float,double>*,<ctype2c>* 2136 <ftype2c> dimension(lx),intent(in,out,copy) :: x 2137 <ftype2c> dimension(ly),intent(in,out,copy) :: y 2138 integer intent(hide),depend(x) :: lx = len(x) 2139 integer intent(hide),depend(y) :: ly = len(y) 2140 <real,double precision> intent(in) :: c 2141 <ftype2c> intent(in) :: s 2142 integer optional, intent(in), check(incx>0||incx<0) :: incx = 1 2143 integer optional, intent(in), check(incy>0||incy<0) :: incy = 1 2144 integer optional, intent(in), depend(lx), check(offx>=0 && offx<lx) :: offx=0 2145 integer optional, intent(in), depend(ly), check(offy>=0 && offy<ly) :: offy=0 2146 integer optional, intent(in), depend(lx,incx,offx,ly,incy,offy) :: n = (lx-1-offx)/abs(incx)+1 2147 check(lx-offx>(n-1)*abs(incx)) :: n 2148 check(ly-offy>(n-1)*abs(incy)) :: n 2149end subroutine <prefix2c>rot 2150 2151subroutine ilaver(major, minor, patch) 2152 integer intent(out) :: major 2153 integer intent(out) :: minor 2154 integer intent(out) :: patch 2155end subroutine ilaver 2156