如何在 Fortran 随机模拟器代码中更改种子数答案

【问题标题】：How to change seed number in Fortran stochastic simulator code如何在 Fortran 随机模拟器代码中更改种子数
【发布时间】：2020-02-06 23:50:25
【问题描述】：

我正在运行一个 Fortran 代码，它对标记的 Poisson 集群过程执行随机模拟。在实践中，事件属性（例如发生时间）是通过反演方法生成的，即通过累积分布函数的随机采样。由于泊松随机性，我希望每个生成的序列都不同，但事实并非如此。我猜原因是伪随机数生成器的种子在每次模拟中都是相同的。我不知道 Fortran，所以我不知道如何解决这个问题。这是与伪随机数生成器有关的部分代码，知道吗？

subroutine pseud0(r)
c     generation of pseudo-random numbers
c     data ir/584287/
      data ir/574289/
      ir=ir*48828125
      if(ir) 10,20,20
   10 ir=(ir+2147483647)+1
   20 r=float(ir)*0.4656613e-9
      return
      end
      subroutine pseudo(random)
c     wichmann+hill (1982) Appl. Statist 31
      data ix,iy,iz /1992,1111,1151/
      ix=171*mod(ix,177)-2*(ix/177)
      iy=172*mod(iy,176)-35*(iy/176)
      iz=170*mod(iz,178)-63*(iz/178)
      if (ix.lt.0) ix=ix+30269
      if (iy.lt.0) iy=iy+30307
      if (iz.lt.0) iz=iz+30323
      random=mod(float(ix)/30269.0+float(iy)/30307.0+
     &            float(iz)/30323.0,1.0)
      return
      end

【问题讨论】：

种子被硬编码到源中。您可以简单地更改那里的值（在data 语句中），但如果您想做一些更奇特的事情，那么这将需要重新设计。
问题是我需要运行这段代码至少 100 次，这样每次手动编辑种子会很麻烦。
一种方法是从命令行传递种子，这样您的代码就不需要更改和重新编译。另一种方法是从/dev/urandom（posix 系统）读取“随机性”。
谢谢@PierredeBuyl，我不应该编辑代码让它接受来自命令行的输入吗？

标签： random fortran simulation random-seed stochastic-process

【解决方案1】：

首先，我将回顾 PRNG 的现代文献并选择一个现代实现。其次，我会用现代 Fortran 重写代码。

您需要遵循@francescalus 的建议并拥有更新种子的方法。在不尝试使您的代码现代化的情况下，这是pseud0 prng 的一种方法

 subroutine init0(i)
    integer, intent(in) :: i
    common /myseed0/iseed
    iseed = i
 end subroutine init0

subroutine pseud0(r)
   common /myseed0/ir
   ir = ir * 48828125
   if (ir) 10,20,20
10 ir = (ir+2147483647)+1
20 r = ir*0.4656613e-9
end subroutine pseud0

program foo
   integer i
   real r1
   call init0(574289)  ! Original seed
   do i = 1, 10
      call pseud0(r1)
      print *, r1
   end do
   print *
   call init0(289574)  ! New seed
   do i = 1, 10
      call pseud0(r1)
      print *, r1
   end do
   print *
end program foo

【讨论】：

谢谢@steve。事实上，这是一个非常古老的代码，这就是为什么对于那些不了解 Fortran 的人来说更加困难的原因。你知道为什么会有pseud0 prng 和pseudo(random) prng 吗？有什么区别？据我所知，在整个代码中，只调用了pseud0 prng。
@francescalus，当然，“现代”是为了避免common。我并没有试图成为现代人，只是向 OP 展示了一种满足 OP 需求的方法。希望 OP 遵循我的第一个建议，并寻找现代 PRNG。
@Angela 对于生成器，您可能比 Wichmann-Hill 做得更差，但您可能会研究 xoroshiro128+ 的 FORTRAN 实现。它在标准测试套件上表现良好，而且速度非常快。
@francescalus，我更新了我的答案，指出我没有尝试对代码进行现代化改造。光看旧代码就眼睛疼。
@steve ，问题是我不明白为什么有 Wichmann-Hill 生成器，程序从不调用伪（随机）子例程。