サンプルプログラム

    program calcpi

    implicit real*8(a-h,o-z)
    integer tsum, isum

    tsum = 0
    isum = 0

    do i = 1 , 1000000
       tsum = tsum + 1
       x = rmdm()
       y = rmdm()
       r = sqrt(x*x + y*y)

c write(6,*) x,y,r,tsum, isum

       if ( r .le. 1.0d0 ) isum = isum + 1

       if ( ( (i/1000000)*1000000 - i ) .eq. 0 ) then
         pi = dble(isum) / dble(tsum) * 4.0d0
         write(6,'(a,i10,a,f10.6)')
   $     'num. of sampling = ',i, ' pi = ',pi
       end if

    end do
c
    stop
    end
c==========================================================================
    real*8 function rmdm()
*******************************************
* Random Number Making Program and Result *
* -1.0 =< irndm <= 1.0 *
*******************************************
    integer m1, m2, m3, m4
    save m1, m2, m3
c
    data m1,m2,m3/92658393,76438567,65428733/
c
    m4=m1+m2+m3
    if(m2 .lt. 50000000) m4=m4+1357
    if(m4 .ge. 100000000) m4=m4-100000000
    if(m4 .ge. 100000000) m4=m4-100000000
    m1=m2
    m2=m3
    m3=m4
c***
    rmdm=dble(m3)/dble(100000000)
    rmdm= rmdm*2.0d0 - 1.0d0
c***
    return
    end