denario-ns2d-levy-v2

Scientist: denario-6 Date: 2026-05-18

Task

Discover the effective physical theory of Lévy flights in 2D turbulence. Specifically: what physical phenomenon in the Navier-Stokes flow produces anomalous diffusion described by the half-Laplacian (fractional diffusion operator)? The goal is mechanistic: identify which physical structures or processes in the velocity field produce power-law tails in displacement PDFs (Lévy flights) and the fractional diffusion equation ∂_t P = -D_α (-Δ)^{α/2} P.

Scientific Context

This is a continuation of prior work that identified the limitations of low-Re 2D turbulence simulations for producing genuine Lévy flights. The critical failure in prior work was:

• T/T_L ≈ 7.9 (below the T/T_L > 10 threshold)
• Vortex condensation (inverse cascade not fully developed)
• Forcing at small k → vortices nearly domain-sized → restorative confinement

The present simulation is specifically designed to overcome these limitations:

• T/T_L ≫ 10 (T=1500 time units, expected T_L ~ 20-40 units → T/T_L ~ 40-75)
• Ekman friction (linear drag α=0.05) prevents condensation into a single vortex
• Forcing at k∈[8,12] creates vortices much smaller than domain → multi-scale hierarchy
• Higher Reynolds number (ν=5×10⁻⁴) for a more developed inverse cascade

Data Generating Process

A 2D incompressible Navier-Stokes simulation on a 512×512 periodic domain:

Equation: ∂_t ω + u·∇ω = ν∇²ω - αω + f(x,t)

where:

• ω = curl(u) is the vorticity (omega_rms ≈ 1.1 at steady state, confirmed non-trivial)
• u = (u,v) is the 2D velocity field derived from ω via the streamfunction ψ: ω = -∇²ψ
• ν = 5×10⁻⁴ (kinematic viscosity)
• α = 0.05 (linear Ekman friction coefficient)
• f(x,t) = stochastic forcing localized in wavenumber band k∈[8,12]

Numerical method:

• Pseudo-spectral solver (FFT-based), 512×512 grid, domain L=2π
• 2/3 dealiasing rule
• Adams-Bashforth 2nd order time integration, dt = 0.01
• Spin-up: 3000 steps (t=0 to t=30) before tracer seeding
• Production: 150000 steps (t=0 to t=1500)

Lagrangian tracers:

• 10,000 passive tracers seeded uniformly at t=0 (after spin-up)
• Advected with RK2 using bilinear interpolation of velocity field
• BOTH wrapped (periodic) and unwrapped positions saved
• Saved every 5 steps (dt_save = 0.05)

Confirmed statistics:

• omega_rms ≈ 1.12 (steady-state, genuine turbulence)
• MSD at tau=1: 0.017, at tau=100: 7.3, at tau=1500: 111 (tracers diffusing ~10 domain sizes)

Files (absolute paths)

All files in /home/node/work/projects/ns2d_levy_v2/data/

• /home/node/work/projects/ns2d_levy_v2/data/traj_x.npy Shape: (30001, 10000), dtype: float32 Tracer x-positions (wrapped to [0, 2π]) at each saved time. Row 0 = t=0.

• /home/node/work/projects/ns2d_levy_v2/data/traj_y.npy Shape: (30001, 10000), dtype: float32 Tracer y-positions (wrapped to [0, 2π]) at each saved time.

• /home/node/work/projects/ns2d_levy_v2/data/traj_x_unwrapped.npy Shape: (30001, 10000), dtype: float64 Tracer x-positions WITHOUT periodic wrapping — essential for displacement statistics. Δx(τ) = traj_x_unwrapped[t+τ] - traj_x_unwrapped[t]

• /home/node/work/projects/ns2d_levy_v2/data/traj_y_unwrapped.npy Shape: (30001, 10000), dtype: float64 Tracer y-positions WITHOUT periodic wrapping.

• /home/node/work/projects/ns2d_levy_v2/data/t_traj.npy Shape: (30001,), dtype: float64 Time values for each saved snapshot. Range: [0, 1500], spacing = 0.05.

• /home/node/work/projects/ns2d_levy_v2/data/vorticity_snapshots.npy Shape: (101, 512, 512), dtype: float32 Eulerian vorticity field ω(x,y,t) at 101 times (Δt=15 apart).

• /home/node/work/projects/ns2d_levy_v2/data/snap_times.npy Shape: (101,), dtype: float64 Times of vorticity snapshots.

• /home/node/work/projects/ns2d_levy_v2/data/energy_spectrum.npy Shape: (256, 2), dtype: float64 Columns: [k, E(k)]. Energy spectrum computed from the last vorticity snapshot.

• /home/node/work/projects/ns2d_levy_v2/data/meta.json Simulation metadata (all parameters, actual shapes, wall clock time).

Physical Scales and Expected Regime

• Forcing scale: L_f = 2π/k_force ≈ 2π/10 ≈ 0.63 (wavenumber 10)
• Domain scale: L = 2π ≈ 6.28 → scale separation ≈ 10×
• Ekman drag scale: L_α = sqrt(ν/α) ≈ 0.1 (smaller than forcing scale)
• Expected Lagrangian integral time: T_L ~ L_f / u_rms (estimate 20-40 time units)
• T/T_L: 1500 / (20-40) = 40-75 → well above the T/T_L > 10 requirement
• Expected energy spectrum: E(k) ~ k^{-5/3} at large k (enstrophy cascade), E(k) ~ k^{-3} or steeper toward small k (inverse cascade)

Suggested Analyses (ordered by scientific priority)

Priority 1: Establish whether Lévy flights actually occur

1. Displacement PDF: Compute P(Δr, τ) for τ = 1, 5, 10, 50, 100, 500 time units.

   • Fit with Lévy stable distribution (characteristic function: φ(k) = exp(-σ|k|^α))
   • Estimate stability index α(τ). If α < 2 and stable across scales, Lévy flights confirmed.
   • Also compute empirical tail exponent via Hill estimator on |Δr| top 5%.

2. MSD and Hurst exponent: Compute ⟨|Δr|²⟩(τ) via time averaging.

   • Running H(τ) = d log⟨|Δr|²⟩ / (2 d log τ) — check if it converges to a stable value.
   • Use only lags τ < T/10 to avoid boundary effects.

3. Verify T/T_L: Compute Lagrangian velocity autocorrelation R_v(τ).

   • T_L = ∫₀^∞ R_v(τ)/R_v(0) dτ
   • Confirm T_simulation / T_L > 10.

Priority 2: Identify the physical mechanism

4. Trapping time distribution:

   • Define trapping: |Δr(t, τ_trap)| < r_trap for connected interval.
   • Measure P(τ_trap). Power law ~ τ^{-μ} with μ < 2 → CTRW (vortex trapping).

5. Jump size distribution:

   • Measure P(ℓ) for displacements between trapping events.
   • Power law P(ℓ) ~ ℓ^{-(1+α)} → Lévy stable.

6. Vortex identification and conditional statistics:

   • Use Okubo-Weiss criterion (Q = s² - ω²/4 < 0) or |ω| > threshold for vortex cores.
   • Tag tracers as "in-vortex" vs "out-of-vortex".
   • Compare H_in vs H_out, trapping time distributions.

7. Eulerian velocity increment statistics:

   • Structure function exponents S_p(r) ~ r^{ζ_p}.
   • Energy spectrum E(k) ~ k^{-β}: extract β.
   • Theoretical prediction: α_Levy = 2/(β-1) (from Kraichnan-Gawędzki theory).

8. Lagrangian velocity increment PDF:

   • δv_L(τ) = v(t+τ) - v(t) for tracers.
   • If these are Lévy-stable distributed, the velocity field drives the mechanism.

Priority 3: Establish the effective theory

9. Fractional diffusion fit:

   • If Lévy flights confirmed: fit ∂_t P = -D_α (-Δ)^{α/2} P to displacement PDF.
   • Verify α is consistent with velocity increment exponent and energy spectrum.

10. CTRW vs Lévy process discrimination:

    • CTRW: P(Δr, τ) has characteristic tent shape, depends on waiting time distribution.
    • Lévy process: P scales as τ^{1/α} P(Δr/τ^{1/α}) (self-similar).
    • Test: does the displacement PDF show self-similar scaling with τ^{1/α}?

Key Scientific Question

After the analysis, the report must answer: what is the effective physical theory of Lévy flights in 2D turbulence?

Candidate answers:

1. CTRW with power-law waiting times: Vortex trapping + ejection. Mechanism: coherent vortex cores.
2. Lévy-stable velocity increments: Power-law energy spectrum → heavy-tailed velocity increments. Mechanism: multi-scale spectral roughness.
3. Fractional Brownian Motion: Long-range VACF memory → anomalous diffusion without fat tails. Mechanism: persistent velocity correlations.
4. Lévy walk: Alternating slow/fast segments with power-law lengths. Mechanism: intermittent strain field.

The effective physical theory is the stochastic model that correctly reproduces P(Δr, τ), ⟨Δr²⟩(τ), and relates the Eulerian velocity spectrum exponent to α.

CRITICAL MEMORY NOTE: The unwrapped position arrays are 2.3 GB each. Do NOT load the full arrays at once. Use memory-mapped loading: np.load(path, mmap_mode='r') and slice only the time indices you need. For MSD computation, subsample tracers (e.g., every 10th tracer) or use a sliding window approach on the time axis. The vorticity snapshots (102 MB) can be loaded fully.