Boundary Files

Boundary Generation

This page explains how solid-boundary geometry is fed into the solver: instead of an explicit mesh, boundaries are rasterized onto a background grid and shipped as a pair of small binary files. It covers the on-disk format, the bit-packing scheme used to tag cells, and a step-by-step recipe (with the existing nozzle generator as a worked example) for writing a generator for a new shape.

Table of contents

• Overview
• The two files
• Cell encoding
• Grid & normal conventions
• Writing a new boundary generator
• Worked example: the nozzle generator
• Validating your output

Overview

A boundary is defined on a 2D grid of Ni x Nj cells. Every cell is either:

• not a boundary (fluid interior or empty space outside the geometry), or
• a solid wall cell, carrying an index into a companion list of outward-facing unit normal vectors. The solver uses the wall cells + their normals to enforce wall boundary conditions (no-slip / no-penetration) without needing an actual surface mesh — the geometry is effectively "painted" onto the grid, voxel by voxel.

Two files are produced per boundary: a grid file (*.bin) and a normals file (*_n.bin). Both are written by two small shared helpers that should be treated as a fixed contract — any new generator should call them rather than writing its own file I/O:

def _write_grid_bin(path, Ni, Nj, norm, array):
    with open(path, "wb") as f:
        np.array([Ni, Nj, norm], dtype=np.int32).tofile(f)
        array.astype(np.uint32).tofile(f)
 
def _write_normals_bin(path, normal):
    if len(normal) > 0:
        np.array(normal, dtype=np.float32).tofile(path)
    else:
        open(path, "wb").close()

The two files

geometry.bin — grid file

field	type	count
Ni	int32	1
Nj	int32	1
norm	int32	1
array	uint32	Ni * Nj

The array is stored row-major in (i, j) order — i.e. for each i (0..Ni-1), all Nj values of j are written contiguously before moving to the next i. This falls out naturally from array.astype(np.uint32).tofile() on a C-order NumPy array shaped (Ni, Nj), so as long as you build array with that shape you get the right layout for free.

norm is a scale factor: generators multiply their geometry dimensions by it, and it's also written into the header so downstream tooling knows what scale the voxel grid was built at.

geometry_n.bin — normals file

A flat stream of float32, laid out as nx0, ny0, nx1, ny1, ... — one 2D unit vector per wall cell, in ascending order of the index value packed into that cell (see below). If there are zero wall cells, the file is still created, just empty (0 bytes) — readers should tolerate that rather than treating a missing/empty file as an error.

Cell encoding

Each grid cell is a single uint32, built with two fixed helpers:

def pack_uint16_to_uint32(high_16, low_16) -> int:
    return (high_16 << 29) | low_16
 
def unpack_uint32_to_uint16(packed_32):
    high_16 = (packed_32 >> 29) & 1
    low_16 = packed_32 & 0xFFFFFF
    return high_16, low_16

Bit layout of a packed cell:

 31        29 28                                              0
+---+---+---+-------------------------------------------------+
|      |  1 |                    index (up to 2^29 - 1)         |
+---+---+---+-------------------------------------------------+
        ^ bit 29 — "this cell is a solid wall" flag

• 0 — not a boundary cell (fluid or empty space). This is the default; you generally never write 0 explicitly, you just leave those cells untouched.
• pack_uint16_to_uint32(1, index) — a solid wall cell. index is the 0-based position of this cell's outward normal inside the normals file.
• Cell count is capped at 2²⁹ − 1 ≈ 536,870,911 wall cells, because index occupies bits 0–28.

Two quirks worth knowing before you write a generator against this format:

1. pack_uint16_to_uint32 does not mask low_16. If index were ever allowed to reach 2^29, it would spill into bit 29 and corrupt the wall flag itself. In practice this is a non-issue as long as your wall-cell count stays under the ~536M cap, but it's why that cap exists and why you shouldn't repurpose the low bits for anything else.
2. unpack_uint32_to_uint16 only ever recovers bit 29, and its own low_16 mask (0xFFFFFF, 24 bits) is narrower than the 29 bits pack actually writes. So a tag written with anything other than 0 or 1 for high_16 (some generators use 2 as an internal-only marker, e.g. for "fluid interior" cells) will not round-trip through unpack_uint32_to_uint16 — it decodes with high_16 == 0, i.e. indistinguishable from background. The visualization code ((array >> 29) & 1) has the same behavior. Treat the format as binary (wall / not-wall) plus an index; don't rely on extra tag values surviving unless you also extend the reader.

Grid & normal conventions

• i is the axial/x dimension (Ni long), j is the transverse/y dimension (Nj long).
• Normals are 2D, (nx, ny), unit length, pointing outward from the solid material into the fluid. The solver assumes they're already normalized — always divide by the vector's length before storing it.

Writing a new boundary generator

The recipe is the same regardless of shape:

1. Pick Ni, Nj, and norm. Leave enough margin around your geometry — nothing stops you from placing geometry flush against the edges, but most existing generators pad the domain so the shape doesn't touch the texture boundary.
2. Allocate array = np.zeros((Ni, Nj), dtype=np.uint32) and normal = [], index = 0.
3. Walk every cell that should be a solid wall. For each one:
   • Compute its outward unit normal (nx, ny).
   • normal.append(np.array([nx, ny]))
   • array[i, j] = pack_uint16_to_uint32(1, index)
   • index += 1
4. Leave every other cell as 0 — fluid interior and exterior/empty space both just stay at the default.
5. Export with the shared helpers:

   _write_grid_bin(destination + "geometry.bin", Ni, Nj, norm, array)
   _write_normals_bin(destination + "geometry_n.bin", normal)

Minimal skeleton for an arbitrary shape, assuming you have (or write) a predicate for "is this cell solid" and a way to compute the outward normal there (analytically, from a signed distance field, from an image mask, etc.):

import numpy as np
 
def create_custom_bc(Ni, Nj, destination, norm=1):
    array = np.zeros((Ni, Nj), dtype=np.uint32)
    normal = []
    index = 0
 
    for i in range(Ni):
        for j in range(Nj):
            if not is_solid(i, j):
                continue  # stays 0
 
            nx, ny = outward_normal(i, j)
            length = np.sqrt(nx**2 + ny**2)
            normal.append(np.array([nx / length, ny / length]))
 
            array[i, j] = pack_uint16_to_uint32(1, index)
            index += 1
 
    _write_grid_bin(destination + "geometry.bin", Ni, Nj, norm, array)
    _write_normals_bin(destination + "geometry_n.bin", normal)

Invariants to keep:

• index in a packed cell must always equal the 0-based position of its normal inside normal — append the normal and pack the cell together, in the same step, so they can never drift out of sync.
• Every wall cell needs exactly one normal, and vice versa — don't append a normal without writing the matching cell, or write a wall cell without appending a normal.
• Total wall-cell count must stay under 2^29.
• Normals must be unit length.
• If your shape is a hollow shell (walls with a gap between them, like a pipe), and the domain boundary cuts through that gap, remember to seal the gap there — otherwise the two shells are structurally disconnected and the enclosed cavity effectively vents into the surrounding fluid.

Worked example: the nozzle generator

create_nozzle_bc(Ni, Nj, destination, norm) builds a converging–diverging (de Laval) nozzle as a hollow double-walled shell:

• Radial profile. A bore radius r(x) is defined piecewise along the axial direction: a straight chamber section, a quadratic converging section down to a throat radius, and a t^1.5 diverging/bell section out to the exit radius. A matching derivative r'(x) is used to tilt each wall's normal so it isn't purely radial on the converging/diverging sections — it follows the local slope of the contour.
• Wall thickness. Measured outward from r(x): an inner wall (d_wall1), a hollow gap (d_gap), then an outer wall (d_wall2).
• Inlet cap. The combustion-chamber end is sealed with a solid disk (normal pointing along -x) so the bore isn't open to the edge of the texture.
• Exit rim. The exit end is intentionally left open (no cap, so exhaust can leave the domain), but the last few axial layers (EXIT_CAP_THICKNESS) fuse the hollow gap into solid wall — this is the "seal the cavity where the domain boundary cuts through it" case described above; without it the gap between inner and outer wall would vent straight into the exhaust flow.
• Two markers, one flag. Wall cells are tagged pack_uint16_to_uint32(1, index); bore-fluid cells are additionally tagged pack_uint16_to_uint32(2, index) purely as an internal marker in the raw array — per the quirk above, this doesn't survive unpack_uint32_to_uint16 or the visualization step as anything other than "not a wall." | knob | meaning | |---|---| | d_wall1, d_gap, d_wall2 | inner wall / gap / outer wall thickness | | R_in, R_t, R_e | inlet / throat / exit radius | | chamber_length, converging_length, exit_length | axial length of each contour section | | margin | empty padding at each end of the texture | | EXIT_CAP_THICKNESS | how many axial layers near the exit fuse inner+outer wall shut |

Validating your output

Before trusting a new generator, check:

• Grid file size should be exactly 12 + 4 * Ni * Nj bytes (3 header int32s + the uint32 array).
• Normals file size should be exactly 8 * N bytes, where N is the final index count (two float32s per normal).
• Visual sanity check — decode the wall flag and plot it before exporting:

  high = (array >> 29) & 1
  plt.imshow(np.where(high == 1, 1.0, 0.0).T, origin="lower", aspect="auto")

• If your geometry doesn't fit inside Ni/Nj (or your last section runs past the domain), grow the grid or shrink the geometry rather than letting it clip — clipped geometry silently produces an incomplete, structurally-broken boundary.

Tips

• Don't modify pack_uint16_to_uint32, unpack_uint32_to_uint16, or the write helpers. Every generator in the project relies on them meaning the same thing; treat them as the shared file-format contract.
• The format is effectively binary (wall / not-wall), not multi-tag, per the bit-29 quirk above. If you need a real multi-category boundary (e.g. distinct inlet/outlet/wall types), that needs a deliberate format change on both the writer and reader side, not just a new high_16 value.
• Index ordering doesn't need to match scan order. You can walk cells in any order you like (row-major, column-major, by feature, etc.) — the only requirement is that each cell's packed index matches its position in the normals list at the moment you write it.