Triangular domino tiling of an almost regular hexagon

Background

An almost regular hexagon is a hexagon where

• all of its internal angles are 120 degrees, and
• pairs of the opposite sides are parallel and have equal lengths (i.e. a zonogon).

The following is an example of an almost regular hexagon, with side lengths 2 (red), 4 (blue), and 3 (yellow). A triangular domino is a domino made of two unit triangles. A triangular domino tiling is a tiling on a shape using triangular dominoes. The following is a possible triangular domino tiling of the above shape (each color represents an orientation of each triangular domino): Challenge

Given the lengths of the three sides of an almost regular hexagon, find the number of distinct triangular domino tilings. The three sides will be always positive integers.

Alternative description

The second image shows that such a tiling can be viewed as an isometric view of stacked unit cubes. Now let's assign three directions to three axes in 3D:

• x = down-right / southeast / SE (blue edges in the first image)
• y = down-left / southwest / SW (red edges)
• z = up / north / N (yellow edges)

Then the stacked unit cubes can be represented as an x-by-y 2D array, whose items represent the height of the stack at that position. So the above challenge is equivalent to the following:

Given three positive integers x, y, and z, find the number of x-by-y arrays whose elements are between 0 and z inclusive and all rows and columns are in decreasing order.

It happens that this is one definition of plane partition in the form of $$\ PL(x,y,z) \$$, and it has a closed-form formula:

$$PL(x,y,z) = \prod_{i=1}^x \prod_{j=1}^y \prod_{k=1}^z \frac{i+j+k-1}{i+j+k-2}$$

Scoring and winning criterion

Standard rules apply. The shortest code in bytes wins.

Note that a submission is valid even if it suffers from integer overflows or floating-point inaccuracies, as long as the underlying algorithm is correct.

Test cases

x,y,z => output
---------------
1,1,1 => 2
1,1,2 => 3
1,2,3 => 10
2,3,1 => 10 (the order of inputs doesn't matter, since it's the same hexagon)
2,3,4 => 490
3,4,2 => 490
3,3,5 => 14112
• @LuisMendo If you're referring to the second figure, each color represents an orientation of each domino. A domino is a rhombus (i.e. diamond) bounded by solid black edges. Nov 18 '19 at 23:46
• Of course, they are rhombi, not triangles. My bad. Still, I'd clarify that the colour just indicates orientation for easier viewing of the figure Nov 18 '19 at 23:48
• Distinct up to what symmetry group? Nov 19 '19 at 8:16
• @PeterTaylor Symmetry is not considered. For example, 180 degrees rotation of the second image is distinct from the original. Nov 19 '19 at 8:23

APL (Dyalog Classic), 14 bytes

×/1+1÷1+∘,1⊥¨⍳

Try it online!

Uses 0 indexing with ⎕IO←0.

Explanation:

⍳     ⍝ Cartesian product of the ranges from 0 to n-1
1⊥¨      ⍝ Sum of each element (using base 1)
,         ⍝ Flattened
∘          ⍝ Composed with
1+1÷1+           ⍝ 1 + 1/(n+1)
×/                 ⍝ And reduced by multiplication

Perl 6, 38 30 bytes

{[*] map 1+1/(*+1),[X+] ^<<@_}

Try it online!

Based on the formula given in the question

Explanation:

{                            }    # Anonymous code block
^<<@_     # Map each input to the range 0 to n-1
[X+]           # Get the cross product of sums
map          ,               # Map these to
1+1/(*+1)                # 1+1/(n+1) = (n+2/n+1)
# Which is the formula compensating for the 0 based range
[*]                              # And reduce by multiplication

Wolfram Language (Mathematica), 28 bytes

Array[1+1/(##-2)&,#,1,1##&]&

Try it online!

05AB1E, 13 11 bytes

L.«â€˜OÍz>P

Explanation:

Uses the derived formula (inspired by @JoKing's answers): $$PL(x,y,z) = \prod_{i=1}^x \prod_{j=1}^y \prod_{k=1}^z \frac{1}{i+j+k-2}+1$$

L            # Map each value in the (implicit) input-list to an inner [1,v]-ranged list
.«          # (Right-)reduce these lists by:
â         #  Taking the cartesian product between two lists
€˜       # Then flatten each inner list
O      # Sum each inner list
Í     # Decrease all by 2
z    # Take 1/v for each value
>   # Increase all by 1
P  # And take the product of that
# (after which the result is output implicitly)
• @Grimmy I thought so as well at first, but from the challenge description: "Note that a submission is valid even if it suffers from integer overflows or floating-point inaccuracies, as long as the underlying algorithm is correct." Nov 19 '19 at 12:59
• Whoops, missed that! Nov 19 '19 at 12:59

foldr(\r k->k+k/(sum r-2))1.mapM(\n->[1..n])

Try it online!

Outputs floats. Thanks to @H.PWiz for 4 bytes using a fold.

• You can get shorter by replacing product.map with a fold Nov 19 '19 at 11:52
• @H.PWiz Nice one, thanks!
– xnor
Nov 20 '19 at 3:25

Jelly, 9 8 bytes

Œp§_2İ‘P

Try it online!

Uses the altered form of the closed form expression. Takes input as a list [x,y,z].

Œp          Cartesian product of the list.
(Each element, being an integer, is implicitly converted to a range.)
§         Sum the items of each triplet,
_2       subtract 2 from each sum,
İ      take the reciprocal of each lowered sum,
‘     increment each reciprocal,
P    and return the product of the increments reciprocals.

K (ngn/k), 17 15 bytes

%/*/'2 1+\:+/!:

Try it online!

last test fails because of an overflow

uses a 0-indexed version of the formula:

$$\PL(x,y,z)=\prod_{i=0}^{x-1}\prod_{j=0}^{y-1}\prod_{k=0}^{z-1}\frac{i+j+k+2}{i+j+k+1}\$$

Bash, 75

Inputs are passed as a comma-separated list.

e={1..${1//,/\}+\{1..}} eval eval echo \\$\[1 *\\$$e-1\\$$ /\\$$e-2\\$$ ]

Try it online! The last test fails due to integer overflow.

JavaScript (ES6),  73 65  64 bytes

(x,y,z)=>(g=n=>!n||g(n-1)/(s=n%z-~(n/z%y)+n/z/y%x|0)*-~s)(x*y*z)

Try it online!

Commented

(x, y, z) => (              // (x, y, z) = input
g = n =>                  // g is a recursive function taking a counter n
!n ||                   //   if n = 0, stop recursion and return 1
g(n - 1) / (            //   otherwise, divide the result of a recursive call by:
s =                   //     s defined as the sum of:
n % z               //       k, 0-indexed: n mod z
- ~(n / z % y)      //       j, 1-indexed: floor((n / z) mod y) + 1
+ n / z / y % x | 0 //       i, 0-indexed: floor((n / z / y) mod x)
)                       //
* -~s                   //   and multiply by s + 1
)(x * y * z)                // initial call to g with n = x * y * z

J, 25 bytes

-4 bytes thanks to Bubbler!

[:*/1+1%1+1#.&>[:,@{i.&.>

Try it online!

J, 29 bytes

[:*/[:(1+1%1++/)@>[:,@{<@i."0

Try it online!

A J port of @JoKing's APL answer (don't forget to upvote it), but twice as long. I'll try to golf it...

• 25 bytes. Mar 31 '20 at 4:00
• @Bubbler Thank you! Mar 31 '20 at 7:00

J, 18 bytes

*/@,@(1+%)1++/&i./

Try it online!

Takes z on its left and x y on the right.

\begin{align} \prod_{i=1}^x \prod_{j=1}^y \frac{i+j+z-1}{i+j-1} &= \prod_{i=1}^x \prod_{j=1}^y \left(1 + \frac{z}{i+j-1}\right) \\ &= \prod_{i=0}^{x-1} \prod_{j=0}^{y-1} \left(1 + \frac{z}{i+j+1}\right) \end{align}

How it works

*/@,@(1+%)1++/&i./  NB. Left =: z, Right =: x y
/  NB. Reduce over x y...
&i.   NB.   Apply range (n -> 0..n-1) to each item and
+/      NB.   Outer product by addition
1+        NB. Increment each, so that each cell has i+j+1
@(1+%)          NB. Compute 1+z/(i+j+1) for each i and j
*/@,                NB. Flatten the matrix and compute the product

Ruby, 85 75 bytes

->x,y,z{w=->t{(0...x*y*z).map{|r|r%x+(r/=x)%y+r/y+t}.reduce &:*};w/w}

Try it online!

Icon, 81, 78 70 bytes

-8 bytes thanks to Peter Taylor

procedure f(x,y,z)
p:=1&p*:=(1+z/(seq()\x+seq()\y-1.))&\u
return p
end

Try it online!

• Does this benefit from telescoping one product to $$\prod_{i=1}^x \prod_{j=1}^y \frac{i+j+z-1}{i+j-1} = \prod_{i=1}^x \prod_{j=1}^y \left(1 + \frac{z}{i+j-1}\right)$$? Nov 19 '19 at 17:12
• Answer: yes, 8-char reduction if the second line is changed to p:=1&p*:=(1+z/(seq()\x+seq()\y-1.))&\u Nov 19 '19 at 17:14
• @Peter Taylor That's great, thanks! Nov 19 '19 at 18:18

Python 3.8 (pre-release), 79 62 bytes

-17 bytes thanks to Jo King and xnor.

lambda x,y,z,t=1:[t:=t+t*z/(i%x-~i//x)for i in range(x*y)][-1]

Try it online!

• 63 bytes, though with a couple of minor rounding errors
– Jo King
Nov 19 '19 at 21:05
• You can also do [-1] rather than and t to extract what t ends p as.
– xnor
Nov 20 '19 at 3:19

Charcoal, 25 bytes

ＮθＮηＦＮＦηＦθ⊞υ⊕⁺λ⁺κιＩ÷Π⊕υΠυ

Try it online! Link is to verbose version of code. Explanation: Uses the 0-indexed version of the given formula.

ＮθＮηＦＮＦηＦθ

Input the three dimensions and loop over the implicit ranges.

⊞υ⊕⁺λ⁺κι

Take each incremented sum and save it in a list, thus flattening the Cartesian product.

Ｉ÷Π⊕υΠυ

Divide the product of the incremented list by the product of the list and cast to string for implicit print.

Japt, 18 bytes

mo rï Ëc rÄ pJ ÄÃ×

Try it

input as array

mo rï               // cartesian product [0,n)
Ëc            // flatten each triplet
rÄ         // reduce with initial value of 1 ( -2 +3=>triplets 0 to 1
pJ      // raised -1
Ä    // +1
Ã×  // reduce by multiplication

Implicit output, -m flag to run the program for each input