CFD Solver to Deterministic LES Accelerator

This repository is the software reference model and future hardware driver for a deterministic LES accelerator.

The current Python codebase is a NumPy-based structured-grid compressible solver that serves as the host-side golden model. The long-term hardware target is not a general-purpose CPU. It is an application-specific, phase-locked streaming tile engine with a small RISC-V control core, wide FP64 datapaths, and large banked on-chip SRAM.

What Is Happening

You are not building a full CPU first. You are building a small, deterministic accelerator slice and proving that it can reproduce a known NumPy stencil result cycle by cycle.

The workflow is:

1. NumPy creates a tiny structured field.
2. The host flattens each cell neighborhood into a packet.
3. Cocotb sends those packets into the RTL one cycle at a time.
4. The RTL returns a fixed-latency result.
5. Cocotb compares the RTL output to the NumPy golden model.
6. Once that is stable, the same style of scheduling grows into RK3 supersteps, halo exchange, and BSP orchestration.

The point of this first stage is not speed alone. The point is to prove deterministic timing, correct data movement, and a host-managed cycle budget before the design gets larger.

Core Goal

Prove that a structured-grid LES solver can run as a fully deterministic, cycle-counted pipeline with no DRAM.

If the problem is larger than on-chip memory, scale by chip count and domain decomposition over an interconnect. For the initial prototype, keep everything on one device, in Verilog simulation, with all working data resident in SRAM.

Numerical Scope

• Structured mesh LES kernels
• Explicit RK3 time integration
• FP64 everywhere
• Flat, contiguous NumPy arrays on the host side
• Static scheduling of stencil movement and compute

Hardware Direction

Prototype architecture

• One small RISC-V control core for sequencing only
• One deterministic FP64 vector/stencil engine
• Multi-banked on-chip SRAM scratchpad
• No cache hierarchy
• No DRAM controller
• No speculative execution
• No dynamic task scheduling

Why this shape

• The solver already has a fixed stencil structure.
• RK3 naturally divides into repeating supersteps.
• The host can compile the mesh into a deterministic cycle plan.
• A fixed-latency pipeline is easier to verify cycle-by-cycle in simulation.

BSP Orchestration Model

The preferred programming model is BSP-style supersteps, but with compile-time scheduling instead of runtime synchronization.

Each RK3 stage is treated as a superstep:

1. Local compute on tile data in SRAM
2. Halo or boundary exchange within SRAM banks
3. Barrier point reached by cycle count, not handshake logic

For the first prototype there is no multi-chip communication. The barrier is a deterministic phase boundary, controlled by the Python host and a fixed cycle budget.

Phase-Locked Streaming Tile Engine

The prototype should behave like a streaming accelerator for a tile of the structured mesh.

Key properties:

• Input neighborhoods are flattened into linear streams
• The datapath consumes a known stencil pattern
• Arithmetic latency is fixed
• Output timing is deterministic
• Every RK3 phase has a known cycle envelope

This is the software-simulation equivalent of a hardware compiler targeting a fixed pipeline.

Simulation Stack

• SystemVerilog for the RTL compute blocks
• Verilator for cycle-accurate simulation
• Cocotb for Python test orchestration
• NumPy for the golden model and host-side scheduling

Why this stack:

• It preserves clock-level behavior.
• It integrates naturally with Python reference code.
• It makes exact cycle counting straightforward.
• It is a realistic bridge from algorithm to hardware.

Verification Strategy

1. Generate a small structured mesh in NumPy.
2. Flatten each stencil neighborhood into a deterministic stream.
3. Drive the RTL cycle by cycle.
4. Capture outputs and exact elapsed cycles.
5. Compare against the NumPy golden model.
6. Confirm throughput after pipeline fill and drain.

The primary proof point is throughput under fixed latency, not just functional correctness.

First RTL Contract

The first prototype uses a single-scalar stencil packet per cell:

• center
• east
• west
• north
• south
• top
• bottom

The first golden operation is a centered six-neighbor derivative on a scalar field:

result = (east - west) + (north - south) + (top - bottom)

The purpose of this kernel is not to model the full LES update yet. It is to prove the deterministic streaming path, pipeline fill/drain timing, and NumPy-to-RTL cycle accounting on a structured-grid stencil.

Formal Packet Contract

1. Field order

Each packet carries one scalar cell value and its six direct neighbors in this exact order:

1. center
2. east
3. west
4. north
5. south
6. top
7. bottom

2. Bit layout

For the first prototype, every field is a 64-bit IEEE-754 floating-point value.

• Total packet width: 7 x 64 = 448 bits
• Encoding: raw fp64 bit pattern, no compression, no packing tricks
• Host side: NumPy float64
• RTL side: logic [63:0] per lane

3. Handshake timing

The contract uses a standard streaming handshake:

• in_valid means the packet on the input lanes is meaningful.
• in_ready means the RTL can accept a new packet on that cycle.
• out_valid means a valid result is present at the output.
• out_ready means the consumer can accept the result.

For the prototype, the target behavior is:

• one packet accepted per cycle when in_valid && in_ready
• fixed pipeline latency after input acceptance
• deterministic output order matching input order
• no dynamic stalls in the ideal case

4. Expected latency

The initial RTL stub uses a fixed pipeline depth parameter:

• PIPE_DEPTH = 4

So the first output should appear after the pipeline fill period, then continue one per cycle if the input stream stays valid and the output remains ready.

In the testbench, total cycles are measured explicitly so we can verify:

$$ ext{Total Cycles} = \text{number of packets} + \text{pipeline depth} - 1 $$

5. Mapping to the NumPy reference

The NumPy reference field is a 3D scalar array. For each cell, the host builds a packet from the local neighborhood:

• center = field[i, j, k]
• east = field[i+1, j, k] with boundary clamping at the edge
• west = field[i-1, j, k] with boundary clamping
• north = field[i, j+1, k] with boundary clamping
• south = field[i, j-1, k] with boundary clamping
• top = field[i, j, k+1] with boundary clamping
• bottom = field[i, j, k-1] with boundary clamping

The first golden operation is:

result = (east - west) + (north - south) + (top - bottom)

That is a simple centered stencil sanity check. It is not the final LES kernel; it is the first proof that the hardware path and the NumPy path agree.

Project Phases

Phase 1

• Define the stream format for stencil neighborhoods
• Stub the SystemVerilog datapath
• Stub the Cocotb testbench
• Establish the NumPy golden reference
• Add a simple superstep controller that can mark compute, exchange, and barrier phases by cycle count

Phase 2

• Implement fixed-latency FP64 operators
• Add tile buffering in banked SRAM
• Model RK3 superstep orchestration

Phase 3

• Add AMReX-style domain decomposition later
• Add multi-chip interconnect later
• Scale beyond single-device SRAM capacity later

Design Rules

• Keep the prototype deterministic.
• Keep memory on chip.
• Keep precision at FP64 unless explicitly justified otherwise.
• Keep the host responsible for orchestration.
• Keep hardware simple, pipelined, and statically scheduled.

Current Repository State

The repository currently contains the Python solver, geometry, boundary, flux, reconstruction, SGS, and time-stepping code needed to define the golden model.

The hardware and Cocotb layers are not yet present and should be added as a separate simulation-oriented subtree.