// Static figures for Chapters 1 to 4.
import { Arrowheads, Figure, Label, Node, arrow } from './svg.tsx'
import { ink, method } from '../plot/palette.ts'

/** The mesh, the three-point stencil, and the band it puts in the matrix. */
export function ThreePointStencil() {
  const L = 70
  const R = 470
  const Y = 96
  const N = 8
  const xs = Array.from({ length: N + 1 }, (_, i) => L + (i * (R - L)) / N)
  const hit = [3, 4, 5]
  return (
    <Figure
      width={716}
      height={200}
      caption="The three-point scheme. Each interior node is coupled to its two neighbours and to nothing else, so the matrix has three non-zero diagonals and the solve is two sweeps."
    >
      <Arrowheads id="tps" colors={[ink.muted]} />

      {/* the mesh */}
      <line x1={L} y1={Y} x2={R} y2={Y} stroke={ink.axis} strokeWidth={1.4} />
      {xs.map((x, i) => (
        <g key={i}>
          <Node x={x} y={Y} r={i === 0 || i === N ? 5 : 3.6} fill={i === 0 || i === N ? method.fd2 : ink.node} />
          <Label x={x} y={Y + 20} size={11} fill={ink.muted} mono>
            {i === N ? (
              <>
                u<tspan fontSize={9} dy={3}>N</tspan>
              </>
            ) : (
              `u${'₀₁₂₃₄₅₆₇₈'[i]}`
            )}
          </Label>
        </g>
      ))}
      <Label x={L} y={Y - 46} size={11.5} fill={method.fd2}>
        u(0) = α
      </Label>
      <Label x={R} y={Y - 46} size={11.5} fill={method.fd2}>
        u(1) = β
      </Label>
      <line x1={L} y1={Y - 36} x2={L} y2={Y - 10} stroke={method.fd2} strokeWidth={1} />
      <line x1={R} y1={Y - 36} x2={R} y2={Y - 10} stroke={method.fd2} strokeWidth={1} />

      {/* h */}
      <line x1={xs[1]} y1={Y + 44} x2={xs[2]} y2={Y + 44} stroke={ink.muted} strokeWidth={1} markerStart={arrow('tps', 0)} markerEnd={arrow('tps', 0)} />
      <Label x={(xs[1] + xs[2]) / 2} y={Y + 60} size={11.5} fill={ink.muted} italic>
        h = 1/N
      </Label>

      {/* the stencil, sitting on one interior node */}
      <path
        d={`M ${xs[hit[0]] - 12} ${Y - 26} H ${xs[hit[2]] + 12} V ${Y + 8} H ${xs[hit[0]] - 12} Z`}
        fill="none"
        stroke={method.fd2}
        strokeWidth={1.4}
        strokeDasharray="4 3"
      />
      {hit.map((i, k) => (
        <Label key={i} x={xs[i]} y={Y - 38} size={12} fill={method.fd2} mono>
          {[' −1', '  2', ' −1'][k]}
        </Label>
      ))}

      {/* the matrix band */}
      <g transform="translate(524 22)">
        <Label x={72} y={-4} size={11.5} fill={ink.muted}>
          the matrix, to scale
        </Label>
        {Array.from({ length: 7 }, (_, r) =>
          Array.from({ length: 7 }, (_, c) => {
            const d = Math.abs(r - c)
            if (d > 1) return null
            return (
              <rect
                key={`${r}-${c}`}
                x={16 + c * 19}
                y={8 + r * 19}
                width={17}
                height={17}
                fill={d === 0 ? method.fd2 : '#1c5cab'}
                opacity={d === 0 ? 0.95 : 0.6}
                rx={2}
              />
            )
          }),
        )}
        <rect x={14} y={6} width={137} height={137} fill="none" stroke={ink.axis} strokeWidth={1} rx={3} />
        <Label x={82} y={160} size={11} fill={ink.muted}>
          three diagonals, O(N) to solve
        </Label>
      </g>
    </Figure>
  )
}

/** Where manufactured data comes from, and what is then compared with what. */
export function ManufacturedLoop() {
  const boxes = [
    { x: 6, w: 126, t: 'pick u(x)', s: 'any function you like' },
    { x: 148, w: 206, t: 'f := −u″,  α := u(0),  β := u(1)', s: 'differentiate by hand' },
    { x: 370, w: 146, t: 'solve for uₕ', s: 'the method under test' },
    { x: 532, w: 178, t: 'error = uₕ − u', s: 'known exactly, everywhere' },
  ]
  const Y = 42
  return (
    <Figure
      width={716}
      height={104}
      caption="The method of manufactured solutions. The data is derived from the answer rather than the other way round, so no reference computation stands between the method and its error."
    >
      <Arrowheads id="mfg" colors={[ink.muted]} />
      {boxes.map((b, i) => (
        <g key={i}>
          <rect x={b.x} y={Y} width={b.w} height={46} rx={6} fill="#201f1e" stroke={i === 3 ? method.fd2 : ink.axis} strokeWidth={1.4} />
          <Label x={b.x + b.w / 2} y={Y + 18} size={12.5} fill={ink.primary}>
            {b.t}
          </Label>
          <Label x={b.x + b.w / 2} y={Y + 35} size={10.5} fill={ink.muted}>
            {b.s}
          </Label>
          {i < boxes.length - 1 ? (
            <line
              x1={b.x + b.w + 4}
              y1={Y + 23}
              x2={boxes[i + 1].x - 4}
              y2={Y + 23}
              stroke={ink.muted}
              strokeWidth={1.3}
              markerEnd={arrow('mfg', 0)}
            />
          ) : null}
        </g>
      ))}
      <path
        d={`M 69 ${Y} C 69 14, 621 14, 621 ${Y}`}
        fill="none"
        stroke={ink.muted}
        strokeWidth={1.1}
        strokeDasharray="4 3"
        markerEnd={arrow('mfg', 0)}
      />
      <Label x={345} y={14} size={11} fill={ink.muted}>
        the same u, kept for comparison
      </Label>
    </Figure>
  )
}

/** The V: discretisation down, roundoff up, and the best attainable error. */
export function RoundoffFloor() {
  const L = 66
  const R = 546
  const T = 24
  const B = 168
  // x is log10 N from 1 to 7; y is log10 error, 1 at the top and 1e-10 at the foot
  const sx = (v: number) => L + ((v - 1) / 6) * (R - L)
  const sy = (v: number) => T + (-v / 10) * (B - T)
  const disc = (n: number) => -2.3 - 2 * (n - 1)
  const round = (n: number) => -16 + 1.925 * n
  const total = (n: number) => Math.log10(10 ** disc(n) + 10 ** round(n))
  const path = (f: (n: number) => number) => {
    const pts: string[] = []
    for (let i = 0; i <= 240; i++) {
      const n = 1 + (6 * i) / 240
      const y = f(n)
      if (y > 0 || y < -10) continue
      pts.push(`${sx(n).toFixed(1)} ${sy(y).toFixed(1)}`)
    }
    return 'M ' + pts.join(' L ')
  }
  const nStar = 4 // where the two contributions are equal
  return (
    <Figure
      width={716}
      height={216}
      caption="Two error sources with opposite trends. Refining reduces the discretisation error and amplifies roundoff, so the total has a minimum: for a second-order scheme in double precision, near N = 10⁴ and an error of 10⁻⁸."
    >
      <line x1={L} y1={T} x2={L} y2={B} stroke={ink.axis} strokeWidth={1.2} />
      <line x1={L} y1={B} x2={R} y2={B} stroke={ink.axis} strokeWidth={1.2} />
      {[0, -2, -4, -6, -8, -10].map((e) => (
        <g key={e}>
          <line x1={L} y1={sy(e)} x2={R} y2={sy(e)} stroke={ink.grid} strokeWidth={1} />
          <Label x={L - 8} y={sy(e)} size={10.5} anchor="end" fill={ink.muted} mono>
            {e === 0 ? '1' : `1e${e}`}
          </Label>
        </g>
      ))}
      {[1, 2, 3, 4, 5, 6, 7].map((n) => (
        <Label key={n} x={sx(n)} y={B + 14} size={10.5} fill={ink.muted} mono baseline="hanging">
          {`1e${n}`}
        </Label>
      ))}
      <Label x={(L + R) / 2} y={B + 34} size={11.5} fill={ink.muted} baseline="hanging">
        intervals N
      </Label>
      <g transform={`translate(${L - 48} ${(T + B) / 2}) rotate(-90)`}>
        <Label x={0} y={0} size={11.5} fill={ink.muted}>
          error
        </Label>
      </g>

      <path d={path(disc)} fill="none" stroke={ink.family} strokeWidth={1.8} strokeDasharray="5 4" />
      <path d={path(round)} fill="none" stroke={ink.family} strokeWidth={1.8} strokeDasharray="5 4" />
      <path d={path(total)} fill="none" stroke={method.fd2} strokeWidth={2.6} strokeLinecap="round" />

      <Label x={sx(2.1)} y={sy(disc(2.1)) - 14} size={11.5} anchor="start" fill={ink.secondary}>
        discretisation, C h²
      </Label>
      <Label x={sx(6.4)} y={sy(round(6.4)) - 14} size={11.5} anchor="end" fill={ink.secondary}>
        roundoff, ε κ(A)
      </Label>

      <line x1={sx(nStar)} y1={sy(total(nStar)) + 8} x2={sx(nStar)} y2={B} stroke={method.fd2} strokeWidth={1} strokeDasharray="3 3" />
      <Node x={sx(nStar)} y={sy(total(nStar))} r={4} fill={method.fd2} />
      <Label x={R + 12} y={sy(total(nStar)) - 8} size={12} anchor="start" fill={method.fd2}>
        best ≈ √ε
      </Label>
      <Label x={R + 12} y={sy(total(nStar)) + 9} size={11} anchor="start" fill={ink.muted}>
        about half the
      </Label>
      <Label x={R + 12} y={sy(total(nStar)) + 24} size={11} anchor="start" fill={ink.muted}>
        digits, and no more
      </Label>
    </Figure>
  )
}

/** Numerov changes the right-hand side and leaves the matrix alone. */
export function CompactScheme() {
  const CX = [268, 336, 404]
  const cell = (x: number, y: number, v: string | null, colour: string) =>
    v === null ? (
      <rect x={x - 30} y={y} width={60} height={34} rx={5} fill="none" stroke={ink.grid} strokeWidth={1.2} strokeDasharray="3 3" />
    ) : (
      <>
        <rect x={x - 30} y={y} width={60} height={34} rx={5} fill="#201f1e" stroke={colour} strokeWidth={1.4} />
        <Label x={x} y={y + 17} size={12.5} fill={colour} mono>
          {v}
        </Label>
      </>
    )
  const row = (y: number, name: string, vals: (string | null)[], colour: string, note: string) => (
    <g>
      <Label x={12} y={y + 17} size={12.5} anchor="start" fill={ink.primary}>
        {name}
      </Label>
      {vals.map((v, i) => (
        <g key={i}>{cell(CX[i], y, v, colour)}</g>
      ))}
      <Label x={452} y={y + 17} size={11.5} anchor="start" fill={ink.muted}>
        {note}
      </Label>
    </g>
  )
  return (
    <Figure
      width={716}
      height={234}
      caption="Numerov's scheme against the three-point scheme. The left-hand side, and therefore the matrix and the cost of the solve, is identical; only the weighting of the data changes."
    >
      {['i−1', 'i', 'i+1'].map((t, i) => (
        <Label key={t} x={CX[i]} y={14} size={11.5} fill={ink.muted}>
          {t}
        </Label>
      ))}

      <Label x={12} y={40} size={11} anchor="start" fill={ink.muted}>
        left-hand side, on the unknowns u
      </Label>
      {row(52, 'both schemes', ['−1', '2', '−1'], method.fd2, '÷ h², tridiagonal')}

      <line x1={12} y1={110} x2={704} y2={110} stroke={ink.grid} strokeWidth={1} />

      <Label x={12} y={130} size={11} anchor="start" fill={ink.muted}>
        right-hand side, on the data f
      </Label>
      {row(142, 'three-point', [null, 'fᵢ', null], ink.secondary, 'order 2')}
      {row(184, 'Numerov', ['f/12', '10f/12', 'f/12'], method.fd4, 'order 4, for one pass over an array')}
    </Figure>
  )
}

/** Equispaced against Chebyshev points, and where they come from. */
export function ChebyshevPoints() {
  const N = 16
  const L = 172
  const R = 562
  const cx = (L + R) / 2
  const rad = (R - L) / 2
  const YARC = 216
  const cheb = Array.from({ length: N + 1 }, (_, j) => 0.5 * (1 - Math.cos((Math.PI * j) / N)))
  const equi = Array.from({ length: N + 1 }, (_, j) => j / N)
  const X = (t: number) => L + t * (R - L)
  return (
    <Figure
      width={716}
      height={300}
      caption="Chebyshev points are equally spaced angles on the semicircle, projected down. The clustering near the ends is what removes the Runge phenomenon, and it is also why the boundary-layer solution turns out to be the easy case for this method."
    >
      <path d={`M ${L} ${YARC} A ${rad} ${rad} 0 0 1 ${R} ${YARC}`} fill="none" stroke={ink.axis} strokeWidth={1.3} />
      {cheb.map((t, j) => {
        const th = Math.PI - (Math.PI * j) / N
        const ax = cx + rad * Math.cos(th)
        const ay = YARC - rad * Math.sin(th)
        return (
          <g key={j}>
            <line x1={ax} y1={ay} x2={X(t)} y2={YARC} stroke={ink.grid} strokeWidth={1} />
            <Node x={ax} y={ay} r={2.6} fill={method.spectral} />
          </g>
        )
      })}
      <Label x={20} y={64} size={11.5} anchor="start" fill={ink.muted}>
        equal angles
      </Label>
      <Label x={20} y={80} size={11.5} anchor="start" fill={ink.muted}>
        on the semicircle
      </Label>

      <line x1={L} y1={YARC} x2={R} y2={YARC} stroke={ink.axis} strokeWidth={1.4} />
      {cheb.map((t, j) => (
        <Node key={j} x={X(t)} y={YARC} r={3.4} fill={method.spectral} />
      ))}
      <Label x={L - 14} y={YARC} size={12} anchor="end" fill={method.spectral}>
        Chebyshev
      </Label>

      <line x1={L} y1={YARC + 44} x2={R} y2={YARC + 44} stroke={ink.axis} strokeWidth={1.4} />
      {equi.map((t, j) => (
        <Node key={j} x={X(t)} y={YARC + 44} r={3.4} fill={ink.family} />
      ))}
      <Label x={L - 14} y={YARC + 44} size={12} anchor="end" fill={ink.muted}>
        equispaced
      </Label>
      <Label x={R + 14} y={YARC + 44} size={11} anchor="start" fill={ink.muted}>
        unusable
      </Label>
      <Label x={R + 14} y={YARC + 59} size={11} anchor="start" fill={ink.muted}>
        past N ≈ 20
      </Label>

      <Label x={cx} y={YARC + 78} size={11.5} fill={ink.muted}>
        N = 16 on [0, 1]
      </Label>
    </Figure>
  )
}