Newton's Method: Root Hunter

// Newton's method root hunter.
// f(x) = x^3 - 2x - 5  (classic Wallis example; real root ≈ 2.09455)
// User clicks anywhere on the plot to set x_0. We then animate the
// iteration x_{n+1} = x_n - f(x_n) / f'(x_n), one step at a configurable
// interval. The current tangent line is drawn from (x_n, f(x_n)) down
// to the x-axis where it intersects at x_{n+1}.

const X_MIN = -3, X_MAX = 4;
const Y_MIN = -8, Y_MAX = 8;
const PAD_L = 46, PAD_R = 14, PAD_T = 38, PAD_B = 26;
const STEP_PERIOD = 1.0;     // seconds per Newton step
const MAX_STEPS = 50;
const CONV_TOL = 1e-12;

let W, H;
let x0 = 1.0;
let xn;                      // current iterate
let steps;                   // number of completed iterations
let timer;                   // seconds since last step
let history;                 // array of {x, fx} for trail
let converged;
let blewUp;
let curvePts;

function f(x) { return x * x * x - 2 * x - 5; }
function fprime(x) { return 3 * x * x - 2; }

function precompute() {
  const M = 320;
  curvePts = new Float32Array(M * 2);
  for (let i = 0; i < M; i++) {
    const xx = X_MIN + (X_MAX - X_MIN) * i / (M - 1);
    curvePts[i * 2] = xx;
    curvePts[i * 2 + 1] = f(xx);
  }
}

function reset(newX0) {
  x0 = newX0;
  xn = newX0;
  steps = 0;
  timer = 0;
  history = [{ x: xn, fx: f(xn) }];
  converged = false;
  blewUp = false;
}

function init({ width, height, ctx }) {
  W = width; H = height;
  precompute();
  reset(1.0);
  ctx.fillStyle = '#06080d';
  ctx.fillRect(0, 0, W, H);
}

function toPx(x, y, plot) {
  return {
    x: plot.px + plot.pw * (x - X_MIN) / (X_MAX - X_MIN),
    y: plot.py + plot.ph - plot.ph * (y - Y_MIN) / (Y_MAX - Y_MIN),
  };
}
function fromPxX(xPx, plot) {
  return X_MIN + (X_MAX - X_MIN) * (xPx - plot.px) / plot.pw;
}

function tick({ ctx, dt, width, height, input }) {
  if (width !== W || height !== H) { W = width; H = height; }
  const plot = {
    px: PAD_L, py: PAD_T,
    pw: W - PAD_L - PAD_R, ph: H - PAD_T - PAD_B,
  };

  // clicks: reset to new x0
  const clicks = input.consumeClicks();
  for (const c of clicks) {
    if (c.x >= plot.px && c.x <= plot.px + plot.pw && c.y >= plot.py && c.y <= plot.py + plot.ph) {
      const nx = fromPxX(c.x, plot);
      reset(nx);
    }
  }

  // advance iteration
  if (!converged && !blewUp && steps < MAX_STEPS) {
    timer += dt;
    while (timer >= STEP_PERIOD) {
      timer -= STEP_PERIOD;
      const fv = f(xn);
      const fp = fprime(xn);
      if (Math.abs(fp) < 1e-10) { blewUp = true; break; }
      const xNext = xn - fv / fp;
      if (!isFinite(xNext) || Math.abs(xNext) > 1e8) { blewUp = true; break; }
      xn = xNext;
      steps += 1;
      history.push({ x: xn, fx: f(xn) });
      if (Math.abs(f(xn)) < CONV_TOL) { converged = true; break; }
    }
  }

  // ---- draw ----
  ctx.fillStyle = '#06080d';
  ctx.fillRect(0, 0, W, H);

  ctx.fillStyle = '#e8ecf4';
  ctx.font = 'bold 13px monospace';
  ctx.textAlign = 'left';
  ctx.fillText("Newton's Method:  f(x) = x³ − 2x − 5", 10, 18);
  ctx.fillStyle = '#7a8398';
  ctx.font = '11px monospace';
  ctx.fillText('click to set initial guess x₀ — real root ≈ 2.09455', 10, 32);

  // panel
  ctx.fillStyle = '#0a0d14';
  ctx.fillRect(plot.px, plot.py, plot.pw, plot.ph);
  ctx.strokeStyle = '#1a2030';
  ctx.lineWidth = 1;
  ctx.strokeRect(plot.px + 0.5, plot.py + 0.5, plot.pw - 1, plot.ph - 1);

  // grid
  ctx.strokeStyle = '#13192a';
  ctx.fillStyle = '#5a6478';
  ctx.font = '10px monospace';
  ctx.textAlign = 'right';
  for (let i = -8; i <= 8; i += 2) {
    const p = toPx(0, i, plot);
    ctx.beginPath(); ctx.moveTo(plot.px, p.y); ctx.lineTo(plot.px + plot.pw, p.y); ctx.stroke();
    ctx.fillText(i.toFixed(0), plot.px - 4, p.y + 3);
  }
  ctx.textAlign = 'center';
  for (let i = -3; i <= 4; i++) {
    const p = toPx(i, 0, plot);
    ctx.beginPath(); ctx.moveTo(p.x, plot.py); ctx.lineTo(p.x, plot.py + plot.ph); ctx.stroke();
    ctx.fillText(i.toFixed(0), p.x, plot.py + plot.ph + 14);
  }
  // axes emphasized
  ctx.strokeStyle = '#2a3450';
  const yZero = toPx(0, 0, plot);
  ctx.beginPath(); ctx.moveTo(plot.px, yZero.y); ctx.lineTo(plot.px + plot.pw, yZero.y); ctx.stroke();
  const xZero = toPx(0, 0, plot);
  ctx.beginPath(); ctx.moveTo(xZero.x, plot.py); ctx.lineTo(xZero.x, plot.py + plot.ph); ctx.stroke();

  // curve
  ctx.strokeStyle = '#5fd3ff';
  ctx.lineWidth = 1.8;
  ctx.beginPath();
  let prev = null;
  for (let i = 0; i < curvePts.length / 2; i++) {
    const xx = curvePts[i * 2];
    const yy = curvePts[i * 2 + 1];
    const yClamped = Math.max(Y_MIN - 1, Math.min(Y_MAX + 1, yy));
    const p = toPx(xx, yClamped, plot);
    if (!prev) ctx.moveTo(p.x, p.y);
    else ctx.lineTo(p.x, p.y);
    prev = p;
  }
  ctx.stroke();

  // analytical root marker at ≈ 2.09455148...
  // (we don't bother re-deriving; numeric high-precision)
  const ROOT = 2.0945514815423265;
  const pr = toPx(ROOT, 0, plot);
  ctx.fillStyle = 'rgba(124,240,138,0.9)';
  ctx.beginPath(); ctx.arc(pr.x, pr.y, 5, 0, Math.PI * 2); ctx.fill();
  ctx.strokeStyle = 'rgba(124,240,138,0.3)';
  ctx.setLineDash([3, 3]);
  ctx.beginPath(); ctx.moveTo(pr.x, plot.py); ctx.lineTo(pr.x, plot.py + plot.ph); ctx.stroke();
  ctx.setLineDash([]);

  // history trail: vertical drops & tangent segments
  for (let i = 0; i < history.length; i++) {
    const h = history[i];
    const isLast = i === history.length - 1;
    const yClamp = Math.max(Y_MIN, Math.min(Y_MAX, h.fx));
    const pPt = toPx(h.x, yClamp, plot);
    const pBase = toPx(h.x, 0, plot);
    const fade = isLast ? 1 : 0.25 + 0.5 * (i / Math.max(1, history.length - 1));

    // vertical drop from (x_i, 0) up to (x_i, f(x_i))
    ctx.strokeStyle = `rgba(255,207,102,${fade.toFixed(3)})`;
    ctx.lineWidth = isLast ? 1.8 : 1;
    ctx.beginPath(); ctx.moveTo(pBase.x, pBase.y); ctx.lineTo(pPt.x, pPt.y); ctx.stroke();

    // tangent segment to next iterate
    if (i < history.length - 1) {
      const next = history[i + 1];
      const pNextBase = toPx(next.x, 0, plot);
      ctx.strokeStyle = `rgba(255,122,209,${(0.4 + 0.5 * (i / Math.max(1, history.length - 1))).toFixed(3)})`;
      ctx.lineWidth = 1.4;
      ctx.beginPath(); ctx.moveTo(pPt.x, pPt.y); ctx.lineTo(pNextBase.x, pNextBase.y); ctx.stroke();
    }

    // point at (x_i, f(x_i))
    ctx.fillStyle = isLast ? '#fff0a1' : `rgba(255,207,102,${fade.toFixed(3)})`;
    ctx.beginPath(); ctx.arc(pPt.x, pPt.y, isLast ? 4.2 : 3, 0, Math.PI * 2); ctx.fill();
  }

  // x-axis label of current iterate
  if (history.length) {
    const last = history[history.length - 1];
    const pBase = toPx(last.x, 0, plot);
    ctx.fillStyle = '#fff0a1';
    ctx.font = 'bold 11px monospace';
    ctx.textAlign = 'center';
    ctx.fillText('x' + steps, pBase.x, pBase.y + 14);
  }

  // HUD
  ctx.textAlign = 'right';
  ctx.font = '11px monospace';
  const hudX = W - 14;
  let hy = plot.py + 14;
  ctx.fillStyle = '#9aa3b8';
  ctx.fillText('x₀   = ' + x0.toFixed(4), hudX, hy); hy += 14;
  ctx.fillStyle = '#fff0a1';
  ctx.fillText('x_' + steps + ' = ' + xn.toFixed(8), hudX, hy); hy += 14;
  ctx.fillStyle = '#5fd3ff';
  ctx.fillText('|f(x_n)| = ' + Math.abs(f(xn)).toExponential(3), hudX, hy); hy += 14;
  ctx.fillStyle = '#7cf08a';
  ctx.fillText('root = ' + ROOT.toFixed(8), hudX, hy); hy += 14;
  ctx.fillStyle = '#9aa3b8';
  ctx.fillText('iter = ' + steps, hudX, hy); hy += 14;
  if (converged) {
    ctx.fillStyle = '#7cf08a';
    ctx.fillText('CONVERGED', hudX, hy);
  } else if (blewUp) {
    ctx.fillStyle = '#ff7a7a';
    ctx.fillText('DIVERGED', hudX, hy);
  } else if (steps >= MAX_STEPS) {
    ctx.fillStyle = '#ff7ad1';
    ctx.fillText('MAX ITER', hudX, hy);
  }

  // bottom formula
  ctx.fillStyle = '#7a8398';
  ctx.font = '10px monospace';
  ctx.textAlign = 'center';
  ctx.fillText('x_{n+1} = x_n − f(x_n) / f\'(x_n)     f\'(x) = 3x² − 2',
               plot.px + plot.pw / 2, H - 6);
}