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 costruito ha estremi
    // elementi consecutivi di F(2,53) oppure f(z) = 0.
    //
    // z: approssimazione suggerita.
    // s,d: estremo sinistro e destro, rispettivamente,
    //      dell'ultimo intervallo calcolato.
    // iter: numero di iterazioni effettuate.
    //
    iter = 0;
    if ( f(a) == 0 | f(b) == 0 | sign(f(a)) == sign(f(b)) ) then
       error('la funzione non assume valori di segno opposto agli estremi');
    else
        s = a; d = b;
        z = (s+d)/2;
        while ( f(z) <> 0 & d > nearfloat('succ',s)),
              if sign(f(s)) <> sign(f(z)) then d = z;
              else s = z; end;
              z = (s+d)/2;
              iter = iter + 1;
        end;
    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
//
clf(); 
a = 0;
b = 3;
[zF,s,d,iter] = Bisezione(F,a,b)
printf('\n[zF,s,d,iter] = Bisezione(F,a,b)')
//
printf('\n\n   zF = %10.9e , F(zF) = %3.2e , iter = %d , ampiezza = %3.2e\n\n', ...
       zF,F(zF),iter,d-s);
printf('\n...immettere <resume> / <abort> per riprendere / interrompere l''esecuzione\n');
pause;
//a = 0;
//b = 3;
[zG,s,d,iter] = Bisezione(G,a,b)
printf('\n[zG,s,d,iter] = Bisezione(G,a,b)')
//
printf('\n\n   zG = %10.9e , G(zG) = %3.2e , iter = %d , ampiezza = %3.2e\n\n', ...
       zG,G(zG),iter,d-s);
printf('\n...immettere <resume> / <abort> per riprendere / interrompere l''esecuzione\n');
pause;
//
x = linspace(-1,1,901)'; x = x + 2;
plot(x,G(x),'-r',x,F(x),'-b'); xgrid();
xlabel('x'); legend('G(x)','F(x)','in_upper_left');
halt('...premere INVIO per continuare ...');
//
clf();
x = linspace(-0.23,0.23,700)'; x = x + 2;
plot(x,G(x),'-r',x,F(x),'-b'); xgrid();
xlabel('x'); legend('G(x)','F(x)','in_upper_left');
assi = gca();
assi.children(2).children(1).thickness=3;
plot(zG,G(zG),'ko');
halt('...premere INVIO per continuare ...');
clf();
x = linspace(-0.23,0.23,700)'; x = x + 2;
plot(x,abs(G(x)-F(x)),'-r'); xgrid();
xlabel('x','fontsize', 2); ylabel('|G(x) - F(x)|','fontsize', 2);
xtitle('Grafico per stimare max|G(x) - F(x)|')
