Simulation methods and reproducibility

New audit implementation, 6 September 2026. The numerical code is in preview/lib/dust-physics.ts; the interactive controls are in preview/components/ResearchLab.tsx; invariant checks are in preview/scripts/test-dust-physics.mjs. Existing simulations and original essays are preserved.

Scope and units

These are deterministic, optically thin, monochromatic dust demonstrations. The amplitudes are polarized specific intensities at 353 GHz, in arbitrary common units. They are not fractions of extinction or total intensity. A layer uses

The code uses h/k = 0.04799243073366221 K/GHz and expm1 to evaluate the ratios. Polarization angles are emission angles in degrees modulo 180°. No actual field map, grain-size distribution, extinction inversion or distance reconstruction is inferred. Layer temperatures and indices are selected inputs.

The lab does not include thermodynamic CMB conversion, bandpass/color corrections, noise, telescope beams, masks, calibration, CMB, synchrotron, lensing or an r likelihood. Common frequency-dependent unit conversion would leave angle diagnostics invariant, but changing the unit of only one amplitude would not. All comparisons here are internally in specific-intensity units.

Sightline experiment

The reference layer has (A,T,β,ψ)=(1,20 K,1.6,0°). The default second layer has (0.7,16 K,1.9,65°). Frequency, second-layer amplitude, temperature, β and angular separation can be varied. Presets expose equal SEDs, aligned layers and exact cancellation.

The chart samples 90–405 GHz at 2.1 GHz spacing. It plots the signed wrapped angle difference from 353 GHz, using a half-angle atan2 of the Stokes cross and dot products. The chart breaks connections at changes larger than 45° instead of implying a continuous angle through an unsampled zero or branch jump. This is a plotting safeguard, not an adaptive physical root finder. The selected-frequency diagnostics evaluate the exact formula at that frequency.

Reported coherence is $C=|P|/\sum_i|P_i|$, not polarization fraction $p=|P|/I$. The angle derivative and its upper bound are derived in the exact diagnostics note. The code returns an undefined slope near numerical cancellation (C < 10⁻¹⁰); reference-angle and scalar-fit checks also suppress divisions by numerically vanishing vectors. These thresholds protect floating-point displays and are not noise-based data cuts.

The local residual minimizes $|P_\nu-aP_{353}|^2$ over an unconstrained real scalar $a$ and divides by $|P_\nu|^2$. Negative coefficients are allowed. Therefore an antiparallel vector has zero local fit residual despite a 90° angle flip. This is a geometry diagnostic, not a fit to a physical positive MBB. It is neither a BB-power residual nor a tensor-to-scalar bias.

Synthetic sky experiment

The domain is a periodic 25 × 25 grid with no assigned angular size. It is deliberately odd-sized to avoid ambiguous Nyquist-frequency orientations. At each integer (x,y), define a=2πx/25, b=2πy/25 and

The layer-1 amplitude is $1+0.35\cos a+0.2\sin b$, with $T_1=20\,K$, $\beta_1=1.6+\delta\beta$ and angle $\psi$. The layer-2 amplitude is $f[1+0.3\cos(a-b)]$, with adjustable $T_2$, $\beta_2=1.6+\Delta\beta+\delta\beta$ and angle $\psi+\Delta\psi$. Both amplitudes are nonnegative over the permitted controls. The spatial-variation slider s is the maximum amplitude of this deterministic β pattern; it is not a Gaussian standard deviation.

Generate Q/U at 150 and 353 GHz. Scale the latter by a single fixed $s_{150}(20 K,1.6)$ and subtract it from the former. For each of target, scaled template and residual, a separable direct discrete Fourier transform uses 1/N² forward normalization and no inverse normalization. Apply

omit the zero mode, and inverse-transform to the displayed real B map. No random-number seed is required: the sky is given entirely by these equations and slider values. The blue/amber maps share a symmetric color scale including the residual, so stronger residuals cannot silently rescale themselves to appear small.

The displayed B correlation uses the summed cross-power of all nonzero modes divided by the geometric mean of the two auto-powers. It is not a multipole-bin Planck decorrelation estimator. The fixed-template residual is total residual B power divided by true 150 GHz dust B power and can exceed 100%. The oracle best global scalar residual is $1-R^2$, using the same modes and unconstrained scalar. The oracle has access to the true synthetic dust B field; it is a diagnostic lower bound within this one-scalar family, not an achievable component-separation method on noisy observations.

The exact-template null gives correlation 1 and numerical-zero residual for a common SED. The spatial SED only preset aligns the clouds and gives them identical local SEDs, while changing β across the grid. Pixel angles stay fixed with frequency, but the B maps can decorrelate. This demonstrates why the repository's local-vector and angular-power statements cannot be interchanged.

Separation experiment

The third panel plots the exact condition number of the two-column known-angle Stokes design matrix. It is singular at 0° and 90° and optimally conditioned at 45°. It assumes dust-only data, exactly known angles, two free amplitudes and equal uncorrelated Q/U noise. The κ > 10 display label is a pedagogical heuristic. It is not an observational requirement or evidence that the amplitudes can be recovered at useful signal-to-noise.

Verification and export

From preview/, using Node 22.6+ with type stripping support (verified here with Node 24):

node --experimental-strip-types scripts/test-dust-physics.mjs
npm run build

Checks cover reference normalization, the MBB slope by finite difference, equal-SED angle invariance, analytic rotation derivative, zero-vector handling, both inversion singularities, the multilayer coherence bound, DFT round-trip and Parseval normalization, pure E/B Fourier modes, exact-template cleanup, and the spatial-SED counterexample. An independent review also checked general-direction pure E/B waves and compared the oracle residual with direct real-map least squares, with discrepancies below 1.5 × 10⁻¹⁴ in the tested cases. These validate the mathematics implemented here; they do not validate an astrophysical foreground population.

Both data experiments export CSV containing numerical outputs and #-prefixed metadata lines with their selected parameters and assumptions. Readers can skip comment lines when loading the table. Undefined ratios export empty cells. The exports do not contain observed sky data.

Next research steps are specified in the predictive-validation protocol: real bandpasses and covariance, uncertain ancillary information, matched empirical/sham priors, held-out regions and frequencies, then an instrument-specific cosmological likelihood.