--- a/src/DLD-FUNCTIONS/lsode.cc
+++ b/src/DLD-FUNCTIONS/lsode.cc
@@ -475,6 +475,78 @@
 }
 
 /*
+
+%% dassl-1.m
+%%
+%% Test lsode() function
+%%
+%% Author: David Billinghurst (David.Billinghurst@riotinto.com.au)
+%%         Comalco Research and Technology
+%%         20 May 1998
+%%
+%% Problem
+%%
+%%    y1' = -y2,   y1(0) = 1
+%%    y2' =  y1,   y2(0) = 0
+%%
+%% Solution
+%%
+%%    y1(t) = cos(t)
+%%    y2(t) = sin(t)
+%!function xdot = f (x, t)
+%!  xdot = [-x(2); x(1)];
+%!test
+%! 
+%! x0 = [1; 0];
+%! xdot0 = [0; 1];
+%! t = (0:1:10)';
+%! 
+%! tol = 500 * lsode_options ("relative tolerance");
+%! 
+%! 
+%! x = lsode ("f", x0, t);
+%! 
+%! y = [cos(t), sin(t)];
+%! 
+%! assert(all (all (abs (x - y) < tol)));
+
+%!function xdotdot = f (x, t)
+%!  xdotdot = [x(2); -x(1)];
+%!test
+%! 
+%! x0 = [1; 0];
+%! t = [0; 2*pi];
+%! tol = 100 * dassl_options ("relative tolerance");
+%! 
+%! x = lsode ("f", x0, t);
+%! 
+%! y = [1, 0; 1, 0];
+%! 
+%! assert(all (all (abs (x - y) < tol)));
+
+%!function xdot = f (x, t)
+%!  xdot = x;
+%!test
+%! 
+%! x0 = 1;
+%! t = [0; 1];
+%! tol = 100 * dassl_options ("relative tolerance");
+%! 
+%! x = lsode ("f", x0, t);
+%! 
+%! y = [1; e];
+%! 
+%! assert(all (all (abs (x - y) < tol)));
+
+%!test
+%! lsode_options ("absolute tolerance", eps);
+%! assert(lsode_options ("absolute tolerance") == eps);
+
+%!error <Invalid call to lsode_options.*> lsode_options ("foo", 1, 2);
+
+*/
+
+/*
 ;;; Local Variables: ***
 ;;; mode: C++ ***
 ;;; End: ***