[SciPy-user] Wrapping TWPBVPLC BVP solver
Pauli Virtanen
pav at iki.fi
Sun Nov 16 18:55:07 EST 2008
Sun, 16 Nov 2008 15:03:02 -0800, John Salvatier wrote:
> I want to wrap the Fortran Boundary Value Problem Solver TWPBVPLC
> (http://
> www.ma.ic.ac.uk/~jcash/BVP_software/readme.php<http://www.ma.ic.ac.uk/%
7Ejcash/BVP_software/readme.php>)
> which I believe I have
> arranged to have released under a BSD license. Has someone else done
> this before? I don't want to repeat work unnecessarily. I want to wrap
> the Fortran Boundary Value Problem Solver TWPBVPLC (http://<br> <a
> href="http://www.ma.ic.ac.uk/%7Ejcash/BVP_software/readme.php"
> target="_blank">www.ma.ic.ac.uk/~jcash/BVP_software/readme.php</a>)
> which I believe I have<br> arranged to have released under a BSD
> license. Has someone else done<br> this before? I don't want to
> repeat work unnecessarily.
Not TWPBVP*, but COLNEW is here:
http://www.iki.fi/pav/software/bvp/index.html
COLNEW is of course non-free, but TWPBVP should go along the same lines
-- the Python wrapper is BSD, so you can lift whatever you want from the
wrapper part. Just avoid reading anything in Fortran from there, to stay
within a clean room environment :)
I think it would be useful if you considered these things when writing
the TWPBVP* wrapper:
* To cut down Python <-> Fortran call overhead, modify TWPBVP* so that
FSUB and DFSUB evaluate the result for all X-points in the mesh in one
call. Also collapse GSUB and DGSUB by vectorizing the `i` index away.
This makes writing efficient vectorised Python code quite a bit easier.
* Make DFSUB and DGSUB optional, by falling back to simple numerical
differentiation if the user omits them. Writing the DGSUB and DFSUB
takes typically much more work than GSUB and FSUB, and it's nice not
to have to do it when it's not necessary. (I think you can lift some
code from this directly from the colnew wrapper.)
I remember Ascher and Bader (?) mentioning in some paper that naive
numerical differentiation typically is enough, and this indeed seems to
be usually the case.
--
Pauli Virtanen
More information about the SciPy-User
mailing list