CelerisLab/docs/validation/duetsch_re100_kc5_protocol.md

Prescribed in-line oscillation: Dütsch Re 100 / KC 5 — PROTOCOL (review draft)

Production under test: claude/moving-bodies-s3a at 653fd997 (unchanged). Primary source: Dütsch, Durst, Becker & Lienhart, J. Fluid Mech. 360, 249–271 (1998), doi:10.1017/S002211209800860X. Source facts below were verified from the paper (Slack C0C78UFKGRM/1793335557.460509, TRON 1793356317–1793356724, page numbers given). Everything marked ENGINEERING is our a-priori choice, fixed before any benchmark result; it is not a paper uncertainty.

0. Source facts (with provenance)

• Motion (p.251): x_c = -A sin(2 pi f t); paper phase phi = 0 is the mean position moving left at -U_max. KC = 5 gives A/D = 5/(2 pi) = .7957747.
• Start (p.254): the computation starts a quarter cycle before phi = 0, at x = +A with zero velocity (phi = -90°).
• Experiment (pp.252, 254–256): D = 10 mm, L = .65 m, tank 1.2 × .7 × .7 m, rho = 998.2 kg/m³, mu = 1.00818 mPa·s, Re = 100, KC = 5; LDA with 1° phase bins and 200 samples per node and degree; no quantified experimental error recovered.
• Profiles (Fig. 8, p.260): LAB-FIXED stations x1 = -6, 0, 6, 12 mm, i.e. x/D = -.6, 0, .6, 1.2, at phi = 180°, 210°, 330°; u1 and u2 versus x2.
• Original numerics (§2): cylinder-fixed accelerated frame (x_i = x~_i + x_i,s, outer velocity -u_c, inertial source -rho du_c/dt), O-grid outer diameter 120 D, illustrated 192 circumferential nodes, 720 steps/cycle; tank influence neglected.
• Force fit (verified): force ON THE BODY per unit span F1 = -0.5 rho D c_d xdot|xdot| - (pi/4) rho D² c_i xddot with x = -A sin phi; normalized by 0.5 rho D U_max²: c_d cos phi |cos phi| - (pi² c_i / KC) sin phi. c_i is the ADDED-MASS coefficient; total inertia c_m = c_i + 1 is a different convention and is not used.
• Table 1 computed Morison coefficients (c_d, c_i) for 96 / 192 / 384 perimeter nodes at dt/T = 1/720: (2.15, 1.42), (2.10, 1.45), (2.09, 1.45). These are fitted computed coefficients, not a measured force and not an error bar. The paper states the coarse direct force is closer to the finest direct force than the Morison fit, but gives no waveform-error norm.

1. Kinematics, time levels and phase bookkeeping

• Our run is in the LAB frame: the circle moves, the fluid is initially at rest (rho = 1, u = 0); no fictitious inertial source is added. Lab-fixed profile stations are compared directly (no frame transform needed).
• Geometry reference: the circle's base centre is the INITIAL position, the right turning point x_mean + A; the declared offset box is x in [-2A, 0], y = 0. The body is never teleported: offset 0 at s = 0.
• Elapsed time s = n dt from the start. Analytic position relative to the mean x(s) = A cos(2 pi s / T), velocity u(s) = -U_max sin(2 pi s / T); paper phase phi = (360 s / T - 90) mod 360.
• Step n moves from s_n to s_{n+1}: target offset x(s_{n+1}) - A; wall velocity = the single-step interval average (x(s_{n+1}) - x(s_n)) / dt (analytic, supplied, not differenced from stored offsets). The per-step force of run(1) / read_force(normalize=True) is attributed to that same interval.
• Paper cycle k = 1, 2, ... spans s in [T/4 + (k-1) T, T/4 + k T); the initial quarter cycle (s < T/4) is cycle 0 and never scored. Every logged sample carries (k, phi). Example: phi = 330° of cycle 1 is at s = 7T/6.
• T is chosen as a multiple of 360 native steps, so phi = 180°, 210°, 330° fall exactly on step boundaries and 1° bins hold whole steps.

2. Resolution ladder (uniform grid, fixed Mach), ENGINEERING

D (cells)	T (steps)	U_max = 5 D / T	nu_lat = U_max D / 100	tau (3 nu + .5)
20	2,520	.039683	.0079365	.5238
40	5,040	.039683	.015873	.5476
80	10,080	.039683	.031746	.5952

• MRT collision. Mach fixed (Ma ≈ .0687); a separate Mach check at D = 40 with T = 10,080 (U_max halved, tau .5238) is labelled MACH and not part of the ladder.
• Per-step displacement ≤ U_max ≈ .04 cells (bound .5).

3. Domain and boundaries (finite-domain reproduction), ENGINEERING

The paper's 120 D O-grid is not reproduced; our rectangular box is a reproduction choice. The native dense channel (zero-velocity inlet, pressure-NEQ outlet, bounce-back y walls) is a finite-domain approximation, not infinite still fluid. The G4 x-periodic / y-slip override is a separately labelled sensitivity family and never mixed into the primary ladder. The runner owner confirms the zero-velocity inlet is supported as configured before any run.

Order (closed):

1. Smoke + ledger at D = 40, native BC, 40 D × 40 D (finite box).
2. Domain sensitivity before any literature statement: same BC, D = 40, L ∈ {20 D, 40 D, 60 D}; record the change in Fx waveform, c_d, c_i and the scored profiles.
3. Optional labelled BC sensitivity: one x-periodic / y-slip case at D = 40, 40 D; reported as a delta against native only.
4. Only after step 2 passes: the D = 20 / 40 / 80 ladder on the chosen BC and domain.

Domain criteria (a-priori, ENGINEERING; separate from the Table 1 bands): the reference is 60 D at D = 40; 20 D and 40 D are compared with it on the common window (same cycle indices) and the common valid point set:

• RMS(Fx_L - Fx_60D) / F_ref ≤ .01 over the window (per-step waveforms);
• |c_d,L - c_d,60D| ≤ .02 and |c_i,L - c_i,60D| ≤ .02 (window fits);
• max(|u_L - u_60D|, |v_L - v_60D|) over common points / U_max ≤ .02. The chosen domain is the smallest of {20 D, 40 D} passing all three; if neither passes, the domain verdict is FAIL and the ladder is not run.
• Memory/runtime: measured by a smoke run before each new size; measurements inform scheduling only, never tolerances.

4. Observables (raw, per step, preserved)

• obs Fx, Fy (CME) and all MovingTables.DIAG columns per step: CME, CME + Wen and CME + cover − uncover reported separately; crossings, eq_only, no_donor (no_donor must stay 0).
• Velocity u1, u2 on lab-fixed lines x/D = -.6, 0, .6, 1.2 (relative to the mean cylinder position) versus x2, at the scored phases. Lab frame throughout; if any body-frame intermediate is used, position and velocity are both transformed.
• Sampling cost: raw force and diagnostics counters are read every step (small). Flags and velocity fields are read only at the steps inside the three prescribed 1° phase bins (φ ∈ [180, 181), [210, 211), [330, 331), T/360 steps each) of the cycles needed for settling and scoring, plus sparse declared audits; never every step.
• Point set: fixed x2/D sample points per station, identical on every grid and declared before running; per point the interpolation weights are recorded; the scored set is the common subset valid on all compared grids at that phase, with its coverage recorded; point sets never change silently.
• Interpolation: bilinear on cell centres using only cells that are fluid at that step (flags read back via get_flags()); weights renormalized over the fluid corners; a sample with no fluid corner, or inside the circle at that step, is MISSING and logged, never set to zero. One-degree bins keep valid fluid weights. The numerical average need not mimic the 200 LDA samples.

5. Periodic-regime gate (every case, before any statistic), ENGINEERING

Ledger rows, same schema for every case: dx, dt (in D and T units); U_max, nu, Re, KC (lattice and nondimensional); T (native steps, nondimensional); duration (native steps, s/T, paper cycles); discard; accepted window; cycle ledger; verdict PASS or FAIL-CLOSED.

Scales: F_ref = 0.5 rho U_max² D, U_ref = U_max; never a near-zero mean force as denominator.

• Settling metrics, cycle k−1 → k, phase-aligned by step within the cycle (error budget: one fifth of the tightest comparison band of the same observable, so residual drift cannot obscure a domain or resolution pass):
  • waveform RMS(dFx) / F_ref ≤ .002;
  • |d mean|, |d RMS|, |d peak| of Fx, each / F_ref ≤ .002;
  • scored profiles (fluid-only) |du| / U_max, |dv| / U_max ≤ .004 at the three phases × four stations;
  • Morison |dc_d| ≤ .004, |dc_i| ≤ .004.
• Conservative comparison: every domain and cross-resolution criterion is applied to |difference| + drift_a + drift_b, where drift is the largest pair metric of that case inside the scored window.
• These constants (protocol.py SETTLE, DOMAIN, CROSS) are frozen before any physical result; no post-result loosening.
• Rule: all metrics for N = 3 consecutive whole paper cycles → settled.

Cycle indexing (worked example, binding). Paper cycle k spans s/T in [0.25 + (k-1), 0.25 + k), i.e. φ0 → φ0; cycle 0 is the startup quarter s/T < 0.25. Pair p_k compares cycle k−1 with cycle k.

• Discarded: cycles 0, 1, 2 (s/T < 2.25). Pairs are formed only between cycles ≥ 3, so the first pair is p_4 (3 → 4).
• Settling proof for a case: three consecutive pairs all pass; K_case = the first cycle AFTER the last of those pairs. Earliest possible: pairs 3→4, 4→5, 5→6 pass ⇒ K_case = 7.
• Common window: K = max over compared cases of K_case; every case scores the same three complete cycles K, K+1, K+2 and must stay settled continuously: every pair p_j for j = its own K_case .. K+2 must pass (including the wait between its local settle and the common K), otherwise that case is FAIL-CLOSED.
• Cap: cycle k ends at s/T = k + 0.25, so the last admissible scored cycle is 11 (ends 11.25 T) and K ≤ 9. Examples: K = 7 → cycles 7–9 (ends 9.25 T); K = 9 → cycles 9–11 (ends 11.25 T); any K_case > 9 → FAIL-CLOSED.
• Run length: cases run sequentially; each runs to the end of cycle 11 (11.25 T: 28,350 / 56,700 / 113,400 steps for D = 20 / 40 / 80; the 12 T cap is 30,240 / 60,480 / 120,960) and keeps its scalar and profile data, then the common K is chosen. An unscorable partial cycle 12 is not run.
• Profiles: an empty common valid point set is FAIL-CLOSED; coverage, missing points, bin time weights (equal per step) and the original fluid corner weight sums are reported.
• The earlier G4 timing check (100 coarse vs 400 diffusive-fine steps) was a matched short transient, not periodic evidence.

6. Metrics and per-observable criteria (ENGINEERING, not paper uncertainty)

Three criteria stay separate; the G4 same-grid 2% target is not reused.

• Morison fit: least squares of the normalized per-step force over the accepted window against the two basis functions, each averaged over its single step interval, fixed paper phase, no free phase shift; c_d and added-mass c_i per force candidate.
• Cross-resolution: every scored observable's error against its reference must decrease monotonically D20 → D40 → D80 (necessary, not sufficient).
• Against Table 1 (computed fit, not measured force): at D ≥ 40, |c_d - 2.09| ≤ .15 and |c_i - 1.45| ≤ .15, and monotone.
• Force waveform against its own Morison reconstruction: residual reported only, no pass/fail (the paper reports fit error above its coarse numerical error and gives no percentage).
• Fig. 8 LDA profiles: QUALITATIVE ONLY. The original figure could not be accessed or rendered; no digitized ordinates and no measurement uncertainty exist. After a digitization with figure number, tool and uncertainty note: RMS(|du|, |dv|) / U_max ≤ .15 at D ≥ 40 and monotone; LDA error stays UNQUANTIFIED and is stated with every comparison.
• Fy sanity (not literature): |mean Fy| / F_ref ≤ .05 in the accepted window.
• No force candidate is selected for an RHS by these results; each is reported.

Superseded (earlier drafts): the velocity-jump start, guessed stations and phases, "C_M" for c_i, and the uncertainty wording; cycles 4–6, the 30-cycle cap and the earlier tolerance numbers.

6b. Amendment: duration-only diagnostic (declared before the run)

• Background: the D40 domain stage (20 D / 40 D / 60 D, 11.25 T each) failed settling in every box (no window). That result stands.
• One declared diagnostic: D = 40, 60 D box, from rest, through cycle 31 (157,500 steps = 31.25 T). Physics, settling tolerances, bands and production unchanged; three proof pairs plus a continuously settled K..K+2 with K ≤ 29.
• Invocation: run.py --D 40 --box 60 --cycles 31 then analyze.py --horizon 31; the JSON records horizon and horizon_mode. Without --horizon the analyzer is the original cap-11 scoring.
• Failure at 31.25 T: report and stop; no automatic extension, matrix, D80 or VIV. A pass shows only that this box settles with time; it is not a domain or literature pass.
• The phase-locked high-frequency force component (residual of the mean waveform after its Morison reconstruction, above 20/T) is reported per case as unqualified; it needs separate qualification before any Morison or literature statement.

6a. Driver unit tests (CPU, before any run)

• Phase/cycle arithmetic: s ↔ (k, φ), cycle boundaries, bin membership, window/cap selection for K = 7, 8, 9, 10.
• Solid-aware interpolation: all-fluid, partial-solid (renormalized weights), no-fluid (MISSING) and inside-circle cases.
• Interval-averaged Morison bases and the two-parameter fit on a synthetic waveform with known c_d, c_i.

7. Order and scope

Protocol review → one-shot Dütsch driver (analytic set_body + run(1), cycle ledger, lab-fixed fluid-only station sampler, two-basis Morison fit; no reusable runner until the first smoke shows shared helpers pay off; no production change) → smoke (memory/runtime) → domain sensitivity → ladder → physical-validation verdict → only then single-cylinder VIV. No AMR, no pinball, no manuscript work.