function [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,LAMBDA,HMIN)
//
// Integra numericamente, sull'intervallo [t0,tf], il Problema
// di Cauchy in R(n):
//
// x' = F(t,x)
// x(t0) = x0
//
// con il metodo TS(1) - Eulero esplicito - a passo variabile.
//
// x0: condizione iniziale (colonna di n elementi)
// t0: istante iniziale (numero reale)
// tf: istante finale (numero reale)
// F: function che definisce l'equazione differenziale - F(t,x) deve
//    essere una colonna di n numeri reali
// G: function che restituisce la derivata seconda in t della soluzione
//    dell'equazione differenziale che all'istante t assume valore x.
// E: valore massimo della stima dell'errore locale (numero reale)
// LAMBDA: numero reale che stabilisce il valore massimo del passo
//         (OPZIONALE - valore predefinito: 1d-5)
// HMIN: valore minimo consentito del passo
//       (OPZIONALE - valore predefinito: (tf - t0) / 1d6)
//
// T = [t(0),...,t(N)]: riga contenente gli istanti di integrazione
// X = [x(0),...,x(N)]: matrice n x (N+1) contenente le approssimazioni
// PASSO = [h(0),...,h(N-1)]: riga contenente i passi di integrazione
//
// Valore degli argomenti opzionali
//
if ~exists('LAMBDA','l') then LAMBDA = 1d-5; end;
if ~exists('HMIN','l') then HMIN = (tf - t0) / 1d6; end;
//
// Inizializzazione delle variabili di uscita
//
T(1,1) = t0;
X(:,1) = x0;
PASSO = [];
//
// ciclo principale
//
while (T(1,$) < tf), // arresta la costruzione se ha raggiunto tf
  //
  // scelta del passo
  //
  Nd2x = norm(G(T(1,$),X(:,$)));
  d = max(LAMBDA, Nd2x);
  PASSO(1,$+1) = min(sqrt(2*E/d), tf - T(1,$));
  //
  // calcolo approssimazione e nuovo istante di integrazione
  //
  X(:,$+1) = X(:,$) + F(T(1,$),X(:,$)) * PASSO(1,$);
  T(1,$+1) = T(1,$) + PASSO(1,$);
  //
  // arresta la costruzione se il passo calcolato risulta troppo
  // piccolo e non ha raggiunto tf
  //
  if (PASSO(1,$) < HMIN) & (T(1,$) < tf) then break; end;
  //
end;
//
// Verifica se l'integrazione ha raggiunto tf
//
if T(1,$) < tf then
   printf("\n\nIntegrazione interrotta a T = %3.2e", T(1,$));
end;
//
endfunction
//
// Nota sui valori di LAMBDA e HMIN
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E);
//
//         => LAMBDA = valore predefinito, HMIN = valore predefinito
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,y)
//
//    equivalente a
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,LAMBDA = y)
//
//         => LAMBDA = y, HMIN = valore predefinito
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,y,z)
//
//    equivalente a
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,LAMBDA = y,HMIN = z)
//
//         => LAMBDA = y, HMIN = z
//
//    [T,X,PASSO] = TS_1_pv(x0,t0,tf,F,G,E,HMIN = z)
//
//         => LAMBDA = valore predefinito, HMIN = z
//
