Week 8

   

Progress this Week

This week I worked on the Runge Kutta Fehlberg numerical integration algorithm. While the integrator is set up and running, some of the test cases fail due to issues with error estimation. 

For example, in model 3 in the core SBML test suite, the simulator calculates the values as follows:

$ cargo run ../testsuites/core-semantic/00003/00003-sbml-l3v2.xml 5 50 1e-4 1e-5
5 seconds with 50 steps.
t               S1                      S2
0.000000        0.0150000000000000      0.0000000000000000
0.100000        0.0135747032487483      0.0028505935025034
0.200000        0.0122848378861052      0.0054303242277897
0.300000        0.0111175352508572      0.0077649294982856
0.400000        0.0100611494591923      0.0098777010816154

Whereas, the expected values are:

time            S1                      S2

0               0.015                   0

0.1             0.01357256127053939     0.002854877458921213

0.2             0.01228096129616973     0.005438077407660545

0.3             0.01111227331022577     0.007775453379548465

0.4             0.01005480069053459     0.009890398618930821

At t = 0.4, for species S2, even though the estimated error is within the tolerance in the simulator, this is not the case when the test is being run. 

The tolerance during the test is calculated as follows:

tol = relative_tol * S2                   + absolute_tol

    = 1e-4         * 0.009890398618930821 + 1e-5
    = 1.0989039861893083e-05

error = abs(expected_values - simulator_result)

      = abs(0.009890398618930821, 0.0098777010816154)

      = 1.2697537315422072e-05

I have gone through the following reference implementations so far, but I have been unable to point out exactly what the issue is. I'll continue working on this in the coming week.

1. Numerical methods using MATLAB/John H. Mathews, Kurtis D. Fink—4th ed.
2. Python implementation on Dr John Burkardt's website (link)
3. Python implementation by amirh0ss3in on Github (link)
4. C++ implementation in Boost (link)
5. SciPy reference (link)

Progress so Far

As per discussion with my project mentors, we had diverged from the initially proposed timeline. The first half was completed as proposed and the next two weeks were spent in building support for more SBML features instead of the Runge-Kutta-Fehlberg integrator. Work on the RKF integrator was started in week 8 and will continue in the coming week.

Plan for the Coming Week

Fix errors in the Runge-Kutta-Fehlberg integrator