function [z,s,d,iter] = Bisezione(f,a,b)
    //
    // Applica il metodo di bisezione alla funzione f a partire dall'intervallo
    // [a,b]. L'iterazione si arresta quando l'intervallo non varia oppure
    // f(z) = 0.
    //
    // z: vettore contenente la parte calcolata della successione numerica.
    // s,d: estremo sinistro e destro, rispettivamente, degll'ultimo
    //      intervallo calcolato.
    // iter: numero di iterazioni effettuate.
    //
    iter = 0;
    if f(a)*f(b) > 0 then
       error('la funzione non assume valori di segno opposto agli estremi');
    else
        z = (a+b)/2;
        while ( z($) <> a & z($) <> b & f(z($)) <> 0 ),
              if sign(f(a)) <> sign(f(z($))) then b = z($);
              else a = z($); end;
              z($+1) = (a+b)/2;
              iter = iter + 1;
        end;
        s = a; d = b;
    end;
endfunction
//
function y = F(x)
    y = (x - 2).^13;
endfunction
//
function y = G(x)
    p = poly(2*ones(1,13),'x');
    y = horner(p,x);
endfunction
//
a = 0;
b = 3; //1;
[z,s,d,iter] = Bisezione(F,a,b)
//
printf('\n   z = %10.9e , F(z) = %3.2e , iter = %d , ampiezza = %3.2e\n', ...
       z($),F(z($)),iter,d-s);
//
pause;
//
a = 0; //+ rand();
b = 3; //+ rand();
[z,s,d,iter] = Bisezione(G,a,b)
//
printf('\n   z = %10.9e , G(z) = %3.2e , iter = %d , ampiezza = %3.2e\n', ...
       z($),G(z($)),iter,d-s);
//
pause
//
x = linspace(-0.23,0.23,700)'; x = x + 2;
clf(); plot2d(x,[G(x),F(x)],style = [5,2]); xgrid();
xlabel('x'); legend('G(x)','F(x)');
plot2d(z($),G(z($)),style = -5);
assi = gca();
assi.children(3).children(1).thickness=2;
//
//  Esercizio:
//
//  Adattare quanto fatto al caso:
//
//   f(x) = (3*x - 2)^13
//
