Παρασκευή 17 Ιουλίου 2020

Finding π

Finding the pi using random point set in a square with a circle inside.










\\ finding π with monte carlo method
\\ Using internal random generator
cls 0, 0
Gradient 1, 5
Double : Report 2, "Finding π" : Normal
dx=(scale.x/3 div twipsX)*twipsX
radius=(scale.x/6 div twipsX)*twipsX
minX=dx
maxX=minX*2
maxY=(scale.y*2/3 div twipsY) *twipsY
minY=maxY-minX
white=#FFFFFF
red=#FF0000
Flush
N=10
MM=3000
profiler
for k=1 to N {
      Refresh 500
      cls 5, -height div 7
      Print format$("Try {0}/{1}", k, N)
      move dx, maxY
      Pen 5 {Polygon 5, dx, 0, 0, -dx, -dx, 0, 0, dx}
      move scale.x/2, minY+radius
      circle fill 15, radius, 1, 15
      countwhite=0 : counttotal=0
      For i=1 to MM {
            move rnd*dx+minX, rnd*dx+minY
            if point=white Then countwhite++ else.if point=red Then continue
            counttotal++ :      pset red
      } : Push 4*countwhite/counttotal
}
Print TimeCount
Refresh 100
N=stack.size
If N=0 then exit
\\ get all data from stack using [] converting stack object to array object using array()
Print  "π=";array([])#sum()/N




Second variation using external random generator


\\ finding π with monte carlo method
\\ Using external random generator: CryptGenRandom from advapi32.dll
Declare CryptAcquireContext Lib "advapi32.CryptAcquireContextW" {Long &hProv, pszContainer$,pszProvider$, long dwProvType, long dwFlags}
Declare CryptReleaseContext Lib "advapi32.CryptReleaseContext" {Long hProv, Long dwFlags}
Declare CryptGenRandom Lib"advapi32.CryptGenRandom" {Long hProv, Long dwLen, Long ByPointer}
Const PROV_RSA_FULL As Long = 1
Const VERIFY_CONTEXT As Long = 0xF0000000&
Long hProv=0
Cls 0, 0
Gradient 1, 5
Double
Report 2, "Finding π"
Normal
dx=(scale.x/3 div twipsX)*twipsX
radius=(scale.x/6 div twipsX)*twipsX
minX=dx
maxX=minX*2
maxY=(scale.y*2/3 div twipsY) *twipsY
minY=maxY-minX
white=#FFFFFF
red=#FF0000
N=20
MM=2000
\\ we need 2 values for MM*N points
buffer clear rdata as integer*MM*N*2
\\  len(rdata) is equal to MM*N*4, because integer is 2 bytes unsigned value
Call void CryptAcquireContext(&hProv, "", "", PROV_RSA_FULL, VERIFY_CONTEXT)
\\ with one call we get all the random numbers
Call Void CryptGenRandom( hProv, len(rdata), rdata(0))
Call Void CryptReleaseContext(hProv, 0&)
group one {
Private:
      \\ a now is a pointer to rdata
      a=rdata, c=0
      d=dx/0x10000
Public:
      value {
            =eval(.a,.c)*.d : .c++
      }
}
\\ flush the current stack
Flush
Profiler
for k=1 to N {
      Refresh 500
      Cls 5, -height div 7
      Print format$("Try {0}/{1}", k, N)
      move dx, maxY
      Pen 5 {Polygon 5, dx, 0, 0, -dx, -dx, 0, 0, dx}
      move scale.x/2, minY+radius
      circle fill 15, radius, 1, 15
      countwhite=0 : counttotal=0
      For i=1 to MM {
            move one+minX, one+minY
            if point=white Then countwhite++ else.if point=red Then continue
            counttotal++ :      pset red
      }
      \\ send result to stack
      Push 4*countwhite/counttotal
}
Print timecount
Refresh 100
N=stack.size
If N=0 then exit
\\ get all data from stack using [] converting stack object to array object using array()
Print  "π=";array([])#sum()/N


Δεν υπάρχουν σχόλια:

Δημοσίευση σχολίου

You can feel free to write any suggestion, or idea on the subject.