To generate N positive numbers that sum to a positive number M at random, where each possible combination is equally likely:
Generate N exponentially-distributed random variates. One way to generate such a number can be written as—
number = -ln(1.0 - RNDU())
where ln(x) is the natural logarithm of x and RNDU() is a method that returns a uniform random variate greater than 0 and less than 1. Note that generating the N variates with a uniform distribution is not ideal because a biased distribution of random variate combinations will result. However, the implementation given above has several problems, such as being ill-conditioned at large values because of the distribution's right-sided tail, especially when the implementation involves floating-point arithmetic. Another implementation is given in another answer.
Divide the numbers generated this way by their sum.
Multiply each number by M.
The result is N numbers whose sum is approximately equal to M (I say "approximately" because of rounding error). See also the Wikipedia article Dirichlet distribution.
This problem is also equivalent to the problem of generating random variates uniformly from an N-dimensional unit simplex.
However, for better accuracy (compared to the alternative of using floating-point numbers, which often occurs in practice), you should consider generating n random integers that sum to an integer m * x, and treating those integers as the numerators to n rational numbers with denominator x (and will thus sum to m assuming m is an integer). You can choose x to be a large number such as 232 or 264 or some other number with the desired precision. If x is 1 and m is an integer, this solves the problem of generating random integers that sum to m.
The following pseudocode shows two methods for generating n uniform random integers with a given positive sum, in random order. (The algorithm for this was presented in Smith and Tromble, "Sampling Uniformly from the Unit Simplex", 2004.) In the pseudocode below—
- the method
PositiveIntegersWithSum returns n integers greater than 0 that sum to m, in random order,
- the method
IntegersWithSum returns n integers 0 or greater that sum to m, in random order, and
Sort(list) sorts the items in list in ascending order (note that sort algorithms are outside the scope of this answer).
METHOD PositiveIntegersWithSum(n, m)
if n <= 0 or m <=0: return error
ls = [0]
ret = NewList()
while size(ls) < n
c = RNDINTEXCRANGE(1, m)
found = false
for j in 1...size(ls)
if ls[j] == c
found = true
break
end
end
if found == false: AddItem(ls, c)
end
Sort(ls)
AddItem(ls, m)
for i in 1...size(ls): AddItem(ret,
ls[i] - ls[i - 1])
return ret
END METHOD
METHOD IntegersWithSum(n, m)
if n <= 0 or m <=0: return error
ret = PositiveIntegersWithSum(n, m + n)
for i in 0...size(ret): ret[i] = ret[i] - 1
return ret
END METHOD
Here, RNDINTEXCRANGE(a, b) returns a uniform random integer in the interval [a, b).