From 799cd8be78736446c9f8b09245882a0fe71c5d14 Mon Sep 17 00:00:00 2001 From: "N. Brouard" Date: Sat, 3 Feb 2024 23:37:11 +0000 Subject: [PATCH] Summary: Random from Brent in Algol W --- src/random.alw | 38 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 38 insertions(+) create mode 100644 src/random.alw diff --git a/src/random.alw b/src/random.alw new file mode 100644 index 0000000..052a308 --- /dev/null +++ b/src/random.alw @@ -0,0 +1,38 @@ +COMMENT: RANDOM NUMBER GENERATOR + *********************** + + PROCEDURE RANDOM RETURNS A LONG REAL RANDOM NUMBER UNIFORMLY + DISTRIBUTED IN (0,1) (INCLUDING 0 BUT NOT 1), + RANINIT(R) WITH R ANY INTEGER MUST BE CALLED FOR + INITIALIZATION BEFORE THE FIRST CALL TO RANDOM, AND THE + DECLARATIONS OF RAN1, RAN2 AND RAN3 MUST BE GLOBAL, + THE ALGORITHM RETURNS X(N)/2**56, WHERE + X(N) = X(N-1) + X(N-127) (MOO 2**56), + SINCE 1 + X + X**127 IS PRIMITIVE (MOD 2), THE PERIOD IS AT + LEAST 2**127 - 1 > 10**38, SEE KNUTH (1969), PP. 26, 34, 464, + X(N) IS STORED IN A LONG REAL WORD AS + RAN3 = X(N)/2**56 - 1/2, AND ALL FLOATING POINT ARITHMETIC + IS EXACT; + + LONG REAL PROCEDURE RANDOM(INTEGER VALUE NAUGHT); + BEGIN + LONG REAL RAN1; INTEGER RAN2; LONG REAL ARRAY RAN3 (0::126); + INTEGER R; LOGICAL INIT; + INIT := FALSE; + IF INIT THEN GO TO L3; + R := ABS(NAUGHT) REM 8190 + 1; + RAN2 := 127; WHILE RAN2 > 0 DO + BEGIN RAN2 := RAN2 - 1; RAN1 := -2L**55; + FOR I := 1 UNTIL 7 DO + BEGIN R := (1756*R) REM 8191; + RAN1 := (RAN1 + (R DIV 32) )*( 1/256) ; + END; + RAN3 (RAN2) := RAN1 + END; + INIT := TRUE; + L3: RAN2 := IF RAN2 = 0 THEN 126 ELSE RAN2 - 1; + RAN1 := RAN1 + RAN3 (RAN2); + RAN3 (RAN2) := RAN1 := IF RAN1 < 0L THEN RAN1 + 0.5L + ELSE RAN1 - 0.5L; + RAN1 + 0.5L + END RANDOM. -- 2.43.0