IOAI Training: Learn AI by solving problems

The problem in brief

A space station is surrounded by $x$ meteorites, and exactly $g$ of them carry gold. Every minute one meteorite, chosen among those still out there, is pulled in, so the meteorites arrive in a uniformly random order and none comes back. You read $x$ and $g$ ( $1 \le g \le x < 100\,000$ ) and print the expected number of minutes until the first gold meteorite arrives. The whole task is to find a closed formula: the program itself is one line.

How you are scored

This is a classic stdin/stdout judge task (1 second, 64 MB): your printed number is compared with the reference value on hidden tests, presumably with a small floating-point tolerance. The only published test is 1 1, answer 1. Since the official tests are hidden, we validate the formula against two independent exact computations on all 820 pairs with $1 \le g \le x \le 40$ , plus a Monte Carlo simulation.

Answer:

\boxed{\;\mathbb{E}[T] = \frac{x+1}{g+1}\;}

First look: what is random here?

Let $T$ be the minute at which the first gold meteorite arrives. $T$ is at least 1 (the very first pull may be gold) and at most $x-g+1$ (in the worst order, all $x-g$ empty rocks come first). Two quick sanity cases:

$g = x$ : every meteorite is gold, $T = 1$ . The formula gives $(x+1)/(x+1) = 1$ .
$g = 1$ : the single gold meteorite is equally likely to be at any of the $x$ positions, so $\mathbb{E}[T] = (1+2+\dots+x)/x = (x+1)/2$ . The formula agrees.

The tempting wrong answer is $x/g$ . That is the expectation of a geometric law, the waiting time when each pull is gold with probability $g/x$ independently, which would mean throwing meteorites back after looking at them. Without replacement, every empty rock you pull makes the next pull more likely to be gold, so the wait is shorter. The right panel of the figure shows the gap: for $x = 20, g = 1$ it is $20$ against the true $10.5$ .

Left: the law of T for x = 20 and g = 1, 3, 8, with its mean as a dashed line. For g = 1 the law is flat (uniform position), for larger g it decays fast. Right: E[T] as a function of g for x = 20: the formula (x+1)/(g+1) passes exactly through the simulated points, while x/g, the answer with replacement, is too large for small g.

The key idea

Look at the gaps, not at the minutes. Lay the $x$ meteorites out in their arrival order. The $g$ gold ones cut the line of $x-g$ empty rocks into $g+1$ gaps:

\underbrace{\circ\circ\circ}_{G_0}\ \bullet\ \underbrace{\circ}_{G_1}\ \bullet\ \underbrace{\ }_{G_2}\ \bullet\ \underbrace{\circ\circ}_{G_3}

(here $x = 9$ , $g = 3$ , gaps of sizes $3, 1, 0, 2$ ). The first gold meteorite arrives right after the first gap, so

T = G_0 + 1 .

Now the symmetry argument. Take one particular empty rock. Together with the $g$ gold meteorites it forms a group of $g+1$ objects whose relative order is uniformly random, so the rock is equally likely to be in any of the $g+1$ slots "before the 1st gold", "between the 1st and 2nd", ..., "after the last". In particular it falls into gap $G_0$ with probability $1/(g+1)$ . By linearity of expectation (no independence needed), summing over the $x-g$ empty rocks:

\mathbb{E}[G_0] = \frac{x-g}{g+1}, \qquad \mathbb{E}[T] = \frac{x-g}{g+1} + 1 = \frac{x+1}{g+1}.

All gaps have the same expected size, so the average gap is the total $x-g$ shared equally among $g+1$ gaps. This "indicator per object + linearity" trick is the single most useful tool for expectation problems: it replaces a messy distribution by one probability per object.

A second derivation, by brute force

If you do not trust symmetry, compute the law of $T$ directly. $T = k$ means that the first $k-1$ pulls are empty, the $k$ -th is gold, and the remaining $g-1$ gold ones are somewhere among the last $x-k$ positions. Counting sets of gold positions:

P(T = k) = \frac{\binom{x-k}{g-1}}{\binom{x}{g}}, \qquad k = 1,\dots, x-g+1 .

Then $\mathbb{E}[T] = \sum_k k\,P(T=k)$ . Using the tail-sum formula $\mathbb{E}[T] = \sum_{k \ge 0} P(T > k)$ with $P(T > k) = \binom{x-g}{k}/\binom{x}{k}$ (the first $k$ pulls are all empty) gives a third route. solution.py evaluates both sums with exact fractions and checks they equal $(x+1)/(g+1)$ in all 820 cases with $x \le 40$ . (Proving the identity algebraically uses the hockey-stick identity $\sum_{m} \binom{m}{g} = \binom{x+1}{g+1}$ ; the gap argument gets there without any algebra.)

The full solution

x, g = map(int, input().split())
print(f"{(x + 1) / (g + 1):.10g}")

It runs in constant time, far inside the limit even for $x = 99\,999$ . Two practical details:

Use floating-point division and print enough digits: $(x+1)/(g+1)$ is usually not an integer (for 10 3 the answer is 2.75).
%.10g prints 1 for the sample, exactly as expected, and up to 10 significant digits otherwise. Simulating the process instead (random shuffles until the first gold) would be both too slow for $x \approx 10^5$ and not precise enough for a numeric checker.

Results

approach	value for `x=10, g=3`	verdict
closed form $(x+1)/(g+1)$	2.75	equals both exact sums on all 820 pairs with $x \le 40$
exact sum $\sum_k k\,P(T=k)$ (fractions)	11/4	agrees
tail sum $\sum_k P(T>k)$ (fractions)	11/4	agrees
Monte Carlo, 200 000 shuffles	2.7507	agrees within noise
geometric law $x/g$ (with replacement)	3.333	wrong: ignores that rocks are not put back
$(x-g)/g$ or $x/(g+1)$ (off-by-one variants)	2.333, 2.5	wrong, fails the sample or `g = x`

What to remember

Check the extreme cases first ( $g = x$ , $g = 1$ ): they kill wrong formulas in seconds and suggest the right one.
Without replacement is not geometric. Draws without replacement are dependent; the answer $(x+1)/(g+1)$ is always smaller than the independent-draws answer $x/g$ .
Indicators plus linearity of expectation: write the quantity as a sum of 0/1 variables, one per object, and add their probabilities. No independence is required.
Symmetry of random orders: in a uniformly random permutation, any object is equally likely to sit in any of the $g+1$ gaps defined by $g$ others.
Verify a formula three ways (exact sum, independent identity, simulation) before submitting: it is cheap and catches off-by-one mistakes.

"""FAIO 2024 qualification, Problem H "Meteorites with Gold".

x meteorites, g of them with gold, pulled one per minute in a uniformly random
order (without replacement). Expected minute of the first gold one?
Answer: (x + 1) / (g + 1).

Run as a judge program:   echo "10 3" | python solution.py
Run the checks + figure:  python solution.py --check
"""
import sys
from fractions import Fraction
from math import comb


def answer(x: int, g: int) -> float:
    # x - g empty rocks are split by the g gold ones into g + 1 gaps of equal
    # expected size (x - g) / (g + 1); the first gold arrives after the first gap.
    return (x + 1) / (g + 1)


def exact_by_definition(x: int, g: int) -> Fraction:
    """E[T] = sum_k k * P(T = k), with P(T = k) = C(x-k, g-1) / C(x, g)
    (first k-1 draws empty, draw k gold, the other g-1 golds among the last x-k)."""
    total = comb(x, g)
    return sum(Fraction(k * comb(x - k, g - 1), total) for k in range(1, x - g + 2))


def exact_by_tail_sum(x: int, g: int) -> Fraction:
    """E[T] = sum_{k>=0} P(T > k), P(T > k) = C(x-g, k) / C(x, k) (first k all empty)."""
    return sum(Fraction(comb(x - g, k), comb(x, k)) for k in range(0, x - g + 1))


def monte_carlo(x, g, trials=200_000, seed=0):
    import numpy as np
    rng = np.random.default_rng(seed)
    # positions of gold = g distinct random positions among 1..x; T = the smallest
    keys = rng.random((trials, x))
    gold_pos = np.argsort(keys, axis=1)[:, :g]          # random g-subset of positions
    return (gold_pos.min(axis=1) + 1).mean()


def check():
    # 1) the closed form equals two independent exact computations
    for x in range(1, 41):
        for g in range(1, x + 1):
            f = Fraction(x + 1, g + 1)
            assert exact_by_definition(x, g) == f == exact_by_tail_sum(x, g), (x, g)
    print("closed form = exact sums for all 1 <= g <= x <= 40 (820 cases)")
    # 2) a simulation agrees
    for x, g in [(1, 1), (5, 1), (10, 3), (30, 7)]:
        print(f"x={x:2d} g={g:2d}  formula={answer(x, g):.4f}  monte-carlo={monte_carlo(x, g):.4f}")
    # 3) extreme input: fast, printed with enough digits
    print("x=99999 g=1 ->", f"{answer(99999, 1):.10f}")
    figure()


def figure():
    import numpy as np
    import matplotlib
    matplotlib.use("Agg")
    import matplotlib.pyplot as plt
    x = 20
    fig, axes = plt.subplots(1, 2, figsize=(9.5, 3.8), dpi=150)
    ax = axes[0]
    for g, col in [(1, "#4c72b0"), (3, "#dd8452"), (8, "#55a868")]:
        ks = np.arange(1, x - g + 2)
        p = [comb(x - k, g - 1) / comb(x, g) for k in ks]
        ax.plot(ks, p, "o-", ms=3.5, color=col, label=f"g = {g}")
        ax.axvline((x + 1) / (g + 1), color=col, ls="--", lw=1)
    ax.set_xlabel("minute k of the first gold meteorite")
    ax.set_ylabel("P(T = k)")
    ax.set_title(f"x = {x}: distribution of T (dashed: mean)")
    ax.legend(frameon=False)
    ax = axes[1]
    gs = np.arange(1, x + 1)
    ax.plot(gs, (x + 1) / (gs + 1), "k-", label="(x+1)/(g+1)")
    ax.plot(gs, x / gs, ":", color="#c44e52", label="x/g (with replacement)")
    sims = [monte_carlo(x, int(g), trials=20_000, seed=int(g)) for g in gs]
    ax.plot(gs, sims, "o", ms=4, color="#4c72b0", label="simulation")
    ax.set_xlabel("number of gold meteorites g")
    ax.set_ylabel("expected minutes E[T]")
    ax.set_title(f"x = {x}")
    ax.legend(frameon=False)
    for a in axes:
        a.spines[["top", "right"]].set_visible(False)
    fig.tight_layout()
    fig.savefig("figures/distribution.png", facecolor="white")
    print("saved figures/distribution.png")


if __name__ == "__main__":
    if "--check" in sys.argv:
        check()
    else:
        x, g = map(int, sys.stdin.read().split()[:2])
        print(f"{answer(x, g):.10g}")

H. Meteorites with Gold