[SciPy-user] RREF - upper triangular
nwagner at mecha.uni-stuttgart.de
Thu Feb 23 02:17:57 CST 2006
arnd.baecker at web.de wrote:
>On Wed, 22 Feb 2006, Alan G Isaac wrote:
>>On Wed, 22 Feb 2006, Nils Wagner apparently wrote:
>>>here is a code. However it doesn't make use of LAPACK. So
>>>it might be to slow.
>>This exchange suggests to me that Robert thinks that
>>wrapping the LAPACK code is approximately trivial but that
>>Nils does not see it that way. Apologies for any
>But in this particular case previous exposure to wrapping lapack
>functions is granted:
>The number of arguments of DGEQP3 is the same as for dgehrd,
>so on this level it is not more difficult.
>>Anyway, coming from GAUSS as a *user*, I would not know
>>where to begin if I needed to wrap a new LAPACK function.
>Thanks to f2py, wrapping Fortran code is (with a bit of effort)
>trivial in many cases. For complicated functions
>requiring many arguments the wrapper can become longish.
>Fortunately, many things can be learnt from
>looking at ``scipy/Lib/linalg/generic_flapack.pyf``
>In particular, the documentation at
>I also found the f2py notes by Fernando Perez very helpful,
>Let me try to give some general remarks on how to start
>(the real authority on all this is of course Pearu, so
>please correct me if I got things wrong here):
>- first find a routine which will do the job you want:
> - If the lapack documentation is installed properly
> on Linux you could do
> apropos keyword
> - www.netlib.org provides a nice decision tree
>- make sure that that it does not exist in scipy:
> from scipy.lib import lapack
> lapack.clapack.<TAB> (assuming Ipython)
> Remark: routines starting with c/z are for double/single complex
> and routines for d/s for double/single real numbers.
> The calling sequence for c/z and d/s are (I think always) the same and
> sometimes they are also the same for the real and complex case.
>- Then one has to download the fortran file for the lapack routine
> of interest.
>- Generate wrapper by calling
> f2py -m wrap_lap -h wrap_lap.pyf lapack_routine.f
>- Generate library
> f2py -c wrap_lap.pyf lapack_routine.f -latlas -llapack -lblas
>- You can use this by
> import wrap_lap
> Note, that this is not yet polished (this is the part on
> which has to spent some effort ;-), i.e. one
> has to tell which variables are input, which are output
> and which are optional. In addition temporary
> storage has to be provided with the right dimensions
> as described in the documentation part of the lapack routine.
>Concrete (and very simple) example (non-lapack):
>Wrapping Hermite polynomials
>Download code (found after hours of googling ;-)::
>and extract ``hermite.f``
>Generate wrapper framework::
> # only run the following line _once_
> # (and never again, otherwise the hand-modified hermite.pyf
> # goes down the drains)
> f2py -m hermite -h hermite.pyf hermite.f
>Then modify ``hermite.pyf``
>Create the module::
> f2py2.3 -c hermite.pyf hermite.f
> # add this if you want:
> -DF2PY_REPORT_ON_ARRAY_COPY=1 -DF2PY_REPORT_ATEXIT
> import hermite
> import scipy
>A more complicated example about
>how to wrap routines for band matrices can be found at
>under "Python and Co - some recent developments".
>It would be nice if someone could add notes
>on all this to the scipy wiki.
>SciPy-user mailing list
>SciPy-user at scipy.net
If it is so trivial to add the rank revealing QR decomposition to scipy,
I would be glad if the several steps on building a wrapper for DGEQP3
documented on the Wiki as an exercise.
The best way to learn something is to have a well-elaborated example. :-)
More information about the SciPy-user