I was reading a thread on Math Forum that inspired todays post.

The question was how to go about generating a set of pseudo-random numbers within a particular range that met some constraint. The obvious solution is to keep generating numbers in that range until the constraint is satisfied and that is where the thread ended. But this got me wondering if this is the best approach. Before diving into another possible approach let’s first explore the tools Mathematica provides for generating pseudo-randomness.

*RandomInteger* and *RandomReal* are the most basic functions for random number generation. They can generate 1 or many numbers in specified range.

In[1]:= RandomInteger[{1, 30}]

Out[1]= 23

In[2]:= RandomInteger[{1, 30}, 5]

Out[2]= {6, 1, 30, 6, 7}

Here I generate a single integer in the range {1,30} and then 5 integers in that range. One of the first constraints you might have in generating a set random numbers is that there be no repetitions. **RandomInteger** is not useful here because it provides no such guarantee. Calling **RandomInteger** repeatedly until there are no repetitions is way to brutish for anyone’s taste and certainly unnecessary because Mathematica provides the proper tools.

*RandomChoice* and *RandomSample* are two functions that select random numbers from a set you provide. RandomSample will never sample an element of the set more than once, so for this application it is exactly what we need. You simple generate the *Range* of numbers you want and **RandomSample** does the rest.

In[3]:= RandomSample[Range[30], 5]

Out[3]= {14, 22, 30, 29, 4}

Given these tools we are close to solving the problem in the method suggested on Math Forum. All we need is a function that specifies the constraints. In the post the constraints were stated to be:

- 5 integers in the Range[1,30]
- At most 2 integers <= 12
- At least one integer >=25
- Each of the integers 3,4 and 7 can’t be a divisor for more than 3 of the random chosen integers
- At most two consecutive integers.

The first constraint is satisfied by simply using **Range**. The remainder can be tested using the following function. This function returns the list if the constraint is satisfied otherwise returns the empty list.

`constraint[l_List] := If[`

Count[l, x_Integer /; x <= 12] <= 2 &&

Count[l, x_Integer /; x >= 25] >= 1 &&

Count[l, x_Integer /; Or @@ Divisible[x, {3, 4, 7}]] <= 3 &&

Not[MatchQ[

Sort[l], {___, x1_, y1_, ___, x2_, y2_, ___} /;

x1 == y1 - 1 && (x2 == y2 - 1 || y1 == x2 - 1)]], l, {}]

We can now generate random numbers that satisfy the constraint as follows.

In[5]:= Module[{s},

While[Length[s = constraint[Evaluate[RandomSample[Range[30], 5]]]] == 0];s]

Out[5]= {19, 27, 6, 13, 20}

This is basically the solution arrived at in Math Forum, except that solution uses Reap and Sow, (and seems to change the problem definition by the end of the post), but in spirit these solutions are the same.

What bothered me about this solution in general is that there could be many calls to RandomSample before the constraint is satisfied. To see how many, we can make a small modification.

In[6]:= Module[{s, n = 1},

While[Length[s = constraint[Evaluate[RandomSample[Range[30], 5]]]] ==

0, n++];

{s, n}]

Out[6]= {{2, 23, 29, 12, 22}, 2}

This code counts how many calls were made before the constraint was met and returns a list with the value and the count. This time we got lucky and met the constraint in 2 calls to **RandomSample**. Let’s wrap this in a table to see how we do over several calls.

` In[8]:= Table[`

Module[{s, n = 1},

While[Length[

s = constraint[Evaluate[RandomSample[Range[30], 5]]]] == 0,

n++];

{s, n}], {30}]

Out[8]= {{{5, 20, 3, 29, 28}, 3}, {{30, 14, 21, 11, 13},

7}, {{18, 1, 29, 22, 28}, 4}, {{25, 29, 12, 28, 6},

1}, {{22, 14, 30, 18, 26}, 1}, {{12, 18, 13, 1, 30},

2}, {{25, 4, 22, 18, 15}, 3}, {{26, 24, 9, 16, 19},

2}, {{5, 20, 6, 26, 18}, 5}, {{16, 8, 22, 1, 30},

6}, {{23, 22, 27, 13, 8}, 6}, {{25, 26, 30, 23, 12},

3}, {{3, 6, 22, 25, 24}, 3}, {{15, 14, 29, 1, 26},

3}, {{27, 24, 22, 7, 26}, 4}, {{2, 29, 26, 23, 19},

1}, {{17, 23, 28, 5, 21}, 1}, {{6, 27, 20, 19, 17},

1}, {{9, 29, 17, 19, 30}, 1}, {{29, 15, 25, 7, 27},

3}, {{10, 28, 29, 27, 18}, 2}, {{28, 20, 17, 19, 26},

1}, {{22, 28, 18, 21, 5}, 3}, {{5, 13, 7, 25, 24},

1}, {{18, 10, 25, 30, 4}, 5}, {{1, 26, 12, 21, 17},

3}, {{27, 22, 8, 26, 1}, 2}, {{28, 11, 26, 21, 10},

1}, {{23, 24, 19, 27, 9}, 5}, {{27, 19, 20, 14, 22}, 4}}

Here you can see that it took up to 7 calls before RandomSample found a match. This is not that bad and for this problem it is probably a satisfactory solution. However, I tend to think in general terms when programming. In general terms, we could have a constraint that was very restrictive such that only a handful of values in the space we are sampling met the constraint. This would be bad because then **RandomSample** could be called many more times before in wandered onto a set that met the constraint.

After mulling over this for a day, another approach hit me. Why not generate every set that meets the constraint and then simply use **RandomSample** to select from that set. This solution has a potential problem of its own, in that the sample space might be huge! But let’s see how we do in this particular case.

constraint[l_List] :=

Count[l, x_Integer /; x <= 12] <= 2 &&

Count[l, x_Integer /; x >= 25] >= 1 &&

Count[l, x_Integer /; Or @@ Divisible[x, {3, 4, 7}]] <= 3 &&

Not[MatchQ[

Sort[l], {___, x1_, y1_, ___, x2_, y2_, ___} /;

x1 == y1 – 1 && (x2 == y2 – 1 || y1 == x2 – 1)]]

In[19]:= sampleSpace = Select[Subsets[Range[30], {5}], constraint];

Length[sampleSpace]

Out[20]= 54047

It took a few seconds to complete, but we see that the sample space is a manageable 54047 elements. Note, I tweaked the constraint function to simply return **True** or **False** so using it with **Select** is easier. Here **Subsets** is used to generate all subsets of the integers 1-30 of exactly length 5 (that is what {5} means – the {} is important without the curly’s it will generate all subsets of 5 items or less).

Now that we captured *sampleSpace* we can simply use **RandomChoice** to generate as many sets as we need and be guaranteed that each answer will meet the constraint.

In[22]:= RandomChoice[sampleSpace, 25]

Out[22]= {{11, 22, 24, 26, 28}, {2, 12, 19, 22, 30}, {2, 15, 19, 23,

25}, {3, 4, 21, 26, 29}, {8, 14, 23, 26, 29}, {6, 18, 23, 28,

29}, {8, 18, 19, 22, 26}, {2, 9, 23, 28, 29}, {9, 11, 16, 22,

28}, {9, 13, 15, 22, 30}, {6, 8, 13, 16, 26}, {11, 13, 15, 27,

28}, {2, 14, 18, 19, 27}, {2, 7, 26, 27, 30}, {3, 15, 17, 22,

29}, {3, 5, 13, 16, 30}, {2, 5, 18, 19, 28}, {1, 7, 21, 22, 26}, {1,

8, 18, 22, 26}, {10, 14, 16, 21, 26}, {3, 15, 21, 22, 29}, {2, 4,

20, 23, 25}, {2, 13, 16, 19, 25}, {11, 16, 19, 20, 28}, {2, 6, 14,

25, 26}}

## Conclusion

Two viable solutions are presented to the problem of generating sets of random numbers subject to constraints. Which one is better depends on the particular parameters of your problem. For this particular problem, either is okay. For problems with a huge sample space, the first solution is better because finding a value should not take too many calls to **RandomSample**. For small sample spaces, the second technique is much better. This is a classic time/space trade-off and you’ll need to decide for yourself which technique suits the problem at hand.

Perhaps you have another approach? If so, I’d love to hear it!

Download Notebook