Question

fortran程序解读

今天遇到这样一个东西，我没接触过Fortran，看不懂。而我的需要是用C++来用同样的算法实现产生随机数。求达人给出下面程序的解释，最主要的是其所用的算法。100分奉上~
程序如下：
subroutine randomz(ia,ib,x) !伪随机数产生器子程序
implicit none
integer ia,ib
real*8 x
integer,save:: initial(15)
real*8 ,save:: iy(15),iz
selectcase(ia)
case(1) ! 初始
iz=100000001.d0
read(5,'(5i8)') initial
iy=dble(initial)
x=iy(ib)*1.0d-8
case(2) ! 产生随机数
iy(ib)=mod(329.d0*iy(ib),iz)
x=iy(ib)*1.0d-8
case(3) ! 复位为初始iy
iy=dble(initial)
endselect
end

Answer 1

double randomz (int ia, int ib)
{

double x; //返回值

static int initial[15]={0};
static double iz,iy[15]={0.0};

//使用static类型，为下次调用保留值，不然每次都要从头开始

switch(ia)
{
case 1: //ia参数为1，从键盘输入种子初始化随机数生成器
iz=100000001.0

read(5,'(5i8)') initial 这句话直接翻成C很难，
//意思就是从键盘(5号设备)连续读入5个整数（难道不是15？
// 我怀疑源码写错了，应该是15i8），存入整型数组 initial ；
//而且Fortran的格式描述符i8强制每个整数都是8位（千万位），
//如果输入不足8位就切换到下一个数进行输入，则Fortran会
//自动在不足8位的数的右侧补上足够的0，以放大到千万位

iy=dble(initial) 这句话直接翻成C很难，
//本句用到Fortran90/95特色的数组整体操作，C/C++要用循环；
//是把整型数组initial逐项复制给double数组iy，转换成double型

x = iy[ib]*1.0E-8 ;

case 2: //ia参数为2，继续使用已经初始化好了的生成器
iy[ib] = 329.0 * iy[ib] % iz ;
x = iy[ib]*1.0E-8 ;

case 3: //ia参数为3，重新初始化，但还使用原先的种子
iy=dble(initial) 这句话直接翻成C很难，用循环完成。
}
//switch结束

return x;
}

=================================

算法的主要思想就是“线性同余法”，
linear-congruential method
其基本迭代公式为
X[n+1] = ( A * X[n] + B ）% C
X的初始值随便取

在你给的源代码里：
A= 329.0
B= 0
C= 100000001.0

源代码中的关键是 iy[ib] = 329.0 * iy[ib] % iz
另外 x = iy[ib]*1.0E-8 是为了将结果归一化到0~1之间再返回

你可以自己手工算几个数，就能看出这个算法的奥妙了。

另外需要指出的是——源代码里用static就是为了每次case 2时候的调用，都是在对上一次的结果进行迭代。而ib参数的用处就是保持有几组不同的独立迭代序列可用，防止不同用途的几处“生成伪随机数”调用互相干扰。