! Signatures for f2py-wrappers of FORTRAN LAPACK Positive Definite Matrix functions. ! subroutine pstrf(n,a,lda,piv,rank_c,tol,work,info,lower) ! c,x,info = pstrf(a,b,lower=0,overwrite_a=0,overwrite_b=0) ! Solve A * X = B. ! A is symmetric/hermitian positive semidefinite ! P**T * A * P = U**T * U , if UPLO = 'U', ! P**T * A * P = L * L**T, if UPLO = 'L', ! where U is an upper triangular matrix and L is lower triangular, and ! P is stored as vector PIV. callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,piv,&rank_c,&tol,work,&info) callprotoargument char*,F_INT*,*,F_INT*,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) dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a integer intent(hide),depend(a) :: lda=max(1,shape(a, 0)) integer intent(out),dimension(n),depend(n) :: piv integer intent(out) :: rank_c optional, intent(in) :: tol = -1.0 intent(hide),dimension(2*n),depend(n) :: work integer intent(out) :: info end subroutine pstrf subroutine pstf2(n,a,lda,piv,rank_c,tol,work,info,lower) ! c,x,info = pstf2(a,b,lower=0,overwrite_a=0,overwrite_b=0) ! Solve A * X = B. ! A is symmetric positive semidefinite ! P**T * A * P = U**T * U , if UPLO = 'U', ! P**T * A * P = L * L**T, if UPLO = 'L', ! where U is an upper triangular matrix and L is lower triangular, and ! P is stored as vector PIV. callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,piv,&rank_c,&tol,work,&info) callprotoargument char*,F_INT*,*,F_INT*,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) dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a integer intent(hide),depend(a) :: lda=max(1,shape(a, 0)) integer intent(out),dimension(n),depend(n) :: piv integer intent(out) :: rank_c optional,intent(in) :: tol = -1.0 intent(hide),dimension(2*n),depend(n) :: work integer intent(out) :: info end subroutine pstf2 subroutine posv(n,nrhs,a,b,info,lower) ! c,x,info = posv(a,b,lower=0,overwrite_a=0,overwrite_b=0) ! Solve A * X = B. ! A is symmetric positive defined ! A = U^T * U, C = U if lower = 0 ! A = L * L^T, C = L if lower = 1 ! C is triangular matrix of the corresponding Cholesky decomposition. callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,a,&n,b,&n,&info) callprotoargument char*,F_INT*,F_INT*,*,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) integer depend(b),intent(hide):: nrhs = shape(b,1) dimension(n,n),intent(in,out,copy,out=c) :: a check(shape(a,0)==shape(a,1)) :: a dimension(n,nrhs),intent(in,out,copy,out=x),depend(n):: b check(shape(a,0)==shape(b,0)) :: b integer intent(out) :: info end subroutine posv subroutine posvx(fact,n,nrhs,a,lda,af,ldaf,equed,s,b,ldb,x,ldx,rcond,ferr,berr,work,irwork,info,lower) ! Solve A * X = B for Symmetric/Hermitian A ! "expert" version of the ?POSV routines threadsafe callstatement (*f2py_func)(fact,(lower?"L":"U"),&n,&nrhs,a,&lda,af,&ldaf,equed,s,b,&ldb,x,&ldx,&rcond,ferr,berr,work,irwork,&info) callprotoargument char*,char*,F_INT*,F_INT*,*,F_INT*,*,F_INT*,char*,*,*,F_INT*,*,F_INT*,*,*,*,*,*,F_INT* character optional,intent(in):: fact = "E" 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) 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) optional,intent(in,out,out=lu),dimension(n,n),depend(n):: af integer depend(af),intent(hide):: ldaf = shape(af,0) character optional,intent(in,out):: equed = "Y" optional,dimension(n),depend(n),intent(in,out):: s dimension(n,nrhs),check(shape(b,0)==n),depend(n),intent(in,copy,out,out=b_s):: b integer depend(b),intent(hide):: ldb = shape(b,0) dimension(n,nrhs),depend(n,nrhs),intent(out):: x integer depend(x),intent(hide):: ldx = shape(x,0) intent(out):: rcond intent(out),dimension(nrhs),depend(nrhs):: ferr intent(out),dimension(nrhs),depend(nrhs):: berr intent(hide),dimension(<3*n,3*n,2*n,2*n>),depend(n):: work intent(hide),dimension(n),depend(n):: irwork integer intent(out):: info end subroutine posvx subroutine pocon(uplo,n,a,lda,anorm,rcond,work,irwork,info) ! Computes the 1- or inf- norm reciprocal condition number estimate ! for a positive definite symmetric/hermitian matrix. threadsafe callstatement (*f2py_func)(uplo,&n,a,&lda,&anorm,&rcond,work,irwork,&info) callprotoargument char*,F_INT*,*,F_INT*,*,*,*,*,F_INT* character optional,intent(in):: uplo = 'U' integer depend(a),intent(hide):: n = shape(a,0) dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in):: a integer depend(a),intent(hide):: lda = shape(a,0) intent(in):: anorm intent(out):: rcond depend(n),dimension(<3*n,3*n,2*n,2*n>),intent(hide,cache):: work depend(n),dimension(n),intent(hide,cache):: irwork integer intent(out):: info end subroutine pocon subroutine potrf(n,a,lda,info,lower,clean) ! c,info = potrf(a,lower=0,clean=1,overwrite_a=0) ! Compute Cholesky decomposition of symmetric positive defined matrix: ! A = U^T * U, C = U if lower = 0 ! A = L * L^T, C = L if lower = 1 ! C is triangular matrix of the corresponding Cholesky decomposition. ! clean==1 zeros strictly lower or upper parts of U or L, respectively callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,&info); if(clean){F_INT i,j;if(lower){for(i=0;i\*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(clean==0||clean==1) :: clean = 1 integer depend(a),intent(hide) :: n = shape(a,0) dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a integer depend(n),intent(hide) :: lda = MAX(1,n) integer intent(out) :: info end subroutine potrf subroutine potrf(n,a,lda,info,lower,clean) ! c,info = potrf(a,lower=0,clean=1,overwrite_a=0) ! Compute Cholesky decomposition of symmetric positive defined matrix: ! A = U^H * U, C = U if lower = 0 ! A = L * L^H, C = L if lower = 1 ! C is triangular matrix of the corresponding Cholesky decomposition. ! clean==1 zeros strictly lower or upper parts of U or L, respectively callstatement (*f2py_func)((lower?"L":"U"),&n,a,&lda,&info); if(clean){F_INT i,j,k;if(lower){for(i=0;i\r=(a+k)->i=0.0;}} else {for(i=0;i\r=(a+k)->i=0.0;}}} callprotoargument char*,F_INT*,*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer optional,intent(in),check(clean==0||clean==1) :: clean = 1 integer depend(a),intent(hide):: n = shape(a,0) dimension(n,n),check(shape(a,0)==shape(a,1)),intent(in,out,copy,out=c) :: a integer depend(n),intent(hide) :: lda = MAX(1,n) integer intent(out) :: info end subroutine potrf subroutine potrs(n,nrhs,c,b,info,lower) ! x,info = potrs(c,b,lower=0=1,overwrite_b=0) ! Solve A * X = B. ! A is symmetric positive defined ! A = U^T * U, C = U if lower = 0 ! A = L * L^T, C = L if lower = 1 ! C is triangular matrix of the corresponding Cholesky decomposition. callstatement (*f2py_func)((lower?"L":"U"),&n,&nrhs,c,&n,b,&n,&info) callprotoargument char*,F_INT*,F_INT*,*,F_INT*,*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(c),intent(hide):: n = shape(c,0) integer depend(b),intent(hide):: nrhs = shape(b,1) dimension(n,n),intent(in) :: c check(shape(c,0)==shape(c,1)) :: c dimension(n,nrhs),intent(in,out,copy,out=x),depend(n):: b check(shape(c,0)==shape(b,0)) :: b integer intent(out) :: info end subroutine potrs subroutine potri(n,c,info,lower) ! inv_a,info = potri(c,lower=0,overwrite_c=0) ! Compute A inverse A^-1. ! A = U^T * U, C = U if lower = 0 ! A = L * L^T, C = L if lower = 1 ! C is triangular matrix of the corresponding Cholesky decomposition. callstatement (*f2py_func)((lower?"L":"U"),&n,c,&n,&info) callprotoargument char*,F_INT*,*,F_INT*,F_INT* integer optional,intent(in),check(lower==0||lower==1) :: lower = 0 integer depend(c),intent(hide):: n = shape(c,0) dimension(n,n),intent(in,out,copy,out=inv_a) :: c check(shape(c,0)==shape(c,1)) :: c integer intent(out) :: info end subroutine potri