# OK to post homework
# Robert Dougherty-Bliss
# 2021-04-11
# Assignment 20

# 1.

read "C20.txt":

AvHD := proc(S)
    add(Nei(S, s), s in S) / nops(S):
end:

# 2.

AvAv := proc(dim, J, K)
    local B:
    B := Box(dim):
    evalf(add(AvHD(randcomb(B, J)), k=1..K) / K):
end:

# AvAv([2,2,2,2,2,2,2,2],63,100);

# 3.

# Let a(n, k) be the average number of neighbors in a subset of size k of
# vertices in Box([2 $ n]).

# My first observation is that a(n, k) is proportional to k. Let's see some
# exact values.

a := proc(n, k)
    local B:
    B := Box([2 $ n]):
    add(AvHD(S), S in choose(B, k)) / binomial(2^n, k):
end:

an := n -> [seq(a(n, k), k=1..2^n)]:

# Denominators are 2^n - 1.
# Let's normalize by that.

bn := n -> [seq((2^n - 1) * a(n, k), k=1..2^n)]:

# Now it's clear! The numerators are n * (k - 1). So:

#   a(n, k) = n(k - 1) / (2^n - 1).

# Or at least that's our conjecture.

# Proof?

(*

Let's take a look at that sum.

    sum(AvHD(S), S in choose(B, k))
        = sum(sum(number of neighbors of x in S, x in S), S in choose(B, k)) / k.

Instead of summing over subsets, let's sum over elements.

Each element will contribute j to the above sum for each subset that it is in
which contains precisely j of its neighbors. So our sum is

    sum(j * number of k-subsets of B that contain x and exactly j neighbors of x, x in B) / k.

How many k-subsets of B contain x and exactly j of its neighbors? There are n
neighbors total, so we can choose them in binomial(n, j) ways. We must contain
x, and that leaves us with 2^n - n - 1 elements to fill out the remaining k - j
- 1 elements. So there are

    binomial(n, j) * binomial(2^n - n - 1, k - j - 1)

such subsets. Note that this is independent of x, so the total sum is

    2^n / k * sum(j * binomial(n, j) * binomial(2^n - n - 1, k - j - 1), j).

This sum looks intimidating, but it is actually just Vandermode--Chu with some
extra steps. Absorb the j into the first term, shift the summation index up by
1, then apply V-C. (Maple can do this for you.) After some simplifications, you
get

    n 2^n binomial(2^n - 2, k - 2) / k.

Dividing by binomial(2^n, k) (to take the average) and simplifying gives the
result that we're looking for.

Nice!

This same argument tells us how to compute other things. Say that we want to
compute the average value of f(number of neighbors) over k-subsets. (The
original problem is f(x) = x). Then we just need to evaluate the sums

    sum(f(j) * binomial(n, j) * binomial(2^n - n - 1, k - j - 1), j).

Note that this is NOT a sum of hypergeometric terms in (j, n), but it IS
hypergeometric in (j, k). (EKHAD says nice things if you pass in f(j) =
binomial(j, c).)

*)

# 4.

# I'm not sure I understand this problem. What value are we sampling? If it's
# MaxHD, then we'll obviously see a max of n at some point. If it's MinHD, then
# we'll see a min of 0 at some point. Thus it makes the most sense to sample
# AvHD and look for...the min, I suppose.

# After all, if the average seems to be > sqrt(n) for all subsets S, then,
# probably, every subset of the appropriate size will have a vertex with at
# least sqrt(n) neighbors.

SimuHao := proc(dim, J, K)
    local B, champ, min, S, k:
    B := Box(dim):
    champ := FAIL:
    min := infinity:

    for k from 1 to K do
        S := randcomb(B, J):
        if AvHD(S) < min then
            min := AvHD(S):
            champ := S:
        fi:
    od:

    min, champ:
end:

SimuHao([2 $ 10], 513, 10)[1];