c モンテカルロ法による円の面積の計算
      write(6,*) 'ntrial '
      read(5,*) ntrial
      write(6,*) 'iseed ?'
      read(5,*) iseed

      r0=1.0               /* 円の半径*/
      r02=r0**2

      nin=0
      do 10 i=1,ntrial
        x=randf(iseed)
        y=randf(iseed)
        r2=x**2+y**2       /* 原点からの距離の2乗 */
        if(r2.le.r02) nin=nin+1   /* 円の中 */
 10   continue
      s=4.0 *(float(nin)/float(ntrial))
      write(6,*) 'ntrial= ',ntrial,' s= ',s
      stop
      end 

      function irandf(iseed)
      iseed=iseed*48828125
      if(iseed.le.0) iseed=(iseed+2147483647)+1
      irand=iseed
      return
      end

      function randf(iseed)
      iseed=iseed*48828125
      if(iseed.le.0) iseed=(iseed+2147483647)+1
      randf=float(iseed)/2147483647
      return
      end