线性同余生成器 - 输出全为 0？

Question

我一直在尝试在 Fortran 77 中制作一个非常基本的 LCG 伪随机数生成器，以将 1000 个随机数打印到一个文件中，但无论出于何种原因，输出只是 1000 个 0。整个代码很短，所以我多次梳理它并尝试改变一些东西，但我终究无法找出问题所在。我有一种预感，这可能是一个范围问题（如果这样的概念在 Fortran 中甚至有用的话），但这确实是没有根据的。

      PROGRAM RANDOM
      COMMON ISEED, RANDOMNUMBER
      ISEED = 123
      OPEN (UNIT=1,FILE='rand.in',STATUS='UNKNOWN')

      J=1

    7 CALL RANDU(ISEED)
      J=J+1
      WRITE(1,*) RANDOMNUMBER
      IF(J<1000)GOTO 7

      STOP
      END

      SUBROUTINE RANDU(ISEED)
      PARAMETER (IMAX = 2147483647, IMAXINV = 1./IMAX)
      ISEED = ISEED * 65539
      IF(ISEED<0) ISEED = ISEED + IMAX + 1
      RANDOMNUMBER = ISEED * IMAXINV
      RETURN
      END

这里有人有什么想法吗？我刚出来。

Answer 1

我已经有几十年没用 Fortran 编程了，但我会尽力提供帮助。

首先，IMAXINV 是一个整数变量，因为名称以 I 开头并且您没有将其声明为浮点数。所以除法结果将被截断为整数值 0，这解释了你的零输出。无论如何，为了正确性和速度，您的随机数生成器应该坚持使用整数运算而不是引入泛洪点运算。

Fortran 77 支持 return 值的函数，是吗？与将子例程的结果存储在全局变量中相比，这将更加简洁和模块化。

IIRC，COMMON 语句用于在模块之间共享全局值，这对随机数生成器的私有状态来说是有风险的。

您有一个名为 ISEED 的 COMMON 全局变量和一个同名的子例程形式参数（除非我记错了 Fortran 子例程声明的工作方式）。这会使事情变得混乱，应该加以解决。让子例程更新其参数 ISEED 而不是全局变量将导致它在每次循环调用它时 return 具有相同的值。也就是说，除非形式参数是实际参数的按引用调用别名——在此代码中具有相同的名称。你看，这很混乱。

你有调试器吗？如果是这样，单步执行程序并观察变量将很快揭示程序偏离您预期的地方。

Answer 2

现在可以补充@Jerry101 的答案，我已经编写了修改后的代码。这里，关键问题是IMAXINV没有显式声明为REAL，所以被解释为INTEGER（导致原代码中IMAXINV = 1.0 / IMAX一直为0 ).另外，我从 COMMON 块中删除了 ISEED（因为它作为参数传递）并在 RANDU 中放置了另一个 COMMON 语句以在例程之间共享变量。通过这些修改，程序似乎可以正常工作。

      PROGRAM RANDOM
      COMMON RANDOMNUMBER    !<--- ISEED is deleted from here

      ISEED = 123
      J=1

    7 CALL RANDU(ISEED)
      J=J+1
      WRITE(*,*) RANDOMNUMBER       !<--- write to STDOUT for test
      IF (J < 100) GOTO 7
      END

      SUBROUTINE RANDU(ISEED)
      real IMAXINV                   !<--- this is necessary
      COMMON RANDOMNUMBER            !<--- this is also necessary to share variables
      PARAMETER (IMAX = 2147483647, IMAXINV = 1./IMAX)

      ISEED = ISEED * 65539
      IF(ISEED<0) ISEED = ISEED + IMAX + 1
      RANDOMNUMBER = ISEED * IMAXINV
      END

正如另一个答案中所建议的，我们也可以直接使用 FUNCTION 到 return 一个变量。那么我们就不需要使用COMMON，这样代码就变得干净一些了。

      PROGRAM RANDOM
      ISEED = 123
      J=1

 7    RANDOMNUMBER = RANDU(ISEED)
      J=J+1
      WRITE(*,*) RANDOMNUMBER
      IF (J < 100) GOTO 7
      END

      FUNCTION RANDU(ISEED)
      real IMAXINV
      PARAMETER (IMAX = 2147483647, IMAXINV = 1./IMAX)

      ISEED = ISEED * 65539
      IF(ISEED<0) ISEED = ISEED + IMAX + 1
      RANDU = ISEED * IMAXINV                !<--- "RANDU" is the return variable
      END

但注意当使用FUNCTION时，如果函数名不符合隐式规则，则return变量的类型应在调用例程中显式声明。（在上面的代码中，RANDU 没有显式声明，因为它被解释为 REAL）。所以无论如何，Fortran77中的隐式类型规则有很多注意事项...

补充说明：

为了避免这些陷阱，我建议使用 Fortran >=90（而不是 Fortran77），因为它提供了许多防止此类错误的功能。例如，经过最少修改的代码可能如下所示：

module mymodule
contains

subroutine randu ( istate, ran )
    implicit none
    integer, parameter :: IMAX = 2147483647
    real, parameter :: IMAXINV = 1.0 / IMAX
    integer, intent(inout) :: istate
    real,    intent(out)   :: ran

    istate = istate * 65539
    if ( istate < 0 ) istate = istate + IMAX + 1
    ran = istate * IMAXINV
end subroutine

end module

program main
    use mymodule, only: randu
    implicit none
    integer :: j, istate
    real    :: randomnumber

    istate = 123    !! seed for RANDU()

    do j = 1, 99
        call randu ( istate, randomnumber )
        write(*,*) randomnumber
    enddo
end program

这里，

implicit none 用于显式强制声明所有变量。这有助于避免错误键入变量（例如问题中的 IMAXINV！）。
子程序RANDU包含在module中，这样编译器就可以提供显式接口和许多有用的检查（简而言之，module类似于C++中的命名空间）。 module 也可用于以比 COMMON.
我使用 do ... enddo 构造来循环 j 而不是手动递增它并使用 goto。前者实际上更易于使用，而且 goto 往往会使代码的可读性降低...
我把程序文件命名为"test.f90"（注意后缀.f90），允许自由格式。另外，变量用小写字母也是可以的。
[此外，因为 iseed 存储了（伪）随机数生成器的当前状态信息，最好使用一些不同的变量名称（如 istate 等？）来提醒它的值需要在通话期间保持。]

因此，如果您有兴趣，请考虑使用更现代的 Fortran 版本（而不是 Fortran77），这样我们可以编写更安全、更健壮的代码:)

线性同余生成器 - 输出全为 0？

Linear congruential generator - output is all 0s?

random

fortran

scope

fortran77