#!/usr/bin/env python
#
# I wrote this to calculate the minimum number of dice rolls required in order to have greater than a 50% chance
# of hitting at least one roll per dice value.
import gmpy
import sys
from multiprocessing import Process, Lock, Value

coefficients = list()

def main_loop(trial, found, lock, sides):
	# Keep going until we've found the break-even point
	while not found.value:
		# Fetch the next value
		lock.acquire()
		trial.value = trial.value + 1
		rnd = trial.value
		lock.release()

		terms = list()
		# Calculate the value of each term in the polynomial
		for x in xrange(len(coefficients)):
			terms.append( pow(sides - x, rnd) * coefficients[x] * pow(-1, x) )
		#print terms
		# Sum each term of the polynomial
		val = sum(terms)
		# Determine the denominator
		den = pow(sides, rnd)
		# Calculate the quotient between
		probability = val / gmpy.mpf(den)
		# Report
		if rnd % 100 == 0:
			print "Probability after %d rolls is %s" % (rnd, gmpy.fdigits(probability, 10, 10))
		if probability > 0.5:
			print "Probability after %d rolls is %s" % (rnd, gmpy.fdigits(probability, 10, 10))
			found.value = True
			
def main():
	if len(sys.argv) != 2:
		print "Please provide exactly one number"
		sys.exit(-1)
	# The number of sides on the die in question
	sides = int(sys.argv[1])
	
	# First we generate the coefficients of the equation
	n = gmpy.mpz(0)
	while n < sides:
		coefficients.append( gmpy.bincoef(sides, n) )
		n += 1
	
	# Next we start at the number of sides and increment, making a calculation at each row.
	trial = Value('i', sides)
	found = Value('b', False)
	lock  = Lock()
	
	for num in range(4):
		Process(target=main_loop, args=(trial, found, lock, sides)).start()

if __name__ == '__main__':
	main()