[Shootout-list] Science-related benchmarks

[Sebastien Loisel]

Jon, I'm answering your email and also describing what I've worked on so fa=
r.

Aside from the integrator you've tweaked, I've proposed the new
spectral norm test and the implicitode tests; my next bench is to be
something multi-dimensional involving either LU, sparse LU or GMRES
and other algorithms.

My goal was to add some numerical benchmarks that used nontrivial sets
of numerical algorithms in an interrelated way, with interesting data
structures and hopefully some polymorphism. The test implicitode is a
(possibly flawed) example of what I was thinking of. It solves several
odes using the trapezoidal method, with a newton solver which uses
automatic differentation(*) to compute f' with a variety of scalar
fields to force polymorphism.

I think the polymorphism in particular has been poorly received
because it makes the implementation longer in most languages and more
importantly more difficult to understand. I will take that into
account in future tests.

On 4/25/05, Jon Harrop <jon@ffconsultancy.com> wrote:
> On Monday 25 April 2005 15:28, Sebastien Loisel wrote:
> > Looks good to me. Why didn't you use the conjugate gradient method thou=
gh?
>=20
> IIRC, Newton's method is actually better here because it is more numerica=
lly

Newton is almost always better but it has a matrix inverse, so people
don't use it for large scale problems, that's why I suggested CG.

> The Mathematica JIT compiler I wrote for Wolfram Research included this
> functionality. Coupled with the unusually powerful pattern matching

Interesting. Did that particular bit end up in the consumer product?

> Yes, I'm not sure it would be FP bound. Perhaps a Newton-Raphson based qu=
adric

I think that's a good idea.

> We could try FP intensive tests (e.g. bounding box) instead and see if th=
ey
> were faster. I'd be surprised if they were though, unless we're also doin=
g

I don't think they'd be faster, but the algorithm would be more
representative of a real algorithm.

> The C++ was virtually identical - just cosmetic improvements. I shaved 10=
% or
> so off the OCaml (mainly by type specialising the heap, IIRC).

In that case, for apples to apples, it would be good to write a
specialized heap for C++ as well, although I suspect there won't be a
gain. Also, I saw you preallocate the heap in the ocaml, I think we
should either disallow that or make everyone preallocate it.

Sebastien Loisel

(*) I know this is not a reasonable place to use AD. My ODE algorithms
typically use AD for sensitivity analysis, but this would've made the
program bigger.