[Haskell-cafe] Numerical Analysis

[Pierre-Etienne Meunier]

Hello Cafe,

Being a complete beginner in the field of numerical analysis, but anyway needing it to solve "real problems", I wrote a few functions recently to solve systems of polynomial equations using the "projected polyhedron" method by Maekawa and Patrikakalis.
This requires solving systems of linear equations precisely, thus the simple Gauss method was not enough, and I had to write also an algorithm for the "SVD decomposition".

Upon discovering the algol / fortran specifications of these :-( , heavily published in important journals, I thought it would be nice to provide the world with fast reliable implementation of these numerical methods (i.e. not simply bindings to lapack). Moreover, writing numerical things in haskell is much more pleasant than I thought at first. Here are a few random thoughts on this :

- The haskell 98 norm does not require enough about IEEE-754 compliance, thus C bindings are still needed to guess for instance the machine epsilons, or manipulating ULPs. Moreover, taking advantage of hardware rounding is not easy, even if the hardware is IEEE-compliant : calling a C function from haskell screws up the speed advantages of hardware rounding, for instance. Maybe the new LLVM backend will make this possible ?

- The current Array library is definitely not adapted to production code. It makes debugging tricky, requires a heavy use of Debug.Trace to actually see what happens, and does not seem as fast as one could expect. It seems that each algebra library on hackage redefines part of it, but a unified version would be nice : a discussion within the haskell community seems to be needed...

- A numerical analysis library should really take advantage of the parallelism in GHC, especially with the arrival of hardware such as fermi (anyway, I do not know how much haskell is compilable to fermi code). The love for loops and side-effects among this community is hard to understand, but that's more of a cultural problem.

Finally, as stated by William Stein, the creator of SAGE, of course it would take thousands of man-years to rewrite these codes in python, but if a language like haskell, and a compiler do 90% of the work, how many man-years are left ?

If anyone here has got the time, the team and the will to start such a project, I'd love to contribute !

Cheers,
PE