\ Numerical solution of first-order ordinary differential equation \ by 4th order Runge-Kutta algorithm following \ Abramowitz & Stegun p. 896, 25.5.10 \ --------------------------------------------------- \ (c) Copyright 1999 Julian V. Noble. \ \ Permission is granted by the author to \ \ use this software for any application pro- \ \ vided this copyright notice is preserved. \ \ --------------------------------------------------- \ This is an ANS Forth program requiring the \ FLOAT, FLOAT EXT, FILE and TOOLS EXT wordsets. \ Environmental dependencies: \ Assumes independent floating point stack \ Solves dx/dt = f(x,t) \ For simplicity and clarity, the program is intended to be compiled \ using a FORmula TRANslator \ Example \ : fnb ( f: x t -- t^2*exp[-x]) \ f^2 FSWAP FNEGATE FEXP F* ; \ use( fnb 0e0 0e0 5e0 0.1e0 )runge MARKER -runge INCLUDE ftran111.f : fvariables 0 DO FVARIABLE LOOP ; 9 fvariables t h h2 x tf xk1 xk2 xk3 xk4 v: fdummy : x' \ integration step f" xk1 = h2*fdummy(x,t)" \ compute k1 f" t = t+h2 " \ t -> t + h/2 f" xk2 = h*fdummy(x+xk1,t)" f" xk3 = h*fdummy(x+xk2/2,t)" f" t = t+h2 " \ t -> t + h/2 f" xk4 = h2*fdummy(x+xk3,t)" f" x = x + (xk1 + xk2 + xk3 + xk4)/3 " ; : display CR t F@ FS. 5 SPACES x F@ FS. 5 SPACES f" ln(1+t^3/3)" FS. ; : done? t F@ tf F@ F> ; : )runge ( xt --) ( f: x0 t0 tf h -- ) defines fdummy \ vector fn_name h F! tf F! t F! x F! \ initialize variables f" h2 = h/2" BEGIN display x' \ indefinite loop done? UNTIL ;