Is this calculation of mining probability using the Bernoulli trial formula accurate?

Question

I'm trying to understand the formulas for successfully mining a block (without a pool) based on hashrate over time (i.e. to chart the probability of finding a block at various historical points compared to network hashrate at the time).

It seems like this post is exactly the formula I'm looking for, but it doesn't actually add up. I'm awful at math so I can't tell if I'm missing a step or if the post is wrong. If there's a better/more accurate way of figuring this out, please share, but I'm specifically interested in learning if this is wrong or not and why.

It says the probability of being the one to successfully mine a block with 13.5 TH/s and a network hashrate of 11,000,000 TH/s = 0.00000122

Then it says:

Now, in a month, there are N=6*24*30 mined blocks. To find the probability of winning one award in a month, we use the Bernoulli trial formula

N*p1*(1-p1)^(N-1) = 0.0026

Now when I put that in Excel with N=4320, I get .005243, I only come up with the same answer when I divide N by two for N=2160, but I can't figure out any reason for that to be the case. Am I missing something or is this just a typo? Is this formula even a good way of approximating this?

Answer 1

The formula is a reasonable approximation. It looks to me like your calculation is right, and the blog author has made an error.

There is a minor inaccuracy in that the number of blocks per month is not fixed at 6*24*30, but is subject to random variation. So a slightly improved solution is to work out the average number of terahashes needed for one block (h = 11000000 * 600, assuming that the network hash rate is stable and the block difficulty has adjusted accordingly). Then your probability of winning a block on each of your hashes is q1 = 1 / h, and you will perform M = 12.5 * 60 * 60 * 24 * 30 hashes in one month. If you treat the number of blocks you win as a Poisson distribution, the probability of winning one block in a month is given by λ * exp(-λ) where λ = M*q1. This yields the slightly larger number 0.005273.