Hi,
do you consider this answer below the right answer?
Maxima branch_5_30_base_103_g0860c83_dirty http://maxima.sourceforge.net
using Lisp SBCL 1.1.1.0.debian
(%i1) fpprintprec:12$
(%i2) r:0$
(%i3) f(x):=(if r<5 then r:r+1, r)$
(%i4) rk('(f(x)),x,0,[t,0,1,0.2]);
(%o4) [[0.0, 0.0], [0.2, 0.5], [0.4, 1.5], [0.6, 2.5], [0.8, 3.5], [1.0,
4.5]]
Notice that I'm using Maxima 5.30post. In Maxima 5.27 rk() was not able
to handle the kind of thing that you want to do. And notice also that
you do not have to load dynamics (it will be loaded automatically the
first time you use rk). Please update your version of Maxima, or you you
cannot do that, locate the diretory share/dynamics in your system and
replace the files dynamics.mac and complex-dynamics.lisp in that
directory with their newest versions from the Maxima Git repository.
In any case, did you expect f(x) to have the values 1, 2, 3, 4, 5 and 5
at the points 0, 0.2, 0.4, 0.6, 0.8 and 1?
It will not work that way because at each step rk will calculate f(x)
for times. Different implementations of the same 4th order Runge-Kutta
algorithm will give different results for your f(x) function. For
example, if for some reason the algorithm calculated f(x) more than 4
times at one step (perhaps to check for the expected error of the
result?) your f(x) function would become equal to 5 before you expected.
Regards,
Jaime
On 18-09-2013 14:49, Ether Jones wrote:
>
>
> What is the proper method in Maxima to integrate a stateful function
> with rk()?
>
> Below is a greatly simplified example to illustrate the problem. %i4
> gives the wrong answer, and %i5 crashes.
>
> Thank you.
>
> ------------------ begin example -----------------------------
>
> Maxima 5.27.0 http://maxima.sourceforge.net
> using Lisp GNU Common Lisp (GCL) GCL 2.6.8 (a.k.a. GCL)
> Distributed under the GNU Public License. See the file COPYING.
> Dedicated to the memory of William Schelter.
> The function bug_report() provides bug reporting information.
>
> (%i1) load("dynamics")$
>
> (%i2) r:0$
>
> (%i3) f(x):=(if r<5 then r:r+1, r)$
>
> (%i4) rk(f(x),x,0,[t,0,1,0.2]);
>
> (%o4) [[0, 0], [0.2, 0.2], [0.4, 0.4], [0.6, 0.6], [0.8, 0.8], [1.0, 1.0]]
>
> (%i5) rk('(f(x)),x,0,[t,0,1,0.2]);
>
> Expecting a number when the initial state is replaced in the
> equations, but instead found:
> r
> #0: rk(odes=f(x),state=x,initial=0,domain=[t,0,1,0.2])
> -- an error. To debug this try: debugmode(true);
>
> -------------------- end --------------------------------------
>