Gaussian RNG?, numerics in general...



>>>>> "Ray" == Raymond Toy <toy@rtp.ericsson.se> writes:

>>>>> "Liam" == Liam M Healy <Liam.Healy@nrl.navy.mil> writes:
    Liam> #+allegro
    Liam> (foreign-functions:def-foreign-call
    Liam> ;;; See http://sources.redhat.com/gsl/ref/gsl-ref_7.html#SEC94
    Liam> (jacobian-elliptic-function-int "gsl_sf_elljac_e")
    Liam> ((u :double)
    Liam> (m :double)
    Liam> (sn (:array :double)) (cn (:array :double)) (dn (:array :double)))
    Liam> :convention :stdcall
    Liam> :returning :int)

    Ray> I assume gsl_sf_elljac_e is declared in C to be

    Ray> (double u, double m, double* sn, double* cn, double* dn)

Exactly.

    Ray> So you really have to have an array for this?  I'm not familar with
    Ray> ACL's FFI.

To tell you the truth, I don't know if you _have_ to.  I a) needed to
get some elliptic function results from my CL b) found this worked in
ACL.  ACL experts can tell you if you there's a better way.


    Ray> CMUCL's version would be something like

    Ray>   (alien:def-alien-routine ("gsl_sf_elljac_e" jacobian-elliptic-function-int) int
    Ray>     (u double :in)
    Ray>     (m double :in)
    Ray>     (sn double :out)
    Ray>     (cn double :out)
    Ray>     (dn double :out))

    Ray> I didn't test this, though.

    Liam> to do this by hand.  Second, if m>1, the elliptic functions are
    Liam> undefined.  GSL signals an error, but it is handled by default by

    Ray> A&S gives a transformation for m > 1 to the case where m < 1, but
    Ray> there are constraints then on the range of u.  This is one deficiency
    Ray> of GSL:  it doesn't handle complex elliptic functions.  maxima should
    Ray> strive for this whenever possible and reasonable.

Though interesting, that's not the issue, really.  The question is,
how do you handle conditions signalled from the foreign code?

    Liam> Having said that, I think using the existing library is the only way
    Liam> to go.  The suggestions to reimplement in CL, while having a certain
    Liam> appeal, are just not practical.  We need to leverage the existing
    Liam> resources, and reimpmelementaion would take a huge amount of time and
    Liam> effort, and on an ongoing basis as well.

    Ray> I still think f2cl is a decent solution, but I guess I'm in the
    Ray> minority here.  There are Fortran implementations for just about
    Ray> everything we would want and I suspect that f2cl could convert all of
    Ray> them correctly.

I guess I won't convince you otherwise, but you will need to duplicate
the amount of effort that they are putting into the C library.  The
amount of work put into GSL is impressive - I don't see any other
project trying to make a coherent, reasonably well-documented, up to
date, tested scientific/math library.  And it appears there are a
dozen or so active (if not full time, significant time) participants.
I prefer, where possible, to build on others' work instead of
duplicating it.

Though obviously there are serious obstacles to my approach as well! 


    Ray> I think there's more than enough for everyone to do with maxima, the
    Ray> symbolic math program, instead of maxima, the amazingly fast number
    Ray> cruncher.  I'm more than happy with slow, correct numbers.  If I want
    Ray> fast, I'll use octave, matlab, scilab, or whatever.

My intersection with maxima is pretty small.  I see this as well
beyond the world of maxima (hence my previous suggestion to take this
to comp.lang.lisp).  I do extensive sci/num computation in CL
and other languages, and only a little in Maxima.  I want good,
robust, tested, documented, actively supported libraries that I can
use from CL.  This has only incidentally to do with Maxima.

I don't understand your statement about "slow, correct numbers".  Why
not fast, correct numbers, from libraries that have widespread use,
and are documented and tested?


Liam