Will the buffer overflow
A router’s buffer holds N packets, and it currently has k. Every microsecond exactly one event happens: with probability p a packet arrives and there is one more, otherwise one packet is sent on and there is one fewer. Events in different microseconds are independent. If the buffer empties, the router goes to sleep; if it holds N packets, it overflows and losses begin.
With what probability does the buffer overflow before it empties?
Input
The first line holds the number of queries q (1 ≤ q ≤ 105). Then q lines N k a: the buffer capacity, the current number of packets and the arrival probability p = a / 109 (1 ≤ N ≤ 109, 0 ≤ k ≤ N, 1 ≤ a ≤ 109 − 1).
Output
For each query, the probability of overflow with an absolute error of at most 10−6.
Examples
Example 1
4 2 1 500000000 10 3 600000000 5 5 300000000 7 0 900000000
0.5000000000 0.7161223611 1.0000000000 0.0000000000
Example 2
4 1000000000 500000000 499999999 1000000000 10 500000001 1000000000 999999990 499999990 1000000000 123456789 500000000
0.1192029220 0.0000000407 0.9999996000 0.1234567890
Limits
- Time · Python
- 4000 ms
- Time · Go
- 1000 ms
- Memory
- 256 MB
- Judging
- 6 tests
Walkthrough
Open it step by step: each step gives away a little more. You can jump straight to the solution, but a problem solved with one hint stays with you longer.
1Hint 1
Call the answer Pk and analyse the first microsecond: Pk = p·Pk+1 + q·Pk−1, where q = 1 − p, and P0 = 0 and PN = 1. As in the article, solve this equation for the differences of neighbouring P.
2Hint 2
The formula will be short, but try computing the second example with it, a buffer of a billion packets. What happens to (q/p)N?
3The idea
Rewrite the equation as p(Pk+1 − Pk) = q(Pk − Pk−1): each next difference is r = q/p times the previous one. The differences form a geometric progression, their sum from 0 to N is 1, and for p ≠ ½
Pk = (rk − 1) / (rN − 1),
while for p = ½ all differences are equal and Pk = k/N, as in the article’s theorem. That is all the mathematics of the problem. The difficulty is in the arithmetic: with N = 109 and p slightly below ½ the number rN fits in no floating-point type, and with p slightly above ½ the numerator and denominator are differences of numbers close to one, and the subtraction eats all the significant digits.
4Walkthrough
Switch to exponents: r = ex, where x = ln(q/p) = ln(1 + (q − p)/p). Compute the difference q − p = (109 − 2a)/109 in integers, and the logarithm with the function log1p, which is accurate for small arguments. If x < 0, then Pk = expm1(kx) / expm1(Nx), where expm1(t) = et − 1 is computed without subtraction. If x > 0, divide numerator and denominator by rN: Pk = e(k − N)x · (1 − e−kx) / (1 − e−Nx); now all the powers are at most one. Handle separately a = 5·108 (the answer is k/N) and k = 0 or N.
In the second example a buffer of a billion packets is half full, and packets arrive with probability 0.499999999, two billionths less often than they leave. It seems an almost fair coin, and the answer should be about ½. But Nx ≈ 4: over a billion steps the tiny bias has time to take effect, and the buffer overflows with probability only 0.119.
4+Solution
5Complexity
O(1) per query.
6Common mistakes
- The “textbook” formula. (1 − rk) / (1 − rN) overflows for large N, and for p = ½ divides zero by zero, on the very first example.
- q = 1 − p in fractions. For p = ½ + 10−9 the difference q − p computed from fractions loses about half of its sixteen significant digits; it is safer to compute it in integers.
- Simulation. The mean duration here reaches 1017 microseconds, and you cannot wait for even one stroll to finish.