# Nome del file: equazioni_non_lineari.mpl6
# Autore: Claudio Marsan
# Ultima revisione: 28 novembre 2001
# Realizzato per Maple V release 6.02
# Testato su MS-Windows 2000
# Contenuto: alcuni metodi per la risoluzione di equazioni non lineari
#  (1) Il metodo di bisezione
#  (2) Il metodo di falsa posizione (regula falsi)
#  (3) Il metodo delle secanti
#  (4) Il metodo delle tangenti (metodo di Newton)
#  (5) Il metodo di iterazione


##############################
# (1) Il metodo di bisezione #
##############################

metodo_di_bisezione := proc(f::algebraic,      # la funzione di cui si cerca uno zero
                            a,                 # estremo sinistro dell'intervallo
                            b,                 # estremo destro dell'intervallo
                            max_passi::posint, # numero massimo di passi da eseguire
                            epsilon)           # precisione richiesta
  local F, incognite, x0, x1, x2, n;
  incognite := op(select(type, indets(f), name));
  F := unapply(f, incognite);
  x0 := evalf(a); x1 := evalf(b); 
  if (x0 >= x1) then
    ERROR(`Intervallo inesistente!`);
  end if;
  if evalf(F(x0)*F(x1)) >= 0.0 then
    ERROR(`Intervallo non valido!`);
  end if;
  n := 0; 
  while (n < max_passi) do 
    n := n + 1;
    x2 := 0.5*(x0 + x1);
    if evalf(F(x0)*F(x2)) < 0.0 then 
      x1 := x2;
    else 
      x0 := x2;
    end if;
    lprint(n, x2, F(x2));
    if abs(F(x2)) < epsilon then
      break;
    end if;
  end do;
end;


###################################################
# (2) Il metodo di falsa posizione (regula falsi) #
###################################################

regula_falsi := proc(f::algebraic,      # la funzione di cui si cerca uno zero
                     a,                 # una ascissa 
                     b,                 # un'altra ascissa
                     max_passi::posint, # numero massimo di passi da eseguire
                     epsilon)           # precisione richiesta
  local F, incognite, x0, x1, x2, y0, y1, y2, n;
  incognite := op(select(type, indets(f), name));
  F := unapply(f, incognite);
  x0 := evalf(a); x1 := evalf(b); 
  if (x0 >= x1) then
    ERROR(`Intervallo inesistente!`);
  end if;
  y0 := evalf(F(x0)); y1 := evalf(F(x1));
  if evalf(F(x0)*F(x1) >= 0.0) then
    ERROR(`Intervallo non valido!`);
  end if;
  n := 0; 
  while (n < max_passi) do
    n := n + 1;
    x2 := x1 - y1*(x1 - x0)/(y1 - y0);
    y2 := evalf(F(x2));
    if y0*y2 < 0.0 then 
      x1 := x2;
      y1 := y2;
    else 
      x0 := x2;
      y0 := y2;
    end if;
    lprint(n, x2, y2);
    if abs(y2) < epsilon then
      break;
    end if;
  end do;
end;



###############################
# (3) Il metodo delle secanti #
###############################

metodo_delle_secanti := proc(f::algebraic,      # la funzione di cui si cerca uno zero
                             a,                 # estremo sinistro dell'intervallo
                             b,                 # estremo destro dell'intervallo
                             max_passi::posint, # numero massimo di passi da eseguire
                             epsilon)           # precisione richiesta
  local F, incognite, x0, x1, x2, y0, y1, y2, n;
  incognite := op(select(type, indets(f), name));
  F := unapply(f, incognite);
  x0 := evalf(a); x1 := evalf(b); 
  y0 := evalf(F(x0)); y1 := evalf(F(x1));
  n := 0; 
  while (n < max_passi) do
    n := n + 1;
    x2 := x1 - y1*(x1 - x0)/(y1 - y0);
    y2 := evalf(F(x2));
    x0 := x1;
    y0 := y1;
    x1 := x2;
    y1 := y2;
    lprint(n, x2, y2);
    if abs(y2) < epsilon then
      break;
    end if;
  end do;
end;



###################################################
# (4) Il metodo delle tangenti (metodo di Newton) #
###################################################

metodo_di_Newton := proc(f::algebraic,      # la funzione di cui si cerca uno zero
                         a,                 # valore iniziale
                         max_passi::posint, # numero massimo di passi da eseguire
                         epsilon)           # precisione richiesta
  local F, incognite, x0, xn, yn, n;
  incognite := op(select(type, indets(f), name));
  F := unapply(f, incognite);
  xn := evalf(a);
  n := 0; 
  while (n < max_passi) do
    n := n + 1;
    x0 := xn;
    xn := x0 - F(x0)/(D(F)(x0)); 
    yn := F(xn);
    lprint(n, xn, yn);
    if abs(yn) < epsilon then
      break;
    end if;
  end do;
end;


###############################
# (5) Il metodo di iterazione #
###############################

metodo_di_iterazione := proc(f::algebraic,      # la funzione di cui si cerca uno zero
                             a,                 # valore iniziale
                             max_passi::posint, # numero massimo di passi da eseguire
                             epsilon)           # precisione richiesta
  local F, incognite, x1, x2, n;
  incognite := op(select(type, indets(f), name));
  F := unapply(f, incognite);
  x1 := evalf(a); x2 := F(x1);
  n := 0; 
  while (n < max_passi) do 
    n := n + 1;
    x1 := x2;
    x2 := F(x1);
    lprint(n, x1, x2);
    if abs(x2 - x1) < epsilon then
      break;
    end if;
  end do;
end;