1 | ! $Id: lmdz_ran1.f90 5117 2024-07-24 14:23:34Z evignon $ |
---|
2 | |
---|
3 | MODULE lmdz_ran1 |
---|
4 | IMPLICIT NONE; PRIVATE |
---|
5 | PUBLIC ran1 |
---|
6 | CONTAINS |
---|
7 | |
---|
8 | FUNCTION RAN1(IDUM) |
---|
9 | IMPLICIT NONE |
---|
10 | REAL :: RAN1 |
---|
11 | REAL, SAVE :: R(97) |
---|
12 | REAL, PARAMETER :: RM1 = 3.8580247E-6, RM2 = 7.4373773E-6 |
---|
13 | INTEGER, SAVE :: IFF = 0 |
---|
14 | INTEGER, save :: ix1, ix2, ix3 |
---|
15 | INTEGER, PARAMETER :: M1 = 259200, IA1 = 7141, IC1 = 54773 |
---|
16 | INTEGER, PARAMETER :: M2 = 134456, IA2 = 8121, IC2 = 28411 |
---|
17 | INTEGER, PARAMETER :: M3 = 243000, IA3 = 4561, IC3 = 51349 |
---|
18 | INTEGER :: IDUM, J |
---|
19 | |
---|
20 | IF (IDUM<0.OR.IFF==0) THEN |
---|
21 | IFF = 1 |
---|
22 | IX1 = MOD(IC1 - IDUM, M1) |
---|
23 | IX1 = MOD(IA1 * IX1 + IC1, M1) |
---|
24 | IX2 = MOD(IX1, M2) |
---|
25 | IX1 = MOD(IA1 * IX1 + IC1, M1) |
---|
26 | IX3 = MOD(IX1, M3) |
---|
27 | DO J = 1, 97 |
---|
28 | IX1 = MOD(IA1 * IX1 + IC1, M1) |
---|
29 | IX2 = MOD(IA2 * IX2 + IC2, M2) |
---|
30 | R(J) = (REAL(IX1) + REAL(IX2) * RM2) * RM1 |
---|
31 | END DO |
---|
32 | IDUM = 1 |
---|
33 | ENDIF |
---|
34 | IX1 = MOD(IA1 * IX1 + IC1, M1) |
---|
35 | IX2 = MOD(IA2 * IX2 + IC2, M2) |
---|
36 | IX3 = MOD(IA3 * IX3 + IC3, M3) |
---|
37 | J = 1 + (97 * IX3) / M3 |
---|
38 | IF(J>97.OR.J<1) stop 1 |
---|
39 | RAN1 = R(J) |
---|
40 | R(J) = (REAL(IX1) + REAL(IX2) * RM2) * RM1 |
---|
41 | |
---|
42 | END FUNCTION RAN1 |
---|
43 | END MODULE lmdz_ran1 |
---|