// Share of linear (transverse) energy above fc over time: nonlinear strike vs the same strike scaled down.
// Usage: tsx tests/cascade.ts script.m duration fc scale...
import { readFileSync, readdirSync } from "fs";
import { numblEig, runScript } from "../src/sim/runScript";
import { simulate } from "../src/sim/solver";
import { computeModes, isPlanar } from "../src/sim/modes";
const helpers = readdirSync("matlab/helpers").map(name => ({ name, source: readFileSync("matlab/helpers/" + name, "utf8") }));
const { sys } = runScript(readFileSync(process.argv[2], "utf8"), helpers, () => {});
sys.duration = Number(process.argv[3]);
const fc = Number(process.argv[4]);
const modes = computeModes(sys, numblEig, isPlanar(sys) ? "transverse" : "full");
const nm = modes.freqs.length;
const scaleHits = (s: number) => sys.hits.map(h => ({ ...h, velocity: h.velocity.map(v => v * s) as [number, number, number], force: h.force.map(v => v * s) as [number, number, number] }));
const share = (r: ReturnType<typeof simulate>, t: number) => {
  const f = r.frameTimes.findIndex(x => x >= t); if (f < 0) return NaN;
  let hi = 0, tot = 0; for (let k = 0; k < nm; k++) { const e = r.modeEnergy![f * nm + k]; tot += e; if (modes.freqs[k] > fc) hi += e; } return hi / tot;
};
const ts = [0.002, 0.01, 0.03, 0.1, 0.3].filter(t => t < sys.duration);
for (const sc of process.argv.slice(5).map(Number)) {
  const big = simulate({ ...sys, hits: scaleHits(sc) }, { modes });
  const small = simulate({ ...sys, hits: scaleHits(sc * 1e-3) }, { modes });
  console.log(`scale ${sc}: share of energy above ${fc} Hz (linear -> actual): ` + ts.map(t => `${t}s ${(100*share(small, t)).toFixed(1)}%->${(100*share(big, t)).toFixed(1)}%`).join("  ") + ` ${big.aborted ?? ""}`);
}