Runge-Kutta method: Difference between revisions

</pre>

=={{header|REXX}}==

<lang rexx>/*REXX program uses the Runge-Kutta method to solve the differential */

/*equation: y'(t)=t²√y(t) which has the exact solution: y(t)=(t²+4)²/16*/

numeric digits 40; d=digits()%2 /*use forty digits, show ½ that. */

x0=0; x1=10; dx=.1; n=1 + (x1-x0) / dx; y.=1

do m=1 to n-1; mm=m-1

y.m=Runge_Kutta(dx, x0+dx*mm, y.mm)

end /*m*/

say center(x,13,'─') center(y,d,'─') ' ' center('relative error',d,'─')

do i=0 to n-1 by 10; x=(x0+dx*i)/1; y2=(x*x/4+1)**2

relE=format(y.i/y2-1,,13)/1; if relE=0 then relE=' 0'

say center(x,13) right(format(y.i,,12),d) ' ' left(relE,d)

end /*i*/

exit /*stick a fork in it, we're done.*/

/*──────────────────────────────────RATE subroutine─────────────────────*/

rate: return arg(1)*sqrt(arg(2))

/*──────────────────────────────────Runge_Kutta subroutine──────────────*/

Runge_Kutta: procedure; parse arg dx,x,y

k1 = dx * rate(x , y )

k2 = dx * rate(x+dx/2 , y+k1/2 )

k3 = dx * rate(x+dx/2 , y+k2/2 )

k4 = dx * rate(x+dx , y+k3 )

return y + (k1+2*k2+2*k3+k4) / 6

/*──────────────────────────────────SQRT subroutine─────────────────────*/

sqrt: procedure;parse arg x; if x=0 then return 0; d=digits()

numeric digits 11; g=.sqrtG()

do j=0 while p>9; m.j=p; p=p%2+1; end; do k=j+5 to 0 by -1

if m.k>11 then numeric digits m.k

g=.5*(g+x/g); end; numeric digits d; return g/1

.sqrtG: numeric form; m.=11; p=d+d%4+2

parse value format(x,2,1,,0) 'E0' with g 'E' _ .; return g*.5'E'_%2</lang>

'''output'''

<pre style="overflow:scroll">

──────X────── ─────────Y────────── ───relative error───

0 1.000000000000 0

1 1.562499854278 -9.3262010934636E-8

2 3.999999080521 -2.2986980018794E-7

3 10.562497090438 -2.7546153355698E-7

4 24.999993765091 -2.4939637458826E-7

5 52.562489180303 -2.058444217357E-7

6 99.999983405404 -1.6594596403363E-7

7 175.562476482271 -1.3395644712797E-7

8 288.999968434799 -1.092221504E-7

9 451.562459276840 -9.018277747572E-8

10 675.999949016709 -7.5419068845528E-8