[SciPy-dev] newscipy, band matrix routines
Nils Wagner
nwagner at mecha.uni-stuttgart.de
Wed Nov 9 04:36:41 CST 2005
Arnd Baecker wrote:
>Hi,
>
>a while ago I created f2py wrappers for
>- dsbev
>- dsbevd
>- dgbtrf
>- dgbtrs
>- dgbrfs
>which deal with eigenvalue/vector computation,
>LU decomposition, ... for band matrices.
>
>The present version of this, including examples/tests and
>a little bit of documentation is at
> http://www.physik.tu-dresden.de/~baecker/python/Band2.zip
>
>I would like to get this into newscipy but we need some assistance:
>- should the wrappers be added to generic_flapack.pyf ?
>
> (+adding the new functions to setup_linalg.py so that the
> corresponding single precision wrappers can be skipped)
>
>- Are the only scipy specific tricks the use of
> <tchar=s,d,c,z> and <type_in>
> to support the different types?
>- For testing/development purpuses, is there
> a way to put this into a separate sandbox
> or is it just easier to use all the tools which
> are provided by Lib/linalg/?
>- it would be great if an expert (Pearu - are you listening? ;-)
> has a look at the wrapper (either already now, or
> when the above steps are done) - see also the questions below.
>
>Any further guidance/comments/suggestions are welcome!
>
>Best, Arnd
>
>
>Some more technical questions
>
>Concerning dsbevx:
> a) abstol parameter:
> from dsbev.f:
> "Eigenvalues will be computed most accurately when ABSTOL is
> set to twice the underflow threshold 2*DLAMCH('S'), not zero."
>
> Maybe the easiest (best?) is to wrap DLAMCH as well (easy)
> and let the user provide the value, if necessary?
> ((Even nicer would be to compute this value internally -
> is this possible with f2py ? If not, I am not sure if adding this
> is worth the effort...))
>
> b) Presently I suppress the array q
> Q (output) DOUBLE PRECISION array, dimension (LDQ, N)
> If JOBZ = 'V', the N-by-N orthogonal matrix used in the
> reduction to tridiagonal form.
> If JOBZ = 'N', the array Q is not referenced.
> Should one better return this as well and therefore replace
> double precision dimension(ldq,ldq),intent(hide),depend(ldq) :: q
> by
> double precision dimension(ldq,ldq),intent(out),depend(ldq) :: q
> ?
>
> c) Should one make vl,vu,il,iu, optional as well ?
> Background: dsbevx allows to compute
> range: 0 = 'A': all eigenvalues will be found;
> 1 = 'V': all eigenvalues in the half-open interval (vl,vu]
> will be found;
> 2 = 'I': the il-th through iu-th eigenvalues will be found.
> so we may have that either vl,vu or il,iu or neither of
> the two pairs are used.
>
> If one makes all vl,vu,il,iu, optional, one could set
> set the defaults to
> il=0 iu=10
> vl=0.0 vu=20.0 # this of course absolutely
> # dependent on the problem,
> # i.e in general nonsense...
>
>
>_______________________________________________
>Scipy-dev mailing list
>Scipy-dev at scipy.net
>http://www.scipy.net/mailman/listinfo/scipy-dev
>
There is also some progress in improving LAPACK in this context.
http://www.cs.berkeley.edu/~demmel/Sca-LAPACK-Proposal.pdf
Nils
More information about the Scipy-dev
mailing list