# Sticky polyhexes

## Background

A polyhex of size $$\n\$$ is a contiguous shape made from joining $$\n\$$ unit regular hexagons side-by-side. As an example, the following image (from Wikipedia) contains all 7 distinct tetrahexes (polyhexes of size 4):

A dihex is a polyhex of size 2.

Now, let's define a sticky polyhex as follows:

• A single dihex is a sticky polyhex.
• If a polyhex X is sticky, X plus a dihex Y is also sticky if X and Y share at least two unit sides.
• Any polyhex that cannot be formed by the above definition is not sticky.

The following is the only size-2 sticky polyhex:

 _
/ \
\ /
/ \
\_/


The following are the only two essentially different size-4 sticky polyhexes:

 _   _
/ \_/ \
\ /  _/
/ \_/
\_/
_
_/ \
/ \ /
\ / \
/ \_/
\_/


The following is a simplified illustration of all 15 size-6 sticky polyhexes:

    * * *    * * *    *        * * * *    * * *
* *      * *    * * * *    * *          * * *
*          *        *

* *       *        * *       * *     * * *
* * *   * * * *    * * *   * * *       * *
*       *      *         *         *

* * *    * * *    * * * *   * * *    * * *
* *    * * *    *   *       *        * *
*                       * *        *


## Challenge

Given a positive integer $$\n\$$, count the number of sticky polyhexes of size $$\2n\$$. Two polyhexes which are equivalent under rotation and reflection are counted only once.

If you want, you may take the value of $$\2n\$$ as input instead (and assume the value is always even).

Standard rules apply. The shortest code wins.

## Test cases

The following are the expected results for $$\1 \le n \le 7\$$. Generated using reference implementation in Python (it times out shortly after calculating $$\n = 6\$$; the result for $$\n = 7\$$ was found on my local machine, and the program segfaulted while calculating $$\n = 8\$$).

1, 2, 15, 110, 1051, 10636, 113290


This sequence is now on OEIS as A342963.

• Add it to the OEIS! Apr 1, 2021 at 0:58
• @PeterKagey Draft submitted. :) Apr 1, 2021 at 1:13
• I think a fastest-code version of this challenge would also be very interesting :) Apr 1, 2021 at 17:27

length<$>iterate(\h->nub[r$x:y:u|u<-h,t<-u,x<-n t\\u,y<-n x\\u,([x,x,y]\$$u>>=n))/=[x,y]])[[l,0<l]] r u=minimum[sort(!!i).p.map(*j).m a<>u|a<-u,j<-l,i<-[0..5]] n=(<>p l).m m=zipWith(-) p=permutations l=[-1..1] import Data.List  Try it online! An infinite list whose $$\n\$$-th element (1-indexed) is the number of sticky polyhexes of size $$\2n\$$. This is probably the longest I've ever spent golfing a single program, and at the same time it's the answer I'm most proud of. After golfing 90 bytes off of the first iteration, it's time explain in some depth how this solution works and share the things I've learned. ## The grid The first thing to decide is how to handle the hexagonal grid. At first I was demotivated by the apparent lack of regularity in any coordinate system I could find, but then I stumbled across this website1. Along with some amazing visualization tools, which helped me later in debugging my code, it describes several reasonable coordinate systems for a hexagonal grid. The one I settled for in the end is called cube coordinates and is very clever. The following animation (taken from the same website) explains the concept beautifully. We start from the standard lattice grid in $$\3\$$-dimensional space, and then we cut along the plane $$\x+y+z=0\$$. On this plane, we get a hexagonal grid, where each hexagon is encoded as a triplet of integers $$\[x:y:z]\$$ with the constraint that $$\x+y+z=0\$$. Coordinate systems such as this one are called homogeneous, since they have one extra degree of freedom (three numbers for two dimensions) but they also have one additional constraint. While using three numbers to represent two dimensions may seem wasteful, the redundancy of cube coordinates is abundantly outweighted by the elegance and symmetry they provide. The main advantage is that, since we are operating in a vector subspace2 of $$\\mathbb{R}^3\$$, all the linear operations we may want to perform (translations, rotations, reflections) become astoundly straightforward (more on this later). For now, the following picture3 should give you an intuition of how cube coordinates work. Before moving on to the actual problem, let's code(golf) a few grid-related functions. As anticipated, each hexagon is encoded as a triplet of integers with sum equal to $$\0\$$. In Haskell4, we will represent it as a list [x,y,z]. • Since, as I said, we are working in a vector subspace, there is a well-defined notion of displacement vector between two points: coordinate-wise, we just take the difference of each coordinate. The relevant function is m=zipWith(-) (mnemonic: minus), and it works as follows: m[x1,y1,z1][x2,y2,z2] ≡ [x1-x2,y1-y2,z1-z2]. • Unsurprisingly, we will later need a function to compute the neighboring hexes of a given hexagon. By looking at the above picture, it's easy to see that all the neighbors of point $$\[x:y:z]\$$ can be obtained by summing (coordinate-wise) $$\[-1:0:1]\$$ or a permutation thereof. After importing Data.List (I know, that's very costly, but totally worth it5) and creating the aliases p=permutations (mnemonic: permutations) and l=[-1..1]≡[-1,0,1] (mnemonic: list6), we can define our function n=(<>p l).m (mnemonic: neighbors). If you're not fluent in Haskell, this is what n would look like if I weren't codegolfing: n x=[m x y|y<-p l]. In other words, we take each permutation y of l and subtract it (coordinate-wise) from x; adding or subtracting is the same thing, since -l is a permutation of l. ## Attaching dihexes Now that we have our grid set up, it's time to move on to the actual problem. A polyhex will be very naturally represented as a list of hexes, i.e. a list of lists of integers. Assume we are given a list h containing all the sticky polyhexes of size $$\2n\$$. Our task for this section will be to compute a list of all the sticky polyhexes of size $$\2n+2\$$, possibly with duplicates. This is the relevant code, where we abuse7 the list difference operator (\$$.

[x:y:u|                         -- the list of  polyhexes obtained by adding x and y to u, where
u<-h,                       -- u is a polyhex in h
t<-u,                       -- t is a hexagon in u
x<-n t\\u,                  -- x is a neighbor of t not in u
y<-n x\\u,                  -- y is a neighbor of x not in u
([x,x,y]\\(u>>=n))/=[x,y]   -- ?? see below
]


Everything should be pretty straightforward. The only thing that's missing is checking that the polyhex x:y:u is actually sticky, i.e. that the dihex [x,y] has at least two sides in common with u. The most concise way of checking this I could find is ([x,x,y]\\(u>>=n))/=[x,y]. How does it work? u>>=n is the list of neighbors of u, counted with multiplicity. Since x appears at least once in this list (by definition x is a neighbor of t), the difference [x,x,y]\\(u>>=n) is equal to [x,y] if and only if the dihex [x,y] shares exactly one (i.e. less than two) side with u.

## Canonical representatives

The only thing that's left is to remove duplicates. We say that two polyhexes are equivalent if they are related by a finite sequence of translations, rotations and reflections. A standard approach for this kind of task is to choose a canonical representative for every polyhex, with two rules in mind:

1. a polyhex is equivalent to its canonical representative;
2. if two polyhexes are equivalent, then they have the same canonical representative (by rule 1., the converse also holds).

If we find a suitable canonical representative for every polyhex, the task is solved: simply replace each polyhex with its representative, and then remove duplicates the standard way (which, in Haskell, is nub).

A reasonable choice for canonical representatives is the following: given a polyhex $$\u\$$, we define its representative $$\r(u)\$$ as the lexicographically minimum polyhex which is equivalent to $$\u\$$ and contains the hex $$\[0:0:0]\$$. We will always assume that the hexes in a polyhex are lexicographically sorted, so that this definition is well-posed.

This choice of representatives may seem somewhat arbitrary, but the idea behind it is really simple: given a polyhex $$\u\$$, we generate the set of all the equivalent polyhexes (with the restriction of containing $$\[0:0:0]\$$, otherwise there would be infinitely many). Clearly, if we apply this procedure to every polyhex equivalent to $$\u\$$ we will end up with the same set. Now we just need to pick a fixed element of this set, and (with bytes in mind) we select the minimum.

Now that we have an objective in mind, how do we implement it? This is where the linearity of cube coordinates comes into play. Equivalent polyhexes are always related by linear8 transformations, which should be relatively easy to implement. Before we look at the code, let's analyze these transformations in detail. Translations are trivial: the translation sending hexagon x to [0,0,0] is just m x.

As far as rotations and reflections are concerned, by looking at the picture above you can convince yourself that they all behave as follows: permute the three coordinates and, possibly, flip all the signs. For instance, the following are rotations or reflections: $$[x:y:z]\longmapsto[y:z:x]\\ [x:y:z]\longmapsto[-x:-z:-y].$$

We are now ready to implement the function r (mnemonic: representative).

r u=minimum[sort$(!!i).p.map(*j).m a<$>u|a<-u,j<-l,i<-[0..5]]


In plain English, the code does the following. Given a polyhex u, we want to compute its canonical representative as defined above, so we generate a list of all the polyhexes containing [0,0,0] and equivalent to u. To this aim, we fix:

• a hexagon a in u (which will be mapped to [0,0,0] with a translation);
• a sign j in [-1,1];
• a permutation index i in [0,1,2,3,4,5].

The code (!!i).p.map(*j).m a<$>u takes u and applies a translation (m a, mapping a to [0,0,0]), possibly a sign flip (map(*j)), and finally a permutation of the coordinates ((!!i).p). The sort ensures that the polyhex we obtain is "normalized" (i.e. lexicographically sorted). Finally, we take the lexicographically minimum polyhex in our list, and this is the canonical representative of u. A small codegolfing note: you may have noticed the j<-l instead of the expected j<-[-1,1]. This saves 5 bytes at the cost of producing some "degenerate" polyhexes, but the output is the same since a degenerate polyhex will never be selected as the minimum. ## Conclusion Now that we have all the ingredients, let's put them together! This is the last bit of code (with ellipses instead of the section already explained above). length<$>                   -- take the length of each list in
iterate(                    -- the list generated by repeatedly applying
\h->                    --     the function that takes the list h of sticky polyhexes of size 2*n
nub[                --     and removes duplicates from
r$-- the list of canonical representatives of x:y:u|...]) -- all the polyhexes obtained by suitably attaching dihexes to elements of h [[l,0<$l]]              -- starting from the list of sticky dihexes (there is just one)


0<\$l is short for [0,0,0].

And that's it! Thanks for reading through this huge post (or just skipping to the end to read the footnotes, most likely). Despite having spent several hours golfing this program, I believe there are still at least a few bytes that can be somehow shaved off: if you have any ideas, please let me know!

1 Ok maybe I should've resorted to Google earlier, it's literally the first result for "hex grid".
2 Don't worry if you don't know what I'm talking about, I have this bad habit of inserting random mathematical terms even when it's completely unnecessary.
3 Guess where it's taken from.
4 Please don't stop reading, I promise there will be no catamorphisms.
5 Trust me, I've spent several hours trying every possible variation of the code.
6 Creativity is my middle name.
7 Can you blame me? It's just two bytes and it's amazing!
8 +10 internet points if you shouted "Translations aren't linear!" with a disgusted voice.

# C++, Not competing

After spending so long on the golfing side of this challenge, I was curious to see what I could do from a perspective. Here is my C++ attempt, which should be able to run up to $$\n=8\$$ on every computer, and $$\n=9\$$ if you have enough RAM and patience. The approach is the same as my other answer, but optimized for speed instead of bytecount. Here are the results I got so far.

$$\n\$$ Number of sticky polihexes of size $$\2n\$$
1 1
2 2
3 15
4 110
5 1051
6 10636
7 113290
8 1234189
9 13674761
10 153289285

It's highly unlikely I'll ever be able to go further, it looks like 256GB of RAM are not enough to compute the list of sticky polyhexes for $$\n=11\$$.

#include <algorithm>
#include <iostream>
#include <iterator>
#include <vector>

using T = char;

struct hex {
T x, y;
constexpr T z() const {
return -x -y;
}
bool operator < (const hex& other) const {
return x < other.x or (x == other.x and y < other.y);
}
bool operator == (const hex& other) const {
return x == other.x and y == other.y;
}
template<typename F>
void on_neighbors(F f) const {
f(hex{ x + 1, y });
f(hex{ x + 1, y - 1});
f(hex{ x, y + 1 });
f(hex{ x, y - 1 });
f(hex{ x - 1, y + 1 });
f(hex{ x - 1, y });
}
};

using polyhex = std::vector<hex>;

std::vector<polyhex> extend(const std::vector<polyhex>& h) {
auto nub = [&](auto& v) {
std::sort(v.begin(), v.end());
v.erase(std::unique(v.begin(), v.end()), v.end());
};
auto in = [&](const auto& v, const hex& a) {
return binary_search(v.begin(), v.end(), a);
};
std::vector<polyhex> ans0;
for(const auto& u : h) {
std::vector<hex> layer;
for(const hex& t : u) {
t.on_neighbors([&](hex&& a) {
if(not in(u, a)) {
layer.push_back(a);
}
});
}
nub(layer);
std::vector<polyhex> ans1;
for(const hex& a : layer) {
int cnt0 = 0;
a.on_neighbors([&](hex&& b) {
if(in(u, b)) {
++cnt0;
}
});
a.on_neighbors([&](hex&& b) {
if(not in(u, b) and (a < b or not in(layer, b))) {
int cnt = cnt0;
b.on_neighbors([&](hex&& c) {
if(in(u, c)) {
++cnt;
}
});
if(cnt < 2) {
return;
}
polyhex v = u;
v.push_back(a);
v.push_back(b);
v.shrink_to_fit();
std::sort(v.begin(), v.end());
ans1.push_back(v);
}
});
}
nub(ans1);
ans0.insert(ans0.end(), std::make_move_iterator(ans1.begin()), std::make_move_iterator(ans1.end()));
}
nub(ans0);
std::vector<polyhex> ans;
for(auto it = std::make_move_iterator(ans0.begin()); it != std::make_move_iterator(ans0.end()); ++it) {
polyhex v = *it, best = v;
auto trans = [&](auto&& f) {
polyhex vv = v;
for(auto& c : vv) {
c = f(c);
}
std::sort(vv.begin(), vv.end());
hex m = vv.back();
for(auto& c : vv) {
c = hex{ c.x - m.x, c.y - m.y };
}
if(vv < best) {
best = std::move(vv);
}
};
for(int s = -1; s <= 1; s += 2) {
trans([&](hex& c) { return hex{ s * c.x, s * c.y }; });
trans([&](hex& c) { return hex{ s * c.x, s * c.z() }; });
trans([&](hex& c) { return hex{ s * c.y, s * c.x }; });
trans([&](hex& c) { return hex{ s * c.y, s * c.z() }; });
trans([&](hex& c) { return hex{ s * c.z(), s * c.x }; });
trans([&](hex& c) { return hex{ s * c.z(), s * c.y }; });
}
ans.push_back(best);
}
nub(ans);
return ans;
}

int main() {
std::vector<polyhex> h = { { hex{ -1, 0 }, hex{ 0, 0 } } };
for(int n = 1; ; ++n) {
std::cout << "n=" << n << " -> " << h.size() << std::endl;
h = extend(h);
}
return 0;
}


Try it online!

• Even in C++, it's nice that you were able to calculate several terms beyond what Bubbler's reference implementation could. +1. Apr 3, 2021 at 0:45

# Ruby, 290 270 bytes

->n,p=[[2,3]],s=2{s<n ?(q=[]
p.map{|a|a.product(d=[w=n*3,1,~w,-w,-1,-~w],d).map{|i|r=[];t=0
c=a+b=[i[0]+i[1],i.sum]
a&b==r&&12.times{|j|t+=a.count b[j%2]+d[j/2];c=c.map!{|k|j%2<1?k%w*(2+w)-k:k/w+k%w*w}.sort!.map{|m|m-c[0]+n};r<<c*1}
t>1&&q<<r.min}}
f[n,q&q,s+2]):p.size}


Try it online!

A function that accepts the total number of hexagons as an argument n (as permitted by the rules) and returns the number of sticky polyhexes. Will run up to 6 dihexes / 12 hexes on TIO in about 7 seconds.

This is about 3 times faster than the old version, because the code now checks if the new hexes overlap with the old hexes before calculating all 12 symmetries and 12 possible shared edges, hence saving a lot of wasted time. It also saves a byte, because r is empty at this stage, so we can use it to compare with a&b to ensure this array intersection is empty, instead of comparing with the literal empty array [].

Explanation

Each hex is saved as a single number; polyhexes are saved as arrays.

The coordinate system is a square grid of infinite height, and width w three times the number of hexes in the last generation to be investigated. The hexagonal grid is mapped onto this. Hexes with coordinates that differ by the following amounts are considered adjacent: 1=horizontal, w= vertical, and ~w = -w-1= diagonal. As it is a hexagonal grid, coordinates on the other diagonal are not considered adjacent.

The function proceeds by adding a pair of hexes to the previous generation and checking if the new hexes share at least 2 edges with the old ones; this requires 6 directions x 2 new hexes = 12 checks. To ensure the newly added hexes do not wrap around the grid, the old polyhex is positioned such that the hex with the lowest coordinate is at n (where n is the total number of hexes in the last generation to be investigated.)

Symmetries

Coincidentally, there are also 12 symmetries to be checked, so the same loop performs both duties. All symmetries are cycled through by alternately performing 2 different reflections. After each reflection, the hexes are sorted and the coordinates are adjusted so that the lowest coordinate is equal to n.

The first reflection is a vertical flip  k = k%w - k/w*w + k%w*w . The first term keeps the x coordinate the same and the second one (with integer division to discard the x coordinate) flips the y coordinate. The third term is a correction to maintain the skew. This formula can be simplified to k = k%w*(2+w)-k. This reflection is done first in case any new hexes have negative y coordinates.

The second is a diagonal flip which swaps the x and y coordinates: k = k/w+k%w*w. Doing this second ensures there are no negative y coordinates to handle (because after every reflection, the lowest hex's coordinate is adjusted to be equal to n.)

These two reflections combined give a 1/6 rotation (see below) hence 12 reflections cycle through all symmetries.

 SYMMETRY EXAMPLE
Start    Flip      Skew          Flip
Pattern  Vertical  Correction    Diagonal
01         43      54            50
5 2   ->  5 2  ->  0 3  ------>  4 1
43       01        12            32


Once the polyhex's coordinates under all 12 symmetry operations have been generated, the lexically lowest set of coordinates is added to the list of polyhexes.

Commented code (original version)

f=->n,p=[[2,3]]{                               #Assume initial hexes are at 2 and 3 if no info provided.
s=p[0].size                                   #Check size of first polyhex in list, others should be the same.
d=[w=n*3,1,~w,-w,-1,-~w]                      #d contains offsets for the six possible directions of travel. w is the width of the grid.
s<n ?(q=[]                                    #If size of polyhex less than required n, make array q to hold new generation.
p.map{|a|                                    #Iterate through list p of previous generation polyhexes a.
(s*36).times{|i|r=[];t=0                    #Try adding a chain of two hexes in 6*6 possible directions to each existing hex.
c=a+b=[g=a[i/36]+d[i/6%6],g+d[i%6]]        #b contains the new hexes chained onto old hex a[i/36]. Add them to existing polyhex a to make new polyhex c.
12.times{|j|t+=a.count(b[j%2]+d[j/2])      #Iterate through 2 new hexes x 6 directions = 12 to find the number of edges shared with the old polyhex.
c=c.map!{|k|j%2<1?k%w*(2+w)-k:k/w+k%w*w}. #Reflect the new polyhex...
sort!.map{|m|m-c[0]+n}                    #...sort the coordinates, and adjust them so that the smallest equals n
r<<c*1}                                   #Add the current symmetry of the new polyhex to r (the c*1 ensures a new copy is added, not just a pointer)
a&b==[]&&t>1&&q<<r.min                      #If new hexes b don't overlap old hexes a && there's more than one common edge add lowest value by symmetry to q
}
}
f[n,q&q]):p(p.size)                          #q&q to remove duplicates. Recursively call to build next generation. If s already equal to n, output the total.
}