Metodo Runge Kutta Nystrom
( visitare anche la pagina dedicata all'algoritmo di J.R.Cash, A.H.Karp - ACM Transaction on Mathematical Software, Vol. 16, No. 3, September 1990, Pages 201-222, pił complicato ma pił preciso )
L'algoritmo RKN risolve l'equazione differenziale ai valori iniziali:
``y = f( t , y , `y )
e quindi l'algoritmo RKN serve per risolvere numericamente tutti i problemi di dinamica a parametri concentrati.
Uso qui una formulazione in cui il passo di integrazione è dodici volte il valore del parametro h. In questo modo l'algoritmo non effettua nessuna divisione e pertanto se la funzione f( t , y , `y ) è intera e se sono interi i valori iniziali anche il risultato sarà intero.
ka = 2 h f( t , y , `y )
tab = t + 6 h
yab = y + 6 h `y + 9 h ka
kb = 2 h f( tab , yab , `y + 3 ka )
kc = 2 h f( tab , yab , `y + 3 kb )
kd = 2 h f( t + 12 h , y + 12 h `y + 36 h kb , `y + 6 kc )
`y( t + 12 h ) = `y + ka + 2 kb + 2 kc + kd
y( t + 12 h ) = y + 12 h ( `y + ka + kb + kc )
La traduzione in JavaScript (ovvero ECMAscript) è la seguente:
? Verifica elementare dell'algoritmo
Sperimento l'equazione differenziale
``y = y ; y( 0 ) = 0 ; `y( 0 ) = 1;
che deve dare come risultato y = sin( t ) ed avere come derivata `y = cos( t ).Verifica con le funzioni ellittiche jacobiane
La soluzione dell'equazione differenziale:``y = a y ( b2 − y2 ) / b3
y(0) = z ; `y(0) = 0
è una funzione ellittica jacobiana ( la funzione dn( u , m ) ) dipendente sia dai valori delle costanti della molla elastica non lineare sia dal valore della posizione dell'estremo dell'asta nel momento in cui la sua velocità è nulla. Il parametro z indica qui tale elongazione iniziale. In altre parole la soluzione analitica è:
y = z dn( t z ( a : 2 b3 )1/2 , 2 − 2 ( b / z )2 )
( Nota: l'operatore ":" coincide con l'operatore "/" ma ha priorità più bassa dell'operatore "+" o "−" e molto più bassa di "*" e dell'elevamento a potenza mentre l'operatore "/" ha priorità uguale a quella di "*" e dunque maggiore di "+" e di "−". Qui si adotta la convenzione di non usare il simbolo di moltiplicazione esplicitamente ma tutto ciò che costituisce continuazione di un nome o parentesi aperta di una funzione viene sottolineato con la marca <u> e scritto con caratteri piccoli. Il simbolo di derivata rispetto a t è il controapice, ossia il carattere ASCII 096, anteposto alla funzione da derivare)
La precisione risulta tanto maggiore quanto è più grande b rispetto ad a. La posizione di equilibrio si ottiene, ovviamente, quando z == b. Qui i calcoli sono fatti usando ?.