Hexagonal Grids: Integer Coordinate Mapping and Barycentric Computation Without Floating-Point Arithmetic

Hexagonal grids see wide use in geometric modeling, spatial indexing, map tiling, and graphics rendering pipelines. This article targets graphics and GIS developers who need to build their own hexagonal grid storage and rendering. Starting from the negative-index pain point of southward traversal under axial coordinates, it introduces a joint design of alternating face arrangement and integer barycentric coordinates, giving you a grid mapping approach that works entirely without floating-point arithmetic.

1. Why Southward Traversal Produces Negative Indices

Hexagonal grids are typically described with axial coordinates (q, r) or cube coordinates (x, y, z), where the latter satisfies x + y + z = 0. Axial coordinates map naturally to 2D array access, but offsets in many directions mix positive and negative signs — moving south increases the row number while the column may shift by either +1 or -1. In large-scale traversals or GPU layout scenarios, this introduces extra branching.

Take triangular faces in the equatorial region as an example. If faces of the same orientation are laid out consecutively, the row/column relationship can flip when crossing face boundaries, so a single traversal loop cannot cover an entire row uniformly. The fix is to organize faces in alternating orientations: the first face is an up-triangle, the second a down-triangle, and so on. With this arrangement, shared edges between adjacent faces show a consistent row/column growth trend in array coordinates, southward traversal no longer produces negative indices, and no extra subtraction branches are needed.

2. How Alternating Arrangement Unifies Row and Column Stepping

Let (row, col) denote a position in the 2D array, and introduce y = row - N_offset, where N_offset is a row offset constant that converts row numbers into offsets relative to the equator. A larger row number means moving south, and the column represents horizontal position. In up-triangle regions the geometric range is col < y; in down-triangle regions it is col >= y.

Moving one step south from any position increases the row by one and the column by at most one, so y and col grow in the same direction. The condition col < y continues to hold inside up-triangle regions, and col >= y no longer flips inside down-triangle regions. Within the entire rectangular block there is no need to switch offset formulas based on face type — the traversal logic reduces to fixed row/column stepping. For storage structures that pack a hexagonal grid into a rectangular array, this constraint significantly simplifies boundary checks and loop control.

3. Deriving Barycentric Coordinates with Integer Arithmetic

Inside a triangular sub-cell you often need barycentric coordinates (u, v, w) to express a point’s position relative to the vertices. The conventional approach uses vector cross products or area ratios, which involve floating-point division. When the grid satisfies row/column monotonicity, you can use integer arithmetic instead.

Let the sum of the barycentric components be N_granularity (a grid granularity parameter), satisfying u + v + w = N_granularity. Comparing col with y distinguishes which region the point falls in:

• Up-triangle (col < y):

  u = N_granularity - col
  v = y - col
  w = col

• Down-triangle (col >= y):

  u = N_granularity - y
  v = col - y
  w = y

You can verify this: in the up-triangle case, u + v + w = N_granularity - col + y - col + col = N_granularity + y - col. When the point lies on the boundary line col == y, we get v = 0, u = N_granularity - col, and w = col, so the sum is N_granularity. The down-triangle case works the same way. On the boundary line, both formulas yield {N_granularity - col, 0, col}, guaranteeing interpolation continuity.

Compared with the traditional floating-point approach — which needs roughly three divisions plus two subtractions — the integer approach needs only two subtractions plus one comparison. The computation is more deterministic and free of rounding error.

4. Reusing the Equatorial Formulas in the Polar Regions

Hexagonal grids are usually divided into an equatorial band and two polar regions. A polar region can be treated as a degenerate triangle fan: one vertex sits at the pole and the other two lie on the adjacent equatorial boundary. Its row/column relationship in the array can still be expressed as a linear condition between y and col. Therefore, you only need a single pole special case before entering the function — treat the pole as a fixed vertex of its triangle and return directly — while every other case reuses the same set of integer formulas as the equatorial band.

This compresses the coordinate conversion for the entire grid into one function, with no need to maintain separate logic for the equatorial band, the north pole, and the south pole.

5. Pseudocode Structure and Input Constraints

// Input constraints: row and col are non-negative integers giving the face's
// position in the 2D array; N_granularity > 0 is the sum of the barycentric
// components (the granularity parameter); N_offset is a row offset constant
// that must match the storage layout.
std::array<int, 3> get_barycentric(int row, int col, int N_granularity, int N_offset) {
    // Pole special case: if this face belongs to a polar triangle fan,
    // return the pole vertex coordinates
    if (is_pole(row, col)) {
        return pole_barycentric(row, col, N_granularity);
    }
    int y = row - N_offset;
    if (col < y) {
        // Up-triangle region
        int u = N_granularity - col;
        int v = y - col;
        int w = col;
        return {u, v, w};
    } else {
        // Down-triangle region
        int u = N_granularity - y;
        int v = col - y;
        int w = y;
        return {u, v, w};
    }
}
// is_pole checks whether (row, col) falls inside a polar triangle fan
// pole_barycentric returns the barycentric coordinates at the pole,
// typically {N_granularity, 0, 0} or a permutation of it

The branches above depend only on integer comparisons between col and y, with no floating-point operations. Note that when N_granularity or col grows large, u, v, and w may exceed the range of int; use an integer type of 32 bits or more, or add range checks up front.

6. Design Constraints and Applicability

This scheme relies on the premise that face column offsets and row offsets grow in the same direction. If your system uses a different row-major or column-major order, normalize the coordinate axes first; if the array index direction is geometrically opposite to north-south, adjust the meaning of the row increment and change the sign of y accordingly.

The scheme mainly optimizes traversal and interpolation within row/column-monotonic regions. If you need frequent jumps across regions (for example, from the equator to the opposite pole), you still need extra region-location logic, but it can be decoupled from the local barycentric computation. In large grids, integer computation is exactly reproducible, which makes unit testing and cross-platform consistency checks easier; compared with floating-point division, it avoids rounding discrepancies at boundaries and makes symmetry easier to guarantee.

7. Practical Recommendations

1. First fix the mapping between array major order and geometric direction. Every optimization here builds on that convention; changing it midway invalidates the formulas.
2. Treat face orientation and region conditions as first-class citizens. Abstract them into the two cases col < y and col >= y instead of writing if-else logic per face.
3. Prefer integer comparisons over floating-point area ratios. When the grid granularity is fixed, the integer formulas are fully deterministic.
4. Reuse a single local-coordinate function wherever possible. Encapsulate the equatorial/polar difference inside the pole special case.
5. Design consistency tests for boundary lines. For example, when col == y == 5 and N_granularity == 10, both the up-triangle and down-triangle formulas should return {5, 0, 5} — this prevents interpolation seams.

The hard part of hexagonal grid design is often not the math itself but coordinating storage layout, traversal logic, and coordinate computation. By constraining the face arrangement order and building integer barycentric coordinates from row/column conditions, you can achieve a more unified and more maintainable implementation across the equatorial band and the polar regions — a solid foundation for a self-built grid module.