import random, math

class calcPi:
    def __init__(self):
        self.times = pow(10,6)
        self.i = 0
        self.isnot = 0

    def IsOnCircle(self,x,y):
        if math.sqrt(x**2+y**2) < 1:
            return True
        else:
            return False
    
    def run(self):
        for x in range(self.times):
            x,y = random.random(),random.random()
            if self.IsOnCircle(x,y):
                self.i+=1
            else:
                self.isnot+=1
    
    def getResults(self):
        return (float(self.i), float(self.isnot))

    def getPi(self):
        self.run()
        r = self.getResults()
        return r[0]/(r[0]+r[1])*4

if __name__ == '__main__':
    print calcPi().getPi()

#!/usr/local/bin/ruby -w

require 'benchmark' 

class Integer    
   def prime?         # cf. http://snippets.dzone.com/posts/show/4636
     n = self.abs()
     return true if n == 2
     return false if n == 1 || n & 1 == 0
     d = n-1
     d >>= 1 while d & 1 == 0
     20.times do                               # 20 = k from above
       a = rand(n-2) + 1
       t = d
       y = ModMath.pow(a,t,n)                  # implemented below
       while t != n-1 && y != 1 && y != n-1
         y = (y * y) % n
         t <<= 1
       end
       return false if y != n-1 && t & 1 == 0
     end
     return true
   end
end
 
module ModMath
   def ModMath.pow(base, power, mod)
     result = 1
     while power > 0
       result = (result * base) % mod if power & 1 == 1
       base = (base * base) % mod
       power >>= 1;
     end
     result
   end
end



class Integer

   def primes   # cf. http://snippets.dzone.com/posts/show/3734

      sieve = []
      3.step(self, 2) { |i| sieve[i] = i }
      sieve[1] = nil
      sieve[2] = 2

      3.step(Math.sqrt(self).floor, 2) do |i| 
         next unless sieve[i]
         (i*i).step(self, i) do |j|
            sieve[j] = nil
         end
      end

      sieve.compact!

   end 


   def primes2       # cf. http://snippets.dzone.com/posts/show/3734

      primes = [nil, nil].concat((2..self).to_a)

      (2 .. Math.sqrt(self)).each do |i|
         next unless primes[i]
            (i*i).step(self, i) do |j|
               primes[j] = nil
            end
      end
	
      primes.compact!

   end

end



class Integer

   @primes_cache ||= 100.primes
   #@primes_cache ||= [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31]
   class << self; attr_accessor :primes_cache; end

   @nthprime_limit ||= 5_000_000
   class << self; attr_reader :nthprime_limit; end
   #class << self; attr_accessor :nthprime_limit; end


   def sieve_size
      num = self.to_i.abs
      logn = Math.log(num)
      if num < 15985     # cf. http://primes.utm.edu/howmany.shtml
         ( num * (logn + Math.log(logn) - 1 + 1.8 * Math.log(logn) / logn) ).floor
      else
         ( num * (logn + Math.log(logn) - 0.9427) ).floor
      end
   end


   def nthprime

      num = self.to_i.abs

      # set a limit; cf. http://primes.utm.edu/nthprime/: The 5,000,000th prime is 86,028,121.
      raise "#{num}: number too big for Integer#nthprime" if num > Integer.nthprime_limit      

      primes_cache_size = Integer.primes_cache.size

      if num > primes_cache_size

         # optional; cf. Integer#nthprime_add and Integer#nthprime_add_mr below
         if num >= 50_000 && num < 100_000 && (num - primes_cache_size) < 5000 then return num.nthprime_add end
         if num >= 100_000 && num < 200_000 && (num - primes_cache_size) < 15000 then return num.nthprime_add end
         if num >= 200_000 && (num - primes_cache_size) < 35000 then return num.nthprime_add end
         
         logn = Math.log(num)

         if num < 15985       # cf. http://primes.utm.edu/howmany.shtml
            limit = ( num * (logn + Math.log(logn) - 1 + 1.8 * Math.log(logn) / logn) ).floor
            Integer.primes_cache = limit.primes
         else
            limit = ( num * (logn + Math.log(logn) - 0.9427) ).floor
            Integer.primes_cache = limit.primes
         end

=begin
         elsif num >= 15985 && num <= 1_000_000
            limit = ( num * (logn + Math.log(logn) - 0.9427) ).floor
            Integer.primes_cache = limit.primes
         elsif num > 1_000_000
            Integer.primes_cache = 20_000_000.primes
         end
=end

      else
         return Integer.primes_cache.at(num-1)
      end


      if num <= Integer.primes_cache.size
         return Integer.primes_cache.at(num-1)
      end

      if Integer.primes_cache.size > 2_500_000
         num.nthprime_add_mr 
      else
         num.nthprime_add   # faster for a (relatively) small prime cache
      end

   end


   def nthprime_add       #  add primes to Integer.primes_cache up to the nth prime

      num = self.to_i.abs

      if num <= Integer.primes_cache.size
         return Integer.primes_cache.at(num-1)
      end

      last_prime = Integer.primes_cache.last
      last_prime_divmod = last_prime.divmod(6)

      if last_prime_divmod.last == 1
         i = last_prime_divmod.first
         Integer.primes_cache.pop      # avoid a duplicate prime
      else
         i = last_prime_divmod.first + 1
      end


      while Integer.primes_cache.size < num

         n1 = 6*i+1       # cf. http://betterexplained.com/articles/another-look-at-prime-numbers/ and
         n2 = n1+4        # http://everything2.com/index.pl?node_id=1176369
         i += 1

         [n1, n2].each do |p| 

            next if p % 5 == 0 || p % 7 == 0

            next_num_sqrt = Math.sqrt(p).floor

            Integer.primes_cache.each do |d| 
               if d > next_num_sqrt 
                  Integer.primes_cache << p
                  #print "\r\e[0K#{Integer.primes_cache.size}"
                  #$stdout.flush
                  break
               elsif p % d == 0
                  break
               end  
            end
         end   # each  

      end  # while

      return Integer.primes_cache.at(num-1)

   end


   def nthprime_add_mr       #  add Miller-Rabin primes to Integer.primes_cache up to the nth prime

      num = self.to_i.abs

      if num <= Integer.primes_cache.size
         return Integer.primes_cache.at(num-1)
      end

      last_prime = Integer.primes_cache.last
      last_prime_divmod = last_prime.divmod(6)

      if last_prime_divmod.last == 1
         i = last_prime_divmod.first
         Integer.primes_cache.pop
      else
         i = last_prime_divmod.first + 1
      end


      while Integer.primes_cache.size < num

         search_next_prime = true

         while search_next_prime

            n1 = 6*i+1       
            n2 = n1+4        
            i += 1

            [n1, n2].each do |p| 

               next if p % 5 == 0 || p % 7 == 0

               if p.prime?
                  Integer.primes_cache << p
                  search_next_prime = false
                  #print "\r\e[0K#{Integer.primes_cache.size}"
                  #$stdout.flush
               end
            end
         end

      end  #  first while loop

      return Integer.primes_cache.at(num-1)

   end


#------------------------------------------- some additional prime methods


   def next_primes_in_cache    # next primes in Integer.primes_cache

      search_next_primes = self.to_i.abs

      last_prime = Integer.primes_cache.last
      last_prime_divmod = last_prime.divmod(6)

      if last_prime_divmod.last == 1
         i = last_prime_divmod.first
         Integer.primes_cache.pop
      else
         i = last_prime_divmod.first + 1
      end

      while search_next_primes > 0

         n1 = 6*i+1       
         n2 = n1+4        
         i += 1

         [n1, n2].each do |p| 

            next if p % 5 == 0 || p % 7 == 0
            next_num_sqrt = Math.sqrt(p).floor

            Integer.primes_cache.each do |d| 
               if d > next_num_sqrt 
                  Integer.primes_cache << p
                  search_next_primes -= 1
                  #print "\r\e[0K#{Integer.primes_cache.size}"
                  #$stdout.flush
                  break
               elsif p % d == 0
                  break
               end  
            end
         end 
      end   # while

      Integer.primes_cache.last(self.to_i.abs)

   end


   def next_mr_prime     # next Miller-Rabin prime
     
      last_num = self.to_i.abs
      last_num_divmod = last_num.divmod(6)

      if last_num_divmod.last == 1
         i = last_num_divmod.first
      else
         i = last_num_divmod.first + 1
      end

      next_prime = nil
      search_next_prime = true

      while search_next_prime

         n1 = 6*i+1       
         n2 = n1+4        
         i += 1

         [n1, n2].each do |p| 

            next if p % 5 == 0 || p % 7 == 0

            if p > last_num && p.prime?
               next_prime = p
               search_next_prime = false 
               break
            end
         end
      end

         next_prime

   end



   def next_mr_primes(n)     # next Miller-Rabin primes

      last_num = self.to_i.abs
      last_num_divmod = last_num.divmod(6)

      if last_num_divmod.last == 1
         i = last_num_divmod.first
      else
         i = last_num_divmod.first + 1
      end

      next_primes = []
      search_next_primes = n.to_i.abs

      while search_next_primes > 0

         n1 = 6*i+1      
         n2 = n1+4       
         i += 1

         [n1, n2].each do |p| 

            next if p % 5 == 0 || p % 7 == 0

            if p > last_num && p.prime?
               next_primes << p
               search_next_primes -= 1 
            end
         end
      end

      next_primes

   end

end



num1 = 10_000
num1 = 1_000
num1 = 5_000

num2 = 210_349
num2 = 100_125
num2 = 55_000

ret1 = nil
ret2 = nil
ret3 = nil
ret4 = nil


Benchmark.bm(16) do |t| 

   t.report("first #{num1} primes: ") do
      ret1 = num1.nthprime_add
   end 

   Integer.primes_cache.clear
   Integer.primes_cache.concat(100.primes)

   t.report("first #{num1} primes: ") do
      ret2 = num1.nthprime_add_mr
   end 

   Integer.primes_cache.clear
   Integer.primes_cache.concat(100.primes)

   t.report("first #{num1} primes: ") do
      ret3 = num1.nthprime
   end 

   Integer.primes_cache.clear
   Integer.primes_cache.concat(100.primes)

   t.report("first #{num2} primes: ") do
      ret4 = num2.nthprime
   end 

end


puts
puts "the #{num1}th prime number: #{ret1}"
puts "the #{num1}th prime number: #{ret2}"
puts "the #{num1}th prime number: #{ret3}"
puts "the #{num2}th prime number: #{ret4}"
puts
puts "the #{num1}th prime: #{Integer.primes_cache.at(num1-1)}"
puts "the #{num2}th prime: #{Integer.primes_cache.at(num2-1)}"
puts
p Integer.primes_cache.first(10)
p Integer.primes_cache.last(10)
p Integer.primes_cache.size
puts


p 15.next_primes_in_cache

puts 594_213.next_mr_prime

p 149_137.next_mr_primes(5)

puts
p 1_000.sieve_size
p 1_000_000.sieve_size
p 2_500_000.sieve_size
p 5_000_000.sieve_size
p 10_000_000.sieve_size
p 100_000_000.sieve_size


=begin

# check with primegen.c, http://cr.yp.to/primegen.html
primes 2 48611 | nl | tail -n 1
primes 2 104729 | nl | tail -n 1
primes 2 679277 | nl | tail -n 1
primes 2 1301423 | nl | tail -n 1
primes 2 2903603 | nl | tail -n 1
primes 1303129 1303283 | nl
primes 594213 $((594213 + 500)) | head -n 1
primes 149137 $((149137 + 1000)) | head -n 5

# further prime tools from http://libtom.org
- LibTomMath
   bn_mp_prime_is_prime.c
   bn_mp_prime_miller_rabin.c
- TomsFastMath
   fp_isprime.c
   fp_prime_miller_rabin.c

=end

class String
  def is_numeric?
    Float self rescue false
  end
end

# cf. http://www.ruby-forum.com/topic/97105

6.times { k=0; p $*.map!{|i|k+k=i} << 1 }

# ... or ...

ar=[]; 6.times { k=0; p ar.map!{|i|k+k=i} << 1 }


# a more general approach (for polynomials)

class Polynomial

   def triangle(nterms, row, pos=nil)

      return nil if nterms < 2 || row < 1
      nterms = nterms - 2
      num_of_rows = row

      var1 = 0 + nterms   
      var2 = 1 + nterms
      var3 = 3 + nterms 
  
      ar1 = [0, 1, 0]   # first row
      var1.times { ar1.push(0) }
      var1.times { ar1.unshift(0) }

      ar2 = []
      ar3 = []
      ar4 = [[1]]

      for num in 0..(num_of_rows - 1)  

         nextnum = ar1.size - var2

         for nextn in 1..nextnum
            sum = 0
            count = 0
            ar1.each do |n|  
               count += 1 
               if count < var3 then t = sum += n; ar2 << t else break end 
            end

            ar3 << ar2.last
            ar2.clear
            ar1.shift

         end   # second for-loop

         ar1.clear
         ar1 << ar3
         ar1.flatten!

         var2.times { ar1.push(0) }
         var2.times { ar1.unshift(0) }

         ar4 << ar3
         ar3 = []

      end  # first for-loop

      if !pos.nil?
         ret = ar4.at(row).at(pos)
         return "No such position: #{pos} in row: #{row}" if ret.nil?
         ret
      else
         ar4.map! { |r| r.join('-') }
         ar4
      end
   end 
end 


puts Polynomial.new.triangle(2, 5)
puts Polynomial.new.triangle(3, 5)
puts Polynomial.new.triangle(4, 5)
puts Polynomial.new.triangle(5, 5)
puts Polynomial.new.triangle(5, 4, 8)
puts Polynomial.new.triangle(4, 9)
puts Polynomial.new.triangle(4, 9, 10)


#------------------------


class Integer
   def fak
      f=1
      (2..self).each { |i| f *= i   }
      f
   end
end

module Enumerable
   def sum
      inject { |n, m| n + m  }
   end
end

# cf. http://blade.nagaokaut.ac.jp/~sinara/ruby/math/combinatorics/array-rep_perm.rb
class Array
  def rep_perm(n)
    if n < 0
    elsif n == 0
      yield([])
    else
      rep_perm(n - 1) do |x|
	each do |y|
	  yield(x + [y])
	end
      end
    end
  end
end


nterms = 2
exponent = 80
exponent = 8

# create the same number of variable names as there are terms 
# example: ['a', 'b'] for (a+b)**3

var_names = ('a'..'z').to_a.slice(0, nterms)
#var_names = (('a'..'z').to_a << ('A'..'Z').to_a << ('aa'..'zz').to_a).flatten!.slice(0, nterms)

ar1 = []
(0..exponent).to_a.rep_perm(nterms) { |x| p x; ar1 << x if x.sum == exponent }   # example: ... if [2,6].sum == 8
ar1.reverse!

#p ar1
#puts ar1

ar2 = []

ar1.each do |term|

   #puts "term: #{term.inspect}"  # example: term: [5, 0, 0]
   count = 0
   var1 = 1
   term.each { |i| var1 *= i.fak }
   var2 = exponent.fak / var1 
   var3 = "#{var2}   (  #{ term.join('-') << '-' } )"  # prepare term for parsing with gsub below
   ar2 << var3

end

#p ar2

result = ar2.collect do |term|
   p term
   count = -1
   term.gsub!(/(\d+)-/)  { count += 1;  "#{var_names.at(count)}" << '**' << $1 << ' ' }
   term.gsub!(/^(\d+)( +)/, '\1\2*\2')
end

puts result

# Pseudo random number generator
# From: http://gnuvince.wordpress.com/2005/11/24/pseudo-random-number-generator
# Author: Vincent Foley-Bourgon

def random_number
  t = Time.now.to_f / (Time.now.to_f % Time.now.to_i)
  random_seed = t * 1103515245 + 12345;
  (random_seed / 65536) % 32768;
end

10.times { puts random_number }

cnt = Array.new(10,0)
20000.times { cnt[(random_number % 10).to_i] += 1 }
p cnt
puts


# Ruby Gaussian Random Number Generator
# Author: Glenn
# http://webhost101.net/rails/typo/articles/2007/07/31/ruby-gaussian-random-number-generator

def gaussian_rand 
   u1 = u2 = w = g1 = g2 = 0  # declare
   begin
     u1 = 2 * rand - 1
     u2 = 2 * rand - 1
     w = u1 * u1 + u2 * u2
   end while w >= 1
   
   w = Math::sqrt( ( -2 * Math::log(w)) / w )
   g2 = u1 * w;
   g1 = u2 * w;
   # g1 is returned  
end

sum = 0
sumsq = 0
n = 1000
0.upto(n) do 
  #r = gaussian_rand
  r = gaussian_rand * 100 + 50   # new_random_number = gaussian_rand * standard_deviation + average
  #puts r
  sum += r
  sumsq += (r*r)
end

ave = sum/n
stddev = Math::sqrt(sumsq/n - ave * ave)
puts "Average = #{ave}"
puts "StdDev = #{stddev}"
puts


# Random Number Generator
# http://scutigena.sourceforge.net/test-random.html
# http://scutigena.sourceforge.net/sources/ruby-1.7.2/random.html
# $Id: random.ruby,v 1.2 2003/12/30 01:25:05 davidw Exp $

IM = 139968
IA = 3877
IC = 29573

$last = 42.0
def gen_random (max) (max * ($last = ($last * IA + IC) % IM)) / IM end

10.times do
    puts gen_random(100000.0)
end

printf "%.5f\n", gen_random(100.0)
puts


# From: http://matt.blogs.it/entries/00002641.html
# Author: Matt Mower
# cf. http://www.taygeta.com/random/gaussian.html
# cf. http://www.bearcave.com/misl/misl_tech/wavelets/hurst/random.html

def box_mueller( mean = 0.0, stddev = 1.0 )
  x1 = 0.0, x2 = 0.0, w = 0.0

  until w > 0.0 && w < 1.0
    x1 = 2.0 * rand - 1.0
    x2 = 2.0 * rand - 1.0
    w = ( x1 * x2 ) + ( x2 * x2 )
  end

  w = Math.sqrt( -2.0 * Math.log( w ) / w )
  r = x1 * w

  mean + r * stddev
end

10.times { puts box_mueller(5.0, 1.0) }
puts


# Generating random numbers with a specified distribution
# From: http://www.cs.nyu.edu/~michaels/blog/?p=24
# Author: Michael Schidlowsky

def gen
  (x=rand)>0.5 ? x : (x < rand/2.0) ? 1.0 - x : x
end

def gen2
   (x = rand) && rand ? 1.0 - x : x
end

def compute_distribution(numSamples, &block)
   samples = []
   values = 10
   numSamples.times{samples << yield}
   dist = Array.new(values, 0)
   samples.each{ |s| dist[(s*values).floor] += 1 }
   dist
end

p compute_distribution(1000){rand}
p compute_distribution(1000){gen}
p compute_distribution(1000) { gen2 }
puts


# Random number generation for a triangular distribution
# From: http://www.brighton-webs.co.uk/distributions/triangular.asp
# cf. http://en.wikipedia.org/wiki/Triangular_distribution

def rng(m, low, high)

   # cf. the parameter info at http://www.brighton-webs.co.uk/distributions/triangular.asp
   return nil unless high > low && m > low && m < high    

   u = rand

   if u <= (m-low)/(high-low)
      r = low+ Math.sqrt(u*(high-low)*(m-low))
   else
      r = high - Math.sqrt((1.0-u)*(high-low)*(high-m))
   end

end

10.times do
  #puts rng(0.0, -25.0, 25.0)
  puts rng(50.0, 25.0, 75.0)
end
puts


# From: http://www.ruby-forum.com/topic/95931
# Author: Christoffer Lernö

class Range
  def rand
    return first if exclude_end? && last == first
    Kernel::rand(last - first + (exclude_end? ? 0 : 1)) + first
  end
end

r1 = -9..9
r2 = -90...100
p r1.rand, r2.rand

r3 = 0..9

num = []
10.times do
   num << r3.rand
   if num.first.zero? then num.shift; redo end
end

p num
puts num.to_s.to_i
puts


# See:
# http://redcorundum.blogspot.com/2008/01/randomizing-array-revisited.html
# http://redcorundum.blogspot.com/2006/12/randomizing-array-and-other-faqs.html
# http://szeryf.wordpress.com/2007/06/19/a-simple-shuffle-that-proved-not-so-simple-after-all/

class Array
  def shuffle
    array = dup
    size.downto 2 do |j|
      r = rand j
      array[j-1], array[r] = array[r], array[j-1]
    end
    array
  end
end

array = (1..50).to_a

10.times { p array.shuffle.first(10) }


#-------------------


class RNG

   def num(min=8,max=min+5,iter=1)

      return nil unless min.is_a?(Integer) && max.is_a?(Integer) && max > min && iter.is_a?(Integer)

      ret = []
      stats = Hash.new(0)     # optional; cf. stats[random_num] += 1
      digits = Array(0..9).sort_by {rand}
      #digits = (Array(0..9) * (rand(4)+1)).sort_by {rand}
      digits_size = digits.size

      iter.times do
         count = 0
         len = min + rand(max-min+1)
         ar = []
         while count < len
            i = rand(digits_size)    # get a random array index
            random_num = digits.at(i)
            stats[random_num] += 1
            ar << random_num
            digits = digits.sort_by {rand}
            count += 1
            if count == 1 && ar.first.zero? then ar.shift; stats[0] -= 1; count = 0 end
         end
         ret << ar.to_s.to_i
      end  # iter
      #ret
      ret << stats   # optional
   end

end

puts

min = 8
max = 13
iter = 10000
result =  RNG.new.num(min, max, iter)
stats = result.pop

t = 0
stats.each_value { |v| t += v }  # get the overall frequency

n = 0
t = t * 1.0
m = t / 10.0 

stats.sort.each do |k,v| 
   i = (v / t) * 100.0
   x = v > m ? v - m : m - v
   y = (x / m) * 100.0
   n += i
   puts "#{k}  ::  #{"%.2f" % m}  ::  #{v}  ::  #{"%.2f" % i} %  ::  #{ "%.2f" % ((m-v) * -1) }  ::  #{"%.2f" % y} %"
end

puts n, m, t

puts
puts RNG.new.num(20, 25, 10)[0..-2]

puts
rn = RNG.new.num(80, 100)
puts rn[0..-2]
p rn.last