The following problems appeared in The Riddler. And it features an interesting combination of an ancient game with the four-color theorem.
The famous four-color theorem states, essentially, that you can color in the regions of any map using at most four colors in such a way that no neighboring regions share a color. A computer-based proof of the theorem was offered in 1976.
Some 2,200 years earlier, the legendary Greek mathematician Archimedes described something called an Ostomachion (shown below). It’s a group of pieces, similar to tangrams, that divides a 12-by-12 square into 14 regions. The object is to rearrange the pieces into interesting shapes, such as a Tyrannosaurus rex. It’s often called the oldest known mathematical puzzle.
Your challenge today: Color in the regions of the Ostomachion square with four colors such that each color shades an equal area. (That is, each color needs to shade 36 square units.) The coloring must also satisfy the constraint that no adjacent regions are the same color.
Extra credit: How many solutions to this challenge are there?

Here are the details of how I solved the problem:
[Show Solution]
I’d like to start by pointing out that Hector Pefo also posted a nice solution to this problem. There were also diagrams of solutions posted on Twitter, such as those from @DimEarly and from @chattinthehatt. The solution I’ll describe here is similar to Hector’s, but I’ll explain how to model it as an integer linear program and solve it that way. Of course both solutions produce the same answers (because they’re both correct!) but it’s nice to see different approaches. So without further ado…
Integer linear program
Let’s suppose there are $N$ shapes and we’d like to use $C$ colors. Define:
\[
x_{ij} = \begin{cases}1&\text{if shape }i\text{ is colored using color }j\\0&\text{otherwise}\end{cases}
\]The first constraint is that each shape must be assigned to exactly one color. We can express this algebraically as follows:
\[
\sum_{j=1}^C x_{ij} = 1\qquad\text{for }i=1,\dots,N
\]Now suppose that shape $j$ has area $a_j$. The constraint that each color has an equal share of the area means that the sum of the areas of the shapes of each color must be $A_\text{tot}/C$, where $A_\text{tot}$ is the total area of the Ostomachion. Algebraically, we can write this as:
\[
\sum_{i=1}^N a_i x_{ij} = A_\text{tot}/C\qquad\text{for }j=1,\dots,C
\]Finally, we must deal with the edge constraints. We say that $(p,q)$ is an “edge” if shape $p$ shares an edge with shape $q$. Shapes that share an edge must be different colors. We can encode this algebraically as:
\[
x_{pj} + x_{qj} \le 1\qquad\text{for all edges }(p,q)\text{ and for }j=1,\dots,C
\]And that’s it! We described the feasible set of solutions as the set of matrices $x \in \{0,1\}^{N\times C}$ that satisfy the three sets of linear constraints above. This sort of optimization model is a modified version of a Knapsack problem. More generally, it’s an Integer Linear Program (ILP). Such problems are known to be NP-hard, which roughly means that in the worst-case, there is no efficient way to solve them short of trying all possible colorings. However, modern ILP solvers are quite sophisticated and can usually solve these types of problems very efficiently.
Finding all solutions
When using an ILP model, most commercial solvers are designed to find a single solution. For this problem, we are interested in finding all solutions. One way to get around the limitation is to use Barvinok’s algorithm. One can also use variants of branch-and-cut. Yet another way, which is not as systematic but much easier to explain and implement, is to use a randomized approach.
Since all our variables are binary (0 or 1), the convex hull of the feasible set cannot contain any interior solution points. All solutions will be on the boundary. This means that any solution can be found by specifying the right objective function for the optimization. As stated, our optimization is just a feasibility problem and there is no objective specified. For this problem, we start by letting $r_{ij}$ be random numbers (we can draw them i.i.d. from a standard normal distribution, for example). Then, we simply add the objective function:
\[
\text{Minimize}\qquad \sum_{i=1}^N\sum_{j=1}^C r_{ij} x_{ij}
\]and solve the problem. If we do this repeatedly with different realizations of the $\{r_{ij}\}$, we will eventually enumerate all possible solutions of the optimization problem.
To jump straight to the solution (and pretty pictures!) see below.
[Show Solution]
Solutions
I implemented the ILP described earlier using Julia and JuMP and I solved everything using the Gurobi solver. Although the randomized approach for finding all solutions works well, Gurobi has built-in functionality that allows one to directly find all solutions. In addition to doing this, I also included logic that removes duplicate solutions arising from different re-colorings of the same solution. e.g. if you take a solution and just swap red for green, I counted that as the same solution. Here are the solutions I found. I included the time it took Gurobi (on my laptop) to find all solutions for each case:
For $C=3$ colors, there are $2$ distinct solutions (took 1.6 ms to solve)

For $C=4$ colors, there are $9$ distinct solutions (took 32 ms to solve):

For $C=6$ colors, there are $33$ distinct solutions (took 700 ms to solve):

Interestingly, using fewer colors makes it easier to satisfy the equal-area constraint, but it makes it harder to satisfy the edge constraints. Ultimately, there are only solutions for the cases above ($C=3,4,6$). There are no solutions for any other values of $C$.
If you are interested in seeing my code, here is my IJulia notebook.
The small triangle in the lower-left corner is the same color (green) in all 44 solutions.
What do you make of that?
Good observation! I counted two colorings as being the same solution if they had the same coloring pattern but with permuted colors. For example, if you had a shape with regions (R1,R2,R3,R4), and they were colored with (red,blue,green,blue), respectively, then this would be equivalent to the coloring (blue,green,red,green), since the defining characteristic of the first coloring is that R2 and R4 are the same color. Every 4-coloring belongs to a family of 24 equivalent colorings (since 4! = 24) which come from just permuting the colors.
In order to make the duplicate removal procedure more efficient so that I only count one of the 24 equivalent colorings for each solution, I constrained one of the triangles to always be green — precisely the triangle you pointed out.
Have you ever considered that the (O)Stomachion might be a 2-D rendering of a 3-D shape?