ATTENZIONE ! QUI CARICO UNA PAGINA DI SORGENTE Javascript MASCHERATA da documento HTML per cui la posso commentare bene, a colori, e posso leggerla con un qualsiasi browser, tipo Internet Explorer o Firefox .. Questa pagina, se il PC non è connesso in rete, avrebbe bisogno che nella sua cartella ci fosse il file: libreria-dormand-prince.js.html

Dormand Prince in Javascript

Trascrivo qui i risultati di quattro casi prova di crescente difficoltà

Il risultato finale del primo caso prova dovrebbe essere [1,7.38905609893065]

...100 ... !Javascript ORA NON ATTIVO!

...101 ... !Javascript ORA NON ATTIVO!

...102 ... !Javascript ORA NON ATTIVO!

Se ha girato come previsto viene : { [true] [1] [ { [1.0000000000000002] [7.389056098930669] } ] [888] }

Il problema da risolvere, scritto per esseri umani, è questo: per chi sa programmare ossia per un quasi marziano:

Ossia per esseri matematici di dubbia umanità :

x' = 1
y' = e2·x + 3·y2

con la condizione iniziale  [0,1]  e il valore finale della variabile indipendente  x = 1  e dunque la funzione  y , dipendente da  x , che vale  e2·x  deve valere  e2  ossia, come dice wxMaxima a cui faccio pubblicità in http://www.elegio.it/max/ : 7.389056098930650227230427460575007813180315570551847324087127823b0 .

Il risultato finale del secondo caso prova dovrebbe essere [-1,0]

...200 ... !Javascript ORA NON ATTIVO!

...201 ... !Javascript ORA NON ATTIVO!

...202 ... !Javascript ORA NON ATTIVO!

Se ha girato come previsto viene : { [true] [4.71238898038469] [ { [-1.0000000000000002] [5.435642170745216e-14] } ] [801] }

Per i marziani informatici questo sistemino di tre variabili, la prima,   x   l'indipendente ossia il primo argomento obbligatorio ma non usato nel sistema, e le due dipedenti indicate da   y   e da   z   per maggiore leggibilità :

ossia, per chi digerisce la matematica, ho il seguente sistema:

y' = z2
z' = −2·y

con le condizioni iniziali per la variabile indipendente ( ossia quella che ha implicitamente derivata prima 1 ) ossia: x = 0   e per le due variabili dipendenti   y = 0 e z = 2   e con il valore imposto finale della variabile indipendente dato da   4.71238898038469   ossia  3·π2 .

Dato che la soluzione è   y = sin( x )   e   z = 2·cos( x )  , quando   x   vale  3·π2  allora   y   vale   −1   e   z   vale   0   ma, per piccoli errori numerici, vale un numero piccolissino ossia un qualcosa nettamente minore di   1e-12.

Il risultato finale del terzo caso prova dovrebbe essere [146079760576,0,0,30500]

... ossia la distanza dal Sole dovrebbe essere circa 146 miliardi di metri e la velocità tangenziale, lungo y dovrebbe essere quella data all'inizio ossia 30500. L'orbita dovrebbe essere completata e la durata dell'orbita ossia il tempo finale è stato determinato usando la soluzione esatta del moto ellittico. Usando le formule si è fissata la lunghezza del transitorio ossia un anno ovvero un intero periodo di rivoluzione che è risultato essere di 31556606 secondi... 31 milioni circa...

Notare che questo è un test più impegnativo dei precedenti dato che richiede la soluzione di 4 equazioni differenziali.

...300 ... !Javascript ORA NON ATTIVO!

...301 ... !Javascript ORA NON ATTIVO!

...302 ... !Javascript ORA NON ATTIVO!

Se ha girato come previsto viene : { [true] [31556606.083602715] [ { [146079760591.51837] [-49.05691683292389] [0.000009993655297080295] [30499.999998315412] } ] [242] }

Per gli informatici il sistema da risolvere è questo:

E per determinare il periodo ossia la durata dell'anno e le quattro condizioni iniziali del moto bidimensionale , ho usato questa funzione che fa partire l'asteroide simil Terra dal perielio e dunque con velocità radiale nulla :

dove   1.327581e20 rappresenta il prodotto della costante di gravitazione universale G per la massa del Sole M,   149.61e9 rappresenta l'asse maggiore   a   e   30500.0   rappresenta la velocità in metri al secondo, da capire se al perielio o all'afelio a seconda che risulti maggiore o minore della velocità media circolare che non dipende dalla eccentricità ma da   G·M   e dal semiasse maggiore solitamente indicato da   a. La formuletta è semplicissima ossia   vm benscritta al quadrato vm2   vale   G·Ma dove la velocità media non è la media sulla vera traiettoria ma sul cerchio di raggio uguale al semiasse maggiore.

Dalla velocitità media e da quella minima all'afelio e massima al perielio ( legate dalla relazione   vp·va = vm·vm ) si ricava l'eccentricità dell'ellisse data da   vp − vavp + va   e dalla eccentricità ricavo la distanza dell'asteroide dal Sole al perielio ovvero a·(1 − ec ).

Quarto caso prova. Un sistema a tre corpi non ammette una soluzione analitica ma non è molto difficile da calcolare con questo algoritmo...

Purtroppo il solo test possibile è il confronto tra calcoli fatti a diversa tolleranza degli errori...

...400 ... !Javascript ORA NON ATTIVO!

...401 ... !Javascript ORA NON ATTIVO!

...402 ... !Javascript ORA NON ATTIVO!

...403 ... !Javascript ORA NON ATTIVO!

Se ha girato come previsto viene : { [true] [31556606.083602715] [ { [146084915168.18964] [-3674611.3418972045] [145667408476.40298] [-80256778.7892583] [10.51550926858721] [30511.232113772032] [261.62261764701765] [29601.405412900927] } ] [276] }

La libreria

ORA sto usando UN FILE ESTERNO col sorgente Javascript CAMUFFATO DA NORMALE DOCUMENTO HTML5 ma qui TRASCRIVO LA LIBRERIA in modo da garantire la conservazione della presente libreria. 
//
// Questa e' la funzione fondamentale.
// Richiede in argomento fder(x,y,costa) ossia la funzione da usare
// nel Runge Kutta e che deve produrre come risultato
// un vettore.
// Il primo argomento di fder e' il vettore delle incognite
// ossia le posizioni e le velocita'
// Il secondo argomento di fder e' il vettore delle variazioni
// Il terzo argomento di fder e' il vettore delle costanti usate
// che non vanno cambiate ma solo usate.
//
// Il significato degli altri argomenti e' questo:
//
// x0 == Valore iniziale della variabile indipendente.
// y0 == Valore iniziale del vettore delle variabili dipendenti.
// x1 == Valore finale desiderato della variabile indipendente.
// costa == Costanti usate per calcolare le derivate prime delle
//       variabili dipendenti che all'inizio valgono y1.
// peso == Vettori usato per equilibrare l'importanza degli errori
//         del calcolo. Spesso l'errore di una qualche variabile
//         ha importanza molto diversa dall'errore di una qualche
//         altra variabile.
// tol == parametro usato per limitare l'errore. Piu' e' piccolo
//        e maggiore e' la precisione ossia minore la tolleranza 
//        dell'errore totale ottenuto come somma degli errori 
//        delle varie variabili pesati sul peso della loro importanza.
//
function odedopri(fder,x0,y0,x1,costa,peso,tol,maxiter){
   var ha21,ha31,ha32,ha41,ha42,ha43,ha51,ha52,ha53,ha54,
       ha61,ha62,ha63,ha64,ha65,ha71,ha73,ha74,ha75,ha76;
   var e1,e3,e4,e5,e6,e7;
   var i,j,ifa,x,ny,h,hmax,hmin,verrore,erromax;
   var y=[0],z=[0],k1=[0],k2=[0],k3=[0],k4=[0],k5=[0],k6=[0],k7=[0];
   x=x0;
   ny=y0.length;
   // alert("ny = "+ny);
   for(j=1;ny>j;j++)y[j]=y0[j]; 
   h=10*(x1-x0)/(maxiter*171362822400);
   hmax=h;
   hmin=h/100;
   erromax=0;
   amen=false;
   e1=211228479;
   e3=-728766720;
   e4=6336854370;
   e5=-8716163841;
   e6=7180918272;
   e7=-4284070560;
   k1=fder(x,y,costa);
   //
   // alert("Stampo k1 \n"+k1);
   //
   for(i=1;maxiter>=i;i++){
       ifa=i;
       ha21=34272564480*h;
       for(j=1;ny>j;j++)z[j]=y[j]+ha21*k1[j]; 
       k2=fder(x+34272564480*h,z,costa);
       ha31=12852211680*h;
       ha32=38556635040*h;
       for(j=1;ny>=j;j++)z[j]=y[j]+ha31*k1[j]+ha32*k2[j]; 
       k3=fder(x+51408846720*h,z,costa);
       ha41=167554759680*h;
       ha42=-639754536960*h;
       ha43=609290035200*h;
       for(j=1;ny>j;j++)z[j]=y[j]+ha41*k1[j]+ha42*k2[j]+
          ha43*k3[j]; 
       k4=fder(x+137090257920*h,z,costa);
       ha51=505965644800*h;
       ha52=-1987087872000*h;
       ha53=1683278643200*h;
       ha54=-49833907200*h;
       for(j=1;ny>j;j++)z[j]=y[j]+ha51*k1[j]+ha52*k2[j]+
          ha53*k3[j]+ha54*k4[j]; 
       k5=fder(x+152322508800*h,z,costa);
       ha61=487745760600*h;
       ha62=-1843448544000*h;
       ha63=1526229734400*h;
       ha64=47708967600*h;
       ha65=-46873096200*h;
       hc6=171362822400*h;
       for(j=1;ny>j;j++)z[j]=y[j]+ha61*k1[j]+ha62*k2[j]+
          ha63*k3[j]+ha64*k4[j]+ha65*k5[j]; 
       k6=fder(x+hc6,z,costa);
       ha71=15619007250*h;
       ha73=76982400000*h;
       ha74=111564337500*h;
       ha75=-55243291950*h;
       ha76=22440369600*h;
       for(j=1;ny>j;j++)z[j]=y[j]+ha71*k1[j]+
          ha73*k3[j]+ha74*k4[j]+ha75*k5[j]+ha76*k6[j]; 
       k7=fder(x+171362822400*h,z,costa);
       //
       erromax=0;
       for(j=1; ny>j; j++){verrore=peso[j]*
              Math.abs(e1*k1[j]+e3*k3[j]+e4*k4[j]+e5*k5[j]+
              e6*k6[j]+e7*k7[j]);
          if(verrore>erromax)erromax=verrore;}
       //
       if(tol>erromax/171362822400){avanti=true;}
       else{ 
           if(hc6>hmin){h=h/2;avanti=false;} 
           else { avanti=true; }
           amen=false;}
       //
       if(avanti){
          x=x+hc6;
          for(j=1; ny>j; j++){
             y[j]=y[j]+ha71*k1[j]+ha73*k3[j]+ha74*k4[j]+
                  ha75*k5[j]+ha76*k6[j];
             k1[j]=k7[j]};
          if(amen){break};
          if(hmax>h)h=1.05*h;
          }
       if(x+h*171362822400>x1){ 
          h=(x1-x)/171362822400;amen=true};
      //
      };
   //
   //  amen vale true se ha fatto tutto, altrimenti false.
   //  ifa indica il numero di iterazioni fatte realmente.
   //
   return [0,amen,x,y,ifa];
   }
//
// __________________________________________________________
//
//  CASI PROVA LIBRERIA DORMAND PRINCE
//
//  Ogni caso e' basato su due funzioni. La funzione da integrare
//  e la funzione che fa l'inizializzazione dei dati del calcolo.
//
//  Le seguenti sono le variabili d'ambiente utilizzate:
//
//
var tol=1.0e-13,maxiter=250;
var x_0,x_1,y_0,costa,peso;
//
//
// Primo caso prova, elementare...
//
//
function primafder(x,y,c){
   var rr;
   rr=Math.exp(2*x)/2;
   return [0,1,rr+1.5*y[2]];
   }
//
function ini_fder1(){
   x_0=0;
   x_1=1;
   y_0=[0,0,1];
   costa=[0,1];
   peso=[0,1,1];
   tol=1.0e-13;    
   maxiter=250;
   return [0,x_0,y_0,x_1,costa,peso,tol,maxiter]
   }
//
//
// Secondo caso prova leggermente piu' complicato.
//
//
function secondafder(x,v,c){
   var y,z;
   y=v[1];
   z=v[2];
   return [0,z/2,-2*y];
   }
//
function ini_fder2(){
   x_0=0;
   x_1=1.5*Math.PI;
   y_0=[0,0,2];
   costa=[0,1];
   peso=[0,1,1];
   tol=1.0e-13;    
   maxiter=2000;
   return [0,x_0,y_0,x_1,costa,peso,tol,maxiter]
   }
//
//
// Terzo caso prova: moto piano ( bidimensionale ) di asteroide 
// attorno al Sole.
//
//
function terzafder(t,v,c){
   var x,y,vx,vy,mu,rq,rr,den,x_,y_,vx_,vy_;
   x=v[1];
   y=v[2];
   vx=v[3];
   vy=v[4];
   mu=c[1];
   x_=vx;
   y_=vy;
   rq=x*x+y*y;
   rr=Math.sqrt(rq);
   den=1/(rr*rq);
   vx_=-mu*x*den;
   vy_=-mu*y*den;
   // 
   // Il risultato e' il vettore delle quattro variabili ossia:
   // la derivata della posizione in ascissa, 
   // la derivata della posizione in ordinata,
   // la derivata della velocita' in ascissa,
   // la derivata della velocita' in ordinata.
   //
   return [0,x_,y_,vx_,vy_];
   }
//
function ini_fder3(){
   // velocita' media, al perielio, all'afelio, eccentricita'... 
   var vm,vp,va,ec;
   costa=[0,1.327581e20,149.61e9,30500.0];
   vm=Math.sqrt(costa[1]/costa[2]);
   if(costa[3]>vm){vp=costa[3];va=vm*vm/vp;}
   else {va=costa[3];vp=vm*vm/va};
   ec=(vp-va)/(vp+va);
   x_0=0;
   x_1=2*Math.PI*costa[2]/vm;
   y_0=[0,costa[2]*(1-ec),0,0,vp];
   peso=[0,1.0e-9,1.0e-9,1.0e-4,1.0e-4];
   tol=5.0e-14;    
   maxiter=500;
   return [0,x_0,y_0,x_1,costa,peso,tol,maxiter];
   }
//
//_________________________________________________________________
//
//  Problema di tre corpi nel piano ( Terra e Luna 
//     attorno al Sole immobile )
//  xt e yt sono l'ascissa e l'ordinata della Terra,
//  xl e yl sono l'ascissa e l'ordinata della Luna.
//  x sarebbe il tempo che pero' non interviene esplicitamente
//  nel calcolo delle derivate prime.
//  vxt e vyt sono la velocita' della Terra in ascissa e in ordinata.
//  vxl e vyl sono la velocita' della Luna in ascissa e in ordinata.
//  dunque le variabili sono in totale 8 e di queste la funzione
//  deve produrre le loro derivate prime.
//  Nel caso delle posizioni ossia le prime 4 variabili, le derivate 
//  prime sono ovviamente le velocita'.
//  Nel caso delle velocita' ossia le variabili dalla quinta alla 
//  ottava, le derivate prime sono ovviamente le forze. 
//
//
function quartafder(x,v,c){
   var mu,ts,xt,yt,vxt,vyt,rqt,rrt,dent,xt_,yt_,vxt_,vyt_;
   var dentl,denlt,ls,xl,yl,vxl,vyl,rql,rrl,denl,xl_yl_,vxl_vyl_;
   //
   // Ricopia il vettore delle coordinate solo per rendere 
   // piu' leggibile questa funzione.
   //
   xt=v[1];
   yt=v[2];
   xl=v[3];
   yl=v[4];
   vxt=v[5];
   vyt=v[6];
   vxl=v[7];
   vyl=v[8];
   //
   // mu e' la azione gravitazionale del Sole ossia la sua massa per G 
   //      ossia la costante di gravitazione universale.
   // ts e' il rapporto tra la massa della Terra e  quella del Sole
   // ls e' il rapporto tra la massa della Luna e quella del Sole...
   //
   mu=c[1];
   ts=c[2];
   ls=c[3];
   //
   // Ovviamente le derivate prime delle velocita' della Terra e della
   // Luna sono le velocita' ( bidimensionali, in ascissa e ordinata ) 
   //
   xt_=vxt;
   yt_=vyt;
   xl_=vxl;
   yl_=vyl;
   //
   // rqt e' il quadrato della distanza della Terra dal Sole.
   //
   rqt=xt*xt+yt*yt;
   //
   // rrt e' la distanza della Terra dal Sole.
   //
   rrt=Math.sqrt(rqt);
   //
   // l'inverso del cubo della distanza tra la Terra e in Sole.
   //
   dent=1/(rrt*rqt);
   //
   // rql e' il quadrato della distanza della Luna dal Sole.
   //
   rql=xl*xl+yl*yl;
   //
   // rrl e' la distanza della Luna dal Sole.
   //
   rrl=Math.sqrt(rql);
   //
   // l'inverso del cubo della distanza tra la Luna e il Sole.
   //
   denl=1/(rrl*rql);
   //
   // rqtl e' il quadrato della distanza della Luna dalla Terra.
   //
   rqtl=(xt-xl)*(xt-xl)+(yt-yl)*(yt-yl);
   //
   // rrtl e' la distanza della Luna dalla Terra.
   //
   rrtl=Math.sqrt(rqtl);
   //
   // per velocizzare i calcoli precalcolo queste due quantita'
   //
   dentl=ls/(rqtl*rrtl);
   denlt=ts/(rqtl*rrtl);
   //
   // le forze sulla Terra in ascissa ed in ordinata. Per
   // ottenere una forza attrattiva la forza deve avere il
   // segno meno. 
   //
   vxt_=-mu*(xt*dent+(xt-xl)*dentl);
   vyt_=-mu*(yt*dent+(yt-yl)*dentl);
   //
   // le forze sulla Luna in ascissa ed in ordinata
   //
   vxl_=-mu*(xl*denl+(xl-xt)*denlt);
   vyl_=-mu*(yl*denl+(yl-yt)*denlt);
   return [0,xt_,yt_,xl_,yl_,vxt_,vyt_,vxl_,vyl_];
   }
// _____________________________________________________________
//
// Produce un vettore di dati di questo quarto test ossia:
// 1== L'istante di tempo assunto come iniziale
// 2== Il vettore delle condizioni iniziali ossia le posizioni
//     delle particelle e le velocita' delle particelle.
// 3== L'istante di tempo che bisogna raggiungere.
// 4== Le costanti da usare per calcolare le forze ossia 
//     le accelerazioni.
// 5== Il vettore di equilibramento delle tolleranze degli errori.
// 6== Il parametro di controllo dell'errore, piccolissimo 
//     se molto preciso.
// 7== Il numero massimo delle iterazioni ammesse per fare 
//     evolvere il sistema fino alla fine.  
//
function ini_fder4(){
   var vm,vp,va,ec;
   //
   // Le costanti usate sono:
   // 1== la costante gravitazionale solare ossia
   //     la massa del sole per la G ossia costante di 
   //     gravitazione universale.
   // 2== Il rapporto tra la massa della Terra e la massa del Sole.
   // 3== Il rapporto tra la massa della Luna e la massa del Sole
   //     e la massa della Luna e' circa un ottantesimo della massa
   //     della terra per cui uso questo dato...
   // 4== La distanza iniziale della Terra dal Sole.
   // 5== La velocita' in metri al secondo della Terra.
   //
   costa=[0,1.327581e20, 1/333000.1, 1/(333000.1*80), 
            149.61e9, 30500.0];
   vm=Math.sqrt(costa[1]/costa[4]);
   if(costa[5]>vm){vp=costa[5];va=vm*vm/vp;}
   else {va=costa[5];vp=vm*vm/va; };
   ec=(vp-va)/(vp+va);
   //
   // x_0 e' l'istante iniziale del transitorio. 
   //
   x_0=0;
   //
   // x_1 e' l'istante finale a cui vorrei arrivare nella simulazione.
   //
   x_1=2*Math.PI*costa[4]/vm;
   //
   // y_0 e' il vettore delle incognite: 4 posizioni e 4 velocita'.
   //
   y_0=[0,costa[4]*(1-ec),0,costa[4]*(1-ec),-360e6,0,vp,1100.0,vp];
   //
   // il vettore peso controlla la tolleranza del calcolo
   // Il peso delle posizioni deve essere molto piccolo
   // rispetto ai pesi delle velocita' perche' le
   // posizioni si misurano in miliardi di metri e le
   // velocita' si misurano in decine di migliaia di metri al secondo.
   //
   peso=[0,1e-9,1e-9,1e-9,1e-9,1e-4,1e-4,1e-4,1e-4];
   //
   // Il criterio di tolleranza piuttosto severo...
   //
   tol=0.5e-14; 
   //
   // il massimo numero delle iterazioni consentite:
   //
   maxiter=2000;
   // 
   // Ecco il vettore risultato usabile per inizializzare i calcoli.
   // Notare che tutti questi dati sono variabili globali
   // e non locali mentre il vettore risultato e' solo una raccolta
   // di queste variabili globali inizializzate da questa funzione.
   // 
   return [0,x_0,y_0,x_1,costa,peso,tol,maxiter];
   }
//
// ____________________________________________________________________
//
//PER PILOTARE I CASI PROVA DELLA LIB. DORMAND PRINCE
//
// Devono esistere le seguenti MARCHE HTML:
//
// Per il primo caso   prova : "x100", "x101", "x102".
// Per il secondo caso prova : "x200", "x201", "x202".
// Per il terzo caso   prova : "x300", "x301", "x302".
// Per il quarto caso  prova : "x400", "x401", "x402", "x403".
//
function test_casi_prova_libreria(){
    var obs,ob,p1,ris3,ris4,ris5,inforis1,inforis2;
    //
    // ======================================================
    // Prima prova.
    // Il risultato finale dovrebbe essere [1,7.38905609893065]
    //
    p1=ini_fder1();
    inforis1=[0,"var.indip.x",
       "x considerato dipendente","y","x finale","costa",
       "xpeso","ypeso","tolleranza","iterazioni massime"];
    inforis2=[0,"Bene se true, male se false","x","xdip","y",
       "iterazioni fatte" ];
    ob=document.getElementById("x100");
    ob.innerHTML="Dati di inizializzazione del primo test "+
       vedoris(p1,inforis1);
    ris3=odedopri(primafder,x_0,y_0,x_1,costa,peso,tol,maxiter);
    ob=document.getElementById("x101");
    ob.innerHTML="Avendo fatto al massimo "+maxiter+" iterazioni "+
       " non dovrebbe riuscire a trovare la soluzione "+
       vedoris(ris3,inforis2);
    if(!ris3[1])ris4=odedopri(primafder,ris3[2],ris3[3],x_1,
       costa,peso,tol,4*maxiter);
    ob=document.getElementById("x102");
    ob.innerHTML="Risultato finale avendo "+
       " quadruplicato il massimo delle iterazioni ammesse "+
        vedoris(ris4,inforis2);
    //
    // ======================================================
    //
    // Seconda prova.
    // Il risultato finale dovrebbe essere [-1,0]
    //
    p1=ini_fder2();
    inforis1=[0,"varindip.iniz.x","y","z","varindip.finale",
       "cost","ypeso","zpeso","tolleranza","maxpassi"];
    inforis2=[0,"Bene se true, male se false","xfinale",
       "yfinale","zfinale","iterazionifatte"];
    ob=document.getElementById("x200");
    ob.innerHTML="Dati di inizializzazione del secondo test "+
       vedoris(p1,inforis1);
    ris3=odedopri(secondafder,x_0,y_0,x_1,costa,peso,tol,maxiter);
    ob=document.getElementById("x201");
    ob.innerHTML="Avendo fatto al massimo "+maxiter+" iterazioni "+
       "non trova la soluzione"+
       vedoris(ris3,inforis2);
    if(!ris3[1])ris4=odedopri(secondafder,ris3[2],ris3[3],x_1,
       costa,peso,tol,4*maxiter);
    ob=document.getElementById("x202");
    ob.innerHTML="Risultato finale con il quadruplo di iterazioni"+
      " che non usa se non quando servono al criterio "+
      " di cambio ampiezza passo "+vedoris(ris4,inforis2);
    //
    // ==============================================================
    //
    // Terza prova.
    // Il risultato finale dovrebbe essere :
    //           [146079760576.14456,0,0,30500]
    //
    p1=ini_fder3();
    inforis1=[0,"t-iniziale","xaster","yaster","vxaster","vyaster",
       "t_finale","costaG*M","costa_amaggiore","costavmedia",
       "xpeso","ypeso","vxpeso","vypeso","tolleranza","Maxiterazioni"];
    inforis2=[0,"Bene se true, male se false",
       "t_finale","xfinale","yfinale","vxfinale","vyfinale",
       "iterazionifatte"];
    ob=document.getElementById("x300");
    ob.innerHTML="Dati di inizializzazione del terzo test "+
       vedoris(p1,inforis1);
    ris3=odedopri(terzafder,x_0,y_0,x_1,costa,peso,tol,maxiter);
    ob=document.getElementById("x301");
    ob.innerHTML="Avendo fatto al massimo "+maxiter+" iterazioni "+
       " in questo caso riesce al primo colpo "+
       vedoris(ris3,inforis2);
    if(!ris3[1])ris4=odedopri(terzafder,ris3[2],ris3[3],x_1,
       costa,peso,tol,maxiter)
    else ris4=ris3;
    ob=document.getElementById("x302");
    ob.innerHTML="Risultato finale "+vedoris(ris4,inforis2);
    //
    // ============================================================
    //
    // Quarta prova che commento in dettaglio...
    // Il risultato finale dovrebbe essere ???
    //
    // Innanzi tutto calcola i dati che servono per impostare i calcoli
    //
    p1=ini_fder4();
    inforis1=[0,"t-iniziale","xterra","yterra","xluna","yluna","vxterra",
       "vyterra","vxluna","vyluna","t-finale","costa G*Msole",
       "Mterrasusole","Mlunasusole","Distadalsole","Vtangen",
       "pesoxt","pesoyt","pesoxl","pesoyl","pesovxt",
       "pesovyt","pesovxl","pesovyl","tollera","Maxpassi"];
    inforis2=[0,"Bene se true ma male se false","t-finale",
       "xtfin","ytfin","xlfin","ylfin","vxtfin","vytfin","vxlfin",
       "vylfin","Iterfatte"];
    ob=document.getElementById("x400");
    ob.innerHTML="Dati di inizializzazione del quarto test "+
       vedoris(p1,inforis1);
    //
    // Poi fa un numero massimo di 20 iterazioni.
    //
    ris3=odedopri(quartafder,x_0,y_0,x_1,costa,peso,tol,20);
    ob=document.getElementById("x401");
    ob.innerHTML="Avendo fatto preliminarmente "+20+" iterazioni "+
       vedoris(ris3,inforis2);
    // 
    // Se le sole 20 iterazioni non sono bastate ad ottenere true 
    // fa realmente un altro gruppo di maxiter iterazioni, ma se sono 
    // bastate le prime 20 ad arrivare al tempo x_1, non ripete il calcolo...
    // Notare che ris3[2] e' il tempo a cui e' riuscito ad arrivare con 20 
    // iterazioni e ris3[3] e' il vettore delle condizioni iniziali.
    //
    if(!ris3[1])ris4=odedopri(quartafder,ris3[2],ris3[3],x_1,
       costa,peso,tol,maxiter)
    else ris4=ris3;
    ob=document.getElementById("x402");
    ob.innerHTML="Avendo fatto preliminarmente "+maxiter+
       " iterazioni non ce la fa "+vedoris(ris4,inforis2);
    //
    // Se non e' riuscito a raggiungere il tempo x_1 rifa' un altro 
    // blocco di maxiter iterazioni sperando che almeno questo terzo tentativo 
    // gli consenta di arrivare alla conclusione ossia al tempo x_1 desiderato.
    //
    if(!ris4[1])ris5=odedopri(quartafder,ris4[2],ris4[3],x_1,
        costa,peso,tol,maxiter)
    else ris5=ris4;
    //
    ob=document.getElementById("x403");
    ob.innerHTML="Risultato finale ottenuto ridando "+
      "i due dati, tempo e valori ottenuti "+
      "e facendo un altro blocco di iterazioni con successo "+
      vedoris(ris5,inforis2);
    //
    }
//
//_________________________________________________________________________
//
//  Una funzione di servizio per trasformare vettori 
//  e vettori di vettori in stringhe.
//
function vedoris(rr,inforr){
    // Visualizza il vettore rr
    var j,pn=[],npn,paragrafo,zz,ta=typeof([0]),
        ss=" {\u2192 ",nn=rr.length;
    if(arguments.length>2)paragrafo=arguments[2];
    else {paragrafo=[],inforr[0]=0};
    npn=paragrafo.length;
    for(j=0;npn>j;j++)pn[j]=paragrafo[j];
    pn[npn]=0;
    for(j=1;nn>j;j++){
       pn[npn]=j;
       if(typeof(rr[j])==ta) { zz=vedoris(rr[j],inforr,pn);
            ss+="[\u00ab"+pn+"\u00bb: "+zz+"] "; }
       else { zz=rr[j];
            inforr[0]=Math.min(inforr[0]+1,inforr.length-1);
            ss+="[\u00ab"+pn+":"+inforr[inforr[0]]+
                 "\u00bb: "+zz+"] ";} };
    return ss+" \u2190} ";
    }
//
//
// Funzione generica per salvare nella libreria voluta ossia lib,
// nell'elemento qui ( se e' stato creato oppure in coda )
// i dati contenuti in rr che puo' essere un dato o una lista o
// una lista contenente altre liste...
//
function salvinlib(lib,qui,rr){
   var nn=Math.min(qui,lib.length);
   lib[nn]=[];
   salvorr(lib[nn],rr);
   }
//
function salvorr(lista,rr){
    var j,mm,nn=rr.length,ta=typeof([0]);
    // Salva tutto quello che sta in rr
    for(j=0;nn>j;j++){
       mm=lista.length;
       if(typeof(rr[j])==ta){           
            lista[mm]=[];
            salvorr(lista[mm],rr[j]);}
       else lista[mm]=rr[j];
       }
    }
//

Se tutto è andato ... non resta che affrontare problemi di effettivo interesse...

Buon divertimento !

Giampaolo Bottoni
gpbottoni@gmail.com