Bridge hand distribution

bridgehand.jl

using StatsBase     # sample function
using Combinatorics # combinations function
using Plots         # plotting infrastructure

# Code calculating distribution of HCP in hands of a pair of players in bridge
# Ace is worth 4 points, King 3, Queen 2, Jack 1, all other cards (36 in total) count as 0

# Get the distribution using simulation
function sim_hand()
    n = 1_000_000
    deck = [repmat(1:4, 4); fill(0, 36)] # whole deck in HCP
    distr = fill(0, 41)
    for i in 1:n
        hcp = sum(sample(deck, 26, replace=false))
        distr[hcp + 1] += 1
    end
    distr / n
end

# Calculate the distribution exactly
function exact_hand()
    high = repmat(1:4, 4) # cards that have non zero value in HCP
    distr = fill(0, 41)
    distr[1] = binomial(16, 0) * binomial(36, 26)
    for i in 1:16
        low = binomial(36, 26 - i) # ways to sample cards that have zero value in HCP
        for c in combinations(high, i)
            hcp = sum(c)
            distr[hcp + 1] += low
        end
    end
    distr / sum(distr)
end

# Run the computations and plot the results
# semicolons suppress printing of results in REPL if you copy-paste the code
srand(1);
hcp = 0:40;
sim = sim_hand();
exact = exact_hand();
plotly();
plt = plot(hcp, exact, label="exact", xlab="HCP", ylab="probability");
scatter!(plt, hcp, sim, label="simulated")