/**
 * GPU validation, run in a browser (scripts/test-gpu.mjs drives it headlessly):
 *
 *  1. near weights against the brute-force double precision reference, for the
 *     self pair and the near pairs of a few nodes;
 *  2. Green's identity S_k[du/dn] - D_k[u] = u/2 at all wall nodes through the
 *     full operator (far field plus near weights);
 *  3. the self-test problem: the field of a source outside the room, solved
 *     with the impedance equation and compared with the exact field inside.
 *
 * The surface is kept small so that the software WebGPU adapter can run it.
 */
import { roundedCube } from '../src/surface.ts';
import { green, greenDn, selftestRhs, selftestPoints, relError, ALPHA, type CVec } from '../src/problem.ts';
import { refSelf, refLeaf, momentsToWeights } from '../src/reference.ts';
import { requestSolverDevice } from '../src/gpu/device.ts';
import { Operator } from '../src/gpu/operator.ts';

declare global {
  interface Window {
    __RESULTS__?: { passed: number; failed: number; fatal?: string; lines: string[] };
  }
}

const logEl = document.getElementById('log')!;
const lines: string[] = [];
let passed = 0;
let failed = 0;
const log = (s: string, cls = '') => {
  lines.push(s);
  const span = document.createElement('span');
  span.textContent = s + '\n';
  if (cls) span.className = cls;
  logEl.append(span);
  console.log(s);
};
const check = (name: string, val: number, tol: number) => {
  const ok = val < tol;
  ok ? passed++ : failed++;
  log(`${ok ? 'ok  ' : 'FAIL'} ${name}: ${val.toExponential(2)} (tol ${tol.toExponential(0)})`, ok ? 'ok' : 'fail');
};

async function main() {
  const t0 = performance.now();
  if (logEl.textContent === 'starting…') logEl.textContent = '';
  const gpu = await requestSolverDevice();
  const device = gpu.device;
  const lim = device.limits;
  log(`device: maxStorageBufferBindingSize ${(lim.maxStorageBufferBindingSize / 2 ** 20).toFixed(0)} MB, workgroup storage ${lim.maxComputeWorkgroupStorageSize} B`);

  const params = new URLSearchParams(location.search);
  const m = Number(params.get('m') ?? 3);
  const norder = Number(params.get('norder') ?? 6);
  const k = Number(params.get('k') ?? 4);
  const a = 1;
  const p = Number(params.get('p') ?? 8);
  const surf = roundedCube(a, p, m, norder);
  const npts = surf.npts;
  const npols = surf.npols;
  log(`surface: p = ${p}, ${surf.npatches} patches of order ${norder}, ${npts} nodes; k = ${k}`);
  const xs = surf.srcvals.subarray(0, 3 * npts);
  const nrm = surf.srcvals.subarray(9 * npts, 12 * npts);

  const lam = -0.2;
  const op = new Operator(gpu, surf, k, ALPHA, lam, 0, { progress: () => {} });
  const info = await op.setup();
  log(`operator: ${info.npairs} near pairs, ${info.nleaves} leaves, ${info.noversampled} oversampled nodes, weights ${(info.weightsBytes / 2 ** 20).toFixed(1)} MB, setup ${info.setupSeconds.toFixed(2)} s`);

  // 1. near weights against the reference, for three nodes of one patch
  {
    const { pairs, W } = await op.readWallWeights();
    const tp = 4;
    const worst = new Float64Array(6);
    for (const loc of [0, Math.floor(npols / 2), npols - 2]) {
      const i = tp * npols + loc;
      const xt = [xs[i], xs[npts + i], xs[2 * npts + i]];
      const nt = [nrm[i], nrm[npts + i], nrm[2 * npts + i]];
      for (let pp = pairs.rowPtr[i]; pp < pairs.rowPtr[i + 1]; pp++) {
        const j = pairs.pairP[pp];
        const self = j === tp;
        const I = self ? refSelf(surf, j, surf.uvNodes[i], surf.uvNodes[npts + i], xt, nt, k) : refLeaf(surf, j, xt, nt, k);
        const Wref = momentsToWeights(surf, I);
        const errs: string[] = [];
        for (let mm = 0; mm < 6; mm++) {
          let num = 0;
          let den = 0;
          for (let l = 0; l < 2 * npols; l++) {
            const g = W[2 * ((pp * 6 + mm) * npols) + l];
            num += (g - Wref[mm][l]) ** 2;
            den += Wref[mm][l] ** 2;
          }
          const e = Math.sqrt(num / den);
          worst[mm] = Math.max(worst[mm], e);
          errs.push(e.toExponential(1));
        }
        const d = Math.hypot(xt[0] - surf.centroids[j], xt[1] - surf.centroids[surf.npatches + j], xt[2] - surf.centroids[2 * surf.npatches + j]) / surf.rads[j];
        log(`     node ${i} patch ${j} (${self ? 'self' : 'near'}, dist/rad ${d.toFixed(2)}): rel err per kernel ${errs.join(' ')}`);
      }
    }
    for (let mm = 0; mm < 6; mm++) check(`near weights, kernel ${mm}, worst relative error`, worst[mm], 1e-4);
  }

  // 2. Green's identity through the operator
  {
    const xe = [1.7, 0.6, -0.5];
    const u = green(xs, xe, k);
    const un = greenDn(xs, nrm, xe, k);
    const S = await op.applySingle(4, un);
    const D = await op.applySingle(5, u);
    const res = new Float64Array(2 * npts);
    for (let i = 0; i < 2 * npts; i++) res[i] = S[i] - D[i] - u[i] / 2;
    let num = 0;
    let den = 0;
    for (let i = 0; i < npts; i++) {
      num += (res[2 * i] ** 2 + res[2 * i + 1] ** 2) * surf.wts[i];
      den += (u[2 * i] ** 2 + u[2 * i + 1] ** 2) * surf.wts[i];
    }
    check("Green's identity S[du/dn] - D[u] - u/2 on the walls (relative, weighted)", Math.sqrt(num / den), 1e-4);
  }

  // 3. the self-test problem
  {
    const xe = [1.7, 0.6, -0.5];
    const rhs = selftestRhs(xs, nrm, xe, k, lam);
    const iters: string[] = [];
    const sol = await op.solve(rhs, 1e-7, 100, (it, r) => iters.push(`${it}:${r.toExponential(1)}`));
    log(`     GMRES: ${sol.niter} iterations, residual ${sol.rres.toExponential(2)}, ${sol.seconds.toFixed(2)} s`);
    const xt = selftestPoints(30, a);
    const ut = await op.eval(xt, sol.x, sol.phi1);
    const uex = green(xt, xe, k);
    check('self-test: relative error at 30 interior points', relError(ut, uex), 2e-4);
    const ub = await op.evalOnSurface(sol.x, sol.phi1);
    const ubex = green(xs, xe, k);
    check('self-test: relative error on the walls', relError(ub, ubex), 2e-4);
  }

  op.destroy();
  log(`${passed} passed, ${failed} failed, ${((performance.now() - t0) / 1000).toFixed(1)} s`, failed ? 'fail' : 'ok');
  window.__RESULTS__ = { passed, failed, lines };
}

main().catch((e) => {
  const msg = e instanceof Error ? `${e.message}\n${e.stack ?? ''}` : String(e);
  log(`FATAL ${msg}`, 'fail');
  window.__RESULTS__ = { passed, failed: failed + 1, fatal: msg, lines };
});