# PRNG for generating numbers with n set bits exactly

I'm currently writing some code to generate binary data. I specifically need to generate 64-bit numbers with a given number of set bits; more precisely, the procedure should take some $0 < n < 64$ and return a pseudo-random 64-bit number with exactly $n$ bits set to $1$, and the rest set to 0.

My current approach involves something like this:

1. Generate a pseudorandom 64-bit number $k$.
2. Count the bits in $k$, storing the result in $b$.
3. If $b = n$, output $k$; otherwise go to 1.

This works, but it seems inelegant. Is there some kind of PRNG algorithm which can generate numbers with $n$ set bits more elegantly than this?

What you need is a random number between 0 and ${ 64 \choose n } - 1$. The problem then is to turn this into the bit pattern.

This is known as enumerative coding, and it's one of the oldest deployed compression algorithms. Probably the simplest algorithm is from Thomas Cover. It's based on the simple observation that if you have a word that is $n$ bits long, where the set bits are $x_k \ldots x_1$ in most-significant bit order, then the position of this word in the lexicographic ordering of all words with this property is:

$$\sum_{1 \le i \le k} { x_i \choose i}$$

So, for example, for a 7-bit word:

$$i(0000111) = { 2 \choose 3 } + {1 \choose 2 } + {0 \choose 1} = 0$$ $$i(0001011) = { 3 \choose 3 } + {1 \choose 2 } + {0 \choose 1} = 1$$ $$i(0001101) = { 3 \choose 3 } + {2 \choose 2 } + {0 \choose 1} = 2$$

...and so on.

To get the bit pattern from the ordinal, you just decode each bit in turn. Something like this, in a C-like language:

uint64_t decode(uint64_t ones, uint64_t ordinal)
{
uint64_t bits = 0;
for (uint64_t bit = 63; ones > 0; --bit)
{
uint64_t nCk = choose(bit, ones);
if (ordinal >= nCk)
{
ordinal -= nCk;
bits |= 1 << bit;
--ones;
}
}
return bits;
}


Note that since you only need binomial coefficients up to 64, you can precompute them.

• Cover, T., Enumerative Source Encoding. IEEE Transactions on Information Theory, Vol IT-19, No 1, Jan 1973.
• Beautiful and elegant! Enumerative coding looks like something that's very useful - are there any good resources on it (preferably in textbook form)? – Koz Ross Dec 20 '16 at 7:52
• Does this actually give better performance in practice? (Of course it depends on the speed of the RNG.) If not then there's no point in using more complex code. – Gilles 'SO- stop being evil' Dec 20 '16 at 10:58
• @Giles I interpreted this as a computer science question, since this is cs.se. I only gave source code because I happened to have it lying around from an RRR array implementation. (See, for example, alexbowe.com/rrr for an explanation of what that means.) – Pseudonym Dec 21 '16 at 3:32
• @KozRoss I don't know of any textbook that covers this specifically. The general problem is to enumerate the elements in a finite set, which lies somewhere in the intersection between combinatorics and information theory, so those are the topics I'd start looking at. – Pseudonym Dec 21 '16 at 3:46
• @Gilles To follow up on your question, I implemented both my naive method and the one provided by Pseudonym in Forth. The naive method, even when using a very simple xorshift PRNG, took something in the order of 20 second per number, while Pseudonym's method was almost instantaneous. I used tables of precomputed binomials for this. – Koz Ross Dec 22 '16 at 6:30

Very similar to Pseudonym's answer, obtained by other means.

The total number of available combinations is approachable by the stars and bars method, so it will have to be $c=\binom{64}{n}$. The total number of 64-bit numbers from which you would be trying to sample your number would be obviously much higher than that.

What you need then is a function that can lead you from a pseudorandom number $k$, ranging from $1$ to $c$, to the corresponding 64-bit combination.

Pascal's triangle can help you with that, because every node's value represents exactly the number of paths from that node to the root of the triangle, and every path can be made to represent one of the strings you are looking for, if all left turns are labeled with a $1$, and every right turn with a $0$.

So let $x$ be the number of bits left to determine, and $y$ be the number of ones left to use.

We know that $\binom{x}{y}=\binom{x-1}{y}+\binom{x-1}{y-1}$, and we can use it to properly determine the next bit of the number at each step:

$\mathtt{while}\;\;\; x>0$

$\quad \mathtt{if}\;\;\; x>y$

$\qquad \mathtt{if}\;\;\;k>\binom{x-1}{y}: \;\;\;s \leftarrow s\; + \mathtt{"1"}, \;k\leftarrow k-\binom{x-1}{y}, \;y \leftarrow y-1$

$\qquad \mathtt{else}:\; \;s \leftarrow s\; + \mathtt{"0"}$

$\quad \mathtt{else}: \;\;s \leftarrow s\; + \mathtt{"1"}, \;y \leftarrow y-1$

$\quad x \leftarrow x-1$

You can do the following :

1) Generate random number,$k$ between $1$ and $64$.

2) Set $k$ th $0$ to $1$.

3)Repeat steps 1 and 2 $n$ times

$A[]$ is $64$ bit array with all $0$s

for(i=1 to n)
{
k=ran(1,65-i) % random number between 1 and 65-i
for(x=1;x<65;x++)
{
if(A[x]==0)k--;
if(k==0)break;
}
A[x]=1;
}

• The prose doesn't seem to match your code? The code never assigns 1s to the array. Also it doesn't seem to generate a uniform distribution (and not even numbers that satisfy the constraints) when multiple ks collide – Bergi Dec 20 '16 at 10:13
• @Bergi Ya forgot the $A[x]=1$ line...added it now. And multiple collision of k is handled. See first number is chosen between 1 and 64, second between 1 and "remaining" 63. So it skips the 1 while counting...see the $if(A[x]==0)k--;$ line. And it is uniform distribution. – User Not Found Dec 20 '16 at 12:25
• Ah, I see now. The prose algorithm didn't mention the skipping. – Bergi Dec 20 '16 at 12:38
• @ArghyaChakraborty Are you using 1-based indexing there? – Koz Ross Dec 20 '16 at 19:29
• @KozRoss Start with what happens if $i=1,k=1$ (of course $A$ will be all zeroes) So, it will check $A==0$ and get $true$ meaning $k--;$ which gives $k=0$. So, sets $A=1$ outside the loop. So yeah it is 1-based indexing. To make it 0 based all you have to do is change the inner $for$ to $(x=0;x<64;x++)$ – User Not Found Dec 21 '16 at 3:01