cik@l.cc.purdue.edu (Herman Rubin) (07/16/90)
There have been request for me to post a program which points out some
of the problems with HLLs. This is written in "sort of" C. There may
be minor errors. This program is rather long, and is rather incomplete.
It would have taken me no longer to write this in a reasonable assembler,
but I am writing it this way for easier readability by the audience.
It is not completely optimized, but would do a good job with an optimizing
assembler. Note also the use of assigns. If fall-through is faster than
goto, the rather weird arrangement of the code is time-saving.
Each time G is used, a new random variable with probability(G = i) = 2^(-i)
is obtained. The if(RANBIT == 0) can be any test with probability 1/2
of success, again a new one at each stage. This program is typical of
programs requiring these; an assessment of the timing shows that these
can constitute a large part of the operating time, so that there is a
point in making these hardware instructions.
The switch statement should be optimized as suggested. There is one place
where some code is to be carried out for all subsequent cases; my knowledge
of C is inadequate for this purpose.
The outcome at each stage is either a negative number, or a bit mask, which
may be 0.
for(i = begarray; i <= endarray; i++){
Preliminaries (includes checking that enough G's are available; the
omitted portions, which are uncommon, have additional such checks.)
start: n = G;
switch(n){ /* Optimize by test and jump if greater */
case 1:
a1: out(i)=0; break;
case 2:
if (RANBIT == 0){
assign a1 to j; goto testm;}
else
ab:
out(i)= NEG; break;
case 3:
a2: out(i) = b >> G; break;
case 4:
out(i)= NEG; break;
case 5:
assign a2 to j;
goto testm;
/* this next line applies for n >= 6 */
mask = 0; q = b;
case 6:
a3: m4 = 1;
w3: ....... /* w code will be written together later, but
the blocks should be allocated as shown */
case 7:
a4: k = G;
m4 = k>>1;
if (ODD(k) { m4 += 2; goto w3;}
w4: ...... /* see above remark about w blocks */
case 8:
assign a3 to j; goto testm;
case 9:
assign a4 to j; goto testm;
case 10:
case 11:
case 13;
goto ab;
default:
if(ODD(n))goto high;
m16 = (n-8) .. 2;
if(ODD(n>1)){
assign w5 to j; goto testm;
}
w5: ...... /* see above remark about w blocks */
high: m16 = (n-11) .. 2;
if(EVEN(n>1)){
assign a6 to j; goto testm;
}
a6: nn = G;
m4 = (nn+1)>>1;
m = 6;
if (EVEN(nn)) goto a7;
wh: ...... /* see above remark about w blocks */
a7: ................ /* code omitted; ends up at wh or ab */
}
testm: n = G;
if (n == 1) goto *j;
............. /* code omitted; ends up at *j or ab */
/* Now for the w blocks */
w1: out(i) = mask; break;
w2: out(i) = mask | q >> G; break;
w3: n + G;
w3a: q >> m4;
out(i) = mask | q | q >> (n-1); break;
w4: n = G;
mask |= q >> n;
if (m4 >= n){
if (RANBIT == 0)goto w1;
q >> 1; goto w3a;
}
out(i)= mask | q>> m4;
w5: q >>= m16;
mask |= q;
n = G;
if (RANBIT == 0){
if(ODD(n))goto w2:
goto w1:
}
m4 = (n+1)>>1;
if(ODD(n))goto w3:
goto w4:
wh: ............. /* code omitted; ends up at some w */
}
--
Herman Rubin, Dept. of Statistics, Purdue Univ., West Lafayette IN47907
Phone: (317)494-6054
hrubin@l.cc.purdue.edu (Internet, bitnet) {purdue,pur-ee}!l.cc!cik(UUCP)cik@l.cc.purdue.edu (Herman Rubin) (07/16/90)
One of the readers has complained that he was unable to figure out what my program was doing. Therefore, I am posting a description of the algorithm and how it was arrived at. In article <2362@l.cc.purdue.edu>, cik@l.cc.purdue.edu (Herman Rubin) writes: > > There have been request for me to post a program which points out some > of the problems with HLLs. This is written in "sort of" C. There may > be minor errors. This program is rather long, and is rather incomplete. For more details of the actual program, see the article. This is a fairly typical example of an "infinite precision" procedure. Only a few bits are examined, and the procedure is complete for its intended use except for the copying of random bits. The code produced was produced at one sitting. The algorithm was worked out the previous day, and only the w code had been done before. ----------- It is desired to produce a random variable which is a bit mask m with probability 2(1 - e^(-.5)) and an impossible result the rest of the time, with the property that a uniform random variable with the bits of the mask cleared will have the distribution of the fractional part of an exponential random variable with mean 2. This procedure to accomplish this is to first produce a random variable such that P(J = j) = e^(-.5) * 2^(1-j / j! Then the mask consists of 1's in those places forced 0 by a bit-by-bit process of minimizing j uniform random variables, or something with that distribution. With the remaining probability, an easily tested impossible result is returned. This is to be accomplished using few random bits. Ignoring the exponential factor, the measure for 1 is 1, for 2 is 2^(-2), for 3 and 4 together is 3*2^(-6, for 5 it is C*2^(-11), where C = 16/15. and an upper bound for 6 or more is C*2^(-14). The actual bit patterns are instructive. The additional bits coming from C, as well as the additional bits needed to separate 3 from 4, are useful in obtaining the mask. Also, the exponential is .5 + .125*(1 - 1/6*(1 - 1/8*(1 - 1/10* ...)))) So the procedure is to first get a random variable which has the desired measure with probability .5, and also with probability .125. In the .125 case, appropriate procedures are used to correct for the additional factor. The "6 or more" is further refined, with rejection possibilities. Then the appropriate output procedures are used to get a mask with the right properties. A geometric random variable with p=.5 costs 2 bits average, and a test on one bit costs 1 bit. The expect number of bits used is < 3. -- Herman Rubin, Dept. of Statistics, Purdue Univ., West Lafayette IN47907 Phone: (317)494-6054 hrubin@l.cc.purdue.edu (Internet, bitnet) {purdue,pur-ee}!l.cc!cik(UUCP)