Sunday, March 25, 2012

Parametric IVP solver

Building on lazy analytic functions Veronica and Jan have now completed an implementation of a generic solver for initial value problems. The ideas for the solver are borrowed from the work on validated numerics by Davis and Tucker. The solver computes an interval enclosure of the solution and takes the following parameters:
a)  the field as a function of type (A,A) => A for any type A that can be used for endpoints of intervals and as a representation of the real numbers,
b) one such type
c) the initial value of type A
d) the degree for the Taylor expansion used to approximate the solution at each step
e) the step size
Our implementations are in Scala where we could make use parametric polymorphism, overloading and a lazy datastructure in the procedure for automatic differentiation. As a result the taylor coefficients for the sought function needed at the points where the solution is calculated and at intervals are expressed as a mutually recursive definition using the field and the initial value in the lazy data structure. We have performed some experiments with both constructive reals (a Java implementation by Boehm) and BigDecimals (a Java library). 

The following graph summarizes one experiment with the IVP given by x'(t) = -x(t)^3 with x(0)=1up to t = 4. In the graph below we show the execution time and the width of the enclosing interval at time = 4 for  constructive reals and big decimals with precision 10. 


Future work includes a) completing the library for using big decimals as endpoints, b) completing the library for intervals and c) using these libraries in an implementation of a Picard solver. We will also start writing about the program, the experiments and how to use the solver in the context of the parametric version of Acumen.