[SciPy-user] Wrapping geqp3 (Rank revealing QR)

Nils Wagner nwagner at mecha.uni-stuttgart.de
Fri Feb 24 07:55:50 CST 2006


Hi all,
 
I was trying to wrap geqp3...(my homework as suggested by Arnd ;-) )
It turns out that it is not straightforward (at least for me).

I have used an example taken from generic_flapack.pyf

   subroutine <tchar=s,d,c,z>geqrf(m,n,a,tau,work,lwork,info)

   ! qr_a,tau,work,info = geqrf(a,lwork=3*n,overwrite_a=0)
   ! Compute a QR factorization of a real M-by-N matrix A:
   !   A = Q * R.

     callstatement (*f2py_func)(&m,&n,a,&m,tau,work,&lwork,&info)
     callprotoargument
int*,int*,<type_in_c>*,int*,<type_in_c>*,<type_in_c>*,int*,int*

     integer intent(hide),depend(a):: m = shape(a,0)
     integer intent(hide),depend(a):: n = shape(a,1)
     <type_in> dimension(m,n),intent(in,out,copy,out=qr) :: a
     <type_in> dimension(MIN(m,n)),intent(out) :: tau

     integer optional,intent(in),depend(n),check(lwork>=n||lwork==-1) ::
lwork=3*n
     <type_in> dimension(MAX(lwork,1)),intent(out),depend(lwork) :: work
     integer intent(out) :: info
   end subroutine <tchar=s,d,c,z>geqrf


Some observations:
The number of arguments differ wr.t. to complex (C,Z) and real matrices
(S,D).
How is this handled with scipy ?

      SUBROUTINE CGEQP3( M, N, A, LDA, JPVT, TAU, WORK, LWORK, RWORK,
     $                   INFO )
      SUBROUTINE DGEQP3( M, N, A, LDA, JPVT, TAU, WORK, LWORK, INFO )



Here is my first (incomplete !!) attempt to adapt  the wrapper .  How
can I proceed ?

!    -*- f90 -*-
! Note: the context of this file is case sensitive.

python module wrap_lap ! in
    interface  ! in :wrap_lap
        subroutine <tchar=s,d>geqp3(m,n,a,lda,jpvt,tau,work,lwork,info)
! in :wrap_lap:dgeqp3.f

        ! rrqr_a,tau,work,info = geqp3(a,lwork=...,overwrite_a=0)
        ! Compute the a QR factorization with column pivoting of a matrix A
        ! A*P = Q*R  using Level 3 BLAS

        callstatement (*f2py_func)(&m,&n,a,&m,tau,work,&lwork,&info)
        callprotoargument int*,int*,<type_in_c>*,int*,<type_in_c>*,int*,int*

            integer intent(hide),depend(a):: m=shape(a,0)
            integer intent(hide),depend(a):: n=shape(a,1)
            <type_in> dimension(m,n), intent(in,out,copy,out=rrqr) :: a
            <type_in> dimension(MIN(m,n)),intent(out) :: tau

            integer dimension(*) :: jpvt
            double precision dimension(*), intent(out) :: tau
            double precision dimension(*), intent(out) :: work
            <type_in>
dimension(lwork),intent(hide,cache),depend(lwork):: work
            integer intent(out):: info
        end subroutine <tchar=s,d>geqp3
    end interface
end python module wrap_lap

! This file was auto-generated with f2py (version:2_2163).
! See http://cens.ioc.ee/projects/f2py2e/
~

Nils

 



More information about the SciPy-user mailing list