combinatorics

#!/usr/bin/env ruby
# frozen_string_literal: true

#
# Print out the basic combinatorial state counts.
#

require 'csv'

BOARD_SIZES = (2..4).to_a
MAX_EXPONENTS = (3..11).to_a

def number_with_comma(number)
  number.to_s.reverse.gsub(/(\d{3})(?=\d)/, '\\1,').reverse
end

#
# Basic Estimate: K^C + 1
# (for K = max exponent, C = number of cells)
#
CSV.open('data/blog/total_basic.csv', 'w') do |csv|
  csv << [''] + BOARD_SIZES

  MAX_EXPONENTS.each do |max_exponent|
    row = [2**max_exponent]
    BOARD_SIZES.each do |board_size|
      cells = board_size**2
      count = max_exponent**cells + 1
      row << number_with_comma(count)
    end
    csv << row
  end
end

#
# Slight Refinement: K^C + 1 - (K-2)^C - 2*C
#
CSV.open('data/blog/total.csv', 'w') do |csv|
  csv << [''] + BOARD_SIZES

  MAX_EXPONENTS.each do |max_exponent|
    row = [2**max_exponent]
    BOARD_SIZES.each do |board_size|
      cells = board_size**2
      count = max_exponent**cells + 1 - (max_exponent - 2)**cells - 2 * cells
      row << number_with_comma(count)
    end
    csv << row
  end
end

#
# Layer sums using compositions.
#

# Check whether n is a power of 2 using the bitwise and trick.
def power_of_2?(n)
  n & (n - 1) == 0
end

# Implement the recurrence from
# Chinn, Niederhausen (unknown). Compositions into Powers of 2
# Add one new thing: let the caller manipulate the bounds on the sum index.
def count_compositions(n_max, k_max, i_min = 0, i_max = nil)
  ns = (1..n_max).to_a
  ks = (2..k_max).to_a
  table = ns.map do |n|
    power_of_2?(n) && 2**i_min <= n && (i_max.nil? || n <= 2**i_max) ? [1] : [0]
  end
  ks.each do |k|
    ns.each do |n|
      next if k > n
      i_max_n = i_max || Math.log2(n).floor.to_i
      table[n - 1][k - 1] = (i_min..i_max_n).map do |i|
        n_prev = n - 2**i
        k_prev = k - 1
        k_prev > n_prev ? 0 : table[n_prev - 1][k_prev - 1]
      end.sum
    end
  end
  table
end

#
# Enumerator over all numbers with the given (mixed) radices.
#
def radix_sequence(size, radix)
  array = [0] * size
  Enumerator.new do |y|
    loop do
      y << array.dup
      break unless spin_array(array, radix)
    end
  end
end

def spin_array(array, radix)
  (0...(array.size)).each do |i|
    array[i] += 1
    return true if array[i] <= radix
    array[i] = 0
  end
  false
end

# Explicitly enumerate and check compositions of powers of two.
# n = the sum, k = the number of parts
def check_compositions(n, k, i_min = 0, i_max = nil)
  i_max ||= Math.log2(n + 1).ceil.to_i
  span = i_max - i_min
  candidates = radix_sequence(k, span)
  count = 0
  candidates.each do |exponents|
    # p(exponents.map { |i| 2**i }) if exponents.map { |i| 2**i }.sum == n
    count += 1 if exponents.map { |i| 2**(i_min + i) }.sum == n
  end
  count
end

# Check the compositions table against the explicitly calculated values.
def check_compositions_table(n_max, k_max, i_min = 0, i_max = nil)
  table = count_compositions(n_max, k_max, i_min, i_max)
  p table
  (1..n_max).each do |n|
    (1..k_max).each do |k|
      next if k > n
      check = check_compositions(n, k, i_min, i_max)
      expected = table[n - 1][k - 1]
      if check == expected
        puts "OK #{n}, #{k}: #{expected}"
      else
        puts "mismatch: #{n}, #{k}: table #{expected} != check #{check}"
      end
    end
  end
end
# check_compositions_table(10, 10)
# check_compositions_table(10, 10, 0, 2)
# check_compositions_table(10, 10, 1, 2)

def factorial(n)
  (1..n).inject(:*) || 1
end

def choose(n, k)
  factorial(n) / (factorial(k) * factorial(n - k))
end

#
# Estimate with only the max exponent constraint (no tile constraints).
#
def estimate_layer_counts(board_size, max_exponent)
  num_cells = board_size**2
  i_max = max_exponent - 2
  max_sum = num_cells * 2**i_max
  table = count_compositions(max_sum, num_cells, 0, i_max)
  table.map.with_index do |row, index|
    sum = (index + 1) * 2
    counts = row.map.with_index do |count_k, i|
      k = i + 1
      count_k * choose(num_cells, k)
    end
    [sum, counts.sum]
  end
end
# estimate_layer_counts(2, 5)

CSV.open('data/blog/layers_basic.csv', 'w') do |layer_csv|
  layer_csv << %w[board_size max_exponent layer_sum num_states]
  CSV.open('data/blog/layers_basic_total.csv', 'w') do |total_csv|
    total_csv << %w[board_size max_exponent total_states]
    BOARD_SIZES.each do |board_size|
      MAX_EXPONENTS.each do |max_exponent|
        counts = estimate_layer_counts(board_size, max_exponent)
        counts.each do |layer_sum, count|
          layer_csv << [
            board_size,
            max_exponent,
            layer_sum,
            count
          ]
        end
        # This method doesn't count the zero state or the special win state,
        # which the other 'basic' method does, so add two for comparability.
        total = counts.map { |layer| layer[1] }.sum + 2
        total_csv << [board_size, max_exponent, total]
      end
    end
  end
end

#
# Estimate with max exponent and tile constraints. The approach is to calculate
# the basic estimate but (1) ignore k = 1 (only one tile) and (2) subtract off
# compositions that do not include a power of 2 larger than 4 (must have at
# least one two or four tile).
#
def count_row(row, num_cells)
  row.map.with_index do |count_k, i|
    k = i + 1
    next 0 if k == 1
    count_k * choose(num_cells, k)
  end.sum
end

def estimate_layer_counts_with_restrictions(board_size, max_exponent)
  num_cells = board_size**2
  i_max = max_exponent - 1
  max_sum = num_cells * 2**i_max
  table = count_compositions(max_sum, num_cells, 1, i_max)
  table_4 = count_compositions(max_sum, num_cells, 3, i_max)
  table.zip(table_4).map.with_index do |(row, row_4), index|
    sum = index + 1
    next if sum.odd?
    count = count_row(row, num_cells) - count_row(row_4, num_cells)
    [sum, count]
  end.compact.drop(1)
end

CSV.open('data/blog/layers.csv', 'w') do |layer_csv|
  layer_csv << %w[board_size max_exponent layer_sum num_states]
  CSV.open('data/blog/layers_total.csv', 'w') do |total_csv|
    total_csv << %w[board_size max_exponent total_states]
    BOARD_SIZES.each do |board_size|
      MAX_EXPONENTS.each do |max_exponent|
        counts = estimate_layer_counts_with_restrictions(
          board_size, max_exponent
        )
        counts.each do |layer_sum, count|
          layer_csv << [
            board_size,
            max_exponent,
            layer_sum,
            count
          ]
        end
        # This method doesn't count the win state, which the other refined
        # method does, so add one for comparability.
        total = counts.map { |layer| layer[1] }.sum + 1
        total_csv << [board_size, max_exponent, total]
      end
    end
  end
end
# p estimate_layer_counts_with_24_tile(2, 5)