! Signatures for f2py-wrappers of FORTRAN LAPACK Symmetric/Hermitian Matrix functions. ! subroutine <prefix2>syev(compute_v,lower,n,w,a,lda,work,lwork,info) ! w,v,info = syev(a,compute_v=1,lower=0,lwork=3*n-1,overwrite_a=0) ! Compute all eigenvalues and, optionally, eigenvectors of a ! real symmetric matrix A. ! ! Performance tip: ! If compute_v=0 then set also overwrite_a=1. callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,a,&lda,w,work,&lwork,&info) callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* integer optional,intent(in):: compute_v = 1 check(compute_v==1||compute_v==0) compute_v integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(a):: n = shape(a,0) integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) <ftype2> dimension(n,n),check(shape(a,0)==shape(a,1)) :: a intent(in,copy,out,out=v) :: a <ftype2> dimension(n),intent(out),depend(n) :: w integer optional,intent(in),depend(n) :: lwork=max(3*n-1,1) check(lwork>=3*n-1) :: lwork <ftype2> dimension(lwork),intent(hide),depend(lwork) :: work integer intent(out) :: info end subroutine <prefix2>syev subroutine <prefix2>syev_lwork(lower,n,w,a,lda,work,lwork,info) ! LWORK routines for syev fortranname <prefix2>syev callstatement (*f2py_func)("N",(lower?"L":"U"),&n,&a,&lda,&w,&work,&lwork,&info) callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(n):: lda = MAX(1, n) <ftype2> intent(hide):: a <ftype2> intent(hide):: w integer intent(hide):: lwork = -1 <ftype2> intent(out):: work integer intent(out):: info end subroutine <prefix2>syev_lwork subroutine <prefix2c>heev(compute_v,lower,n,w,a,lda,work,lwork,rwork,info) ! w,v,info = syev(a,compute_v=1,lower=0,lwork=3*n-1,overwrite_a=0) ! Compute all eigenvalues and, optionally, eigenvectors of a ! complex Hermitian matrix A. ! ! Warning: ! If compute_v=0 and overwrite_a=1, the contents of a is destroyed. callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,a,&lda,w,work,&lwork,rwork,&info) callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* integer optional,intent(in),check(compute_v==1||compute_v==0) :: compute_v = 1 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(a):: n = shape(a,0) integer intent(hide),depend(a):: lda = MAX(1,shape(a,0)) <ftype2c> dimension(n,n),check(shape(a,0)==shape(a,1)) :: a intent(in,copy,out,out=v) :: a <ftype2> dimension(n),intent(out),depend(n) :: w integer optional,intent(in),depend(n) :: lwork=max(2*n-1,1) check(lwork>=2*n-1) :: lwork <ftype2c> dimension(lwork),intent(hide),depend(lwork) :: work <ftype2> dimension(3*n-1),intent(hide),depend(n) :: rwork integer intent(out) :: info end subroutine <prefix2c>heev subroutine <prefix2c>heev_lwork(lower,n,w,a,lda,work,lwork,rwork,info) ! LWORK routines for heev fortranname <prefix2c>heev callstatement (*f2py_func)("N",(lower?"L":"U"),&n,&a,&lda,&w,&work,&lwork,&rwork,&info) callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(n):: lda = MAX(1, n) <ftype2c> intent(hide):: a <ftype2> intent(hide):: w <ftype2> intent(hide):: rwork integer intent(hide):: lwork = -1 <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>heev_lwork subroutine <prefix2>syevd(compute_v,lower,n,w,a,lda,work,lwork,iwork,liwork,info) ! w,v,info = syevd(a,compute_v=1,lower=0,lwork=min_lwork,overwrite_a=0) ! Compute all eigenvalues and, optionally, eigenvectors of a ! real symmetric matrix A using D&C. ! ! Performance tip: ! If compute_v=0 then set also overwrite_a=1. callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,a,&lda,w,work,&lwork,iwork,&liwork,&info) callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* integer optional,intent(in),check(compute_v==1||compute_v==0) :: compute_v = 1 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(a):: n = shape(a,0) <ftype2> dimension(n,n),check(shape(a,0)==shape(a,1)) :: a intent(in,copy,out,out=v) :: a integer intent(hide),depend(n) :: lda = MAX(1,n) <ftype2> dimension(n),intent(out),depend(n) :: w integer optional,intent(in),depend(n,compute_v) :: lwork=max((compute_v?1+6*n+2*n*n:2*n+1),1) check(lwork>=(compute_v?1+6*n+2*n*n:2*n+1)) :: lwork <ftype2> dimension(lwork),intent(hide,cache),depend(lwork) :: work integer optional,intent(in),depend(n,compute_v) :: liwork = (compute_v?3+5*n:1) integer dimension(liwork),intent(hide,cache),depend(liwork) :: iwork integer intent(out) :: info end subroutine <prefix2>syevd subroutine <prefix2>syevd_lwork(compute_v,lower,n,w,a,lda,work,lwork,iwork,liwork,info) ! LWORK routines for syevd fortranname <prefix2>syevd callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,&a,&lda,&w,&work,&lwork,&iwork,&liwork,&info) callprotoargument char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(compute_v==0||compute_v==1) :: compute_v = 1 integer intent(hide):: lda = MAX(1,n) <ftype2> intent(hide):: a <ftype2> intent(hide):: w integer intent(hide):: lwork = -1 integer intent(hide):: liwork = -1 <ftype2> intent(out):: work integer intent(out):: iwork integer intent(out):: info end subroutine <prefix2>syevd_lwork subroutine <prefix2c>heevd(compute_v,lower,n,w,a,lda,work,lwork,iwork,liwork,rwork,lrwork,info) ! w,v,info = heevd(a,compute_v=1,lower=0,lwork=min_lwork,overwrite_a=0) ! Compute all eigenvalues and, optionally, eigenvectors of a ! complex Hermitian matrix A using D&C. ! ! Warning: ! If compute_v=0 and overwrite_a=1, the contents of a is destroyed. callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,a,&lda,w,work,&lwork,rwork,&lrwork,iwork,&liwork,&info) callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* integer optional,intent(in):: compute_v = 1 check(compute_v==1||compute_v==0) compute_v integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(a):: n = shape(a,0) <ftype2c> dimension(n,n),check(shape(a,0)==shape(a,1)) :: a intent(in,copy,out,out=v) :: a integer intent(hide),depend(n) :: lda = MAX(1,n) <ftype2> dimension(n),intent(out),depend(n) :: w integer optional,intent(in),depend(n,compute_v) :: lwork=max((compute_v?2*n+n*n:n+1),1) check(lwork>=(compute_v?2*n+n*n:n+1)) :: lwork <ftype2c> dimension(lwork),intent(hide,cache),depend(lwork) :: work integer optional,intent(in),depend(n,compute_v) :: liwork = (compute_v?3+5*n:1) integer dimension(liwork),intent(hide,cache),depend(liwork) :: iwork integer optional,intent(in),depend(n,compute_v) :: lrwork = (compute_v?1+5*n+2*n*n:n) <ftype2> dimension(lrwork),intent(hide,cache),depend(n,lrwork) :: rwork integer intent(out) :: info end subroutine <prefix2c>heevd subroutine <prefix2c>heevd_lwork(compute_v,lower,n,w,a,lda,work,lwork,iwork,liwork,rwork,lrwork,info) ! LWORK routines for heevd fortranname <prefix2c>heevd callstatement (*f2py_func)((compute_v?"V":"N"),(lower?"L":"U"),&n,&a,&lda,&w,&work,&lwork,&rwork,&lrwork,&iwork,&liwork,&info) callprotoargument char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(compute_v==0||compute_v==1) :: compute_v = 1 integer intent(hide):: lda=MAX(1,n) <ftype2c> intent(hide):: a <ftype2> intent(hide):: w integer intent(hide):: lwork=-1 integer intent(hide):: liwork=-1 integer intent(hide):: lrwork=-1 <ftype2c> intent(out):: work <ftype2> intent(out):: rwork integer intent(out):: iwork integer intent(out):: info end subroutine <prefix2c>heevd_lwork subroutine <prefix>sytf2(lower,n,a,lda,ipiv,info) ! Compute the factorization of a symmetric matrix such that ! A = L * D * L^T if lower = 1 ! A = U * D * U^T if lower = 0 callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,ipiv,&info) callprotoargument char*,F_INT*,<ctype>*,F_INT*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1):: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype> dimension(n,n),intent(in,out,copy,out=ldu):: a integer depend(a),intent(hide):: lda = max(shape(a,0),1) integer dimension(n),depend(n),intent(out):: ipiv integer intent(out):: info end subroutine<prefix>sytf2 subroutine <prefix2>sygst(n,a,lda,b,ldb,info,itype,lower) ! c, info = sygst(a,b) ! Transforms the generalized symmetric eigenvalue problem to standard. ! A = inv(U^T) * A * inv(U), if itype == 1 ! A = U^T * A * U or L^T * A * L, if itype == 2 or 3, respectively ! B must contain the factorized U and L from potrf callstatement (*f2py_func)(&itype,(lower?"L":"U"),&n,a,&lda,b,&ldb,&info) callprotoargument F_INT*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,F_INT* integer optional,intent(in),check(itype==1||itype==2||itype==3):: itype = 1 integer optional,intent(in),check(lower==0||lower==1):: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype2> dimension(n,n),intent(in,out,copy,out=c):: a integer depend(a),intent(hide):: lda = max(shape(a,0),1) <ftype2> dimension(n,n),intent(in):: b integer depend(b),intent(hide):: ldb = max(shape(b,0),1) integer intent(out):: info end subroutine <prefix2>sygst subroutine <prefix>sytrf(lower,n,a,lda,ipiv,work,lwork,info) ! Compute the factorization of a symmetric matrix such that ! A = L * D * L^T if lower = 1 ! A = U * D * U^T if lower = 0 ! This is similar to ?SYTF2 but uses BLAS3 blocked calls callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,ipiv,work,&lwork,&info) callprotoargument char*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1):: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype> dimension(n,n),intent(in,out,copy,out=ldu):: a integer depend(a),intent(hide):: lda = max(shape(a,0),1) integer dimension(n),depend(n),intent(out):: ipiv integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1):: lwork = max(n,1) <ftype> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work integer intent(out):: info end subroutine <prefix>sytrf subroutine <prefix>sytrf_lwork(lower,n,a,lda,ipiv,work,lwork,info) ! lwork computation for ?SYTRF fortranname <prefix>sytrf callstatement (*f2py_func)((lower?"L":"U"),&n,&a,&lda,&ipiv,&work,&lwork,&info) callprotoargument char*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype> intent(hide):: a integer depend(n),intent(hide):: lda = max(n,1) integer intent(hide):: ipiv integer intent(hide):: lwork = -1 <ftype> intent(out):: work integer intent(out):: info end subroutine <prefix>sytrf_lwork subroutine <prefix>sysv(n,nrhs,a,lda,ipiv,b,ldb,work,lwork,info,lower) ! Solve A * X = B for symmetric A matrix callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,a,&lda,ipiv,b,&ldb,work,&lwork,&info) callprotoargument char*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) integer depend(b),intent(hide):: nrhs = shape(b,1) <ftype> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=udut):: a integer depend(a),intent(hide):: lda = shape(a,0) integer dimension(n),depend(n),intent(out):: ipiv <ftype> dimension(n,nrhs),check(shape(b,0)==n),depend(n),intent(in,out,copy,out=x):: b integer depend(b),intent(hide):: ldb = shape(b,0) integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1):: lwork = max(n,1) <ftype> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work integer intent(out):: info end subroutine <prefix>sysv subroutine <prefix>sysv_lwork(n,nrhs,a,lda,ipiv,b,ldb,work,lwork,info,lower) ! lwork computation for ?SYSV fortranname <prefix>sysv callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,&a,&lda,&ipiv,&b,&ldb,&work,&lwork,&info) callprotoargument char*,F_INT*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide):: nrhs = 1 <ftype> intent(hide):: a integer depend(n),intent(hide):: lda = n integer intent(hide):: ipiv <ftype> intent(hide):: b integer depend(n),intent(hide):: ldb = n integer intent(hide):: lwork = -1 <ftype> intent(out):: work integer intent(out):: info end subroutine <prefix>sysv_lwork subroutine <prefix>sysvx(n,nrhs,a,lda,af,ldaf,ipiv,b,ldb,x,ldx,rcond,ferr,berr,work,lwork,irwork,info,factored,lower) ! Solve A * X = B for symmetric A matrix ! The expert driver of ?SYSV with condition number, backward,forward error estimates and iterative refinement ! The (c,z) versions assume only symmetric complex matrices. For Hermitian matrices, routine (c,z)HESVX is used threadsafe callstatement (*f2py_func)((factored?"F":"N"),(lower?"L":"U"),&n,&nrhs,a,&lda,af,&ldaf,ipiv,b,&ldb,x,&ldx,&rcond,ferr,berr,work,&lwork,irwork,&info) callprotoargument char*,char*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctype>*,F_INT*,<F_INT,F_INT,float,double>*,F_INT* integer optional,intent(in),check(factored==0||factored==1) :: factored = 0 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) integer depend(b),intent(hide):: nrhs = shape(b,1) <ftype> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,copy,out,out=a_s):: a integer depend(a),intent(hide):: lda = shape(a,0) <ftype> optional,dimension(n,n),depend(n),intent(in,out,out=udut):: af integer optional,depend(af),intent(hide):: ldaf = shape(af,0) integer optional,dimension(n),depend(n),intent(in,out):: ipiv <ftype> depend(n),dimension(n,nrhs),intent(in,copy,out,out=b_s):: b integer depend(b),intent(hide):: ldb = shape(b,0) <ftype> dimension(n,nrhs),intent(out):: x integer depend(b),intent(hide):: ldx = n <ftypereal> intent(out):: rcond <ftypereal> intent(out),dimension(nrhs),depend(nrhs):: ferr <ftypereal> intent(out),dimension(nrhs),depend(nrhs):: berr integer optional,intent(in),check(lwork>=<3*n,3*n,2*n,2*n>||lwork==-1):: lwork = max(3*n,1) <ftype> dimension(MAX(lwork,1)),intent(hide,cache),depend(lwork) :: work <integer,integer,real,double precision> intent(hide,cache),dimension(n),depend(n) :: irwork integer intent(out):: info end subroutine <prefix>sysvx subroutine <prefix>sysvx_lwork(n,nrhs,a,lda,af,ldaf,ipiv,b,ldb,x,ldx,rcond,ferr,berr,work,lwork,irwork,info,factored,lower) ! lwork computation for ?SYSVX fortranname <prefix>sysvx callstatement (*f2py_func)((factored?"F":"N"),(lower?"L":"U"),&n,&nrhs,&a,&lda,&af,&ldaf,&ipiv,&b,&ldb,&x,&ldx,&rcond,&ferr,&berr,&work,&lwork,&irwork,&info) callprotoargument char*,char*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctype>*,F_INT*,<F_INT,F_INT,float,double>*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide) :: factored = 0 integer intent(hide):: nrhs = 1 <ftype> intent(hide):: a integer depend(n),intent(hide):: lda = n <ftype> intent(hide):: af integer depend(n),intent(hide):: ldaf = n integer intent(hide):: ipiv <ftype> intent(hide):: b integer depend(n),intent(hide):: ldb = n <ftype> intent(hide):: x integer depend(n),intent(hide):: ldx = n <ftypereal> intent(hide):: rcond <ftypereal> intent(hide):: ferr <ftypereal> intent(hide):: berr integer intent(hide):: lwork = -1 <integer,integer,real,double precision> intent(hide):: irwork <ftype> intent(out) :: work integer intent(out):: info end subroutine <prefix>sysvx_lwork subroutine <prefix2>sycon(n,a,lda,ipiv,anorm,rcond,work,iwork,info,lower) ! Estimates the reciprocal of the condition number (in the ! 1-norm) of a real symmetric matrix A using the factorization ! A = U*D*U**T or A = L*D*L**T computed by (S/D)SYTRF. ! ! An estimate is obtained for norm(inv(A)), and the reciprocal of the ! condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))) callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,ipiv,&anorm,&rcond,work,iwork,&info) callprotoargument char*,F_INT*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype2> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in):: a integer intent(hide),depend(a) :: lda = max(shape(a,0),1) integer intent(in),dimension(n),depend(n) :: ipiv <ftype2> intent(in) :: anorm <ftype2> intent(out) :: rcond <ftype2> intent(hide),dimension(2*n),depend(n) :: work integer intent(hide),dimension(n),depend(n) :: iwork integer intent(out) :: info end subroutine <prefix2>sycon subroutine <c,z,c,z><sy,\0,he,\2>con(n,a,lda,ipiv,anorm,rcond,work,info,lower) ! Estimates the reciprocal of the condition number (in the ! 1-norm) of a complex symmetric matrix A using the factorization ! A = U*D*U**T or A = L*D*L**T computed by (C/Z)SYTRF. ! ! An estimate is obtained for norm(inv(A)), and the reciprocal of the ! condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))) callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,ipiv,&anorm,&rcond,work,&info) callprotoargument char*,F_INT*,<ctypecomplex>*,F_INT*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypecomplex>*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftypecomplex> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in):: a integer intent(hide),depend(a) :: lda = max(shape(a,0),1) integer intent(in),dimension(n),depend(n) :: ipiv <ftypereal> intent(in) :: anorm <ftypereal> intent(out) :: rcond <ftypecomplex> intent(hide),dimension(2*n),depend(n) :: work integer intent(out) :: info end subroutine <c,z,c,z><sy,\0,he,\2>con subroutine <prefix>syconv(lower,way,n,a,lda,ipiv,e,info) ! ?SYCONV converts A given by ???TRF into L and D and vice-versa. ! Get Non-diag elements of D (returned in workspace) and apply or reverse permutation done in TRF. callstatement (*f2py_func)((lower?"L":"U"),(way?"R":"C"),&n,a,&lda,ipiv,e,&info) callprotoargument char*,char*,F_INT*,<ctype>*,F_INT*,F_INT*,<ctype>*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(way==0||way==1) :: way = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=a):: a integer intent(hide),depend(a) :: lda = max(shape(a,0),1) integer intent(in),dimension(n),depend(n) :: ipiv <ftype> intent(out),dimension(n),depend(n):: e integer intent(out) :: info end subroutine <prefix>syconv subroutine <prefix2c>hegst(n,a,lda,b,ldb,info,itype,lower) ! c, info = hegst(a,b) ! Transforms the generalized Hermitian eigenvalue problem to standard. ! A = inv(U^H) * A * inv(U), if itype == 1 ! A = U^H * A * U or L^H * A * L, if itype == 2 or 3, respectively ! B must contain the factorized U and L from potrf callstatement (*f2py_func)(&itype,(lower?"L":"U"),&n,a,&lda,b,&ldb,&info) callprotoargument F_INT*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* integer optional,intent(in),check(itype==1||itype==2||itype==3):: itype = 1 integer optional,intent(in),check(lower==0||lower==1):: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype2c> dimension(n,n),intent(in,out,copy,out=c):: a integer depend(a),intent(hide):: lda = max(shape(a,0),1) <ftype2c> dimension(n,n),intent(in):: b integer depend(b),intent(hide):: ldb = max(shape(b,0),1) integer intent(out):: info end subroutine <prefix2c>hegst subroutine <prefix2c>hetrf(lower,n,a,lda,ipiv,work,lwork,info) ! Compute the factorization of a hermitian matrix such that ! A = L * D * L^T if lower = 1 ! A = U * D * U^T if lower = 0 ! This is similar to ?HETF2 but uses BLAS3 blocked calls callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,ipiv,work,&lwork,&info) callprotoargument char*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1):: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) <ftype2c> dimension(n,n),intent(in,out,copy,out=ldu):: a integer depend(a),intent(hide):: lda = max(shape(a,0),1) integer dimension(n),depend(n),intent(out):: ipiv integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1):: lwork = max(n,1) <ftype2c> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work integer intent(out):: info end subroutine <prefix2c>hetrf subroutine <prefix2c>hetrf_lwork(lower,n,a,lda,ipiv,work,lwork,info) ! lwork computation for ?HETRF fortranname <prefix2c>hetrf callstatement (*f2py_func)((lower?"L":"U"),&n,&a,&lda,&ipiv,&work,&lwork,&info) callprotoargument char*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2c> intent(hide):: a integer depend(n),intent(hide):: lda = max(n,1) integer intent(hide):: ipiv integer intent(hide):: lwork = -1 <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>hetrf_lwork subroutine <prefix2c>hesv(n,nrhs,a,lda,ipiv,b,ldb,work,lwork,info,lower) ! Solves A * X = B for X ! A is hermitian. For symmetric A see ?SYSV ! A = U * D * U**H if lower = 0 ! A = L * D * L**H if lower = 1 threadsafe callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,a,&lda,ipiv,b,&ldb,work,&lwork,&info) callprotoargument char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) integer depend(b),intent(hide):: nrhs = shape(b,1) <ftype2c> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=uduh):: a integer depend(a),intent(hide):: lda = shape(a,0) integer dimension(n),depend(n),intent(out):: ipiv <ftype2c> dimension(n,nrhs),check(shape(b,0)==n),depend(n),intent(in,out,copy,out=x):: b integer depend(b),intent(hide):: ldb = shape(b,0) integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1):: lwork = max(n,1) <ftype2c> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work integer intent(out):: info end subroutine <prefix2c>hesv subroutine <prefix2c>hesv_lwork(n,nrhs,a,lda,ipiv,b,ldb,work,lwork,info,lower) ! lwork computation for C/ZHESV fortranname <prefix2c>hesv callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,&a,&lda,&ipiv,&b,&ldb,&work,&lwork,&info) callprotoargument char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(in):: n integer intent(hide):: nrhs = 1 <ftype2c> intent(hide):: a integer intent(hide),depend(n):: lda = n integer intent(hide):: ipiv <ftype2c> intent(hide):: b integer intent(hide),depend(n):: ldb = n integer intent(hide):: lwork = -1 <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>hesv_lwork subroutine <prefix2c>hesvx(n,nrhs,a,lda,af,ldaf,ipiv,b,ldb,x,ldx,rcond,ferr,berr,work,lwork,rwork,info,factored,lower) ! Solves A * X = B for X ! Expert driver for ?HESV ! A is hermitian. For symmetric A see ?SYSVX ! A = U * D * U**H if lower = 0 ! A = L * D * L**H if lower = 1 threadsafe callstatement (*f2py_func)((factored?"F":"N"),(lower?"L":"U"),&n,&nrhs,a,&lda,af,&ldaf,ipiv,b,&ldb,x,&ldx,&rcond,ferr,berr,work,&lwork,rwork,&info) callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* integer optional,intent(in),check(factored==0||factored==1) :: factored = 0 integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(a),intent(hide):: n = shape(a,0) integer depend(b),intent(hide):: nrhs = shape(b,1) <ftype2c> dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,copy) :: a integer depend(a),intent(hide):: lda = shape(a,0) <ftype2c> optional,dimension(n,n),depend(n),intent(in,out,out=uduh) :: af integer optional,depend(af),intent(hide):: ldaf = shape(af,0) integer optional,depend(n),dimension(n),intent(in,out):: ipiv <ftype2c> depend(n),dimension(n,nrhs),intent(in,copy) :: b integer depend(b),intent(hide):: ldb = shape(b,0) <ftype2c> depend(n,nrhs),dimension(n,nrhs),intent(out) :: x integer depend(x),intent(hide):: ldx = shape(x,0) <ftype2> intent(out):: rcond <ftype2> intent(out),dimension(nrhs),depend(nrhs):: ferr <ftype2> intent(out),dimension(nrhs),depend(nrhs):: berr <ftype2c> dimension(MAX(1,lwork)),depend(lwork),intent(hide,cache):: work integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1):: lwork = max(2*n,1) <ftype2> intent(hide,cache),dimension(n),depend(n) :: rwork integer intent(out):: info end subroutine <prefix2c>hesvx subroutine <prefix2c>hesvx_lwork(n,nrhs,a,lda,af,ldaf,ipiv,b,ldb,x,ldx,rcond,ferr,berr,work,lwork,rwork,info,factored,lower) ! lwork computation for ?HESVX fortranname <prefix2c>hesvx callstatement (*f2py_func)((factored?"F":"N"),(lower?"L":"U"),&n,&nrhs,&a,&lda,&af,&ldaf,&ipiv,&b,&ldb,&x,&ldx,&rcond,&ferr,&berr,&work,&lwork,&rwork,&info) callprotoargument char*,char*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(in):: n integer intent(hide) :: factored = 0 integer depend(b),intent(hide):: nrhs = 1 <ftype2c> intent(hide) :: a integer depend(n),intent(hide):: lda = n <ftype2c> intent(hide) :: af integer depend(n),intent(hide):: ldaf = n integer intent(hide):: ipiv <ftype2c> intent(hide) :: b integer depend(n),intent(hide):: ldb = n <ftype2c> intent(hide) :: x integer depend(n),intent(hide):: ldx = n <ftype2> intent(hide):: rcond <ftype2> intent(hide):: ferr <ftype2> intent(hide):: berr integer intent(hide):: lwork = -1 <ftype2> intent(hide):: rwork <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>hesvx_lwork subroutine <prefix2>sytrd(lower,n,a,lda,d,e,tau,work,lwork,info) ! Reduce a real symmetric matrix A to real symmetric ! tridiagonal form T by an orthogonal similarity transformation: ! Q**T * A * Q = T. callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,d,e,tau,work,&lwork,&info); callprotoargument char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT* integer intent(hide),depend(a) :: lda=MAX(shape(a,0),1) integer intent(hide),depend(a) :: n=shape(a,1) integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2> dimension(lda,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a <ftype2> dimension(n),intent(out),depend(n) :: d <ftype2> dimension(n-1),intent(out),depend(n) :: e <ftype2> dimension(n-1),intent(out),depend(n) :: tau integer intent(in),optional,depend(n),check(lwork>=1||lwork==-1) :: lwork = MAX(n,1) <ftype2> dimension(lwork),intent(cache,hide),depend(lwork) :: work integer intent(out) :: info end subroutine <prefix2>sytrd subroutine <prefix2>sytrd_lwork(lower,n,a,lda,d,e,tau,work,lwork,info) ! lwork computation for sytrd fortranname <prefix2>sytrd callstatement (*f2py_func)((lower?"L":"U"),&n,&a,&lda,&d,&e,&tau,&work,&lwork,&info); callprotoargument char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2>*,F_INT*,F_INT* integer intent(in) :: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(n) :: lda=MAX(n,1) <ftype2> intent(hide) :: a <ftype2> intent(hide) :: d <ftype2> intent(hide) :: e <ftype2> intent(hide) :: tau <ftype2> intent(out) :: work integer intent(hide) :: lwork = -1 integer intent(out) :: info end subroutine <prefix2>sytrd_lwork subroutine <prefix2c>hetrd(lower,n,a,lda,d,e,tau,work,lwork,info) ! Reduce a complex hermitian matrix A to real symmetric ! tridiagonal form T by an orthogonal similarity transformation: ! Q**H * A * Q = T. callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,d,e,tau,work,&lwork,&info); callprotoargument char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* integer intent(hide),depend(a) :: lda=MAX(shape(a,0),1) integer intent(hide),depend(a) :: n=shape(a,1) integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2c> dimension(lda,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a <ftype2> dimension(n),intent(out),depend(n) :: d <ftype2> dimension(n-1),intent(out),depend(n) :: e <ftype2c> dimension(n-1),intent(out),depend(n) :: tau integer intent(in),optional,depend(n),check(lwork>=1||lwork==-1) :: lwork = MAX(n,1) <ftype2c> dimension(lwork),intent(cache,hide),depend(lwork) :: work integer intent(out) :: info end subroutine <prefix2c>hetrd subroutine <prefix2c>hetrd_lwork(lower,n,a,lda,d,e,tau,work,lwork,info) ! lwork computation for hetrd fortranname <prefix2c>hetrd callstatement (*f2py_func)((lower?"L":"U"),&n,&a,&lda,&d,&e,&tau,&work,&lwork,&info); callprotoargument char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2c>*,<ctype2c>*,F_INT*,F_INT* integer intent(in) :: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(n) :: lda=MAX(n,1) <ftype2c> intent(hide) :: a <ftype2> intent(hide) :: d <ftype2> intent(hide) :: e <ftype2c> intent(hide) :: tau <ftype2c> intent(out) :: work integer intent(hide) :: lwork = -1 integer intent(out) :: info end subroutine <prefix2c>hetrd_lwork subroutine <prefix2>syevr(compute_v,range,lower,n,a,lda,vl,vu,il,iu,abstol,w,z,m,ldz,isuppz,work,lwork,iwork,liwork,info) ! Standard Symmetric/HermitianEigenvalue Problem ! Real - Single precision / Double precision ! ! if jobz = 'N' there are no eigvecs hence 0x0 'z' returned ! if jobz = 'V' and range = 'A', z is (nxn) ! if jobz = 'V' and range = 'V', z is (nxn) since returned number of eigs is unknown beforehand ! if jobz = 'V' and range = 'I', z is (nx(iu-il+1)) callstatement (*f2py_func)((compute_v?"V":"N"),range,(lower?"L":"U"),&n,a,&lda,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,isuppz,work,&lwork,iwork,&liwork,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* <ftype2> intent(in,copy,aligned8),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a integer optional,intent(in),check(compute_v==1||compute_v==0):: compute_v = 1 character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range='A' integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(il>=1&&il<=n),depend(n) :: il=1 integer optional,intent(in),check(n>=iu&&iu>=il),depend(n,il) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>=vl),depend(vl) :: vu=1.0 <ftype2> intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=max(1,26*n)||lwork==-1) :: lwork=max(26*n,1) integer optional,intent(in),depend(n),check(liwork>=max(1,10*n)||liwork==-1):: liwork= max(1,10*n) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(z) :: ldz=max(1,shape(z,0)) <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(liwork),depend(liwork) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2> intent(out),dimension((compute_v?MAX(0,n):0),(compute_v?((*range=='I')?(iu-il+1):MAX(1,n)):0)),depend(n,compute_v,range,iu,il) :: z integer intent(out) :: m ! Only returned if range=='A' or range=='I' and il, iu = 1, n integer intent(out),dimension((compute_v?(2*((*range=='A')||((*range=='I') && (iu-il+1==n))?n:0)):0)),depend(n,iu,il,compute_v,range) :: isuppz integer intent(out) :: info end subroutine <prefix2>syevr subroutine <prefix2>syevr_lwork(lower,n,a,lda,vl,vu,il,iu,abstol,m,w,z,ldz,isuppz,work,lwork,iwork,liwork,info) ! LWORK routines for (s/d)syevr fortranname <prefix2>syevr callstatement (*f2py_func)("N","A",(lower?"L":"U"),&n,&a,&lda,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&isuppz,&work,&lwork,&iwork,&liwork,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* ! Inputs integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 ! Not referenced <ftype2> intent(hide) :: a integer intent(hide),depend(n) :: lda = max(1,n) <ftype2> intent(hide) :: vl=0. <ftype2> intent(hide) :: vu=1. integer intent(hide) :: il=1 integer intent(hide) :: iu=2 <ftype2> intent(hide) :: abstol=0. integer intent(hide) :: m=1 <ftype2> intent(hide) :: w <ftype2> intent(hide) :: z integer intent(hide),depend(n):: ldz = max(1,n) integer intent(hide) :: isuppz integer intent(hide) :: lwork = -1 integer intent(hide) :: liwork = -1 ! Outputs <ftype2> intent(out) :: work integer intent(out) :: iwork integer intent(out) :: info end subroutine <prefix2>syevr_lwork subroutine <prefix2c>heevr(compute_v,range,lower,n,a,lda,vl,vu,il,iu,abstol,w,z,m,ldz,isuppz,work,lwork,rwork,lrwork,iwork,liwork,info) ! Standard Symmetric/HermitianEigenvalue Problem ! Complex - Single precision / Double precision ! ! if jobz = 'N' there are no eigvecs hence 0x0 'z' returned ! if jobz = 'V' and range = 'A', z is (nxn) ! if jobz = 'V' and range = 'V', z is (nxn) since returned number of eigs is unknown beforehand ! if jobz = 'V' and range = 'I', z is (nx(iu-il+1)) callstatement (*f2py_func)((compute_v?"V":"N"),range,(lower?"L":"U"),&n,a,&lda,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,isuppz,work,&lwork,rwork,&lrwork,iwork,&liwork,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* <ftype2c> intent(in,copy,aligned8),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a integer optional,intent(in),check(compute_v==1||compute_v==0):: compute_v = 1 character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range='A' integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(il>=1&&il<=n),depend(n) :: il=1 integer optional,intent(in),check(n>=iu&&iu>=il),depend(n,il) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>vl),depend(vl) :: vu=1.0 <ftype2> intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=max(2*n,1)||lwork==-1) :: lwork=max(2*n,1) integer optional,intent(in),depend(n),check(lrwork>=max(24*n,1)||lrwork==-1) :: lrwork=max(24*n,1) integer optional,intent(in),depend(n),check(liwork>=max(1,10*n)||liwork==-1):: liwork= max(1,10*n) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(z) :: ldz=max(1,shape(z,0)) <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work <ftype2> intent(hide),dimension(lrwork),depend(lrwork) :: rwork integer intent(hide),dimension(liwork),depend(liwork) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2c> intent(out),dimension((compute_v?MAX(0,n):0),(compute_v?((*range=='I')?(iu-il+1):MAX(1,n)):0)),depend(n,compute_v,range,iu,il) :: z integer intent(out) :: m ! MKL implementation has a bug that still accesses isuppz array even if ! range=='A' or range=='I' and il, iu = 1, n which is not the case for ! the reference implementation. Hence here isuppz is allocated regardless ! of the settings. It is wasteful but necessary. The bug is fixed in ! mkl 2020 update 2 and when time comes change this line with ! ! integer intent(out),dimension((compute_v?(2*(*range=='A'||((*range=='I') && (iu-il+1==n))?n:0)):0)),depend(n,iu,il,range,compute_v) :: isuppz ! integer intent(out),dimension(2*max(1,n)),depend(n) :: isuppz integer intent(out) :: info end subroutine <prefix2c>heevr subroutine <prefix2c>heevr_lwork(n,lower,a,lda,vl,vu,il,iu,abstol,m,w,z,ldz,isuppz,work,lwork,rwork,lrwork,iwork,liwork,info) ! LWORK routines for (c/z)heevr fortranname <prefix2c>heevr callstatement (*f2py_func)("N","A",(lower?"L":"U"),&n,&a,&lda,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&isuppz,&work,&lwork,&rwork,&lrwork,&iwork,&liwork,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* ! Inputs integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 ! Not referenced <ftype2c> intent(hide) :: a integer intent(hide),depend(n) :: lda = max(1,n) <ftype2> intent(hide) :: vl=0. <ftype2> intent(hide) :: vu=1. integer intent(hide) :: il=1 integer intent(hide) :: iu=2 <ftype2> intent(hide) :: abstol=0. integer intent(hide) :: m=1 <ftype2> intent(hide) :: w <ftype2c> intent(hide) :: z integer intent(hide),depend(n):: ldz = max(1,n) integer intent(hide) :: isuppz integer intent(hide) :: lwork = -1 integer intent(hide) :: lrwork = -1 integer intent(hide) :: liwork = -1 ! Outputs <ftype2c> intent(out) :: work <ftype2> intent(out) :: rwork integer intent(out) :: iwork integer intent(out) :: info end subroutine <prefix2c>heevr_work subroutine <prefix2>syevx(compute_v,range,lower,n,a,lda,vl,vu,il,iu,abstol,w,z,m,ldz,work,lwork,iwork,ifail,info) ! DSYEVX computes selected eigenvalues and, optionally, eigenvectors ! of a real symmetric matrix A. Eigenvalues and eigenvectors can be ! selected by specifying either a range of values or a range of indices ! for the desired eigenvalues. callstatement (*f2py_func)((compute_v?"V":"N"),range,(lower?"L":"U"),&n,a,&lda,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,&lwork,iwork,ifail,&info) callprotoargument char*,char*,char*,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*,F_INT* integer optional,intent(in),check(compute_v==1||compute_v==0):: compute_v = 1 character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range='A' integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer intent(hide),depend(a) :: n=shape(a,0) integer optional,intent(in),check(il>=1&&il<=n),depend(n) :: il=1 integer optional,intent(in),check(n>=iu&&iu>=il),depend(n,il) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>vl),depend(vl) :: vu=1.0 <ftype2> optional,intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1) :: lwork=max(8*n,1) <ftype2> intent(in,copy),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a integer intent(hide),depend(a) :: lda=max(1,shape(a,0)) integer intent(hide),depend(z) :: ldz=max(1,shape(z,0)) <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(5*n),depend(n) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2> intent(out),dimension((compute_v?MAX(0,n):0),(compute_v?((*range=='I')?(iu-il+1):MAX(1,n)):0)),depend(n,compute_v,range,iu,il) :: z integer intent(out) :: m integer intent(out),dimension((compute_v?n:0)),depend(compute_v,n):: ifail integer intent(out) :: info end subroutine <prefix2>syevx subroutine <prefix2>syevx_lwork(lower,n,a,lda,vl,vu,il,iu,abstol,m,w,z,ldz,work,lwork,iwork,ifail,info) ! LWORK routines for (d/s)syevx fortranname <prefix2>syevx callstatement (*f2py_func)("N","A",(lower?"L":"U"),&n,&a,&lda,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&work,&lwork,&iwork,&ifail,&info) callprotoargument char*,char*,char*,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*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2> intent(hide):: a integer intent(hide),depend(n):: lda = MAX(1,n) integer intent(hide):: il = 1 integer intent(hide):: iu = 0 <ftype2> intent(hide):: vl = 0.0 <ftype2> intent(hide):: vu = 1.0 <ftype2> intent(hide):: abstol = 0.0 integer intent(hide):: m <ftype2> intent(hide):: w <ftype2> intent(hide):: z integer intent(hide),depend(n):: ldz = MAX(1,n) integer intent(hide):: lwork = -1 integer intent(hide):: iwork integer intent(hide):: ifail <ftype2> intent(out):: work integer intent(out):: info end subroutine <prefix2>syevx_lwork subroutine <prefix2c>heevx(compute_v,range,lower,n,a,lda,vl,vu,il,iu,abstol,w,z,m,ldz,work,lwork,rwork,iwork,ifail,info) ! (C/Z)HEEVX computes selected eigenvalues and, optionally, eigenvectors ! of a complex Hermitian matrix A. Eigenvalues and eigenvectors can ! be selected by specifying either a range of values or a range of ! indices for the desired eigenvalues. callstatement (*f2py_func)((compute_v?"V":"N"),range,(lower?"L":"U"),&n,a,&lda,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,&lwork,rwork,iwork,ifail,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* integer optional,intent(in),check(compute_v==1||compute_v==0):: compute_v = 1 character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range='A' integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(il>=1&&il<=n),depend(n) :: il=1 integer optional,intent(in),check(n>=iu&&iu>=il),depend(n,il) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>vl),depend(vl) :: vu=1.0 <ftype2> optional,intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1) :: lwork=max(2*n,1) <ftype2c> intent(in,copy),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(a) :: lda=max(1,shape(a,0)) integer intent(hide),depend(z) :: ldz=max(1,shape(z,0)) <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(5*n),depend(n) :: iwork <ftype2> intent(hide),dimension(7*n),depend(n) :: rwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2c> intent(out),dimension((compute_v*n),(compute_v?((*range=='I')?(iu-il+1):MAX(1,n)):0)),depend(compute_v,range,n,iu,il) :: z integer intent(out) :: m integer intent(out),dimension(compute_v*n),depend(compute_v,n):: ifail integer intent(out) :: info end subroutine <prefix2c>heevx subroutine <prefix2c>heevx_lwork(lower,n,a,lda,vl,vu,il,iu,abstol,m,w,z,ldz,work,lwork,rwork,iwork,ifail,info) ! LWORK routines for (c/z)heevx fortranname <prefix2c>heevx callstatement (*f2py_func)("N","A",(lower?"L":"U"),&n,&a,&lda,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&work,&lwork,&rwork,&iwork,&ifail,&info) callprotoargument char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* integer intent(in):: n integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2c> intent(hide):: a integer intent(hide),depend(n):: lda = MAX(1,n) integer intent(hide):: il = 1 integer intent(hide):: iu = 0 <ftype2> intent(hide):: vl = 0.0 <ftype2> intent(hide):: vu = 1.0 <ftype2> intent(hide):: abstol = 0.0 integer intent(hide):: m <ftype2> intent(hide):: w <ftype2c> intent(hide):: z integer intent(hide),depend(n):: ldz = MAX(1,n) integer intent(hide):: lwork = -1 <ftype2> intent(hide):: rwork integer intent(hide):: iwork integer intent(hide):: ifail <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>heevx_lwork subroutine <prefix2>sygv(itype,jobz,uplo,n,lda,ldb,w,a,b,work,lwork,info) ! Generalized Symmetric Eigenvalue Problem ! Real - Single precision / Double precision ! ! if jobz = 'N' there are no eigvecs returned ! if jobz = 'V' 'a' contains eigvecs callstatement (*f2py_func)(&itype,jobz,uplo,&n,a,&lda,b,&ldb,w,work,&lwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* ! itype,jobz ,uplo , n , a ,lda , b , ldb, w , work ,lwork, info <ftype2> intent(in,copy,aligned8,out,out=v),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer optional,intent(in),check(itype > 0 || itype < 4) :: itype = 1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz='V' character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer optional,intent(in),depend(n),check(lwork>0||lwork==-1) :: lwork=max(3*n-1,1) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(b) :: ldb=max(1,shape(b,0)) <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work <ftype2> intent(out),dimension(n),depend(n) :: w integer intent(out) :: info end subroutine <prefix2>sygv subroutine <prefix2>sygv_lwork(itype,jobz,uplo,n,a,lda,b,ldb,w,work,lwork,info) ! LWORK routine for (S,D)SYGV fortranname <prefix2>sygv callstatement (*f2py_func)(&itype,jobz,uplo,&n,&a,&lda,&b,&ldb,&w,&work,&lwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT* ! Inputs integer intent(in):: n character optional,intent(in),check(*uplo=='U' || *uplo=='L') :: uplo='L' ! Not referenced integer intent(hide) :: itype=1 character intent(hide) :: jobz="N" <ftype2> intent(hide) :: a <ftype2> intent(hide) :: b integer intent(hide) :: lda=max(1,n) integer intent(hide) :: ldb=max(1,n) integer intent(hide) :: lwork=-1 <ftype2> intent(hide) :: w ! Outputs <ftype2> intent(out) :: work integer intent(out) :: info end subroutine <prefix2>sygv_lwork subroutine <prefix2c>hegv(itype,jobz,uplo,n,lda,ldb,w,a,b,work,lwork,rwork,info) ! Generalized Symmetric Eigenvalue Problem ! Complex - Single precision / Double precision ! ! if jobz = 'N' there are no eigvecs returned ! if jobz = 'V' 'a' contains eigvecs callstatement (*f2py_func)(&itype,jobz,uplo,&n,a,&lda,b,&ldb,w,work,&lwork,rwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* ! itype,jobz ,uplo , n , a ,lda , b , ldb, w , work ,lwork,rwork ,info <ftype2c> intent(in,copy,aligned8,out,out=v),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2c> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer optional,intent(in),check(itype > 0 || itype < 4) :: itype = 1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz='V' character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer optional,intent(in),depend(n),check(lwork>0||lwork==-1) :: lwork=max(2*n-1,1) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(b) :: ldb=max(1,shape(b,0)) <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work <ftype2> intent(hide),dimension(MAX(1,3*n-2)),depend(n) :: rwork <ftype2> intent(out),dimension(n),depend(n) :: w integer intent(out) :: info end subroutine <prefix2c>hegv subroutine <prefix2c>hegv_lwork(itype,jobz,uplo,n,a,lda,b,ldb,w,work,lwork,rwork,info) ! LWORK Query routine for (c/z)hegv fortranname <prefix2c>hegv callstatement (*f2py_func)(&itype,jobz,uplo,&n,&a,&lda,&b,&ldb,&w,&work,&lwork,&rwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT* ! Inputs integer intent(in):: n character optional,intent(in),check(*uplo=='U' || *uplo=='L') :: uplo='L' ! Not referenced integer intent(hide) :: itype=1 character intent(hide) :: jobz="N" <ftype2c> intent(hide) :: a <ftype2c> intent(hide) :: b integer intent(hide),depend(n) :: lda=MAX(1,n) integer intent(hide),depend(n) :: ldb=MAX(1,n) integer intent(hide) :: lwork=-1 <ftype2> intent(hide) :: w <ftype2> intent(hide) :: rwork ! Outputs <ftype2c> intent(out) :: work integer intent(out) :: info end subroutine <prefix2c>hegv_lwork ! ! Divide and conquer routines for generalized eigenvalue problem ! subroutine <prefix2>sygvd(itype,jobz,uplo,n,lda,ldb,w,a,b,work,lwork,iwork,liwork,info) ! No call to ILAENV is performed. Hence no need for (d/s)sygvd_lwork. Default sizes are optimal. callstatement (*f2py_func)(&itype,jobz,uplo,&n,a,&lda,b,&ldb,w,work,&lwork,iwork,&liwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* <ftype2> intent(in,copy,aligned8,out,out=v),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer optional,intent(in),check(itype > 0 || itype < 4) :: itype = 1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz='V' character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer optional,intent(in),depend(n),check(lwork>0||lwork==-1) :: lwork=(*jobz=='N'?2*n+1:1+6*n+2*n*n) integer optional,intent(in),depend(n),check(liwork>0||liwork==-1) :: liwork=(*jobz=='N'?1:5*n+3) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(b) :: ldb=max(1,shape(b,0)) <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(liwork),depend(liwork) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w integer intent(out) :: info end subroutine <prefix2>sygvd subroutine <prefix2c>hegvd(itype,jobz,uplo,n,lda,ldb,w,a,b,work,lwork,rwork,lrwork,iwork,liwork,info) ! No call to ILAENV is performed. Hence no need for (c/z)hegvd_lwork. Default sizes are optimal. callstatement (*f2py_func)(&itype,jobz,uplo,&n,a,&lda,b,&ldb,w,work,&lwork,rwork,&lrwork,iwork,&liwork,&info) callprotoargument F_INT*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT*,F_INT* <ftype2c> intent(in,copy,aligned8,out,out=v),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2c> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer optional,intent(in),check(itype > 0 || itype < 4) :: itype = 1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz='V' character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer optional,intent(in),depend(n),check(lwork>0||lwork==-1) :: lwork=(*jobz=='N'?n+1:n*(n+2)) integer optional,intent(in),depend(n),check(lrwork>0||lrwork==-1) :: lrwork=max((*jobz=='N'?n:2*n*n+5*n+1),1) integer optional,intent(in),depend(n),check(liwork>0||liwork==-1) :: liwork=(*jobz=='N'?1:5*n+3) integer intent(hide),depend(a) :: n=shape(a,0) integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(b) :: ldb=max(1,shape(b,0)) <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work <ftype2> intent(hide),dimension(lrwork),depend(lrwork) :: rwork integer intent(hide),dimension(liwork),depend(liwork) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w integer intent(out) :: info end subroutine <prefix2c>hegvd ! ! Expert routines for generalized eigenvalue problem ! subroutine <prefix2>sygvx(itype,jobz,range,uplo,n,a,lda,b,ldb,vl,vu,il,iu,abstol,w,z,m,ldz,work,lwork,iwork,ifail,info) ! (S/D)SYGVX computes selected eigenvalues, and optionally, eigenvectors ! of a real generalized symmetric-definite eigenproblem, of the form ! A*x=(lambda)*B*x, A*Bx=(lambda)*x, or B*A*x=(lambda)*x. Here A ! and B are assumed to be symmetric and B is also positive definite. ! Eigenvalues and eigenvectors can be selected by specifying either a ! range of values or a range of indices for the desired eigenvalues. callstatement (*f2py_func)(&itype,jobz,range,uplo,&n,a,&lda,b,&ldb,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,&lwork,iwork,ifail,&info) callprotoargument F_INT*,char*,char*,char*,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*,F_INT* integer optional,intent(in),check(itype>0||itype<4) :: itype=1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz="V" character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo="L" character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range="A" integer intent(hide),depend(a) :: n=shape(a,0) integer optional,intent(in) :: il=1 integer optional,intent(in),depend(n) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>vl),depend(vl) :: vu=1.0 <ftype2> intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1) :: lwork=max(8*n,1) <ftype2> intent(in,copy),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer intent(hide),depend(a) :: lda=MAX(1,shape(a,0)) integer intent(hide),depend(b) :: ldb=MAX(1,shape(b,0)) integer intent(hide),depend(z) :: ldz=MAX(1,shape(z,0)) <ftype2> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(5*n),depend(n) :: iwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2> intent(out),dimension((jobz[0]=='V'?MAX(0,n):0),(jobz[0]=='V'?((range[0]=='I')?(iu-il+1):MAX(1,n)):0)),depend(n,jobz,range,iu,il) :: z integer intent(out) :: m integer intent(out),dimension((jobz[0]=='N'?0:n)),depend(jobz,n):: ifail integer intent(out) :: info end subroutine <prefix2>sygvx subroutine <prefix2>sygvx_lwork(itype,uplo,n,a,lda,b,ldb,vl,vu,il,iu,abstol,m,w,z,ldz,work,lwork,iwork,ifail,info) ! LWORK Query routine for (s/d)sygvx fortranname <prefix2>sygvx callstatement (*f2py_func)(&itype,"N","A",uplo,&n,&a,&lda,&b,&ldb,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&work,&lwork,&iwork,&ifail,&info) callprotoargument F_INT*,char*,char*,char*,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*,F_INT* integer intent(in):: n character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer intent(hide):: itype = 1 <ftype2> intent(hide):: a integer intent(hide),depend(n):: lda = MAX(1,n) <ftype2> intent(hide):: b integer intent(hide):: ldb = MAX(1,n) integer intent(hide):: il = 1 integer intent(hide):: iu = 0 <ftype2> intent(hide):: vl = 0.0 <ftype2> intent(hide):: vu = 1.0 <ftype2> intent(hide):: abstol = 0.0 integer intent(hide):: m <ftype2> intent(hide):: w <ftype2> intent(hide):: z integer intent(hide),depend(n):: ldz = MAX(1,n) integer intent(hide):: lwork = -1 integer intent(hide):: iwork integer intent(hide):: ifail <ftype2> intent(out):: work integer intent(out):: info end subroutine <prefix2>sygvx_lwork subroutine <prefix2c>hegvx(itype,jobz,range,uplo,n,a,lda,b,ldb,vl,vu,il,iu,abstol,w,z,m,ldz,work,lwork,rwork,iwork,ifail,info) ! (C/Z)HEGVX computes selected eigenvalues, and optionally, eigenvectors ! of a complex generalized Hermitian-definite eigenproblem, of the form ! A*x=(lambda)*B*x, A*Bx=(lambda)*x, or B*A*x=(lambda)*x. Here A and ! B are assumed to be Hermitian and B is also positive definite. ! Eigenvalues and eigenvectors can be selected by specifying either a ! range of values or a range of indices for the desired eigenvalues. callstatement (*f2py_func)(&itype,jobz,range,uplo,&n,a,&lda,b,&ldb,&vl,&vu,&il,&iu,&abstol,&m,w,z,&ldz,work,&lwork,rwork,iwork,ifail,&info) callprotoargument F_INT*,char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* integer optional,intent(in),check(itype>0||itype<4) :: itype=1 character optional,intent(in),check(*jobz=='N'||*jobz=='V') :: jobz='V' character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' character optional,intent(in),check(*range=='A'||*range=='V' ||*range=='I') :: range='A' integer intent(hide),depend(a) :: n=shape(a,0) integer optional,intent(in) :: il=1 integer optional,intent(in),depend(n) :: iu=n <ftype2> optional,intent(in) :: vl=0.0 <ftype2> optional,intent(in),check(vu>vl),depend(vl) :: vu=1.0 <ftype2> intent(in) :: abstol=0.0 integer optional,intent(in),depend(n),check(lwork>=1||lwork==-1) :: lwork=max(2*n,1) <ftype2c> intent(in,copy),check(shape(a,0)==shape(a,1)),dimension(n,n) :: a <ftype2c> intent(in,copy,aligned8),check(shape(b,0)==shape(b,1)),check(shape(b,0)==n),dimension(n,n),depend(n) :: b integer intent(hide),depend(n) :: lda=max(1,n) integer intent(hide),depend(n) :: ldb=max(1,n) integer intent(hide),depend(z) :: ldz=max(1,shape(z,0)) <ftype2c> intent(hide),dimension(lwork),depend(lwork) :: work integer intent(hide),dimension(5*n),depend(n) :: iwork <ftype2> intent(hide),dimension(7*n),depend(n) :: rwork <ftype2> intent(out),dimension(n),depend(n) :: w <ftype2c> intent(out),dimension((jobz[0]=='V'?MAX(0,n):0),(jobz[0]=='V'?((range[0]=='I')?(iu-il+1):MAX(1,n)):0)),depend(n,jobz,range,iu,il) :: z integer intent(out) :: m integer intent(out),dimension((jobz[0]=='N'?0:n)),depend(jobz,n):: ifail integer intent(out) :: info end subroutine <prefix2c>hegvx subroutine <prefix2c>hegvx_lwork(itype,uplo,n,a,lda,b,ldb,vl,vu,il,iu,abstol,m,w,z,ldz,work,lwork,rwork,iwork,ifail,info) ! LWORK Query routine for (c/z)hegvx fortranname <prefix2c>hegvx callstatement (*f2py_func)(&itype,"N","A",uplo,&n,&a,&lda,&b,&ldb,&vl,&vu,&il,&iu,&abstol,&m,&w,&z,&ldz,&work,&lwork,&rwork,&iwork,&ifail,&info) callprotoargument F_INT*,char*,char*,char*,F_INT*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,F_INT*,F_INT*,<ctype2>*,F_INT*,<ctype2>*,<ctype2c>*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,F_INT*,F_INT*,F_INT* integer intent(in):: n character optional,intent(in),check(*uplo=='U'||*uplo=='L') :: uplo='L' integer intent(hide):: itype = 1 <ftype2c> intent(hide):: a integer intent(hide),depend(n):: lda = MAX(1,n) <ftype2c> intent(hide):: b integer intent(hide):: ldb = MAX(1,n) integer intent(hide):: il = 1 integer intent(hide):: iu = 0 <ftype2> intent(hide):: vl = 0.0 <ftype2> intent(hide):: vu = 1.0 <ftype2> intent(hide):: abstol = 0.0 integer intent(hide):: m <ftype2> intent(hide):: w <ftype2c> intent(hide):: z integer intent(hide),depend(n):: ldz = MAX(1,n) integer intent(hide):: lwork = -1 <ftype2> intent(hide):: rwork integer intent(hide):: iwork integer intent(hide):: ifail <ftype2c> intent(out):: work integer intent(out):: info end subroutine <prefix2c>hegvx_lwork subroutine <prefix>syequb(lower,n,a,lda,s,scond,amax,work,info) callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,s,&scond,&amax,work,&info) callprotoargument char*,F_INT*,<ctype>*,F_INT*,<ctypereal>*,<ctypereal>*,<ctypereal>*,<ctype>*,F_INT* integer intent(hide),depend(a) :: n = shape(a,1) integer intent(hide),depend(a) :: lda = MAX(1,shape(a,0)) ! See https://github.com/Reference-LAPACK/lapack/pull/61 for the reason for 3*n <ftype> intent(hide,cache),depend(n),dimension(3*n) :: work integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype> intent(in),dimension(lda,n),check(shape(a,0)==shape(a,1)) :: a <ftypereal> intent(out),dimension(n),depend(n) :: s <ftypereal> intent(out) :: scond <ftypereal> intent(out) :: amax integer intent(out) :: info end subroutine <prefix>syequb subroutine <prefix2c>heequb(lower,n,a,lda,s,scond,amax,work,info) callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,s,&scond,&amax,work,&info) callprotoargument char*,F_INT*,<ctype2c>*,F_INT*,<ctype2>*,<ctype2>*,<ctype2>*,<ctype2c>*,F_INT* integer intent(hide),depend(a) :: n = shape(a,1) integer intent(hide),depend(a) :: lda = MAX(1,shape(a,0)) ! See https://github.com/Reference-LAPACK/lapack/pull/61 for the reason for 3*n <ftype2c> intent(hide,cache),depend(n),dimension(3*n) :: work integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 <ftype2c> intent(in),dimension(lda,n),check(shape(a,0)==shape(a,1)) :: a <ftype2> intent(out),dimension(n),depend(n) :: s <ftype2> intent(out) :: scond <ftype2> intent(out) :: amax integer intent(out) :: info end subroutine <prefix2c>heequb