ported ulmConclusions, ulmRandomGenerators

Former-commit-id: d5a6f185c6
This commit is contained in:
Norayr Chilingarian 2014-01-22 00:43:16 +04:00
parent 1f46569693
commit 39c37311ee
9 changed files with 602 additions and 0 deletions

View file

@ -0,0 +1,419 @@
(* Ulm's Oberon Library
Copyright (C) 1989-1994 by University of Ulm, SAI, D-89069 Ulm, Germany
----------------------------------------------------------------------------
Ulm's Oberon Library is free software; you can redistribute it
and/or modify it under the terms of the GNU Library General Public
License as published by the Free Software Foundation; either version
2 of the License, or (at your option) any later version.
Ulm's Oberon Library is distributed in the hope that it will be
useful, but WITHOUT ANY WARRANTY; without even the implied warranty
of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Library General Public License for more details.
You should have received a copy of the GNU Library General Public
License along with this library; if not, write to the Free Software
Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
----------------------------------------------------------------------------
E-mail contact: oberon@mathematik.uni-ulm.de
----------------------------------------------------------------------------
$Id: RandomGener.om,v 1.9 2004/03/09 21:44:12 borchert Exp $
----------------------------------------------------------------------------
$Log: RandomGener.om,v $
Revision 1.9 2004/03/09 21:44:12 borchert
unpredictable added to the standard set of PRNGs
Revision 1.8 2004/03/06 07:22:09 borchert
Init asserts that the sequence has been registered at Services
Revision 1.7 1998/02/14 22:04:09 martin
Missing calls of Services.Init and Services.CreateType added.
Revision 1.6 1997/10/11 21:22:03 martin
assertion in ValS added, obsolete variable removed
Revision 1.5 1997/10/10 16:26:49 martin
RestartSequence added, range conversions improved,
default implementation replaced.
Revision 1.4 1997/04/01 16:33:41 borchert
major revision of Random:
- module renamed to RandomGenerators
- abstraction instead of simple implementation (work by Frank Fischer)
Revision 1.3 1994/09/01 18:15:41 borchert
bug fix: avoid arithmetic overflow in ValS
Revision 1.2 1994/08/30 09:48:00 borchert
sequences added
Revision 1.1 1994/02/23 07:25:30 borchert
Initial revision
----------------------------------------------------------------------------
original implementation by AFB 2/90
conversion to abstraction by Frank B.J. Fischer 3/97
----------------------------------------------------------------------------
*)
MODULE ulmRandomGenerators;
(* Anyone who considers arithmetical
methods of producing random digits
is, of course, in a state of sin.
- John von Neumann (1951)
*)
IMPORT
Clocks := ulmClocks, Disciplines := ulmDisciplines, Objects := ulmObjects, Operations := ulmOperations, Process := ulmProcess, Services := ulmServices, Times := ulmTimes,
Types := ulmTypes, S := SYSTEM;
TYPE
Sequence* = POINTER TO SequenceRec;
Int32ValSProc* = PROCEDURE (sequence: Sequence): Types.Int32;
LongRealValSProc* = PROCEDURE (sequence: Sequence): LONGREAL;
RewindSequenceProc* = PROCEDURE (sequence: Sequence);
RestartSequenceProc* = PROCEDURE (sequence, seed: Sequence);
SetValSProc* = PROCEDURE (sequence: Sequence; value: Operations.Operand);
CONST
int32ValS* = 0; longRealValS* = 1; rewindSequence* = 2; restartSequence* = 3;
TYPE
CapabilitySet* = SET; (* of [int32ValS..restartSequence] *)
Interface* = POINTER TO InterfaceRec;
InterfaceRec* =
RECORD
(Objects.ObjectRec)
int32ValS* : Int32ValSProc; (* at least one of ... *)
longRealValS* : LongRealValSProc; (* ... these required *)
rewindSequence* : RewindSequenceProc; (* optional *)
restartSequence*: RestartSequenceProc; (* optional *)
END;
SequenceRec* =
RECORD
(Services.ObjectRec)
(* private components *)
if : Interface;
caps: CapabilitySet;
END;
VAR
std* : Sequence; (* default sequence *)
seed*: Sequence; (* sequence of seed values *)
unpredictable*: Sequence;
(* reasonably fast sequence of unpredictable values;
is initially NIL
*)
(* ----- private definitions ----- *)
CONST
modulus1 = 2147483647; (* a Mersenne prime *)
factor1 = 48271; (* passes spectral test *)
quotient1 = modulus1 DIV factor1; (* 44488 *)
remainder1 = modulus1 MOD factor1; (* 3399; must be < quotient1 *)
modulus2 = 2147483399; (* a non-Mersenne prime *)
factor2 = 40692; (* also passes spectral test *)
quotient2 = modulus2 DIV factor2; (* 52774 *)
remainder2 = modulus2 MOD factor2; (* 3791; must be < quotient2 *)
TYPE
DefaultSequence = POINTER TO DefaultSequenceRec;
DefaultSequenceRec =
RECORD
(SequenceRec)
seed1, seed2: LONGINT;
value1, value2: LONGINT;
END;
ServiceDiscipline = POINTER TO ServiceDisciplineRec;
ServiceDisciplineRec =
RECORD
(Disciplines.DisciplineRec)
setValS: SetValSProc;
END;
VAR
service : Services.Service;
serviceDiscID: Disciplines.Identifier;
sequenceType,
defaultSequenceType: Services.Type;
(* ----- bug workaround ----- *)
PROCEDURE Entier(value: LONGREAL): LONGINT;
VAR
result: LONGINT;
BEGIN
result := ENTIER(value);
IF result > value THEN
DEC(result);
END;
RETURN result
END Entier;
(* ----- exported procedures ----- *)
PROCEDURE Init*(sequence: Sequence; if: Interface; caps: CapabilitySet);
(* initialize sequence *)
VAR
type: Services.Type;
BEGIN
ASSERT((if.int32ValS # NIL) OR (if.longRealValS # NIL));
ASSERT(~(int32ValS IN caps) OR (if.int32ValS # NIL));
ASSERT(~(longRealValS IN caps) OR (if.longRealValS # NIL));
ASSERT(~(rewindSequence IN caps) OR (if.rewindSequence # NIL));
Services.GetType(sequence, type); ASSERT(type # NIL);
sequence.if := if;
sequence.caps := caps;
END Init;
PROCEDURE Capabilities*(sequence: Sequence): CapabilitySet;
(* tell which procedures are implemented *)
BEGIN
RETURN sequence.caps
END Capabilities;
PROCEDURE RewindSequence*(sequence: Sequence);
(* re-examine sequence *)
BEGIN
ASSERT(rewindSequence IN sequence.caps);
sequence.if.rewindSequence(sequence);
END RewindSequence;
PROCEDURE RestartSequence*(sequence, seed: Sequence);
(* restart sequence with new seed values *)
BEGIN
ASSERT(restartSequence IN sequence.caps);
sequence.if.restartSequence(sequence, seed);
END RestartSequence;
PROCEDURE ^ LongRealValS*(sequence: Sequence): LONGREAL;
PROCEDURE Int32ValS*(sequence: Sequence): Types.Int32;
(* get random 32-bit value from sequence *)
VAR
real: LONGREAL;
BEGIN
IF int32ValS IN sequence.caps THEN
RETURN sequence.if.int32ValS(sequence)
ELSE
real := LongRealValS(sequence);
RETURN SHORT(Entier( (1. - real - real) * MIN(Types.Int32) ))
END;
END Int32ValS;
PROCEDURE Int32Val*(): Types.Int32;
(* get random 32-bit value from std sequence *)
BEGIN
RETURN Int32ValS(std);
END Int32Val;
PROCEDURE LongRealValS*(sequence: Sequence): LONGREAL;
(* get a uniformly distributed longreal value in [0..1) *)
BEGIN
IF longRealValS IN sequence.caps THEN
RETURN sequence.if.longRealValS(sequence)
ELSE
RETURN 0.5 +
Int32ValS(sequence) / (0. - MIN(Types.Int32) - MIN(Types.Int32))
END;
END LongRealValS;
PROCEDURE LongRealVal*(): LONGREAL;
(* get a uniformly distributed longreal value in [0..1) *)
BEGIN
RETURN LongRealValS(std)
END LongRealVal;
PROCEDURE RealValS*(sequence: Sequence): REAL;
(* get a uniformly distributed real value in [0..1) *)
BEGIN
RETURN SHORT(LongRealValS(sequence))
END RealValS;
PROCEDURE RealVal*(): REAL;
(* get a uniformly distributed real value in [0..1) *)
BEGIN
RETURN SHORT(LongRealValS(std))
END RealVal;
PROCEDURE ValS*(sequence: Sequence; low, high: LONGINT): LONGINT;
(* get a uniformly distributed integer in [low..high] *)
BEGIN
ASSERT(low <= high);
RETURN Entier( low + LongRealValS(sequence) * (1. + high - low) )
END ValS;
PROCEDURE Val*(low, high: LONGINT): LONGINT;
(* get a uniformly distributed integer in [low..high] *)
BEGIN
RETURN ValS(std, low, high)
END Val;
PROCEDURE FlipS*(sequence: Sequence): BOOLEAN;
(* return TRUE or FALSE *)
BEGIN
IF int32ValS IN sequence.caps THEN
RETURN sequence.if.int32ValS(sequence) >= 0
ELSE
RETURN sequence.if.longRealValS(sequence) >= 0.5
END;
END FlipS;
PROCEDURE Flip*(): BOOLEAN;
(* return TRUE or FALSE *)
BEGIN
RETURN FlipS(std)
END Flip;
PROCEDURE Support*(type: Services.Type; setValS: SetValSProc);
(* support service for type *)
VAR
serviceDisc: ServiceDiscipline;
BEGIN
NEW(serviceDisc);
serviceDisc.id := serviceDiscID;
serviceDisc.setValS := setValS;
Disciplines.Add(type, serviceDisc);
Services.Define(type, service, NIL);
END Support;
PROCEDURE SetValS*(sequence: Sequence; value: Operations.Operand);
(* store random value from sequence into already initialized value *)
VAR
baseType : Services.Type;
serviceDisc: ServiceDiscipline;
ok : BOOLEAN;
BEGIN
Services.GetSupportedBaseType(value, service, baseType);
ok := Disciplines.Seek(baseType, serviceDiscID, S.VAL(Disciplines.Discipline, serviceDisc));
ASSERT(ok);
serviceDisc.setValS(sequence, value);
END SetValS;
PROCEDURE SetVal*(value: Operations.Operand);
(* store random value from std sequence into already initialized value *)
BEGIN
SetValS(std, value);
END SetVal;
(* ----- DefaultSequence ----- *)
PROCEDURE CongruentialStep(VAR value1, value2: LONGINT);
BEGIN
value1 :=
factor1 * (value1 MOD quotient1) - remainder1 * (value1 DIV quotient1);
IF value1 < 0 THEN
INC(value1, modulus1);
END;
value2 :=
factor2 * (value2 MOD quotient2) - remainder2 * (value2 DIV quotient2);
IF value2 < 0 THEN
INC(value2, modulus2);
END;
END CongruentialStep;
PROCEDURE DefaultSequenceValue(sequence: Sequence): LONGREAL;
VAR
value: LONGINT;
BEGIN
WITH sequence: DefaultSequence DO
CongruentialStep(sequence.value1, sequence.value2);
value := sequence.value1 - sequence.value2;
IF value <= 0 THEN
INC(value, modulus1);
END;
RETURN (value - 1.) / (modulus1 - 1.)
END;
END DefaultSequenceValue;
PROCEDURE DefaultSequenceRewind(sequence: Sequence);
BEGIN
WITH sequence: DefaultSequence DO
sequence.value1 := sequence.seed1;
sequence.value2 := sequence.seed2;
END;
END DefaultSequenceRewind;
PROCEDURE DefaultSequenceRestart(sequence, seed: Sequence);
BEGIN
WITH sequence: DefaultSequence DO
sequence.seed1 := ValS(seed, 1, modulus1-1);
sequence.seed2 := ValS(seed, 1, modulus2-1);
sequence.value1 := sequence.seed1;
sequence.value2 := sequence.seed2;
END;
END DefaultSequenceRestart;
PROCEDURE CreateDefaultSequences;
VAR
mySeed, myStd: DefaultSequence;
if: Interface;
daytime: Times.Time;
timeval: Times.TimeValueRec;
count: LONGINT;
PROCEDURE Hash(str: ARRAY OF CHAR): LONGINT;
VAR
index,
val: LONGINT;
BEGIN
val := 27567352;
index := 0;
WHILE str[index] # 0X DO
val := (val MOD 16777216) * 128 +
(val DIV 16777216 + ORD(str[index])) MOD 128;
INC(index);
END; (*WHILE*)
RETURN val
END Hash;
BEGIN
(* define interface for all default sequences *)
NEW(if);
if.longRealValS := DefaultSequenceValue;
if.rewindSequence := DefaultSequenceRewind;
if.restartSequence := DefaultSequenceRestart;
(* fake initial randomness using some portably accessible sources *)
NEW(mySeed);
Services.Init(mySeed, defaultSequenceType);
Init(mySeed, if, {longRealValS});
Clocks.GetTime(Clocks.system, daytime);
Times.GetValue(daytime, timeval);
(* extract those 31 bits from daytime that are most likely to vary *)
mySeed.value1 := timeval.usec * 2048 + timeval.second MOD 65536 + 1;
(* generate 31 more bits from the process name *)
mySeed.value2 := Hash(Process.name) MOD (modulus2 - 1) + 1;
(* scramble these values *)
count := 0;
WHILE count < 4 DO
CongruentialStep(mySeed.value1, mySeed.value2);
INC(count);
END;
(* mix them together *)
DefaultSequenceRestart(mySeed, mySeed);
seed := mySeed;
(* now use our seed to initialize std sequence *)
NEW(myStd);
Services.Init(myStd, defaultSequenceType);
Init(myStd, if, {longRealValS, rewindSequence, restartSequence});
DefaultSequenceRestart(myStd, mySeed);
std := myStd;
unpredictable := NIL;
END CreateDefaultSequences;
BEGIN
serviceDiscID := Disciplines.Unique();
Services.Create(service, "RandomGenerators");
Services.CreateType(sequenceType, "RandomGenerators.Sequence", "");
Services.CreateType(defaultSequenceType, "RandomGenerators.DefaultSequence",
"RandomGenerators.Sequence");
CreateDefaultSequences;
END ulmRandomGenerators.