遠い地震で高層ビルだけが大きく揺れるのは、周期が合うから

const G0 = 9.80665; // 重力加速度(m/s²)const STORY = 3.5; // 1 階分の高さ(m)const PLAN = 24; // 平面の 1 辺(m)const FLOOR_M = PLAN * PLAN * 1000; // 1 階分の重さ(kg)。床 1 m² あたり 1 トンとするconst T_PER_H = 0.03; // 「固さ 1 倍」の 1 次の固有周期 ÷ 高さ(秒/m)。T = 高さ × (0.02 + 0.01α) で α = 1 のときconst TAPER = 0.5; // 最上階の層の固さ ÷ 1 階の層の固さconst H = 1 / 240; // 計算の刻み(秒)const ISO_T = 4, ISO_H = 0.2, ISO_M = 1.5; // 免震層の周期(秒)・減衰定数・免震層の上の床の重さ(1 階分の何倍か)const TMD_RATIO = 0.02; // 制振の重りの重さ ÷ 建物の重さconst CY = 0.55, HARDEN = 0.05; // 層が降伏する力の大きさ(ベースシア係数)と、降伏した後の固さの割合const SOFT_K = 0.25, SOFT_Q = 0.6; // 弱い 1 階の固さと降伏する力の割合(ふつうの階と比べて)const T_MIN = 0.1, T_MAX = 10, NP = 240; // 共振曲線の周期の範囲(秒)と点の数const BIN = 0.25, HN = 240; // 揺れの記録の間隔(秒)と個数(1 分分)const FOV = 30; // 視野角(度)

/** * ホルツァー法: ある λ(= 角振動数の 2 乗)で揺れていると仮定し、1 階から順に「層の力 → 次の階のずれ」を求める。 * 最上階の上に余る力が 0 になる λ が固有値。φ と余る力の符号が変わる回数が、λ より小さい固有値の数になる */function holzer(n: number, m: Float64Array, k: Float64Array, lam: number, phi: Float64Array) {  let u = 1, V = k[0], count = 0;  for (let i = 0; i < n; i++) {    phi[i] = u; V -= lam * m[i] * u;    if (i + 1 < n) { const w = u + V / k[i + 1]; if (u * w <= 0) count++; u = w; } else if (u * V < 0) count++;  }  return count;}/** 小さい方から cnt 個の固有周期と振動の形(いちばん大きいずれを 1 にそろえる)を、区間を半分にしていく方法で求める */function modes(n: number, m: Float64Array, k: Float64Array, cnt: number) {  let up0 = 0; for (let i = 0; i < n; i++) up0 = Math.max(up0, (2 * (k[i] + (i + 1 < n ? k[i + 1] : 0))) / m[i]);  const T: number[] = [], phi: Float64Array[] = [], tmp = FA(n);  for (let j = 0; j < cnt; j++) {    let lo = 0, hi = up0;    for (let it = 0; it < 70; it++) { const mid = (lo + hi) / 2; if (holzer(n, m, k, mid, tmp) > j) hi = mid; else lo = mid; }    const lam = (lo + hi) / 2, p = FA(n); holzer(n, m, k, lam, p);    let big = 1e-30; for (let i = 0; i < n; i++) big = Math.max(big, Math.abs(p[i]));    for (let i = 0; i < n; i++) p[i] /= big;    T.push((2 * Math.PI) / Math.sqrt(lam)); phi.push(p);  }  return { T, phi };}/** 重さ・固さ・減衰を決める。揺れている途中のずれと速さはそのまま残す */function tune(b: Bld, stiff: number, hd: number) {  const N = b.N, f = b.f0, m = FA(N).fill(FLOOR_M), k = FA(N);  for (let i = 0; i < N; i++) k[i] = 1 - ((1 - TAPER) * i) / Math.max(1, N - 1);  // 固さの形を決めてから、1 次の固有周期が T = 0.03 × 高さ になるように全体を拡大する。周期は固さの平方根に反比例するので、固さ 4 倍で周期は半分  const w0 = (2 * Math.PI) / modes(N, m, k, 1).T[0], want = (2 * Math.PI) / (T_PER_H * N * STORY), ks = stiff * (want / w0) ** 2;  for (let i = 0; i < N; i++) k[i] *= ks;  if (b.kind === "soft") k[0] *= SOFT_K;  const bare = modes(N, m, k, Math.min(3, N)); b.T = bare.T; b.phi = bare.phi;  const w1 = (2 * Math.PI) / b.T[0], w2 = (2 * Math.PI) / b.T[Math.min(1, b.T.length - 1)];  // レイリー減衰 C = aM + bK。どちらも隣の階としかつながらないので、計算の形(三重対角)が保たれる  const al = b.kind === "iso" ? 0 : (2 * hd * w1 * w2) / (w1 + w2), be = b.kind === "iso" ? (2 * hd) / w1 : (2 * hd) / (w1 + w2);  // 層が降伏する力: 上に載る重さ × 係数 × Ai の分布(上の階ほど大きな力を受けるとする、建築基準法の分布)  let W = 0; const W1 = N * FLOOR_M * G0, Tq = T_PER_H * N * STORY;  for (let i = N - 1; i >= 0; i--) {    W += FLOOR_M * G0;    const a = W / W1, Ai = 1 + ((1 / Math.sqrt(a) - a) * 2 * Tq) / (1 + 3 * Tq);    b.qy[f + i] = b.yields ? CY * Ai * W : Infinity;    b.m[f + i] = m[i]; b.k[f + i] = k[i]; b.c[f + i] = be * k[i]; b.cm[f + i] = al * m[i];  }  if (b.kind === "soft") b.qy[0] *= SOFT_Q;  const M = N * FLOOR_M;  if (b.kind === "iso") {    // 免震層: 上の全部の重さが周期 ISO_T で揺れる固さ    const Mt = M + ISO_M * FLOOR_M, wb = (2 * Math.PI) / ISO_T;    b.m[0] = ISO_M * FLOOR_M; b.k[0] = Mt * wb * wb; b.c[0] = 2 * ISO_H * Mt * wb; b.cm[0] = 0; b.qy[0] = Infinity;    b.Tf = modes(b.n, b.m, b.k, 1).T[0];  } else if (b.kind === "tmd") {    // 制振の重り: デン・ハルトッグの式。1 次の振動の重さ(最上階を基準)に対する比 μ から決める。周期は建物の (1 + μ) 倍、重り自身の周期で測った減衰定数は √(3μ / (8(1 + μ)))    const p = b.phi[0], top = p[N - 1]; let Ms = 0; for (let i = 0; i < N; i++) Ms += (m[i] * p[i] * p[i]) / (top * top);    const md = TMD_RATIO * M, mu = md / Ms, wd = w1 / (1 + mu), z = Math.sqrt((3 * mu) / (8 * (1 + mu))), j = b.n - 1;    b.m[j] = md; b.k[j] = md * wd * wd; b.c[j] = 2 * z * md * wd; b.cm[j] = 0; b.qy[j] = Infinity; b.Td = (2 * Math.PI) / wd;    b.Tf = b.T[0];  } else b.Tf = b.T[0];  factor(b);  for (let i = 0; i < b.n; i++) { const dl = b.x[i] - (i > 0 ? b.x[i - 1] : 0); b.fs0[i] = clamp(b.k[i] * dl, -b.qy[i], b.qy[i]); b.dl0[i] = dl; }  curveOf(b);}function makeBld(N: number, kind: Kind, name: string, x0: number, stiff: number, hd: number): Bld {  const f0 = kind === "iso" ? 1 : 0, n = N + f0 + (kind === "tmd" ? 1 : 0);  const b: Bld = {    kind, name, N, x0, n, f0, top: f0 + N - 1,    m: FA(n), k: FA(n), c: FA(n), cm: FA(n), qy: FA(n), off: FA(n), dd: FA(n), ll: FA(n), x: FA(n), v: FA(n), a: FA(n), xp: FA(n), vp: FA(n), r: FA(n), sh: FA(n), vv: FA(n), fs: FA(n), fs0: FA(n), dl0: FA(n),    T: [], phi: [], Tf: 0, Td: 0, yields: kind === "soft" || kind === "frame", pull: null,    curve: new Float32Array(NP), hA: new Float32Array(HN), hB: new Float32Array(HN), lo: Infinity, hi: -Infinity, loT: 0, hiT: 0,    amp: 0, peakTop: 0, peakAbs: 0, peakDev: 0, ampDev: 0, w0: 0, w1: 0, d0: 0, d1: 0, wt: 0, peak: FA(N), env: FA(N), maxD: 0, maxAt: 1, last: 0, tUp: -1, Tm: 0, noted: false, yielded: false,  };  tune(b, stiff, hd);  return b;}const periodAt = (j: number) => T_MIN * Math.pow(T_MAX / T_MIN, j / (NP - 1));/** 共振曲線: 周期ごとに正弦波の地面の揺れを入れたときの、最上階の揺れ ÷ 地面の揺れ。複素数の三重対角の式を周期の数だけ解く(降伏は考えない) */function curveOf(b: Bld) {  const n = b.n, { m, k, c, cm } = b, cr = FA(n), ci = FA(n), dr = FA(n), di = FA(n);  for (let j = 0; j < NP; j++) {    const w = (2 * Math.PI) / periodAt(j);    let pr = 0, pi = 0, qr = 0, qi = 0; // 1 つ下の行の消去の結果 c'・d'    for (let i = 0; i < n; i++) {      const k1 = i + 1 < n ? k[i + 1] : 0, c1 = i + 1 < n ? c[i + 1] : 0;      const zr = k[i] + k1 - w * w * m[i], zi = w * (c[i] + c1 + cm[i]);      // 下の階とのつながり o = −(k + iωc)      const or = -k[i], oi = -w * c[i], ur = -k1, ui = -w * c1;      const nr = i ? zr - (or * pr - oi * pi) : zr, ni = i ? zi - (or * pi + oi * pr) : zi, den = nr * nr + ni * ni;      const br = w * w * m[i] - (i ? or * qr - oi * qi : 0), bi = -(i ? or * qi + oi * qr : 0);      pr = (ur * nr + ui * ni) / den; pi = (ui * nr - ur * ni) / den;      qr = (br * nr + bi * ni) / den; qi = (bi * nr - br * ni) / den;      cr[i] = pr; ci[i] = pi; dr[i] = qr; di[i] = qi;    }    let xr = dr[n - 1], xi = di[n - 1];    for (let i = n - 2; i >= b.top; i--) { const r = dr[i] - (cr[i] * xr - ci[i] * xi), q = di[i] - (cr[i] * xi + ci[i] * xr); xr = r; xi = q; }    b.curve[j] = Math.hypot(xr, xi);  }}

interface Wave { dt: number; n: number; a: Float32Array; v: Float32Array; d: Float32Array; dur: number; Tc: number; pa: number }/** 合成波の決まり: 濾過の中心の周期 Tc と減衰 zc、ゆっくりした揺れを除く周期 Th、強さの時間変化(t1 まで増え、t2 まで一定、その後は率 dk で減る)、長さ dur、乱数の種 */const WAVES: Record<string, { Tc: number; zc: number; Th: number; t1: number; t2: number; dk: number; dur: number; seed: number }> = {  short: { Tc: 0.5, zc: 0.6, Th: 3, t1: 1.5, t2: 9, dk: 0.18, dur: 40, seed: 7 },  long: { Tc: 5.5, zc: 0.25, Th: 40, t1: 12, t2: 50, dk: 0.12, dur: 90, seed: 11 },};const waveCache: Record<string, Wave> = {};/** * 乱数(白色雑音)を 2 段のばねとダッシュポットに通して、地面の揺れを作る。1 段目は周期 Tc のあたりを残し、2 段目はとてもゆっくりした成分を除く。 * 2 段目の出力を地面のずれとして使うので、速度と加速度はそのずれの微分になり、最後は元の位置に戻る。最大速度を 1 m/s にそろえておき、つまみの値を掛けて使う */function makeWave(api: Api, key: string): Wave {  if (waveCache[key]) return waveCache[key];  const c = WAVES[key], dt = 0.01, n = Math.ceil(c.dur / dt), rnd = api.rng(c.seed);  const a = new Float32Array(n), v = new Float32Array(n), d = new Float32Array(n);  const wf = (2 * Math.PI) / c.Tc, wh = (2 * Math.PI) / c.Th, zh = 0.7;  let xf = 0, vf = 0, y = 0, vy = 0, big = 1e-9;  for (let i = 0; i < n; i++) {    const t = i * dt, e = t < c.t1 ? (t / c.t1) ** 2 : t < c.t2 ? 1 : Math.exp(-c.dk * (t - c.t2));    const w = (e * Math.sqrt(-2 * Math.log(1 - rnd())) * Math.cos(2 * Math.PI * rnd())) / Math.sqrt(dt);    vf += dt * (-2 * c.zc * wf * vf - wf * wf * xf - w); xf += dt * vf;    const g = -(2 * c.zc * wf * vf + wf * wf * xf);    const ay = g - 2 * zh * wh * vy - wh * wh * y;    vy += dt * ay; y += dt * vy;    a[i] = ay; v[i] = vy; d[i] = y; big = Math.max(big, Math.abs(vy));  }  let pa = 0;  for (let i = 0; i < n; i++) { a[i] /= big; v[i] /= big; d[i] /= big; pa = Math.max(pa, Math.abs(a[i])); }  return (waveCache[key] = { dt, n, a, v, d, dur: c.dur, Tc: c.Tc, pa });}interface Ground { t: number; ph: number; ag: number; vg: number; dg: number; T: number; wait: number; ended: boolean; hmin: number; hmax: number; hA: Float32Array; hB: Float32Array }const SWEEP = { T0: 0.3, T1: 5, dur: 90 };const sweepT = (t: number) => SWEEP.T0 * Math.pow(SWEEP.T1 / SWEEP.T0, clamp(t / SWEEP.dur, 0, 1));/** 地面の加速度・速度・ずれを刻み 1 つ分進める。正弦波と周期を上げていく揺れは位相を足して作る(周期のつまみを動かしても波がつながる) */function stepGround(api: Api, s: S, h: number) {  const g = s.g, key = api.str("wave"), V = api.num("amp") / 100;  g.t += h;  if (key === "sine" || key === "sweep") {    const T = key === "sine" ? api.num("period") : sweepT(g.t), w = (2 * Math.PI) / T;    const env = key === "sine" ? smooth(g.t / (2.5 * T)) : smooth(g.t / 3) * smooth((SWEEP.dur - g.t) / 3);    g.T = T; g.ph += w * h;    g.vg = V * env * Math.sin(g.ph); g.ag = V * env * w * Math.cos(g.ph);    g.dg += h * (g.vg - g.dg / 90);    if (key === "sweep" && g.t > SWEEP.dur) endWave(api, s, h);  } else if (key === "short" || key === "long") {    const W = makeWave(api, key), q = g.t / W.dt, i = Math.floor(q), f = q - i;    if (i + 1 < W.n) {      g.ag = V * (W.a[i] + f * (W.a[i + 1] - W.a[i])); g.vg = V * (W.v[i] + f * (W.v[i + 1] - W.v[i])); g.dg = V * (W.d[i] + f * (W.d[i + 1] - W.d[i]));    } else { g.ag = 0; g.vg = 0; g.dg *= Math.exp(-h / 2); endWave(api, s, h); }    g.T = W.Tc;  } else { g.ag = 0; g.vg = 0; g.dg *= Math.exp(-h / 2); g.T = 0; }}/** 地面の揺れが終わったら知らせ、しばらく最上階の揺れを見せてから、同じ揺れをもう一度与える */function endWave(api: Api, s: S, h: number) {  const g = s.g;  if (!g.ended) { g.ended = true; g.wait = 30; api.note("地面の揺れが止まった。建物は固有周期で揺れ続ける", 3200); }  g.wait -= h;  if (g.wait <= 0) { restartGround(s); api.note("同じ揺れをもう一度与える", 2200); }}function restartGround(s: S) {  const g = s.g; g.t = 0; g.ph = 0; g.ended = false; g.wait = 0;  resetPeaks(s); s.sw.fill(0);}function resetPeaks(s: S) {  for (const b of s.bs) { b.peakTop = 0; b.peakAbs = 0; b.peakDev = 0; b.peak.fill(0); b.maxD = 0; b.maxAt = 1; b.noted = false; b.yielded = false; }}

/** 内側の式を解くための下準備(三重対角の行列を前から消しておく)。k̂ = K + (2/h)C + (4/h²)M で、毎刻み同じものを使える */function factor(b: Bld) {  const n = b.n, a0 = 4 / (H * H), a1 = 2 / H;  for (let i = 0; i < n; i++) {    const k1 = i + 1 < n ? b.k[i + 1] : 0, c1 = i + 1 < n ? b.c[i + 1] : 0;    const d = b.k[i] + k1 + a1 * (b.c[i] + c1 + b.cm[i]) + a0 * b.m[i] + (b.pull && b.pull.i === i ? b.pull.kp : 0);    b.off[i] = -(k1 + a1 * c1);    b.ll[i] = i ? b.off[i - 1] / b.dd[i - 1] : 0;    b.dd[i] = i ? d - b.ll[i] * b.off[i - 1] : d;  }}/** * ニューマーク法(β = 1/4、γ = 1/2)で 1 刻み進める。周期の短い振動があっても発散しない。 * 新しい時刻の加速度と速度をずれ x で表し、M a + C v + (層の力) = −M a_g を x について解く。 * 層の力は「傾き k の直線を上限と下限で切った形」(降伏)。降伏するときは残りの力が 0 になるまで解き直す */function stepBld(b: Bld, ag: number) {  const n = b.n, { m, k, c, cm, qy, x, v, a, xp, vp, r, sh, vv, fs, fs0, dl0, dd, ll, off } = b, h = H, a0 = 4 / (h * h), a1 = 2 / h;  for (let i = 0; i < n; i++) { xp[i] = x[i] + h * v[i] + 0.25 * h * h * a[i]; vp[i] = v[i] + 0.5 * h * a[i]; x[i] = xp[i]; }  const iters = b.yields ? 10 : 1, P = b.pull;  for (let it = 0; it < iters; it++) {    let xl = 0, vl = 0;    for (let i = 0; i < n; i++) {      const vi = vp[i] + a1 * (x[i] - xp[i]), dl = x[i] - xl;      let f = fs0[i] + k[i] * (dl - dl0[i]);      if (qy[i] < Infinity) { const hard = HARDEN * k[i] * dl, top = qy[i] * (1 - HARDEN); f = clamp(f, hard - top, hard + top); }      fs[i] = f; vv[i] = vi; sh[i] = f + c[i] * (vi - vl);      xl = x[i]; vl = vi;    }    for (let i = 0; i < n; i++) r[i] = -m[i] * ag - m[i] * a0 * (x[i] - xp[i]) - cm[i] * vv[i] - sh[i] + (i + 1 < n ? sh[i + 1] : 0);    if (P) r[P.i] += P.kp * (P.target - x[P.i]);    for (let i = 1; i < n; i++) r[i] -= ll[i] * r[i - 1];    r[n - 1] /= dd[n - 1];    for (let i = n - 2; i >= 0; i--) r[i] = (r[i] - off[i] * r[i + 1]) / dd[i];    let big = 0;    for (let i = 0; i < n; i++) { x[i] += r[i]; big = Math.max(big, Math.abs(r[i])); }    if (big < 1e-6) break;  }  let xl = 0;  for (let i = 0; i < n; i++) {    const d = x[i] - xp[i]; a[i] = a0 * d; v[i] = vp[i] + a1 * d;    const dl = x[i] - xl;    let f = fs0[i] + k[i] * (dl - dl0[i]);    if (qy[i] < Infinity) { const hard = HARDEN * k[i] * dl, top = qy[i] * (1 - HARDEN); f = clamp(f, hard - top, hard + top); }    fs0[i] = f; dl0[i] = dl; xl = x[i];  }}/** 測る: 層間変形角の最大・最上階の揺れの幅と加速度・上向きに 0 を横切る間隔(測った周期) */function track(api: Api, s: S, b: Bld) {  const j = b.top, xt = b.x[j], T = Math.max(b.Tf, s.g.T || 0), fd = Math.exp(-H / 1.5);  // 揺れの幅: 今の周期のあいだの最大(半周期ごとに 2 つの枠を入れ替えるので、揺れの途中でも値が波打たない)  const dev = b.kind === "iso" ? Math.abs(b.x[0]) : b.kind === "tmd" ? Math.abs(b.x[b.n - 1] - b.x[j]) : 0;  b.wt += H;  if (b.wt >= 0.55 * T) { b.wt = 0; b.w0 = b.w1; b.w1 = 0; b.d0 = b.d1; b.d1 = 0; }  b.w1 = Math.max(b.w1, Math.abs(xt)); b.d1 = Math.max(b.d1, dev);  b.amp = Math.max(b.w0, b.w1); b.ampDev = Math.max(b.d0, b.d1);  b.peakTop = Math.max(b.peakTop, Math.abs(xt));  b.peakAbs = Math.max(b.peakAbs, Math.abs(b.a[j] + s.g.ag));  if (b.last < 0 && xt >= 0 && b.amp > 1e-3) { if (b.tUp >= 0) b.Tm = s.clock - b.tUp; b.tUp = s.clock; }  b.last = xt;  for (let i = 0; i < b.N; i++) {    const q = b.f0 + i, d = Math.abs(b.x[q] - (q > 0 ? b.x[q - 1] : 0)) / STORY;    b.env[i] = Math.max(d, b.env[i] * fd);    if (d > b.peak[i]) b.peak[i] = d;    if (d > b.maxD) { b.maxD = d; b.maxAt = i + 1; }  }  b.peakDev = Math.max(b.peakDev, dev);  if (!b.noted && b.maxD > 0.01 && !b.pull) { b.noted = true; api.note(`${b.name}: 層間変形角が 1/100 を超えた(${b.maxAt} 階)`, 2600); }  if (b.kind === "soft" && !b.yielded && Math.abs(b.fs0[0]) >= b.qy[0] * 0.98) { b.yielded = true; api.note("弱い 1 階が降伏した。ここから変形が 1 階に集まる", 3000); }}function tick(api: Api, s: S, h: number) {  s.clock += h;  stepGround(api, s, h);  const g = s.g;  // 指で引く位置はなめらかに追いかける(急に引くと、波が建物を伝わって一部の階だけが大きくずれる)  if (s.grab) { const P = s.grab.b.pull!; P.target += (s.grab.want - P.target) * (1 - Math.exp(-h / 0.3)); }  for (const b of s.bs) {    stepBld(b, g.ag); track(api, s, b);    const xt = b.x[b.top];    if (xt < b.lo) { b.lo = xt; b.loT = g.t; }    if (xt > b.hi) { b.hi = xt; b.hiT = g.t; }  }  g.hmin = Math.min(g.hmin, g.vg); g.hmax = Math.max(g.hmax, g.vg);  // 周期を上げていく揺れでは、周期ごとの最上階の揺れの比を記録する(共振曲線の上に点で描く)  if (api.str("wave") === "sweep" && !g.ended) {    const j = Math.round(((NP - 1) * Math.log(g.T / T_MIN)) / Math.log(T_MAX / T_MIN)), D = (api.num("amp") / 100 * g.T) / (2 * Math.PI), b = s.bs[s.bs.length - 1];    if (j >= 0 && j < NP) s.sw[j] = Math.max(s.sw[j], Math.abs(b.x[b.top]) / Math.max(1e-6, D));  }  s.bin += h;  if (s.bin >= BIN) {    s.bin -= BIN;    // 区間の中の最小と最大を、起きた順に記録する(ゆっくりした揺れはなめらかな線に、速い揺れは帯になる)    for (const b of s.bs) { const first = b.loT <= b.hiT; b.hA[s.hi] = first ? b.lo : b.hi; b.hB[s.hi] = first ? b.hi : b.lo; b.lo = Infinity; b.hi = -Infinity; }    g.hA[s.hi] = g.hmin; g.hB[s.hi] = g.hmax; g.hmin = Infinity; g.hmax = -Infinity;    s.hi = (s.hi + 1) % HN; s.hn = Math.min(HN, s.hn + 1);  }}

const MAX_COL = 4700, MAX_SLAB = 170, MAX_BEAR = 24;/** 平面の上の柱の位置: 外周の 12 本(画質が低いときは四隅の 4 本) */const COL_P: [number, number][] = [], COL_C: [number, number][] = [];{  const e = PLAN / 2 - 1, q = e / 3;  for (const x of [-e, -q, q, e]) for (const z of [-e, -q, q, e]) { if (Math.abs(x) === e || Math.abs(z) === e) COL_P.push([x, z]); if (Math.abs(x) === e && Math.abs(z) === e) COL_C.push([x, z]); }}interface Meshes { col: THREE.InstancedMesh; tone: Float32Array; toneA: THREE.InstancedBufferAttribute; slab: THREE.InstancedMesh; bear: THREE.InstancedMesh; ball: THREE.Mesh; cable: THREE.LineSegments; cableP: Float32Array; ground: THREE.Group }function makeMeshes(api: Api, v: View3): Meshes {  const box = new THREE.BoxGeometry(1, 1, 1);  const col = new THREE.InstancedMesh(box, v.mat(1), MAX_COL);  const tone = new Float32Array(MAX_COL * 3), toneA = new THREE.InstancedBufferAttribute(tone, 3);  toneA.setUsage(THREE.DynamicDrawUsage); col.instanceColor = toneA; col.instanceMatrix.setUsage(THREE.DynamicDrawUsage);  const slab = new THREE.InstancedMesh(new THREE.BoxGeometry(1, 1, 1), v.mat(0.82), MAX_SLAB);  slab.instanceMatrix.setUsage(THREE.DynamicDrawUsage);  const bear = new THREE.InstancedMesh(new THREE.CylinderGeometry(0.5, 0.5, 1, 12), v.mat(1), MAX_BEAR);  bear.instanceMatrix.setUsage(THREE.DynamicDrawUsage);  for (const m of [col, slab, bear]) { m.castShadow = true; m.receiveShadow = true; m.frustumCulled = false; m.count = 0; }  const ball = new THREE.Mesh(new THREE.SphereGeometry(1, 24, 16), v.mat(0.95, { shine: 0.4 }));  ball.castShadow = true; ball.visible = false;  const cableP = new Float32Array(12), cg = new THREE.BufferGeometry();  cg.setAttribute("position", new THREE.BufferAttribute(cableP, 3));  const cable = new THREE.LineSegments(cg, v.lineMat(1)); cable.frustumCulled = false; cable.visible = false;  // 地面と町並み。建物の土台と同じだけ動く(カメラは揺れない)  const ground = new THREE.Group();  // 地面は奥の面より十分に広くしておく(端が見えると、そこに輪郭線が出る)  const plane = new THREE.Mesh(new THREE.PlaneGeometry(60000, 60000).rotateX(-Math.PI / 2), v.mat(0.14));  plane.receiveShadow = true; ground.add(plane);  const rnd = api.rng(5), CN = 60, city = new THREE.InstancedMesh(new THREE.BoxGeometry(1, 1, 1), v.mat(0.3, { flat: true }), CN), mt = new THREE.Matrix4();  for (let i = 0; i < CN; i++) {    const w = 14 + rnd() * 18, d = 14 + rnd() * 18, h = 5 + rnd() * rnd() * 34, x = -560 + rnd() * 1120, z = -130 - rnd() * 300;    mt.makeScale(w, h, d); mt.setPosition(x, h / 2, z); city.setMatrixAt(i, mt);  }  city.castShadow = true; city.receiveShadow = true; city.frustumCulled = false; ground.add(city);  v.scene.add(col, slab, bear, ball, cable, ground);  return { col, tone, toneA, slab, bear, ball, cable, cableP, ground };}/** 箱 1 個の並べ方の行列を書き込む。縦の軸を shx だけ横へ傾ける(平行四辺形の柱)ので、上下の階にちょうどつながる */function putBox(a: ArrayLike<number> & Float32Array, i: number, sx: number, sy: number, sz: number, shx: number, ox: number, oy: number, oz: number) {  const o = i * 16;  a[o] = sx; a[o + 1] = 0; a[o + 2] = 0; a[o + 3] = 0;  a[o + 4] = shx; a[o + 5] = sy; a[o + 6] = 0; a[o + 7] = 0;  a[o + 8] = 0; a[o + 9] = 0; a[o + 10] = sz; a[o + 11] = 0;  a[o + 12] = ox; a[o + 13] = oy; a[o + 14] = oz; a[o + 15] = 1;}const S_CURVE = (q: number) => q * q * (3 - 2 * q); // 上下の床に直角につながる柱の曲がり方/** 計算の結果を、柱と床の位置に写す。three.js が受け持つのはここから先の描画だけで、毎コマ行列の配列を書き換えるだけにする */function place(api: Api, s: S) {  const K = api.num("scale"), gx = K * s.g.dg, m = s.m, ca = m.col.instanceMatrix.array as Float32Array, sa = m.slab.instanceMatrix.array as Float32Array, ba = m.bear.instanceMatrix.array as Float32Array, ton = m.tone;  const q = api.quality, seg = q >= 2 ? 3 : q === 1 ? 2 : 1, pts = q === 0 ? COL_C : COL_P, CW = 1.1, ST = 0.6;  let ic = 0, is = 0, ib = 0;  m.ground.position.x = gx; m.ball.visible = false; m.cable.visible = false;  for (const b of s.bs) {    const X = b.x0 + gx;    putBox(sa, is++, PLAN + 10, 1, PLAN + 10, 0, X, 0, 0); // 土台    let y = 0.5, xl = 0;    if (b.kind === "iso") {      // 免震層: 積層ゴムの円柱を、土台と床のあいだで傾ける      const xm = K * b.x[0], hb = 2.4;      for (const [cx, cz] of pts) putBox(ba, ib++, 1.9, hb, 1.9, xm, X + cx + xm / 2, y + hb / 2, cz);      y += hb; putBox(sa, is++, PLAN + 2, 0.9, PLAN + 2, 0, X + xm, y + 0.45, 0); y += 0.9; xl = xm;    }    for (let i = 0; i < b.N; i++) {      const xi = K * b.x[b.f0 + i], dl = xi - xl, t = 0.42 + 0.58 * Math.min(1, b.env[i] * 150);      for (const [cx, cz] of pts) {        for (let k = 0; k < seg; k++) {          const a0 = k / seg, a1 = (k + 1) / seg, s0 = S_CURVE(a0), s1 = S_CURVE(a1);          putBox(ca, ic, CW, STORY / seg, CW, dl * (s1 - s0), X + cx + xl + (dl * (s0 + s1)) / 2, y + (STORY * (a0 + a1)) / 2, cz);          ton[ic * 3] = t; ton[ic * 3 + 1] = t; ton[ic * 3 + 2] = t; ic++;        }      }      y += STORY;      putBox(sa, is++, PLAN + 0.8, ST, PLAN + 0.8, 0, X + xi, y - ST / 2, 0);      xl = xi;    }    if (b.kind === "tmd") {      // 屋上の枠から吊った重り。重りのずれは最上階から測り、同じ倍率で描く      const hf = 18, xr = xl, xb = X + xr + K * (b.x[b.n - 1] - b.x[b.top]), yb = y + 8, R = 5.5, W = 9;      for (const cx of [-W, W]) for (const cz of [-W, W]) { putBox(ca, ic, 0.9, hf, 0.9, 0, X + xr + cx, y + hf / 2, cz); ton[ic * 3] = ton[ic * 3 + 1] = ton[ic * 3 + 2] = 0.8; ic++; }      putBox(sa, is++, 2 * W + 2, 0.8, 2 * W + 2, 0, X + xr, y + hf, 0);      m.ball.visible = true; m.ball.scale.setScalar(R); m.ball.position.set(xb, yb, 0);      const c = m.cableP, top = y + hf - 0.4;      c[0] = X + xr; c[1] = top; c[2] = 5; c[3] = xb; c[4] = yb + R * 0.6; c[5] = 2;      c[6] = X + xr; c[7] = top; c[8] = -5; c[9] = xb; c[10] = yb + R * 0.6; c[11] = -2;      m.cable.visible = true; (m.cable.geometry.attributes.position as THREE.BufferAttribute).needsUpdate = true;    }  }  m.col.count = ic; m.slab.count = is; m.bear.count = ib;  m.col.instanceMatrix.needsUpdate = true; m.slab.instanceMatrix.needsUpdate = true; m.bear.instanceMatrix.needsUpdate = true; m.toneA.needsUpdate = true;  // 輪郭線は奥行きの精度に頼るので、見ている距離に合わせて手前と奥の面を決める  const d = s.v.orbit.dist; s.v.camera.near = Math.max(1, d * 0.12); s.v.camera.far = d * 3.2;  if (Math.abs(d - s.fogD) > 0.05 * s.fogD) { s.fogD = d; s.v.set({ fog: [d * 1.2, d * 2.5] }); }}

interface Grab { b: Bld; i: number; px: number; py: number; ex: number; ey: number; x0: number; lim: number; want: number }/** 指の下の建物と階を探す(画面の上で、土台と屋上を結ぶ線からの距離で決める) */function pick(api: Api, s: S) {  const v = s.v, P = api.pointer, K = api.num("scale"), gx = K * s.g.dg, R = api.u(P.kind === "touch" ? 14 : 6);  let best: { b: Bld; f: number } | null = null, bd = Infinity;  for (const b of s.bs) {    const yTop = (b.kind === "iso" ? 3.8 : 0.5) + b.N * STORY, xr = b.x0 + gx + K * b.x[b.top];    const a = v.project([b.x0 + gx, 0, 0]), c = v.project([xr, yTop, 0]), e = v.project([xr + PLAN / 2, yTop, 0]);    const ux = c.x - a.x, uy = c.y - a.y, L2 = Math.max(1, ux * ux + uy * uy), t = clamp(((P.x - a.x) * ux + (P.y - a.y) * uy) / L2, 0, 1);    const dist = Math.hypot(P.x - (a.x + t * ux), P.y - (a.y + t * uy)), half = Math.abs(e.x - c.x) + R;    if (dist < half && dist < bd) { bd = dist; best = { b, f: clamp(Math.round(t * b.N), 1, b.N) }; }  }  return best;}function startPull(api: Api, s: S, hit: { b: Bld; f: number }) {  const b = hit.b, i = b.f0 + hit.f - 1, K = api.num("scale"), gx = K * s.g.dg, X = b.x0 + gx + K * b.x[i], Y = hit.f * STORY;  const a = s.v.project([X, Y, 0]), c = s.v.project([X + 1, Y, 0]);  // 引くばね: 層の固さより十分に固くして、指で引いた階が指についてくるようにする  b.pull = { i, target: b.x[i], kp: 6 * b.k[i] }; factor(b); b.noted = true;  s.grab = { b, i, px: api.pointer.x, py: api.pointer.y, ex: c.x - a.x, ey: c.y - a.y, x0: b.x[i], lim: (hit.f * STORY) / 20, want: b.x[i] };}function movePull(api: Api, s: S) {  const g = s.grab!, p = api.pointer, K = api.num("scale");  const off = ((p.x - g.px) * g.ex + (p.y - g.py) * g.ey) / Math.max(0.04, g.ex * g.ex + g.ey * g.ey);  g.want = clamp(g.x0 + off / K, -g.lim, g.lim);}function endPull(api: Api, s: S) {  const g = s.grab!, b = g.b, moved = Math.abs(b.x[b.top]) > 1e-3;  b.pull = null; factor(b); s.grab = null;  // 引いている間の変形は測った値に入れない。離した後の揺れだけを測り直す  b.peak.fill(0); b.maxD = 0; b.maxAt = 1; b.peakTop = 0; b.peakAbs = 0;  if (moved) { b.tUp = -1; b.Tm = 0; api.note(`離した: ${b.name} は周期 ${fmt(b.Tf, 2)} 秒で揺れ続ける`, 2800); }}

function panels(api: Api) {  const u = api.u, pad = u(10);  if (api.w < 760) {    const H1 = clamp(api.h * 0.1, u(62), u(82)), H2 = clamp(api.h * 0.15, u(84), u(120)), w2 = (api.w - 3 * pad) / 2, y2 = api.h - pad - H2, y1 = y2 - pad * 0.7 - H1;    return { hist: { x: pad, y: y1, w: api.w - 2 * pad, h: H1 }, res: { x: pad, y: y2, w: w2, h: H2 }, dr: { x: 2 * pad + w2, y: y2, w: w2, h: H2 }, frac: (api.h - y1 + pad) / api.h, narrow: true };  }  const ph = clamp(api.h * 0.24, u(100), u(176)), pw = Math.min((api.w - 4 * pad) / 3, u(350)), y = api.h - pad - ph;  return { hist: { x: pad, y, w: pw, h: ph }, res: { x: 2 * pad + pw, y, w: pw, h: ph }, dr: { x: 3 * pad + 2 * pw, y, w: pw, h: ph }, frac: (ph + 2 * pad) / api.h, narrow: false };}/** 上へ切り上げた、きりのよい数 */function nice(x: number) {  const e = Math.pow(10, Math.floor(Math.log10(Math.max(1e-9, x)))), q = x / e;  return (q <= 1 ? 1 : q <= 1.5 ? 1.5 : q <= 2 ? 2 : q <= 3 ? 3 : q <= 5 ? 5 : q <= 7 ? 7 : 10) * e;}const lineStyle = (api: Api, i: number, n: number): { w: number; dash?: number[] } => (i === n - 1 ? { w: api.u(2.2) } : i === n - 2 ? { w: api.u(1) } : { w: api.u(1), dash: [api.u(3), api.u(3)] });const lenText = (m: number) => (m >= 1 ? `${fmt(m, 2)} m` : m >= 0.1 ? `${fmt(m * 100, 0)} cm` : m >= 0.01 ? `${fmt(m * 100, 1)} cm` : `${fmt(m * 1000, 1)} mm`);const roundText = (m: number) => (m >= 1 ? `${m} m` : m >= 0.01 ? `${Math.round(m * 100)} cm` : `${Math.round(m * 1000)} mm`);const angleText = (d: number) => (d < 1e-5 ? "0" : `1/${Math.round(1 / d)}`);function legend(api: Api, s: S, R: Rect, fs: number) {  let lx = R.x + R.w - api.u(6);  for (let i = 0; i < s.bs.length; i++) {    const b = s.bs[i];    api.text(b.name, lx, R.y + fs * 1.05, { size: fs * 0.8, align: "right" });    lx -= api.ctx.measureText(b.name).width + api.u(3);    api.line(lx - api.u(12), R.y + fs * 0.72, lx, R.y + fs * 0.72, lineStyle(api, i, s.bs.length));    lx -= api.u(18);  }}function drawHist(api: Api, s: S, R: Rect, fs: number, narrow: boolean) {  const u = api.u; api.rect(R.x, R.y, R.w, R.h, { fill: "paper", w: u(1) });  api.text(R.w < u(345) && s.bs.length > 2 ? "最上階の揺れ" : "最上階の揺れ(過去 1 分)", R.x + u(6), R.y + fs * 1.05, { size: fs * 0.88 });  legend(api, s, R, fs);  const x0 = R.x + Math.max(u(38), fs * 3.9), w = R.x + R.w - u(6) - x0, y0 = R.y + fs * 1.5, h = R.h - fs * 1.5 - u(5), gh = h * 0.24, ih = h - gh - u(3), n = s.hn, pts = s.pts;  // 地面の速さ(上の帯)  const g = s.g; let vm = 1e-6; for (let k = 0; k < n; k++) vm = Math.max(vm, Math.abs(g.hA[k]), Math.abs(g.hB[k]));  const gy = y0 + gh / 2;  api.line(x0, gy, x0 + w, gy, { w: u(0.6), dash: [u(2), u(4)] });  api.text("地面", x0 - u(5), gy + fs * 0.3, { size: fs * 0.78, align: "right" });  const put = (A: Float32Array, B: Float32Array, y: (q: number) => number) => {    let c = 0;    const dx = w / (HN - 1) / 4;    for (let k = 0; k < n; k++) { const i = (s.hi - n + k + HN) % HN, px = x0 + w - ((n - 1 - k) / (HN - 1)) * w; pts[c++] = px - dx; pts[c++] = y(A[i]); pts[c++] = px + dx; pts[c++] = y(B[i]); }    return c;  };  if (n > 1) api.poly(pts.subarray(0, put(g.hA, g.hB, (q) => gy - (q / vm) * gh * 0.46)), { w: u(0.8) });  // 最上階  let m = 0.001; for (const b of s.bs) for (let k = 0; k < n; k++) m = Math.max(m, Math.abs(b.hA[k]), Math.abs(b.hB[k]));  const top = nice(m), cy = y0 + gh + u(3) + ih / 2, Y = (q: number) => cy - (clamp(q, -top, top) / top) * (ih / 2);  api.line(x0, cy, x0 + w, cy, { w: u(0.7) });  api.text(`±${roundText(top)}`, x0 - u(4), Y(top) + fs * 0.75, { size: fs * 0.72, align: "right" });  if (n > 1) s.bs.forEach((b, i) => { api.poly(pts.subarray(0, put(b.hA, b.hB, Y)), lineStyle(api, i, s.bs.length)); });}function drawRes(api: Api, s: S, R: Rect, fs: number, narrow: boolean) {  const u = api.u; api.rect(R.x, R.y, R.w, R.h, { fill: "paper", w: u(1) });  api.text(narrow ? "共振曲線" : "共振曲線(最上階のずれ ÷ 地面の揺れ幅)", R.x + u(6), R.y + fs * 1.05, { size: fs * 0.88 });  const x0 = R.x + u(narrow ? 22 : 28), w = R.w - u(narrow ? 28 : 36), y0 = R.y + fs * 1.5, h = R.h - fs * 1.5 - fs * 1.35;  const X = (T: number) => x0 + (w * Math.log(clamp(T, T_MIN, T_MAX) / T_MIN)) / Math.log(T_MAX / T_MIN);  let m = 3; for (const b of s.bs) for (let j = 0; j < NP; j++) m = Math.max(m, b.curve[j]);  const top = nice(m), Y = (q: number) => y0 + h - (Math.min(q, top) / top) * h;  api.line(x0, y0 + h, x0 + w, y0 + h, { w: u(0.8) }); api.line(x0, y0, x0, y0 + h, { w: u(0.8) });  for (const T of narrow ? [0.2, 1] : [0.2, 0.5, 1, 2, 5]) { api.line(X(T), y0 + h, X(T), y0 + h + u(3), { w: u(0.8) }); api.text(`${T}`, X(T), y0 + h + fs * 1.1, { size: fs * 0.76, align: "center", font: "mono" }); }  api.text("10 秒", x0 + w, y0 + h + fs * 1.1, { size: fs * 0.76, align: "right", font: "mono" });  api.text(`${top}`, x0 - u(3), y0 + fs * 0.5, { size: fs * 0.76, align: "right", font: "mono" });  api.text("1", x0 - u(3), Y(1) + fs * 0.3, { size: fs * 0.76, align: "right", font: "mono" });  // 入力の周期の印  const key = api.str("wave"), g = s.g;  if (key === "short" || key === "long") {    const Tc = WAVES[key].Tc, a = X(Tc / 1.6), b = X(Tc * 1.6);    api.hatch([[a, y0], [b, y0], [b, y0 + h], [a, y0 + h]], { gap: u(5), w: u(0.7) });  } else if (key === "sine" || key === "sweep") {    const T = g.T || api.num("period");    api.line(X(T), y0, X(T), y0 + h, { w: u(1.1), dash: [u(3), u(3)] });    const lab = `${fmt(T, 2)} 秒`, xL = clamp(X(T), x0 + u(26), x0 + w - u(26));    api.label(lab, xL, y0 + fs * 0.75, { size: fs * 0.78, align: "center", font: "mono", pad: 1 });  }  const pts = s.pts;  s.bs.forEach((b, i) => {    for (let j = 0; j < NP; j++) { pts[2 * j] = X(periodAt(j)); pts[2 * j + 1] = Y(b.curve[j]); }    api.poly(pts.subarray(0, 2 * NP), lineStyle(api, i, s.bs.length));  });  // 周期を上げていく揺れで、実際に測った比を点で重ねる  if (key === "sweep") for (let j = 0; j < NP; j++) if (s.sw[j] > 0) api.dot(X(periodAt(j)), Y(s.sw[j]), u(1.7));}function drawDrift(api: Api, s: S, R: Rect, fs: number, narrow: boolean) {  const u = api.u; api.rect(R.x, R.y, R.w, R.h, { fill: "paper", w: u(1) });  api.text(narrow ? "層間変形角の最大" : "層間変形角の最大(階の高さごと)", R.x + u(6), R.y + fs * 1.05, { size: fs * 0.88 });  const x0 = R.x + u(narrow ? 15 : 22), w = R.w - u(narrow ? 24 : 34), y0 = R.y + fs * 1.5, h = R.h - fs * 1.5 - fs * 1.35;  let m = 0; for (const b of s.bs) m = Math.max(m, b.maxD);  const top = [1 / 400, 1 / 200, 1 / 100, 1 / 50, 1 / 25, 1 / 10, 1 / 5, 1 / 2, 1].find((q) => q >= m * 1.02) ?? 1;  const X = (q: number) => x0 + (Math.min(q, top) / top) * w;  api.line(x0, y0 + h, x0 + w, y0 + h, { w: u(0.8) }); api.line(x0, y0, x0, y0 + h, { w: u(0.8) });  api.text("上", x0 - u(3), y0 + fs * 0.6, { size: fs * 0.76, align: "right" }); api.text("下", x0 - u(3), y0 + h, { size: fs * 0.76, align: "right" });  for (const q of [top, top / 2]) { api.line(X(q), y0, X(q), y0 + h, { w: u(0.6), dash: [u(2), u(4)] }); api.text(`1/${Math.round(1 / q)}`, q === top ? x0 + w : X(q), y0 + h + fs * 1.1, { size: fs * 0.76, align: q === top ? "right" : "center", font: "mono" }); }  const pts = s.pts;  s.bs.forEach((b, i) => {    let c = 0;    for (let f = 0; f < b.N; f++) { const px = X(b.peak[f]); pts[c++] = px; pts[c++] = y0 + h - (f / b.N) * h; pts[c++] = px; pts[c++] = y0 + h - ((f + 1) / b.N) * h; }    api.poly(pts.subarray(0, c), lineStyle(api, i, s.bs.length));  });}function overlay(api: Api, s: S) {  const u = api.u, v = s.v, K = api.num("scale"), gx = K * s.g.dg, fs = clamp(u(12.5), 11, 15), P = panels(api), pad = u(10);  for (const b of s.bs) {    const xt = b.x0 + gx + K * b.x[b.top], yt = (b.kind === "iso" ? 3.8 : 0.5) + b.N * STORY, a = v.project([xt, yt, 0]);    if (!a.on) continue;    const r = u(5), g = s.grab && s.grab.b === b && s.grab.i === b.top;    api.rect(a.x - r, a.y - r, 2 * r, 2 * r, { fill: g ? "ink" : "paper", w: u(1.5) });    const ly = a.y - (b.kind === "tmd" ? u(62) : u(14)), lx = clamp(a.x, u(58), api.w - u(58));    api.label(`${b.name}  T₁ ${fmt(b.Tf, 2)} 秒`, lx, ly - fs * 1.25, { size: fs, align: "center" });    api.label(`±${lenText(b.amp)}`, lx, ly, { size: fs, align: "center", font: "mono" });  }  if (s.grab) {    const g = s.grab, b = g.b, F = b.pull!.kp * Math.abs(b.pull!.target - b.x[g.i]), floor = g.i - b.f0 + 1, y = (b.kind === "iso" ? 3.8 : 0.5) + floor * STORY, a = v.project([b.x0 + gx + K * b.x[g.i], y, 0]);    if (a.on) { api.circle(a.x, a.y, u(7), { w: u(1.6), fill: "paper" }); api.label(`${floor} 階を引く  ${lenText(Math.abs(b.x[g.i]))}  ${fmt(F / 1e6, 2)} MN`, clamp(a.x + u(14), pad, api.w - u(170)), a.y - u(12), { size: fs }); }  }  // 倍率の見本: 1 本の線の長さが、実際のどれだけのずれに見えるか  const main = s.bs[s.bs.length - 1], yt = main.N * STORY, c0 = v.project([main.x0, yt, 0]), c1 = v.project([main.x0 + K, yt, 0]), per = Math.abs(c1.x - c0.x);  let len = 1; for (const q of [0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1, 2, 5]) if (per * q >= u(20)) { len = q; break; }  api.label(`変形は ${K} 倍で表示${P.narrow ? "" : "(地面のずれも同じ)"}`, pad, pad + fs, { size: fs, weight: 600 });  if (per > 1) {    const bx = pad, by = pad + fs * 2, lp = Math.min(per * len, u(170));    api.line(bx, by, bx + lp, by, { w: u(1.8) }); api.line(bx, by - u(3), bx, by + u(3), { w: u(1.2) }); api.line(bx + lp, by - u(3), bx + lp, by + u(3), { w: u(1.2) });    api.label(`= 本当は ${roundText(len)}`, bx + lp + u(5), by + fs * 0.34, { size: fs * 0.88 });  }  drawHist(api, s, P.hist, fs, P.narrow); drawRes(api, s, P.res, fs, P.narrow); drawDrift(api, s, P.dr, fs, P.narrow);}function readout(api: Api, s: S) {  const g = s.g, key = api.str("wave"), bs = s.bs, b = bs[bs.length - 1], a = bs[0], V = api.num("amp");  const tt = (T: number) => `${fmt(T, 2)} 秒`;  let wave = "地面の揺れ  なし";  if (key === "sine" || key === "sweep") {    const T = key === "sine" ? api.num("period") : sweepT(g.t);    wave = `地面の揺れ  ${key === "sine" ? "正弦波" : "周期を上げていく"}  周期 ${tt(T)}    ±${lenText((V / 100 * T) / (2 * Math.PI))}    最大加速度 ${fmt((V / 100 * 2 * Math.PI) / T, 2)} m/s²`;    if (key === "sweep" && g.ended) wave = `地面の揺れ  止まった  あと ${fmt(g.wait, 0)} 秒でもう一度`;  } else if (key === "short" || key === "long") {    const W = makeWave(api, key), name = key === "short" ? "短周期の多い合成波" : "長周期の多い合成波";    wave = g.ended ? `地面の揺れ  止まった  あと ${fmt(g.wait, 0)} 秒でもう一度` : `地面の揺れ  ${name}  残り ${fmt(W.dur - g.t, 0)} 秒    最大加速度 ${fmt((W.pa * V) / 100, 2)} m/s²`;  }  const lines: string[] = [];  if (s.layout === "iso") {    lines.push(`固有周期  そのまま ${tt(a.T[0])}    免震 ${tt(b.Tf)}`, wave,      `最上階の加速度の最大  そのまま ${fmt(a.peakAbs, 2)} m/s²    免震 ${fmt(b.peakAbs, 2)} m/s²`,      `層間変形角の最大  そのまま ${angleText(a.maxD)}    免震 ${angleText(b.maxD)}`,      `免震層の変形の最大  ${lenText(b.peakDev)}`);  } else if (s.layout === "tmd") {    lines.push(`固有周期 T₁  ${tt(b.T[0])}    重りの周期  ${tt(b.Td)}`, wave,      `最上階の揺れ  そのまま ±${lenText(a.amp)}(最大 ${lenText(a.peakTop)})    制振 ±${lenText(b.amp)}(最大 ${lenText(b.peakTop)})`,      `重りの振れ(最上階から)  ±${lenText(b.ampDev)}`);  } else if (s.layout === "soft") {    const res = (q: Bld) => lenText(Math.abs(q.dl0[0]));    lines.push(`固有周期  そのまま ${tt(a.T[0])}    弱い 1 階 ${tt(b.T[0])}`, wave,      `層間変形角の最大  そのまま ${angleText(a.maxD)}(${a.maxAt} 階)    弱い 1 階 ${angleText(b.maxD)}(${b.maxAt} 階)`,      `今の 1 階のずれ  そのまま ${res(a)}    弱い 1 階 ${res(b)}`);  } else {    const Ts = b.T.slice(0, 3).map((T, i) => `T${"₁₂₃"[i]}  ${tt(T)}`).join("    ");    lines.push(`固有周期(${b.N} 階)  ${Ts}`, wave,      `最上階の揺れ  ±${lenText(b.amp)}(最大 ${lenText(b.peakTop)})    測った周期  ${b.Tm > 0 && b.amp > 1e-3 ? tt(b.Tm) : "-"}`,      `層間変形角の最大  ${bs.length > 1 ? bs.map((q) => `${q.name} ${angleText(q.maxD)}`).join("    ") : `${angleText(b.maxD)}(${b.maxAt} 階)`}`);  }  api.readout(lines);}

/** 指の位置を追いかける地面(なめらかにつなぐ 2 次の追従)。その加速度を建物に与え、速度ベルレ法で 1 刻み進める */function advance(api: Api, s: S, h: number) {  const W = 26, M = s.M, n = M.N, hd = api.num("damp") / 100, w1 = (2 * Math.PI) / M.T[0], cc = (2 * hd) / w1;  const ag = W * W * (s.want - s.g) - 2 * W * s.gv;  s.gv += ag * h; s.g += s.gv * h; s.t += h;  for (let i = 0; i < n; i++) { s.v[i] += 0.5 * h * s.a[i]; s.x[i] += h * s.v[i]; }  for (let i = 0; i < n; i++) {    const dl = s.x[i] - (i ? s.x[i - 1] : 0), du = i + 1 < n ? s.x[i + 1] - s.x[i] : 0, vl = s.v[i] - (i ? s.v[i - 1] : 0), vu = i + 1 < n ? s.v[i + 1] - s.v[i] : 0;    const ku = i + 1 < n ? M.k[i + 1] : 0;    // 層の力と、固さに比例する減衰。地面の加速度は全部の階に同じ力として働く    s.a[i] = (-M.k[i] * (dl + cc * vl) + ku * (du + cc * vu)) / M.m - ag;  }  for (let i = 0; i < n; i++) s.v[i] += 0.5 * h * s.a[i];  const xt = s.x[n - 1], f = Math.exp(-h / (1.2 * M.T[0]));  s.amp = Math.max(Math.abs(xt), s.amp * f); s.gAmp = Math.max(Math.abs(s.g), s.gAmp * f);  // 指の周期: 地面の速さが上向きに 0 を横切る間隔  if (s.last < 0 && s.gv >= 0 && s.gAmp > 0.01) { if (s.tUp >= 0) s.Tf = s.t - s.tUp; s.tUp = s.t; }  s.last = s.gv;}