Dominoes

From Rosetta Code
Task
Dominoes
You are encouraged to solve this task according to the task description, using any language you may know.

Warning: If you are looking for pizza delivery you have come to the wrong page. Bloody Google (other search engines are available).

Take a box of dominoes give them a good shuffle and then arrange them in diverse orientations such that they form a rectangle with 7 rows of 8 columns. Make a tableau of the face values e.g.

       05132231
       05505246
       43036620
       06235126
       11300245
       21433466
       64515414

Now torment your computer by making it identify where each domino is.

Do this for the above tableau and one of your own construction.

Extra credit:

How many ways are there to arrange dominoes in an 8x7 rectangle, first ignoring their values, then considering their values, and finally considering values but ignoring value symmetry, i.e. transposing 5 and 4.

F#

// Dominoes: Nigel Galloway. November 17th., 2021.
let cP (n:seq<uint64 list * uint64>) g=seq{for y,n in n do for g in g do let l=n^^^g in if n+g=l then yield (g::y,l)}
let rec fG n g=match g with h::t->fG(cP n h)t |_->fst(Seq.head n)
let solve(N:int[])=let fG=let y=fG [([],0UL)]([for g in 0..47->((N.[g],N.[g+8]),(1UL<<<g)+(1UL<<<g+8))]@[for n in 0..6 do for g in n*8..n*8+6->((N.[g],N.[g+1]),(1UL<<<g)+(1UL<<<g+1))]
                                |>List.groupBy(fun((n,g),_)->(min n g,max n g))|>List.sort|>List.map(fun(_,n)->n|>List.map(fun(n,g)->g))) in (fun n g->if List.contains((1UL<<<n)+(1UL<<<g)) y then "+" else " ")
                   N|>Array.chunkBySize 8|>Array.iteri(fun n g->let n=n*8 in [0..6]|>List.iter(fun y->printf $"%d{g.[y]}%s{fG(n+y)(n+y+1)}"); printfn $"%d{g.[7]}"; [0..7]|>List.iter(fun g->printf $"%s{fG(n+g)(n+g+8)} "); printfn "")

solve [|0;5;1;3;2;2;3;1;
        0;5;5;0;5;2;4;6;
        4;3;0;3;6;6;2;0;
        0;6;2;3;5;1;2;6;
        1;1;3;0;0;2;4;5;
        2;1;4;3;3;4;6;6;
        6;4;5;1;5;4;1;4|]
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +
6 4+5 1+5 4 1+4
solve [|5;6;2;0;0;4;1;4;
        3;6;1;3;0;4;2;2;
        3;5;6;4;3;2;1;1;
        3;5;1;1;3;0;0;5;
        6;0;5;4;3;5;5;2;
        4;4;1;3;6;6;0;2;
        1;2;6;2;6;5;0;4|]
Output:
5+6 2+0 0 4 1+4
        + +
3+6 1 3 0 4 2+2
    + +
3 5 6 4 3+2 1 1
+ +         + +
3 5 1+1 3+0 0 5

6 0 5+4 3+5 5+2
+ +
4 4 1+3 6 6+0 2
        +     +
1+2 6+2 6 5+0 4

Julia

const tableau = [
    0 5 1 3 2 2 3 1;
    0 5 5 0 5 2 4 6;
    4 3 0 3 6 6 2 0;
    0 6 2 3 5 1 2 6;
    1 1 3 0 0 2 4 5;
    2 1 4 3 3 4 6 6;
    6 4 5 1 5 4 1 4
]

const dominoes = [(i, j) for i in 0:size(tableau)[1]-1, j in 0:size(tableau)[2]-1 if i <= j]
sorted(dom) = first(dom) > last(dom) ? reverse(dom) : dom

"""
`patterns` contains solution(s), each containing a partially completed grid, the
dominos used, and steps taken to get to that point in the grid. Proceed via iterating
through possible tile placements from upper left to lower right, adding horizontal and
vertical tile placements, dropping those that require more than one of the same domino.
Consolidate in `patterns`` the newly lengthened layouts each step as moves are added.
"""
function findlayouts(tab = tableau, doms = dominoes)
    nrows, ncols = size(tab)
    patterns = [(zero(tab) .- 1, Tuple{Int, Int}[], Int[])]
    while true
        newpat = empty(patterns)
        for (ut, ud, up) in patterns
            pos = findfirst(x -> x == -1, ut)
            pos == nothing && continue
            row, col = Tuple(pos)
            if row < nrows && ut[row + 1, col] == -1 &&
               !(sorted((tab[row, col], tab[row + 1, col])) in ud)
                newut = copy(ut)
                newut[row:row+1, col] .= tab[row:row+1, col]
                push!(newpat, (newut, [ud; sorted((tab[row, col], tab[row + 1, col]))],
                   [up; [row, col, row + 1, col]]))
            end
            if col < ncols && ut[row, col + 1] == -1 &&
               !(sorted((tab[row, col], tab[row, col + 1])) in ud)
                newut = copy(ut)
                newut[row, col:col+1] .= tab[row, col:col+1]
                push!(newpat, (newut, [ud; sorted((tab[row, col], tab[row, col + 1]))],
                   [up; [row, col, row, col + 1]]))
            end
        end
        isempty(newpat) && break
        patterns = newpat
        length(last(first(patterns))) == length(doms) && break
    end
    return patterns
end

function printlayout(pattern)
    tab, dom, pos = pattern
    bytes = [[UInt8(' ') for _ in 1:(size(tab)[2] * 2 - 1)] for _ in 1:size(tab)[1]*2]
    for idx in 1:4:length(pos)-1
        x1, y1, x2, y2 = pos[idx:idx+3]
        n1, n2 = tab[x1, y1], tab[x2, y2]
        bytes[x1 * 2 - 1][y1 * 2 - 1] = Char(n1 + '0')
        bytes[x2 * 2 - 1][y2 * 2 - 1] = Char(n2 + '0')
        if x1 == x2 # horizontal
            bytes[x1 * 2 - 1][y1 * 2] = Char('+')
        elseif y1 == y2 # vertical
            bytes[x1 * 2][y1 * 2 - 1] = Char('+')
        end
    end
    println(join(String.(bytes), "\n"))
end


for pat in findlayouts()
    printlayout(pat)
end
@time findlayouts()

const t2 = [
    6 4 2 2 0 6 5 0;
    1 6 2 3 4 1 4 3;
    2 1 0 2 3 5 5 1;
    1 3 5 0 5 6 1 0;
    4 2 6 0 4 0 1 1;
    4 4 2 0 5 3 6 3;
    6 6 5 2 5 3 3 4
]
@time lays = findlayouts(t2, dominoes)
printlayout(first(lays))
println(length(lays), " layouts found.")
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +
6 4+5 1+5 4 1+4

  0.000507 seconds (6.06 k allocations: 1.715 MiB)

  0.023503 seconds (92.66 k allocations: 35.817 MiB)
6 4 2 2 0 6+5 0
+ + + + +     +
1 6 2 3 4 1+4 3

2 1 0 2 3+5 5 1
+ + + +     + +
1 3 5 0 5 6 1 0
        + +
4 2 6 0 4 0 1+1
+ + + +
4 4 2 0 5 3 6 3
        + + + +
6+6 5+2 5 3 3 4

2025 layouts found.

Extra credit task

""" From https://en.wikipedia.org/wiki/Domino_tiling#Counting_tilings_of_regions
The number of ways to cover an m X n rectangle with m * n / 2 dominoes, calculated
independently by Temperley & Fisher (1961) and Kasteleyn (1961), is given by
"""
function dominotilingcount(m, n)
    return BigInt(
        floor(
            prod([
                prod([
                    big"4.0" * (cospi(j / (m + 1)))^2 + 4 * (cospi(k / (n + 1)))^2 for
                    k in 1:(n+1)÷2
                ]) for j in 1:(m+1)÷2
            ]),
        ),
    )
end

arrang = dominotilingcount(7, 8)
perms = factorial(big"28")
flips = 2^28

println("Arrangements ignoring values: $arrang")
println("Permutations of 28 dominos: ", perms)
println("Permuted arrangements ignoring flipping dominos: ", arrang * perms)
println("Possible flip configurations: $flips")
println("Possible permuted arrangements with flips: ", flips * arrang * perms)
Output:
Arrangements ignoring values: 1292697
Permutations of 28 dominos: 304888344611713860501504000000
Permuted arrangements ignoring flipping dominos: 394128248414528672328712716288000000   
Possible flip configurations: 268435456
Possible permuted arrangements with flips: 105797996085635281181632579889767907328000000

Mathematica/Wolfram Language

ClearAll[VisualizeState]
VisualizeState[sol_List, tab_List] := Module[{rects},
  rects = Apply[Rectangle[#1 - {1, 1} 0.5, #2 + {1, 1} 0.5] &, sol[[All, 2]], {1}];
  Graphics[{FaceForm[], EdgeForm[Black], rects, MapIndexed[Text[Style[#1, 14, Black], #2] &, tab, {2}]}]
]
ClearAll[FindSolutions]
FindSolutions[tab_] := Module[{poss, possshort, possshorti, posssets, posss, sols},
  poss = Catenate[MapIndexed[#2 &, tab, {2}]];
  possshort = Thread[poss -> Range[Length[poss]]];
  possshorti = Thread[Range[Length[poss]] -> poss];
  posssets = Select[Subsets[poss, {2}], Apply[ManhattanDistance]/*EqualTo[1]];
  posssets = {# /. possshort, Sort[Extract[tab, #]]} & /@ posssets;
  posss = GatherBy[posssets, Last];
  posss = #[[1, 2]] -> #[[All, 1]] & /@ posss;
  posss //= SortBy[Last/*Length];
  
  sols = {};
  ClearAll[RecursePlaceDomino];
  RecursePlaceDomino[placed_List, left_List] := Module[{newplaced, sortedleft, newleft, next},
    If[Length[left] == 0,
     AppendTo[sols, placed];
     ,
     sortedleft = SortBy[left, Last/*Length];
     next = sortedleft[[1]];
     Do[
      newplaced = Append[placed, next[[1]] -> n];
      newleft = Drop[sortedleft, 1];
      newleft[[All, 2]] = newleft[[All, 2]] /. {___, Alternatives @@ n, ___} :> Sequence[];
      If[AnyTrue[newleft[[All, 2]], Length/*EqualTo[0]], Continue[]];
      RecursePlaceDomino[newplaced, newleft]
      ,
      {n, next[[2]]}
      ];
     ]
    ];
  RecursePlaceDomino[{}, posss];
  sols[[All, All, 2]] = sols[[All, All, 2]] /. possshorti;
  sols
]

tab = {{6, 2, 1, 0, 4, 0, 0}, {4, 1, 1, 6, 3, 5, 5}, {5, 4, 3, 2, 0, 
    5, 1}, {1, 3, 0, 3, 3, 0, 3}, {5, 3, 0, 5, 6, 5, 2}, {4, 4, 2, 1, 
    6, 2, 2}, {1, 6, 4, 2, 2, 4, 3}, {4, 6, 5, 6, 0, 6, 1}};
sols = FindSolutions[tab];
Length[sols]
VisualizeState[sols[[1]], tab]
tab = {{6, 4, 4, 1, 2, 1, 6}, {6, 4, 2, 3, 1, 6, 4}, {5, 2, 6, 5, 0, 
    2, 2}, {2, 0, 0, 0, 2, 3, 2}, {5, 5, 4, 5, 3, 4, 0}, {3, 3, 0, 6, 
    5, 1, 6}, {3, 6, 1, 1, 5, 4, 5}, {4, 3, 1, 0, 1, 3, 0}};
sols = FindSolutions[tab];
Length[sols]
VisualizeState[sols[[1]], tab]
Output:
1
[Graphical output of the solution]
2025
[Graphical output of the first solution]

Extra credit task

ClearAll[DominoTilingCount]
DominoTilingCount[m_, n_] := Module[{},
  Round[Product[
    4 Cos[(Pi j)/(m + 1)]^2 + 4 Cos[(Pi k)/(n + 1)]^2,
     {j, Ceiling[m/2]},
     {k, Ceiling[n/2]}
   ]]
]
arrangements = DominoTilingCount[7, 8] // Round;
permutations = 28!;
flips = 2^28;
Print["Arrangements ignoring values: ", arrangements]
Print["Permutations of 28 dominos: ", permutations]
Print["Permuted arrangements ignoring flipping dominos: ", arrangements*permutations]
Print["Possible flip configurations: ", flips]
Print["Possible permuted arrangements with flips: ", flips*arrangements*permutations]
Output:
Arrangements ignoring values: 1292697
Permutations of 28 dominos: 304888344611713860501504000000
Permuted arrangements ignoring flipping dominos: 394128248414528672328712716288000000
Possible flip configurations: 268435456
Possible permuted arrangements with flips: 105797996085635281181632579889767907328000000

Nim

Translation of: Wren
import std/[monotimes, sequtils, strformat, strutils, times]

const
  None = -1
  NoPosition = (None, None)

type
  Value = None..6
  Domino = (Value, Value)
  Position = (int, int)
  Tableau = array[7, array[8, Value]]
  Pattern = (seq[seq[Value]], seq[Domino], seq[Position])

const

  Tableau1 = [[Value 0, 5, 1, 3, 2, 2, 3, 1],
              [Value 0, 5, 5, 0, 5, 2, 4, 6],
              [Value 4, 3, 0, 3, 6, 6, 2, 0],
              [Value 0, 6, 2, 3, 5, 1, 2, 6],
              [Value 1, 1, 3, 0, 0, 2, 4, 5],
              [Value 2, 1, 4, 3, 3, 4, 6, 6],
              [Value 6, 4, 5, 1, 5, 4, 1, 4]]

  Tableau2 = [[Value 6, 4, 2, 2, 0, 6, 5, 0],
              [Value 1, 6, 2, 3, 4, 1, 4, 3],
              [Value 2, 1, 0, 2, 3, 5, 5, 1],
              [Value 1, 3, 5, 0, 5, 6, 1, 0],
              [Value 4, 2, 6, 0, 4, 0, 1, 1],
              [Value 4, 4, 2, 0, 5, 3, 6, 3],
              [Value 6, 6, 5, 2, 5, 3, 3, 4]]


func domino(v1, v2: Value): Domino =
  if v1 > v2: (v2, v1) else: (v1, v2)

func findLayouts(tab: Tableau; domCount: Positive): seq[Pattern] =
  let nrows = tab.len
  let ncols = tab[0].len
  var m = newSeqWith(nrows, repeat(Value None, ncols))
  result = @[(m, newSeq[Domino](), newSeq[Position]())]
  while true:
    var newpats: seq[Pattern]
    for pat in result:
      let (ut, ud, up) = pat
      var pos = NoPosition
      block Outer:
        for j in 0..<ncols:
          for i in 0..<nrows:
            if ut[i][j] == None:
              pos = (i, j)
              break Outer
      if pos == NoPosition:
        continue
      let (row , col) = pos
      if row < nrows - 1 and ut[row+1][col] == None:
        let dom = domino(tab[row][col], tab[row+1][col])
        if dom notin ud:
          var newUt = ut
          newut[row][col] = tab[row][col]
          newut[row+1][col] = tab[row+1][col]
          newpats.add (newut,
                       ud & domino(tab[row][col], tab[row+1][col]),
                       up & @[(row, col), (row+1, col)])
      if col < ncols - 1 and ut[row][col+1] == None:
        let dom = domino(tab[row][col], tab[row][col+1])
        if dom notin ud:
          var newUt = ut
          newut[row][col] = tab[row][col]
          newut[row][col+1] = tab[row][col+1]
          newpats.add (newut,
                       ud & domino(tab[row][col], tab[row][col+1]),
                       up & @[(row, col), (row, col+1)])
    if newPats.len == 0: break
    result = move(newPats)
    if result[0][2].len == domCount:
      break

proc printLayout(pattern: Pattern) =
  let (tab, _, pos) = pattern
  var output = newSeqWith(2 * tab.len, repeat(' ', tab[0].len * 2 - 1))
  var idx = 0
  while idx < pos.len - 1:
    let
      (x1, y1) = pos[idx]
      (x2, y2) = pos[idx+1]
    let
      n1 = tab[x1][y1]
      n2 = tab[x2][y2]
    output[x1*2][y1*2] = chr(48 + n1)
    output[x2*2][y2*2] = chr(48 + n2)
    if x1 == x2:
      output[x1*2][y1*2+1] = '+'
    elif y1 == y2:
      output[x1*2+1][y1*2] = '+'
    inc idx, 2
  for i in 0..output.high:
    echo output[i]

var domCount: Positive
for j in 0..Tableau1[0].high:
  for i in 0..Tableau1.high:
    if i <= j:
      inc domCount

for t in [Tableau1, Tableau2]:
  let start = getMonoTime()
  let lays = t.findLayouts(domCount)
  lays[0].printLayout()
  let lc = lays.len
  let pl = if lc > 1: "s" else: ""
  let fo = if lc > 1: " (first one shown)" else: ""
  echo &"{lc} layout{pl} found{fo}."
  echo &"Took {(getMonoTime() - start).inMilliseconds} ms.\n"
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +        
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +    
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +    
6 4+5 1+5 4 1+4
               
1 layout found.
Took 1 ms.

6 4 2 2 0 6+5 0
+ + + + +     +
1 6 2 3 4 1+4 3
               
2 1 0 2 3+5 5 1
+ + + +     + +
1 3 5 0 5 6 1 0
        + +    
4 2 6 0 4 0 1+1
+ + + +        
4 4 2 0 5 3 6 3
        + + + +
6+6 5+2 5 3 3 4
               
2025 layouts found (first one shown).
Took 42 ms.

Extra credit task

Translation of: Wren
Library: Nim-Integers
import std/[math, monotimes, strutils, times]
import integers

func dominoTilingCount(m, n: Positive): int =
  var prod = 1.0
  for j in 1..((m + 1) div 2):
    for k in 1..((n + 1) div 2):
      let cm = cos(PI * (j / (m + 1)))
      let cn = cos(PI * (k / (n + 1)))
      prod *= (cm * cm + cn * cn) * 4
  result = int(prod)

let
  start = getMonoTime()
  arrang = dominoTilingCount(7, 8)
  perms = factorial(28)
  flips = 1 shl 28

echo "Arrangements ignoring values: ", insertSep($arrang)
echo "Permutations of 28 dominos: ", insertSep($perms)
echo "Permuted arrangements ignoring flipping dominos: ", insertSep($(perms * arrang))
echo "Possible flip configurations: ", insertSep($flips)
echo "Possible permuted arrangements with flips: ", insertSep($(perms * flips * arrang))
echo "\nTook $# µs.".format((getMonoTime() - start).inMicroseconds)
Output:
Arrangements ignoring values: 1_292_697
Permutations of 28 dominos: 304_888_344_611_713_860_501_504_000_000
Permuted arrangements ignoring flipping dominos: 394_128_248_414_528_672_328_712_716_288_000_000
Possible flip configurations: 268_435_456
Possible permuted arrangements with flips: 105_797_996_085_635_281_181_632_579_889_767_907_328_000_000

Took 107 µs.

Perl

Library: ntheory
use v5.36;

# NB: not 'blank' lines, need to be full-width white-space
my $grid1 = <<END;
0 5 1 3 2 2 3 1
               
0 5 5 0 5 2 4 6
               
4 3 0 3 6 6 2 0
               
0 6 2 3 5 1 2 6
               
1 1 3 0 0 2 4 5
               
2 1 4 3 3 4 6 6
               
6 4 5 1 5 4 1 4
END

$grid2 = <<END;
0 0 0 1 1 1 0 2
               
1 2 2 2 0 3 1 3
               
2 3 3 3 0 4 1 4
               
2 4 3 4 4 4 0 5
               
1 5 2 5 3 5 4 5
               
5 5 0 6 1 6 2 6
               
3 6 4 6 5 6 6 6
END

my $grid3 = <<END;
0 0 1 1
       
0 2 2 2
       
1 2 0 1
END

sub find ($rows, $cols, $x, $y, $try, $orig) {
    state $solution;
    my $sum = $rows + $cols;
    my $gap = qr/(.{$sum}) (.{$sum})/s;

    if( $x > $y ) {
        $x = 0;
        $solution = $orig |. $try and return if ++$y == $rows;  # solved!
    }

    while ( $try =~ /(?=(?|$x $gap $y|$y $gap $x))/gx ) {       # vertical
        my $new = $try;
        substr $new, $-[0], 2*($rows+$cols)+3, " $1+$2 ";
        find($rows, $cols, $x + 1, $y, $new, $orig );
    }

    while ( $try =~ /(?=$x $y|$y $x)/g ) {                      # horizontal
        my $new = $try;
        substr $new, $-[0], 3, ' + ';
        find($rows, $cols, $x + 1, $y, $new, $orig );
    }

    $solution
}

say find(7, 8, 0, 0, $grid1, $grid1 ) . "\n=======\n\n";
say find(7, 8, 0, 0, $grid2, $grid2 ) . "\n=======\n\n";
say find(3, 4, 0, 0, $grid3, $grid3 ) . "\n=======\n\n";

use constant PI => 2*atan2(1,0);
use ntheory 'factorial';

sub comma { reverse ((reverse shift) =~ s/(.{3})/$1,/gr) =~ s/^,//r }

my($m,$n, $arrangements)  = (7,8, 1);
for my $j (1 .. $m/2) {
  for my $k (1 .. $n/2) {
    $arrangements *= 4*cos(PI*$j/($m+1))**2 + 4*cos(PI*$k/($n+1))**2
  }
}

printf "%32s:%60s\n", 'Arrangements ignoring values',     comma    $arrangements;
printf "%32s:%60s\n", 'Permutations of 28 dominos',       comma my $permutations = factorial 28;
printf "%32s:%60s\n", 'Flip configurations',              comma my $flips = 2**28;
printf "%32s:%60s\n", 'Permuted arrangements with flips', comma    $flips * $permutations * $arrangements;
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +        
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +    
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +    
6 4+5 1+5 4 1+4

=======

0+0 0+1 1+1 0+2
               
1 2+2 2 0+3 1+3
+     +        
2 3+3 3 0 4 1+4
        + +    
2+4 3+4 4 4 0+5
               
1 5 2+5 3+5 4+5
+ +            
5 5 0+6 1+6 2 6
            + +
3+6 4+6 5+6 6 6

=======

0+0 1+1

0+2 2+2

1+2 0+1

=======

    Arrangements ignoring values:                                                   1,292,697
      Permutations of 28 dominos:                     304,888,344,611,713,860,501,504,000,000
             Flip configurations:                                                 268,435,456
Permuted arrangements with flips: 105,797,996,085,635,281,181,632,579,889,767,907,328,000,000

Phix

with javascript_semantics

function domino_set()
    sequence set = {}
    for i=0 to 6 do
        for j=i to 6 do
            set = append(set,{i,j})
        end for
    end for
    return set
end function

sequence set = domino_set(),
        used = repeat(repeat(false,7),7),
        tags = shuffle(tagset(length(set))),
        grid

function unpack(sequence s)
    s = split(s,' ')
    s = apply(true,join,{s,'?'})
    s = join(s,"\n? ? ? ? ? ? ? ?\n")
    s = split(s,'\n')
    return s
end function

function clear(integer r, c, sequence s)
    if grid[r][c]!='?' then ?9/0 end if
    grid[r][c]='+'
    sequence res = {{r,c}}
    for i=1 to length(s) do
        {r,c} = s[i]
        if r>=1 and r<=13
        and c>=1 and c<=15 then
            integer prev = grid[r][c]
            if prev='?' then
                grid[r][c] = ' '
                res = append(res,{r,c})
            elsif prev!=' ' then
                ?9/0
            end if
        end if
    end for
    return res
end function

procedure restore(sequence s)
    for i=1 to length(s) do
        integer {r,c} = s[i]
        grid[r][c] = '?'
    end for
end procedure

function rand_grid(integer rem=28)
    if rem=0 then
        for r=1 to 13 by 2 do
            for c=2 to 14 by 2 do
                grid[r][c] = '?'
            end for
        end for
        for r=2 to 12 by 2 do
            for c=1 to 15 by 2 do
                grid[r][c] = '?'
            end for
        end for
        return grid
    end if
    for r=1 to 13 by 2 do
        for c=1 to 15 by 2 do
            bool flat = (c<15 and grid[r,c+1]='?'),
                 vert = (r<13 and grid[r+1,c]='?')
            sequence res = {},
                     opt = {}
            if flat then
                opt = append(opt,{{r,c+2,r,c+1},{{r,c-1},{r,c+3},{r-1,c},{r-1,c+2},{r+1,c},{r+1,c+2}}})
            end if
            if vert then
                opt = append(opt,{{r+2,c,r+1,c},{{r-1,c},{r,c-1},{r+2,c-1},{r,c+1},{r+2,c+1},{r+3,c}}})
            end if
            opt = shuffle(opt)
            for i=1 to length(opt) do
                integer {r2,c2,r3,c3} = opt[i][1]
                sequence tile = shuffle(set[tags[rem]])
                grid[r][c] = tile[1]+'0'
                grid[r2][c2] = tile[2]+'0'
                sequence reset = clear(r3,c3,opt[i][2])
                res = rand_grid(rem-1)
                if length(res) then return res end if
                restore(reset)
            end for
            if flat or vert then
                return {}
            end if
        end for
    end for
    return {}
end function

string soln1 = ""
string solnn = ""

function solve(integer rem=28)
    if rem=0 then
        solnn = join(grid,'\n')&"\n\n\n"
        if soln1 = "" then soln1 = solnn end if
        return 1
    end if
    for r=1 to 13 by 2 do
        for c=1 to 15 by 2 do
            bool flat = (c<15 and grid[r,c+1]='?'),
                 vert = (r<13 and grid[r+1,c]='?')
            integer count = 0
            if flat then
                integer {R,C} = sort({grid[r][c]-'0'+1,grid[r][c+2]-'0'+1})
                if not used[R][C] then
                    used[R][C] = true
                    sequence reset = clear(r,c+1,{{r,c-1},{r,c+3},{r-1,c},{r-1,c+2},{r+1,c},{r+1,c+2}})
                    count += solve(rem-1)
                    restore(reset)
                    used[R][C] = false
                end if
            end if
            if vert then
                integer {R,C} = sort({grid[r][c]-'0'+1,grid[r+2][c]-'0'+1})
                if not used[R][C] then
                    used[R][C] = true
                    sequence reset = clear(r+1,c,{{r-1,c},{r,c-1},{r+2,c-1},{r,c+1},{r+2,c+1},{r+3,c}})
                    count += solve(rem-1)
                    restore(reset)
                    used[R][C] = false
                end if
            end if
            if flat or vert then
                return count -- (may still be 0)
            end if
        end for
    end for
    return 0
end function

procedure test(sequence g)
    grid = g
    g = {}
    soln1 = ""
    solnn = ""
    atom t0 = time()
    integer n = solve()
    puts(1,soln1)
    if n>1 then
        if n>2 then puts(1,"...\n\n\n") end if
        puts(1,solnn)
    end if
    printf(1,"%d solution%s found (%s)\n\n\n",{n,iff(n=1?"":"s"),elapsed(time()-t0)})
end procedure
test(unpack("05132231 05505246 43036620 06235126 11300245 21433466 64515414"))
test(rand_grid())
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +
6 4+5 1+5 4 1+4


1 solution found (0.2s)


6+4 2+2 0+6 5 0
            + +
1+6 2+3 4+1 4 3

2+1 0+2 3+5 5+1

1+3 5+0 5+6 1+0

4+2 6 0 4+0 1+1
    + +
4+4 2 0 5 3+6 3
        +     +
6+6 5+2 5 3+3 4


...


6 4 2 2 0 6+5 0
+ + + + +     +
1 6 2 3 4 1+4 3

2 1 0 2 3+5 5 1
+ + + +     + +
1 3 5 0 5 6 1 0
        + +
4 2 6 0 4 0 1+1
+ + + +
4 4 2 0 5 3 6 3
        + + + +
6+6 5+2 5 3 3 4


2025 solutions found (0.1s)

Note that 2025 is not the maximum number of solutions or anything like that, just a higher than average result.

Extra credit

Pretty dumb brute force approach, dreadfully slow.

without js -- too slow
enum IGNORE, CONSIDER, NOSYM
function count(integer what, rem=28, doubles=6)
    if rem=0 then
        return 1
    end if
    atom total = 0
    for r=1 to 13 by 2 do
        for c=1 to 15 by 2 do
            bool flat = (c<15 and grid[r,c+1]='?'),
                 vert = (r<13 and grid[r+1,c]='?')
            sequence res = {},
                     opt = {}
            if flat then
                opt = append(opt,{{r,c+2,r,c+1},{{r,c-1},{r,c+3},{r-1,c},{r-1,c+2},{r+1,c},{r+1,c+2}}})
            end if
            if vert then
                opt = append(opt,{{r+2,c,r+1,c},{{r-1,c},{r,c-1},{r+2,c-1},{r,c+1},{r+2,c+1},{r+3,c}}})
            end if
            for i=1 to length(opt) do
                integer {r2,c2,r3,c3} = opt[i][1]
                sequence reset = clear(r3,c3,opt[i][2])
                if what=IGNORE then
                    total += count(what,rem-1)
                elsif what=CONSIDER then
                    if doubles then
                        total += doubles*count(what,rem-1,doubles-1)
                    end if
                    if rem>doubles then
                        total += 2*(rem-doubles)*count(what,rem-1,doubles)
                    end if
                else -- NOSYM
                    total += 2*rem*count(what,rem-1)
                end if
                restore(reset)
            end for
            if flat or vert then
                return total
            end if
        end for
    end for
    return total
end function

atom t0 = time()
printf(1,"Arrangements ignoring values: %,d\n",count(IGNORE))
--printf(1,"Arrangements considering values: %d\n",count(CONSIDER)) -- too slow
printf(1,"Arrangements ignoring symmetry: %g\n",count(NOSYM))
?elapsed(time()-t0)
Output:
Arrangements ignoring values: 1,292,697
Arrangements ignoring symmetry: 1.05798e+44
"2 minutes and 37s"

Much faster

Translation of: Wren
with javascript_semantics
include mpfr.e

function domino_tiling_count(integer m=7, n=8)
    atom prod = 1
    for j=1 to ceil(m/2) do
        for k=1 to ceil(n/2) do
            atom cm = cos(PI * (j / (m + 1))),
                 cn = cos(PI * (k / (n + 1)))
            prod *= (cm*cm + cn*cn) * 4
        end for
    end for
    return floor(prod)
end function
 
atom start = time()
integer arrang = domino_tiling_count(),
        flips  = power(2,28)
mpz perms = mpz_init()
mpz_fac_ui(perms,28)
 
printf(1,"Arrangements ignoring values: %,d\n", arrang)
printf(1,"Permutations of 28 dominos: %s\n", mpz_get_str(perms,10,true))
mpz_mul_si(perms,perms,arrang)
printf(1,"Permuted arrangements ignoring flipping dominos: %s\n", mpz_get_str(perms,10,true))
printf(1,"Possible flip configurations: %,d\n", flips)
mpz_mul_si(perms,perms,flips)
printf(1,"Possible permuted arrangements with flips: %s\n", mpz_get_str(perms,10,true))
printf(1,"Took %s\n",elapsed(time()-start))
Output:
Arrangements ignoring values: 1,292,697
Permutations of 28 dominos: 304,888,344,611,713,860,501,504,000,000
Permuted arrangements ignoring flipping dominos: 394,128,248,414,528,672,328,712,716,288,000,000
Possible flip configurations: 268,435,456
Possible permuted arrangements with flips: 105,797,996,085,635,281,181,632,579,889,767,907,328,000,000
Took 0.1s

Wren

Translation of: Julia

Basic task

var tableau = [
    [0, 5, 1, 3, 2, 2, 3, 1],
    [0, 5, 5, 0, 5, 2, 4, 6],
    [4, 3, 0, 3, 6, 6, 2, 0],
    [0, 6, 2, 3, 5, 1, 2, 6],
    [1, 1, 3, 0, 0, 2, 4, 5],
    [2, 1, 4, 3, 3, 4, 6, 6],
    [6, 4, 5, 1, 5, 4, 1, 4]
]

var tableau2 = [
    [6, 4, 2, 2, 0, 6, 5, 0],
    [1, 6, 2, 3, 4, 1, 4, 3],
    [2, 1, 0, 2, 3, 5, 5, 1],
    [1, 3, 5, 0, 5, 6, 1, 0],
    [4, 2, 6, 0, 4, 0, 1, 1],
    [4, 4, 2, 0, 5, 3, 6, 3],
    [6, 6, 5, 2, 5, 3, 3, 4]
]

var dominoes = []
for (j in 0...tableau[0].count) {
    for (i in 0...tableau.count) if (i <= j) dominoes.add([i, j])
}

var containsDom = Fn.new { |l, m, n|  // assumes m <= n
    for (i in 0...l.count) {
        var d = l[i]
        if (d[0] == m && d[1] == n) return true
    }
    return false
}

var copyTab = Fn.new { |t|
    var c = List.filled(t.count, null)
    for (r in 0...t.count) c[r] = t[r].toList
    return c
}

var sorted = Fn.new { |dom| (dom[0] > dom[1]) ? [dom[1], dom[0]] : dom }

var findLayouts = Fn.new { |tab, doms|
    var nrows = tab.count
    var ncols = tab[0].count
    var m = List.filled(nrows, null)
    for (i in 0...nrows) m[i] = List.filled(ncols, -1)
    var patterns = [ [m, [], []] ]
    var count = 0
    while (true) {
        var newpat = []
        for (pat in patterns) {
            var ut = pat[0]
            var ud = pat[1]
            var up = pat[2]
            var pos = null
            for (j in 0...ncols) {
                var breakOuter = false
                for (i in 0...nrows) {
                   if (ut[i][j] == -1) {
                       pos = [i, j]
                       breakOuter = true
                       break
                   }
                }
                if (breakOuter) break
            }
            if (!pos) continue
            var row = pos[0]
            var col = pos[1]
            if (row < nrows - 1 && ut[row+1][col] == -1) {
                var dom = sorted.call([tab[row][col], tab[row+1][col]])
                if (!containsDom.call(ud, dom[0], dom[1])) {
                    var newut = copyTab.call(ut)
                    newut[row][col] = tab[row][col]
                    newut[row+1][col] = tab[row+1][col]
                    newpat.add([newut, ud + [sorted.call( [tab[row][col], tab[row+1][col]])],
                        up + [row, col, row+1, col]])
                }
            }
            if (col < ncols - 1  && ut[row][col+1] == -1) {
                var dom = sorted.call([tab[row][col], tab[row][col+1]])
                if (!containsDom.call(ud, dom[0], dom[1])) {
                    var newut = copyTab.call(ut)
                    newut[row][col] = tab[row][col]
                    newut[row][col+1] = tab[row][col+1]
                    newpat.add([newut, ud + [sorted.call([tab[row][col], tab[row][col+1]])],
                        up + [row, col, row, col+1]])
                }
            }
        }
        if (newpat.count == 0) break
        patterns = newpat
        if (patterns[0][-1].count == doms.count) break
    }
    return patterns
}

var printLayout = Fn.new { |pattern|
    var tab = pattern[0]
    var dom = pattern[1]
    var pos = pattern[2]
    var bytes = List.filled(tab.count*2, null)
    for (i in 0...bytes.count) bytes[i] = List.filled(tab[0].count*2 - 1, " ")
    var idx = 0
    while (idx < pos.count-1) {
        var p = pos[idx..idx+3]
        var x1 = p[0]
        var y1 = p[1]
        var x2 = p[2]
        var y2 = p[3]
        var n1 = tab[x1][y1]
        var n2 = tab[x2][y2]
        bytes[x1*2][y1*2] = String.fromByte(48+n1)
        bytes[x2*2][y2*2] = String.fromByte(48+n2)
        if (x1 == x2) { // horizontal
            bytes[x1*2][y1*2 + 1] = "+"
        } else if (y1 == y2) { // vertical
            bytes[x1*2 + 1][y1*2] = "+"
        }
        idx = idx + 4
    }

    for (i in 0...bytes.count) {
        System.print(bytes[i].join())
    }
}

for (t in [tableau, tableau2]) {
    var start = System.clock
    var lays = findLayouts.call(t, dominoes)
    printLayout.call(lays[0])
    var lc = lays.count
    var pl = (lc > 1) ? "s" : ""
    var fo = (lc > 1) ? " (first one shown)" : ""
    System.print("%(lays.count) layout%(pl) found%(fo).")
    System.print("Took %(System.clock - start) seconds.\n")
}
Output:
0+5 1+3 2 2+3 1
        +     +
0 5+5 0 5 2+4 6
+     +        
4 3 0 3 6+6 2 0
  + +       + +
0 6 2 3+5 1 2 6
+         +    
1 1 3 0+0 2 4 5
  + +       + +
2 1 4 3+3 4 6 6
+         +    
6 4+5 1+5 4 1+4
               
1 layout found.
Took 0.014597 seconds.

6 4 2 2 0 6+5 0
+ + + + +     +
1 6 2 3 4 1+4 3
               
2 1 0 2 3+5 5 1
+ + + +     + +
1 3 5 0 5 6 1 0
        + +    
4 2 6 0 4 0 1+1
+ + + +        
4 4 2 0 5 3 6 3
        + + + +
6+6 5+2 5 3 3 4
               
2025 layouts found (first one shown).
Took 0.217176 seconds.

Extra credit (Cli)

Library: Wren-big
Library: Wren-fmt
import "./big" for BigInt
import "./fmt" for Fmt

var dominoTilingCount = Fn.new { |m, n|
    var prod = 1
    for (j in 1..(m/2).ceil) {
        for (k in 1..(n/2).ceil) {
            var cm = (Num.pi * (j / (m + 1))).cos
            var cn = (Num.pi * (k / (n + 1))).cos
            prod = prod * ((cm*cm + cn*cn) * 4)
        }
    }
    return prod.floor
}

var start  = System.clock
var arrang = dominoTilingCount.call(7, 8)
var perms  = BigInt.factorial(28)
var flips  = 2.pow(28)

Fmt.print("Arrangements ignoring values: $,i", arrang)
Fmt.print("Permutations of 28 dominos: $,i", perms)
Fmt.print("Permuted arrangements ignoring flipping dominos: $,i", perms * arrang)
Fmt.print("Possible flip configurations: $,i", flips)
Fmt.print("Possible permuted arrangements with flips: $,i", perms * flips * arrang)
System.print("\nTook %(System.clock - start) seconds.")
Output:
Arrangements ignoring values: 1,292,697
Permutations of 28 dominos: 304,888,344,611,713,860,501,504,000,000
Permuted arrangements ignoring flipping dominos: 394,128,248,414,528,672,328,712,716,288,000,000
Possible flip configurations: 268,435,456
Possible permuted arrangements with flips: 105,797,996,085,635,281,181,632,579,889,767,907,328,000,000

Took 0.00046 seconds.

Extra credit (Embedded)

Library: Wren-gmp

This is just to give what will probably be a rare outing to the Mpf class though (despite their usage in the Julia example) we don't need 'big floats' here, just 'big ints'. Slightly slower than the Wren-cli example as a result.

import "./gmp" for Mpz, Mpf
import "./fmt" for Fmt

var dominoTilingCount = Fn.new { |m, n|
    var prec = 128
    var prod = Mpf.from(1, prec)
    for (j in 1..(m/2).ceil) {
        for (k in 1..(n/2).ceil) {
            var cm = Mpf.pi(prec).mul(Mpf.from(j / (m + 1))).cos.square
            var cn = Mpf.pi(prec).mul(Mpf.from(k / (n + 1))).cos.square
            prod.mul(cm.add(cn).mul(4))
        }
    }
    return Mpz.from(prod.floor)
}

var start  = System.clock 
var arrang = dominoTilingCount.call(7, 8)
var perms  = Mpz.new().factorial(28)
var flips  = 2.pow(28)
 
Fmt.print("Arrangements ignoring values: $,i", arrang)
Fmt.print("Permutations of 28 dominos: $,i", perms)
Fmt.print("Permuted arrangements ignoring flipping dominos: $,i", perms * arrang)
Fmt.print("Possible flip configurations: $,i", flips)
Fmt.print("Possible permuted arrangements with flips: $,i", perms * flips * arrang)
System.print("\nTook %(System.clock - start) seconds.")
Output:
Arrangements ignoring values: 1,292,697
Permutations of 28 dominos: 304,888,344,611,713,860,501,504,000,000
Permuted arrangements ignoring flipping dominos: 394,128,248,414,528,672,328,712,716,288,000,000
Possible flip configurations: 268,435,456
Possible permuted arrangements with flips: 105,797,996,085,635,281,181,632,579,889,767,907,328,000,000

Took 0.00058 seconds.