638.sagews

#some random number examples
points = [[random(),random()] for i in range (2000)]
p = list_plot(points)
p.show(aspect_ratio=1)

#simulation of a double integral approximation
#domain is a circular disk radius = .5 centered at x=.5, y=.5


x,y = var('x y')
f(x,y) = sin(sqrt(log(x+y+1)))
sum = 0
j=0
for i in range(5000):
    x = random()
    y = random()
    if (x-.5)^2 + (y-.5)^2 <= .25:  #if in 2D circular domain
        j = j+1
        sum = sum+f(x,y)

#approximation to integral is 
#(area of domain) * (average height of f over j random points)
val = ((pi)/4)*sum/j
print j, n(val)

3901 0.567296852072786

#rolling a loaded die
count = 0

for i in range(10000):
    x = random()
    if x < .2:
        dice1=1
    elif x<.34:
        dice1=2
    elif x<.56:
        dice1=3
    elif x<.72:
        dice1=4
    elif x<.89:
        dice1=5
    else:
        dice1=6
    
    y = random()
    if y < .2:
        dice2=1
    elif y<.34:
        dice2=2
    elif y<.56:
        dice2=3
    elif y<.72:
        dice2=4
    elif y<.89:
        dice2=5
    else:
        dice2=6
    
    if (dice1 + dice2 == 7):
        count = count+1
        
print "count : ", count, "out of ",(i+1)

count : 1605 out of 10000