XProfan

Window Title "Demo: Runge-Kutta RK4 for systems of Differentialgleichungen höherer Ordnung"
' Angelehnt on: B.Brand: Algorithms to virtually maths, s.240, Oldenbourg
' Testhalber migriert to XProfan-11 in 2014-09 by P.woodpecker, Wien.
' demonstration for private use, no however geartete Gewähr!
Window Style 24:Window 0,0-%maxx,%maxy
declare k1!,k2!,k3!,k4!,l1!,l2!,l3!,l4!
':: example: 'Nichtlinares Differentialgleichungssystem'
':: two Testfunktionen, a gekrümmte expanse form (X: Abszisse, Y: Ordinate, Z: Höhe)
declare x!,y!,z!,f!,g!
declare i&,h!,x![],y![],z![],xend!

proc F :parameters x!,y!,z!

    return (y!-z!)/x!

endproc

proc G :parameters x!,y!,z!

    return (y!+z!)/x!

endproc

' Info: Lines, The with :: begin, stammen of Übersetzter to XProfan-11
print "\n in the example becomes one system from '2 nichtlinearen Differentialgleichungen  "
print " in 3 variables' resolved. The Gleichungen lauten F=(y+z)/x and G=(y-z)/x.   "
print " The Anfangswerte for x, y and z as well as The Schrittweite h are vorzugeben. "
print " moreover becomes one Abbruchkriterium, z.B. the Maximale x-Endwert association. \n"
print " from Found of my comfort go with input 'x=0' The Beispiels-  "
print " values the Originalvorlage taken (Results einsehbar in the Programmtext).\n"
' The Input-values the Original-Beispiels:
'  x!=1 : y!=1 : z!=0 : h!=0.2 : xend!=3
' make for following spending:
'step  X           Y                 Z
' 0       1     1                 0
' 1       1.2   1.18011937557392  0.217584940312213
' 2       1.4   1.32150759540037  0.462240502741668
' 3       1.6   1.42651937035306  0.724648391007437
' 4       1.8   1.49791768072086  0.998168848790414
' 5       2     1.53848713605723  1.27796059664382
' 6       2.2   1.55087497139788  1.56043022121442
' 7       2.4   1.5375315540692   1.84287395296026
' 8       2.6   1.50069575905948  2.12323677414932
' 9       2.8   1.44240050556552  2.39994581694282
' 10      3     1.36448684585809  2.67179234356735
' ---
':: ACHTUNG! too the Endkriterium is each on The Gleichungsstruktur anzupassen,
':: otherwise rechnet the thing, To the memory fully is!!!
print " X = ";:input x! : if x!=0:x!=1:y!=1:z!=0:h!=0.2:xend!=3.0:goto "skip":endif
print " Y = ";:input y!
print " Z = ";:input z!
print " Schrittweite h = ";:input h!
print " Endwert (here: X_end) = ";:input xend!
skip:
set("decimals",17)
clear x![],y![],z![]:clearclip
x![0]=x!:y![0]=y!:z![0]=z!
spending:
print "\n Schritt";tab(25);"    X";tab(50);"    Y";tab(75);"    Z":print mkstr$("-",90)

Repeat

    print "  ";i&,tab(25);stature$("%g",X![i&]);tab(50);stature$("%g",Y![i&]);tab(75);stature$("%g",Z![i&])
    ::putclip st$(i&)+" "+stature$("%g",X![i&])+" "+stature$("%g",Y![i&])+" "+stature$("%g",Z![i&])+"\n"
    k1!=h!*F(x![i&],y![i&],z![i&])
    l1!=h!*G(x![i&],y![i&],z![i&])
    k2!=h!*F(x![i&]+h!/2,y![i&]+k1!/2,z![i&]+l1!/2)
    l2!=h!*G(x![i&]+h!/2,y![i&]+k1!/2,z![i&]+l1!/2)
    k3!=h!*F(x![i&]+h!/2,y![i&]+k2!/2,z![i&]+l2!/2)
    l3!=h!*G(x![i&]+h!/2,y![i&]+k2!/2,z![i&]+l2!/2)
    k4!=h!*F(x![i&]+h!,  y![i&]+k3!,z![i&]+l3!)
    l4!=h!*G(x![i&]+h!,  y![i&]+k3!,z![i&]+l3!)
    y![i&+1]=y![i&]+1/6*(k1!+2*k2!+2*k3!+k4!)
    z![i&+1]=z![i&]+1/6*(l1!+2*l2!+2*l3!+l4!)
    x![i&+1]=x![i&]+h!
    inc i&
    ::if %csrlin>55:waitinput' 55 and 85 on actually Schirmdimensionen adjust!
    ::cls:print:print " Schritt";tab(25);"    X";tab(50);"    Y";tab(75);"    Z":print mkstr$("-",90)
    ::endif

until x![i&]>(xend!*1.000000000000001)'(because of höherer accuracy of Intel-FPUs)

::print "OK.":beep:print " spending standing too in Clipboard!":waitinput
End