better numerical checks

headmyshoulder / odeint-v2

odeint - solving ordinary differential equations in c++ v2

http://headmyshoulder.github.com/odeint-v2/

Other

341 stars 101 forks source link

better numerical checks #145

Closed mariomulansky closed 9 years ago

mariomulansky commented 9 years ago

in light of ee3ea15d5810b6b8d8e34cfada20cb58f01ad303 one could try to check the order of the numerical algorithms by using an ODE that has a polynomial of the same degree of the order of the stepper. The approximate solution should then be exact.

GregorDeCillia commented 9 years ago

As far as i can tell, it should be possible to test all your solvers like this. Lately, I created a script to demonstrate this idea on my git page. Here's the output with the newest odeint version. The elements of the martix show the error of solving x'(t)=t^p from 0 to 1 with 10 steps of size 0.1

If the solvers are implemented correctly, the errors will be just due to machine precision if p+1 <= ord, wich happens to be around 10^-16 for dobles, see wiki.

mariomulansky commented 9 years ago

Yes, this is a very nice idea. If you want, we would very much appreciate if you could convert this code into a unit test with boost.test. Using https://github.com/headmyshoulder/odeint-v2/blob/master/test/numeric/runge_kutta.cpp as a starting point this seems possible without too much effort I think. I'd happily merge a pull request :)

Otherwise I will go ahead and implement this myself at some point.

Thanks again for your help!

GregorDeCillia commented 9 years ago

Ok, I'd like to do that myself especially because of some numerical details (e.g. choosing tolerances, appropriate steppers for initialization etc). It will take me a week or so tough because of the upcoming holidays.

mariomulansky commented 9 years ago

great! no need to rush at all though, enjoy the holidays!

GregorDeCillia commented 9 years ago

thanks. Have a nice time too:)