Create a gist now

Instantly share code, notes, and snippets.

anonymous /
Created Dec 17, 2013

What would you like to do?
For Reddit Math Problem
import random;
def magsq(v):
return sum(x**2 for x in v)
def sub(a,b):
return (x[0]-x[1] for x in zip(a,b))
def in_sph(r0,r,pt):
return magsq(sub(pt,r0)) <= r**2
def rand_pt(domain):
return (random.uniform(domain[0], domain[1]),
random.uniform(domain[2], domain[3]),
random.uniform(domain[4], domain[5]))
def dom_vol(domain):
return (domain[1]-domain[0])*(domain[3]-domain[2])*(domain[5]-domain[4])
def inside(pt):
return in_sph((0,0,0),1.0,pt) or \
in_sph((1,0,0),1.0,pt) or \
in_sph((0,1,0),1.0,pt) or \
in_sph((1,1,0),1.0,pt) or \
in_sph((0,0,1),1.0,pt) or \
in_sph((1,0,1),1.0,pt) or \
in_sph((0,1,1),1.0,pt) or \
Npts = 0
Ninside = 0
D = (-1, 2, -1, 2, -1, 2)
for N in xrange(1000000):
if inside(rand_pt(D)): Ninside += 1
Npts += 1
print 1.0*dom_vol(D)*Ninside/Npts
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment