import { describe, expect, it } from "vitest";
import {
  audibleFreqs,
  isolatedAudible,
  maxAbs,
  reference,
  referenceWav,
  solve,
  solveSingle,
  testSpec,
} from "./helpers";

/**
 * The Rayleigh-Ritz method against the diagnostics of aluminum_plate_sim.py
 * (sim-docs): the same basis, quadrature, hole correction, thread spring,
 * damping, pulse and radiation model, so the numbers should agree to
 * rounding, and the first 50 ms of the peak-normalized sound to the 24-bit
 * quantization of the reference WAV.
 */
describe("Rayleigh-Ritz against the reference simulator", () => {
  it("reproduces the centre hit: modes, contact, sound", async () => {
    const ref = reference("rr_center");
    const { result } = await solve(testSpec("center-hit"), "rr");
    const f = audibleFreqs(result);
    expect(f.length).toBe(ref.modes.length);
    ref.modes.forEach((m, k) => {
      expect(Math.abs(f[k] - m.frequency_hz) / m.frequency_hz).toBeLessThan(1e-7);
    });
    // damping and participations, on modes that are not part of a degenerate pair
    for (const k of isolatedAudible(result)) {
      const m = ref.modes[k];
      expect(result.modes.q[k]).toBeCloseTo(m.q_total, 3);
      expect(Math.abs(result.modes.hit[k])).toBeCloseTo(Math.abs(m.hit_participation), 5);
      expect(Math.abs(result.modes.radiation[k])).toBeCloseTo(
        Math.abs(m.radiation_participation),
        8,
      );
    }
    const c = result.contact;
    expect(c.tc).toBeCloseTo(ref.derived.contact_duration_s, 12);
    expect(c.impulse).toBeCloseTo(ref.derived.impact_impulse_Ns, 10);
    expect(c.peak_force).toBeCloseTo(ref.derived.estimated_peak_contact_force_N, 6);
    expect(result.info.critical_frequency).toBeCloseTo(ref.derived.critical_frequency_hz, 6);
    expect(result.info.plate_mass).toBeCloseTo(ref.derived.plate_mass_kg, 12);
    // the sound: the peak of the reference lies in the first 50 ms, so a 100 ms render normalizes alike
    const wav = referenceWav("rr_center");
    expect(result.audio.fs).toBe(wav.fs);
    const a = result.audio.p;
    const g = result.audio.gain;
    let worst = 0;
    for (let i = 0; i < wav.samples.length; i++)
      worst = Math.max(worst, Math.abs(a[i] * g - wav.samples[i]));
    expect(worst).toBeLessThan(1e-4);
    expect(result.audio.peak_pa).toBeCloseTo(ref.derived.unnormalized_peak_pressure_pa_approx, 6);
  });

  it("reproduces the off-centre hit, which wakes the asymmetric modes", async () => {
    const ref = reference("rr_off");
    const { result } = await solve(testSpec("off-center-hit"), "rr");
    // the modes are the same (the strike does not change them) ...
    const f = audibleFreqs(result);
    ref.modes.forEach((m, k) => {
      expect(Math.abs(f[k] - m.frequency_hz) / m.frequency_hz).toBeLessThan(1e-7);
    });
    // ... but the lowest mode, silent at the centre, now takes part
    expect(Math.abs(result.modes.hit[0])).toBeGreaterThan(1);
    for (const k of isolatedAudible(result))
      expect(Math.abs(result.modes.hit[k])).toBeCloseTo(
        Math.abs(ref.modes[k].hit_participation),
        5,
      );
    const wav = referenceWav("rr_off");
    const a = result.audio.p;
    const g = result.audio.gain;
    let worst = 0;
    for (let i = 0; i < wav.samples.length; i++)
      worst = Math.max(worst, Math.abs(a[i] * g - wav.samples[i]));
    expect(worst).toBeLessThan(1e-4);
  });

  it("is stable in the basis order: order 8 moves the first ten modes by less than 0.05 Hz", async () => {
    const ref = reference("rr_center_o8");
    const spec = testSpec();
    spec.rr.basis_order = 8;
    spec.rr.quadrature_order = 20;
    const { result } = await solve(spec, "rr");
    const f = audibleFreqs(result);
    ref.modes.slice(0, 10).forEach((m, k) => {
      expect(Math.abs(f[k] - m.frequency_hz)).toBeLessThan(0.05);
    });
  });

  it("gives the same result as two stages and as one script", async () => {
    const spec = testSpec();
    spec.render.duration = 0.02;
    const [staged, single] = await Promise.all([solve(spec, "rr"), solveSingle(spec, "rr")]);
    expect(maxAbs(staged.result.audio.p)).toBeGreaterThan(0);
    let worst = 0;
    for (let i = 0; i < single.audio.p.length; i++)
      worst = Math.max(worst, Math.abs(single.audio.p[i] - staged.result.audio.p[i]));
    expect(worst).toBeLessThan(1e-12 * maxAbs(single.audio.p));
    expect(single.frames.w.data.length).toBe(staged.result.frames.w.data.length);
  });

  it("carries the mode shapes, curvatures and frames the views draw", async () => {
    const spec = testSpec();
    spec.render.duration = 0.02;
    const { result } = await solve(spec, "rr");
    const g = spec.rr.grid;
    expect(result.grid.x.length).toBe(g);
    expect(result.modes.shapes.nx).toBe(g);
    expect(result.modes.shapes.nf).toBe(result.modes.n);
    expect(result.frames.w.nf).toBe(spec.render.nframes);
    expect(result.frames.vm.nf).toBe(spec.render.nframes);
    expect(maxAbs(result.frames.w.data)).toBeGreaterThan(0);
    expect(maxAbs(result.envelope.vm.data)).toBeGreaterThan(0);
    expect(result.probe.w.length).toBe(result.t.length);
    expect(result.info.n_rigid).toBe(0);
    expect(result.contact.hammer_v.length).toBe(0);
  });
});