// ==========================================================================
//  Newton.sci: ニュートン法(テキストp.72の式(4.16))で非線形方程式を解く関数
// --------------------------------------------------------------------------
//  f: 関数f(x),  df: 導関数f'(x),  x0: 解の初期値
//  tol: 絶対誤差の許容範囲,  itr: 最大繰り返し回数
// ==========================================================================
function x = Newton(f, df, x0, tol, itr)
x(1)  = x0;                            // xの初期値x0をx(1)とする
y(1)  = f( x(1) );                     // x0に対するf(x0)の値をy(1)とする
dy(1) = df( x(1) );                    // x0に対するf'(x0)の値をdy(1)とする

for k = 1:itr                          // 最大でitr回繰り返す
   x(k+1)  = x(k) - y(k) / dy(k);      // 式(4.16)でxの値を更新
   y(k+1)  = f( x(k+1) );              // 更新したxに対するf(x)の値
   dy(k+1) = df( x(k+1) );             // 更新したxに対するf'(x)の値
   if ( abs( y(k+1) ) < tol ) then     // 式(4.17)で収束判定
      break;                           // 収束したらfor文から抜ける
   end
end
format(12)                             // 値を10桁まで表示する
disp([ [0:k]', x, y ])                 // ニュートン法の繰り返し回数とx, y の値を表示
endfunction