! Calculate pi using the Monte Carlo method
! estimate the area of circle of radius 1
! Program to accompany "Computational Physics" by N. Giordano
! Copyright Prentice Hall 1997
program pi 
   option nolet
   randomize
   n = 100
   do
      call calculate(n,value)
      print n, 4 * value      ! multiply value by 4 to get the area of the entire circle
      n = n * 3
      if key input then exit do
   loop
end
! calculate area of circle with radius = 1    
! consider only portion in first quadrant
! n_total = number of points (darts) to generate   value = integral
sub calculate(n_total,value)
   n_sum = 0
   for i = 1 to n_total
      x = rnd     ! generate the coordinates of a random point
      y = rnd
      if y <= sqr(1 - x^2) then n_sum = n_sum + 1    ! check if point is
                                                     ! the circle
   next i
   value = n_sum / n_total
end sub