Mercurial > octave
changeset 23407:e265ae9e7a6c
Fix hermite cubic interpolation in ode23
* scripts/ode/private/runge_kutta_interpolate.m : fix indices of derivatives being
passed to the cubic interpolation algorithm.
* scripts/ode/ode23.m : add a test to check for accuracy of interpolated values.
| author | Andreas Stahel <Andreas.Stahel@bfh.ch> |
|---|---|
| date | Tue, 18 Apr 2017 12:04:30 +0200 |
| parents | aab07d7e98be |
| children | 0af9a1ae0912 |
| files | scripts/ode/ode23.m scripts/ode/private/runge_kutta_interpolate.m |
| diffstat | 2 files changed, 6 insertions(+), 2 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ode/ode23.m Mon Apr 17 18:14:40 2017 -0400 +++ b/scripts/ode/ode23.m Tue Apr 18 12:04:30 2017 +0200 @@ -413,6 +413,10 @@ %! opt = odeset ("AbsTol", 1e-8, "RelTol", 1e-8); %! sol = ode23 (@fpol, [0 2], [2 0], opt); %! assert ([sol.x(end); sol.y(:,end)], [2; fref'], 1e-3); +%!test # hermite_cubic_interpolation +%! opt = odeset ("RelTol", 1e-8, "NormControl", "on"); +%! [t,sol] = ode23(@(t,x)[x(2);x(1)],linspace(0,1),[1;0],opt); +%! assert(max(abs(sol(:,1)-cosh(t))),0,1e-6) %!test # RelTol and NormControl option -- higher accuracy %! opt = odeset ("RelTol", 1e-8, "NormControl", "on"); %! sol = ode23 (@fpol, [0 2], [2 0], opt);
--- a/scripts/ode/private/runge_kutta_interpolate.m Mon Apr 17 18:14:40 2017 -0400 +++ b/scripts/ode/private/runge_kutta_interpolate.m Tue Apr 18 12:04:30 2017 +0200 @@ -72,8 +72,8 @@ s = (t_out - t(1)) / dt; x_out = ((1 + 2*s) .* (1-s).^2) .* x(:,1) + ... (s .* (1-s).^2 * dt ) .* der(:,1) + ... - ((3-2*s) .* s.^2 ) .* x(:,2) + ... - ((s-1) .* s.^2 * dt ) .* der(:,2); + ((3-2*s) .* s.^2 ) .* x(:,end) + ... + ((s-1) .* s.^2 * dt ) .* der(:,end); endfunction