ken@aiai.ed.ac.uk (Ken Johnson) (12/22/89)
In article <83943@linus.UUCP> rtidd@mwsun.mitre.org writes:
# I was trying to generate a random number bewteen 0 and maxlongint (i.e.
# 2^32) on my machine (80386 with Interactive). The rand() function that
# I have returns a psuedo-random number between 0 and MAXINT (which in
# this case is 32768). Does anyone have any ideas or any code that
# generates random numbers, or anyplace that I can check to find such
# code?
Suppose you have two random numbers, int r1 and int r2, each
chosen so as to lie randomly between 0x0000 and 0xFFFF. Then
long r;
r = ((long) r1 << 16) | ((long) r2) & 0x0000FFFFL);
-- Ken
--
Ken Johnson, AI Applications Institute, 80 South Bridge, Edinburgh EH1 1HN
E-mail ken@aiai.ed.ac.uk, phone 031-225 4464 extension 212
`I have read your article, Mr Johnson, and I am no wiser now than when I
started'. -- `Possibly not, sir, but far better informed.'henseler@uniol.UUCP (Herwig Henseler) (12/22/89)
Hello World,
There seems to be no regular posting 'Answers to frequently asked questions'
in comp.lang.c. Here is the _best_ (!) random number generator dealing with
(long). (You can find the proof in 'Communications of ACM' 1988; I don't
know the exact number, sorry)
-------------------------------- cut here -----------------------------------
#define a_random (long) 16807
#define m_random (long) 2147483647
#define q_random (long) 127773 /* m div a */
#define r_random (long) 2836 /* m mod a */
long seed;
/* Return random value [0..1) */
double
random()
{
seed = a_random * (seed % q_random) - r_random * (seed / q_random);
if( seed <= 0 )
seed += m_random;
return( (double) seed / (double) m_random );
}
main()
{
/* Initialize seed */
(void) time(&seed);
seed += getpid();
random();
while(1) {
printf( "%lf\n", random() );
}
}
-------------------------------- cut here -----------------------------------
I hope this will prevent further questions for a few months.
bye, Herwig
--
## Herwig Henseler (CS-Student) D-2930 Varel, Tweehoernweg 69 | Brain fault- ##
## EMail: henseler@uniol.UUCP (..!uunet!unido!uniol!henseler) | core dumped ##psm@manta.NOSC.MIL (Scot Mcintosh) (12/25/89)
Here's a random number generator that I like. It comes from
Knuth's Seminumerical algorithms, and has a very long period.
(2^55). This one was written for Microsoft C. It should
work equally well if you change the 'int's to 'long's, and
make a few adjustments in the seeding.
static unsigned xa[55];
static int j,k;
static void RandSeed() /* Seed our random number generator */
{
register unsigned i;
srand(1); /* Seed the c library random number generator */
for (i=0; i<55; i++)xa[i] = rand();
j = 23;
k = 54;
}
unsigned Rand() /* Generate a pseudorandom number */
{
register unsigned x;
x = xa[k] += xa[j];
if (j == 0)j = 55; --j;
if (k == 0)k = 55; --k;
return (x);
}
/****************************************************************************************************
* Procedure: GS_frmq *
* *
* Remove a far item from the Queue *
* *
* INPUTS: fitmq *q - Pointer to the Queue to remove an item from *
* *
* RETURNS: fitm far * - Pointer to the item removed from the queue *
* *
* LOCAL VARIABLES: fitmq far *fqtmp - Temporary queue item pointer *
* fitmq far *fqtmp2- 2nd temporary queue item pointer *
* unsigned int intstat - Current processor interrupt status *
* *
* GLOBAL VARIABLES: *
* *
* This routine removes an item from the specified queue. The item removed is the one at the head *
* of the queue. *
* *
****************************************************************************************************/
tMsg far *GS_frmq (q)
fitmq *q;
{
fitm far *fqtmp;
fitm far *fqtmp2;
register unsigned intstat = disable();
/* GS_QueueCheck(q); */ /*DEBUG*/
if( (q)->q_hd != (fitm far *)NULL)
if( ((q)->q_hd)->myq != q) /* Error if top item on q does not belong to this q */
FERROR("Q rem err");
fqtmp2 =
( ((fqtmp = (q)->q_hd) == ((fitm far *) NULL)) ? ((fitm far *) NULL) :\
( ( (((q)->q_hd = fqtmp->i_lnk) == ((fitm far *) NULL)) ? \
( ((q)->q_tl = &(q)->q_hd), (q)->q_nq-- ) : \
(q)->q_nq-- ), fqtmp ) ) ;
enable(intstat);
(fqtmp2->myq) = (fitmq *)NULL; /* Indicate this item is no longer part of this q */
return (tMsg far *)fqtmp2;
}
/****************************************************************************************************
* Procedure: GS_faddq *
* *
* Add a far item to the Queue *
* *
* INPUTS: tMsg far *i - Pointer to the item to be added onto the queue *
* fitmq *q - Pointer to the Queue to add an item to *
* *
* RETURNS: *
* *
* LOCAL VARIABLES: char far *chrptr - Item pointer converted to a character pointer *
* unsigned int intstat - Current processor interrupt status *
* *
* GLOBAL VARIABLES: *
* *
* This routine adds an item onto the specified queue. The item is added at the tail of the queue *
* *
****************************************************************************************************/
void GS_faddq(q,i)
fitmq *q;
tMsg far *i;
{
register unsigned intstat = disable();
char far *chrptr;
/* GS_QueueCheck(q);*/ /*DEBUG*/
/* Error if this item is null or belongs to a q already */
if(( i == (tMsg far *)NULL )|| ( (i->iob).myq != (fitmq *)NULL) )
FERROR("Q add err");
/* DEBUG ONLY-Error if item doesnt start at a valid addr*/
chrptr = (char far *) i;
if( !((unsigned int)(((uns long)chrptr-((uns long)(&g_MNHbufs1[0][0])))%sizeof(tMsg) == (unsigned int)NULL))/* ||
((unsigned int)((uns long)chrptr-((uns long)(&g_MNHbufs2[0][0])))%sizeof(tMsg) == (unsigned int)NULL)) */ )
FERROR("Illegal item addr");
(i->iob).myq = q; /* Mark item as belonging to this queue */
( (((fitm far*) i)->i_lnk = ((fitm far *) NULL)), \
(*(q)->q_tl = ((fitm far*) i)), \
((q)->q_tl = &((fitm far*) i)->i_lnk), ((q)->q_nq++) );
enable(intstat);
}
/* The following is a DEBUG routine that checks a queue for internal
consistency
*/
void GS_QueueCheck(qp)
fitmq *qp;
{
fitm far *ptr = (fitm far *)((qp)->q_hd);
register unsigned i = 0;
while( ptr != (fitm far *)NULL)
{
/* if ((unsigned int)(((uns long)ptr-((uns long)(&g_MNHbufs1[0][0])))%sizeof(tMsg) == (unsigned int)NULL)) ||
((unsigned int)((uns long)ptr-((uns long)(&g_MNHbufs2[0][0])))%sizeof(tMsg) == (unsigned int)NULL))
{
*/
if(++i > 80)
FERROR("Endless Q");
if(ptr->myq != qp)
FERROR("Q has illegal member");
ptr = ptr->i_lnk;
/* }else{
FERROR("Illegal fiorb addr");
}
*/
}
if( i != ((qp)->q_nq) )
FERROR("Q err");
}
/****************************************************************************************************
* *
* PROCEDURE: GS_PostMortem *
* *
* Routine to account for all MNH buffers upon an error condition *
* *
* INPUTS: *
* *
* RETURNS: *
* *
* LOCAL VARIABLES: *
* *
* GLOBAL VARIABLES: *
* *
****************************************************************************************************/
struct tbufacct
{ unsigned char g_nMbq, g_qMbq, g_nBadQ, g_qBadQ, g_nNCHistoryQueue, g_qNCHistoryQueue;
unsigned char g_nOpq, g_qOpq, g_nRrq, g_qRrq, g_nRfq, g_qRfq;
unsigned char g_nmRmq, g_qmRmq, g_nrRmq, g_qrRmq;
unsigned char g_nGWq, g_qGWq, g_nInMsgQueue, g_qInMsgQueue;
unsigned char g_nRlyMsgQueue, g_qRlyMsgQueue, g_nSelfSendQueue, g_qSelfSendQueue;
unsigned char g_nTcq, g_qTcq;
unsigned char g_nTmq, g_qTmq, g_nOcq, g_qOcq, g_nCcq, g_qCcq, g_nCiq, g_qCiq;
}ba;
fitmq g_lost;
fitmq g_invalid;
GS_PostMortem()
{
fitmq *temp;
unsigned char i;
fclrq (&g_lost);
fclrq (&g_invalid);
ba.g_nMbq = 0;
ba.g_nBadQ = 0;
ba.g_nNCHistoryQueue = 0;
ba.g_nOpq = 0;
ba.g_nRrq = 0;
ba.g_nRfq = 0;
ba.g_nmRmq = 0;
ba.g_nrRmq = 0;
ba.g_nGWq = 0;
ba.g_nInMsgQueue = 0;
ba.g_nRlyMsgQueue = 0;
ba.g_nSelfSendQueue = 0;
ba.g_nTcq = 0;
ba.g_nTmq = 0;
ba.g_nOcq = 0;
ba.g_nCcq = 0;
ba.g_nCiq = 0;
ba.g_qMbq = g_Mbq.q_nq;
ba.g_qBadQ = g_BadQ.q_nq;
ba.g_qNCHistoryQueue = g_NCHistoryQueue.q_nq;
ba.g_qOpq = m_Opq.q_nq;
ba.g_qRrq = m_Rrq.q_nq;
ba.g_qRfq = m_Rfq.q_nq;
ba.g_qmRmq = m_Rmq.q_nq;
ba.g_qrRmq = r_Rmq.q_nq;
ba.g_qInMsgQueue = n_InMsgQueue.q_nq;
ba.g_qRlyMsgQueue = n_RlyMsgQueue.q_nq;
ba.g_qSelfSendQueue = n_SelfSendQueue.q_nq;
ba.g_qTcq = t_Tcq.q_nq;
ba.g_qTmq = t_Tmq.q_nq;
ba.g_qOcq = t_Ocq.q_nq;
ba.g_qCcq = t_Ccq.q_nq;
ba.g_qCiq = t_Ciq.q_nq;
for (i=0; i<NUMBUFS/2; i++)
{
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &g_Mbq)
ba.g_nMbq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &g_BadQ)
ba.g_nBadQ++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &g_NCHistoryQueue)
ba.g_nNCHistoryQueue++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &m_Opq)
ba.g_nOpq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &m_Rrq)
ba.g_nRrq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &m_Rfq)
ba.g_nRfq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &m_Rmq)
ba.g_nmRmq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &r_Rmq)
ba.g_nrRmq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &n_InMsgQueue)
ba.g_nInMsgQueue++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &n_RlyMsgQueue)
ba.g_nRlyMsgQueue++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &n_SelfSendQueue)
ba.g_nSelfSendQueue++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &t_Tcq)
ba.g_nTcq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &t_Tmq)
ba.g_nTmq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &t_Ocq)
ba.g_nOcq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &t_Ccq)
ba.g_nCcq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *) &t_Ciq)
ba.g_nCiq++;
else
if ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq) == (fitmq *)NULL)
{ faddq (&g_lost, &g_MNHbufs1[i][0]);}
else
{ temp = ((fitmq *)(((fiorb far *)&g_MNHbufs1[i][0])->myq));
((fiorb far *)&g_MNHbufs1[i][0])->myq = (fitmq *)NULL;
faddq (&g_invalid, &g_MNHbufs1[i][0]);
((fiorb far *)&g_MNHbufs1[i][0])->myq = (fitmq *)temp;
} ;
}
}
--
----
Scot McIntosh
Internet: psm@helios.nosc.mil
UUCP: {ihnp4,akgua,decvax,decwest,ucbvax}!sdscvax!nosc!psmpem@zyx.SE (Per-Erik Martin) (02/01/90)
In article <1549@mbf.UUCP> mark@mbf.UUCP (MHR {who?}) writes: > >Nope, no such luck. You see, what I am usually looking for is a good >(best would be _fine_!) RNG which uses integer arithmetic only. >[...] >So, how about it? > The article mentioned previously is "Random number generators: Good ones are hard to find" by Stephen K. Park and Keith W. Miller, and was published in Communications of the ACM, October 1988 Volume 31 Number 10. The generator is based on the function f(z) = 16807z mod 2147483647 which has a full period. The article do have an integer arithmetic version of the generator which could look like this in C: ----------------------------------------------------------------------------- /* Seed may be any number between 1 and 2147483646 [1, 2^31-2]. */ static unsigned long Seed = 4711; #define RAND_A 16807L #define RAND_M 2147483647L #define RAND_Q 127773L /* m div a */ #define RAND_R 2836L /* m mod a */ /* Returns a pseudo random number between 1 and 2147483646 [1, 2^31-2]. */ unsigned long random() { register long lo, hi, test; hi = Seed/RAND_Q; lo = Seed - hi*RAND_Q; /* Seed mod RAND_Q */ test = RAND_A*lo - RAND_R*hi; if (test < 0) test += RAND_M; return Seed = test; } ----------------------------------------------------------------------------- It obviously depends on 32-bit longs. Note that it does not have a full long period (but it's close enough) and never generates 0. All seeds within the given interval are equally good, no folklore magic is needed (such as: "at least five digits with the last one odd"). To check that the implementation of the generator is correct, Seed should have the value 1043618065 after running this code: Seed = 1; i = 10000; while (i--) random(); Finally, I should say that I'm no expert at all on random number generators, I just happened to need a good one and didn't trust home made bit fiddlers. If you want to know anything about the theories behind the generator above, you should read the article. -- ------------------------------------------------------------------------------- - Per-Erik Martin, ZYX Sweden AB, Bangardsgatan 13, S-753 20 Uppsala, Sweden - - Email: pem@zyx.SE - -------------------------------------------------------------------------------
caprio@uhccux.uhcc.hawaii.edu (Mike Caprio) (02/05/90)
>Nope, no such luck. You see, what I am usually looking for is a good >(best would be _fine_!) RNG which uses integer arithmetic only. Why? >Because I'm a cheapskate who refuses to invest in a floating point >processor for my PC, and also if I ever get around to writing a program >good enough to be worth distributing or (choke) selling, I want it to be >able to run well on as many systems as possible, so FP is out - too many >emulators are TOO SLOW! > >So, how about it? > The following random number generator might be of interest to you. It produces integers in the range 0-2^31, with a period of approx. 2^31. It was published in Modeling and Simulation on Micro- computers (see ref. below), and was tested extensively by the authors. The original code was fortran, but I've run this on Turbo C and the output matched the test output listed in the article. I would check it against that listing if you use a different compiler to be sure. If you can't get ahold of the journal, E-MAIL me and I'll send along at least part of the listing. For the curious, the bit about seed<0 and then adding 2x31 is apparently a quick way to do modulo 2^31 with signed longs. Incidently, this generator is quite fast, quicker than the integer random generator built into Turbo C. long seed; long lrand () { seed = seed * 1220703125L; if (seed < 0) seed += 2147483648L; } From: Dudewicz, E.J., Z. A. Karian and R. J. Marshall III. 1985. Random number generation on microcomputers. pp. 9-14. In R. G. Lavery [ed.], Modeling and simulation on microcomputers:1985. Society for Computer Simulation, La Jolla. _Mike_ caprio@uhccux Opinons? I had one last week...
chris@mimsy.umd.edu (Chris Torek) (02/05/90)
In article <6482@uhccux.uhcc.hawaii.edu> caprio@uhccux.uhcc.hawaii.edu (Mike Caprio) quotes an algorithm from >Dudewicz, E.J., Z. A. Karian and R. J. Marshall III. 1985. >Random number generation on microcomputers. pp. 9-14. In R. G. >Lavery [ed.], Modeling and simulation on microcomputers:1985. >Society for Computer Simulation, La Jolla. which goes as follows: >long seed; > >long lrand () >{ > > seed = seed * 1220703125L; > if (seed < 0) > seed += 2147483648L; >} This code is not portable, because it depends on integer overflow being ignored. (It also fails to return a value from lrand---no doubt simply an oversight.) It can be made portable quite simply: replace the `long's with `unsigned long's, add a mask (in case ULONG_MAX > 0xffffffff, e.g., if longs are 64 or 256 bits or something [yes, 256 bit longs are legal, if unlikely]), and use a cast or two: unsigned long seed; long lrand() { seed *= 1220703125; /* the L suffix is unnecessary */ seed &= 0xffffffff; /* mask to 32 bits */ if ((long)seed < 0) seed += 0x80000000; return ((long)seed); } Next, the test for ((long)seed < 0) can be eliminated by observing the following: 1. (long)seed will be < 0 in both 1s and 2s complement arithmetic only when (a) long is 32 bits and (b) bit 31 is set. 2. In both cases, adding 0x80000000 will (since long is 32 bits and unsigned arithmetic is modular) simply clear bit 31. Hence, to obtain the same behaviour as before, we can simply mask bit 31 while masking other possible bits: long lrand() { seed *= 1220703125; seed &= 0x7fffffff; return ((long)seed); } Finally, although I am no numerical analyst nor statistician and get my random numbers from Knuth vol. 2 (or more simply, from the 4BSD `random()' library), one observation: >[lrand()] produces integers in the range 0-2^31 .... As written, lrand() must never produce the value 0, because if it does it will continue to produce 0 indefinitely. Most multiplicative random number generators involve an addition step as well, so that the algorith is seed = (A * seed) + B; for two (fixed) nonzero values of A and B. A zero value for either one results in a generator with at least one fixed point (when seed is 0). -- In-Real-Life: Chris Torek, Univ of MD Comp Sci Dept (+1 301 454 7163) Domain: chris@cs.umd.edu Path: uunet!mimsy!chris