patterncsharpModerate
Calculate all possible distributions of n objects to k containers
Viewed 0 times
objectsallpossibledistributionscalculatecontainers
Problem
Let's say I have \$5\$ apples, thus \$n = 5\$. I have three bags (\$k = 3\$), each having the capacities of \$4\$, \$4\$ and \$2\$:
$$c = { \{4, 4, 2 \}}$$
I'd like to calculate the number of ways the \$5\$ apples can be distributed into those bags; each bag can be empty. The results for this example would be following output:
$$\{4,1,0\}\\
\{4,0,1\}\\
\{1,4,0\}\\
\{0,4,1\}\\
\{3,2,0\}\\
\{3,0,2\}\\
\{2,3,0\}\\
\{0,3,2\}\\
\{3,1,1\}\\
\{1,3,1\}\\
\{2,2,1\}\\
\{2,1,2\}\\
\{1,2,2\}\\$$
Thus there are \$13\$ possible distributions for \$n = 5\$, \$k = 3\$ and \$c = \{4, 4, 2\}\$.
I have written a simple program which can calculate this, the basic algorithm is like following:
$$\{4,1,0\}\\
\{3,2,0\}\\
\{3,1,1\}\\
\{2,2,1\}\\$$
$$\\\{3,1,1\}\\
\{1,3,1\}\\
\{1,1,3\}\\$$
I'm storing the partitions in a
My code works well for small numbers. Another example would yield \$48\$ possible distributions:
$$\\n = 10\\
k = 3\\
c = \{10,8,5\}$$
The problem of my approach is that it isn't efficient for large numbers at all. Following input could take days, if not weeks to compute:
$$c = { \{4, 4, 2 \}}$$
I'd like to calculate the number of ways the \$5\$ apples can be distributed into those bags; each bag can be empty. The results for this example would be following output:
$$\{4,1,0\}\\
\{4,0,1\}\\
\{1,4,0\}\\
\{0,4,1\}\\
\{3,2,0\}\\
\{3,0,2\}\\
\{2,3,0\}\\
\{0,3,2\}\\
\{3,1,1\}\\
\{1,3,1\}\\
\{2,2,1\}\\
\{2,1,2\}\\
\{1,2,2\}\\$$
Thus there are \$13\$ possible distributions for \$n = 5\$, \$k = 3\$ and \$c = \{4, 4, 2\}\$.
I have written a simple program which can calculate this, the basic algorithm is like following:
- Generate all integer partitions for \$n\$. I'm optimizing this by limiting each partition to each number in \$c\$; as such, the partition \$\{5, 0, 0\}\$ is impossible. For the example above, this would generate following partitions:
$$\{4,1,0\}\\
\{3,2,0\}\\
\{3,1,1\}\\
\{2,2,1\}\\$$
- Permute each partition. Each partition will have \$(\frac{k!}{i!})\$ possible permutations, whereas \$i\$ describes the number of identical integers. The set \$\{3, 1, 1\}\$ would have precisely \$3\$ permutations, as \$(\frac{3!}{2!} = 3)\$:
$$\\\{3,1,1\}\\
\{1,3,1\}\\
\{1,1,3\}\\$$
- After each permutation is calculated, validate them based on \$c\$. The permutation \$\{1, 1, 3\}\$ would be invalid, as it wouldn't fit into \$c\$.
- Finally, if the permutation is valid, increase the counter.
I'm storing the partitions in a
BlockingCollection to avoid running out of memory in case the numbers get too large (permutations will always be slower than partitions): \$n = 1000\$ could potentially yield \$2.4 \cdot 10^{31}\$ partitions!My code works well for small numbers. Another example would yield \$48\$ possible distributions:
$$\\n = 10\\
k = 3\\
c = \{10,8,5\}$$
The problem of my approach is that it isn't efficient for large numbers at all. Following input could take days, if not weeks to compute:
Solution
The structure of the code makes it difficult to test. Suppose you had two different implementations, and wanted to test that they returned the same results? The code doesn't lend itself well to this situation.
My first attempt at speeding up the code was to use the answer to TopinFrassi's question on Mathematics.SE.
As the answer explains, the solution is given by the coefficient of \$x^n\$ in
$$\prod_{c \in C} f_c(x),$$
where \$f_c(x) = 1 + x + x^2 + \cdots + x^c\$.
First, a bare-bones implementation of
Some helper methods:
And putting it all together:
This works fine for your example of \$n = 30\$, \$C = \left\{ 1, 2, 3, \ldots, 20 \right\}\$, calculating the solution \$6209623185136\$ in about 0.03s (solution cross-checked against Joe Wallis's answer).
I then translated Joe's answer to C#, only to find that it performed faster still, taking about 0.01s for the same problem.
The difference got worse as I tried bigger numbers. For example, Joe's solution takes just 0.8s to solve \$n = 30\$ and \$C = \left\{ 1, 2, 3, \ldots, 60 \right\}\$, while my solution above takes about 4.4s.
My next attempt was based on the fact that we only need one coefficient of the product. If we lazily calculate the coefficients of the polynomials, can we cut down on the time the method takes?
To delay evaluation, I changed the coefficients from
This turned out to be considerably faster than both my first attempt, and my translation of Joe's code. Some timings for \$n = 30\$ are below:
$$
\begin{array}{lrr}
& C = \left\{ 1, \ldots, 20 \right\} & C = \left\{ 1, \ldots, 60 \right\} \\
\hline
\text{Polynomial} & 0.032s & 4.39s \\
\text{Translated} & 0.010s & 0.80s \\
\text{LazyPolynomial} & 0.002s & 0.05s \\
\end{array}$$
I've now got some more accurate timings using BenchmarkDotNet:
BenchmarkDotNet=v0.7.8.0
OS=Microsoft Windows NT 6.2.9200.0
Processor=Intel(R) Core(TM) i7-4712HQ CPU @ 2.30GHz, ProcessorCount=8
HostCLR=MS.NET 4.0.30319.42000, Arch=32-bit
Type=ThirtyTwentyCompetition Mode=Throughput Platform=HostPlatform Jit=HostJit .NET=HostFramework
```
Method | AvrTime | StdDev | op/s |
--------------------------- |----------- |---------- |------- |
PolynomialLazyThirtyTwenty | 1.2705 ms | 0.0114 ms | 787.08 |
PolynomialThirtyTwenty | 13.4990 ms | 0.2331 ms | 74.08 |
TranslatedThirtyTwenty | 9.0769 ms | 0.4108 ms | 110.17 |
Method | AvrTime | StdDev | op/s |
-------------------------- |-------------- |------------ |------ |
PolynomialLazyThirtySixty | 12.3559 ms | 0.2563 ms | 80.93 |
PolynomialThirtySixty | 2,610.2568 ms | 138.6060 ms | 0.38 |
TranslatedThirtySixty | 712.4170 m
My first attempt at speeding up the code was to use the answer to TopinFrassi's question on Mathematics.SE.
As the answer explains, the solution is given by the coefficient of \$x^n\$ in
$$\prod_{c \in C} f_c(x),$$
where \$f_c(x) = 1 + x + x^2 + \cdots + x^c\$.
First, a bare-bones implementation of
Polynomial:public struct Polynomial
{
private readonly IReadOnlyList _coefficients;
public Polynomial(IEnumerable coefficients)
{
if (coefficients == null)
{
throw new ArgumentNullException(nameof(coefficients));
}
_coefficients = coefficients.ToArray();
}
public static readonly Polynomial One = new Polynomial(new[] { BigInteger.One });
public static Polynomial operator *(Polynomial a, Polynomial b)
{
var degree = a.Degree + b.Degree;
var coefficients = new BigInteger[degree + 1];
for (var k = 0; k a[i] * b[k - i])
.Sum();
}
return new Polynomial(coefficients);
}
public BigInteger this[int n] => n >= 0 && n _coefficients.Count - 1;
}Some helper methods:
internal static class EnumerableExtensions
{
public static BigInteger Sum(this IEnumerable source)
{
return source.Aggregate(BigInteger.Zero, (x, y) => x + y);
}
public static Polynomial Product(this IEnumerable source)
{
return source.Aggregate(Polynomial.One, (x, y) => x * y);
}
}And putting it all together:
public static BigInteger Solve(int n, IEnumerable bags)
{
return bags.Select(bag => new Polynomial(Enumerable.Repeat(BigInteger.One, bag + 1)))
.Product()[n];
}This works fine for your example of \$n = 30\$, \$C = \left\{ 1, 2, 3, \ldots, 20 \right\}\$, calculating the solution \$6209623185136\$ in about 0.03s (solution cross-checked against Joe Wallis's answer).
I then translated Joe's answer to C#, only to find that it performed faster still, taking about 0.01s for the same problem.
The difference got worse as I tried bigger numbers. For example, Joe's solution takes just 0.8s to solve \$n = 30\$ and \$C = \left\{ 1, 2, 3, \ldots, 60 \right\}\$, while my solution above takes about 4.4s.
My next attempt was based on the fact that we only need one coefficient of the product. If we lazily calculate the coefficients of the polynomials, can we cut down on the time the method takes?
To delay evaluation, I changed the coefficients from
BigIntegers to Lazys:public struct LazyPolynomial
{
private readonly IReadOnlyList> _coefficients;
public LazyPolynomial(IEnumerable> coefficients)
{
if (coefficients == null)
{
throw new ArgumentNullException(nameof(coefficients));
}
_coefficients = coefficients.ToArray();
}
public static readonly LazyPolynomial One = new LazyPolynomial(new[] { new Lazy(() => BigInteger.One) });
public static LazyPolynomial operator *(LazyPolynomial a, LazyPolynomial b)
{
var degree = a.Degree + b.Degree;
var coefficients = new Lazy[degree + 1];
for (var k = 0; k (() => Enumerable.Range(0, tmp + 1)
.Select(i => a[i] * b[tmp - i])
.Sum());
}
return new LazyPolynomial(coefficients);
}
public BigInteger this[int n] => n >= 0 && n _coefficients.Count - 1;
}This turned out to be considerably faster than both my first attempt, and my translation of Joe's code. Some timings for \$n = 30\$ are below:
$$
\begin{array}{lrr}
& C = \left\{ 1, \ldots, 20 \right\} & C = \left\{ 1, \ldots, 60 \right\} \\
\hline
\text{Polynomial} & 0.032s & 4.39s \\
\text{Translated} & 0.010s & 0.80s \\
\text{LazyPolynomial} & 0.002s & 0.05s \\
\end{array}$$
I've now got some more accurate timings using BenchmarkDotNet:
iniBenchmarkDotNet=v0.7.8.0
OS=Microsoft Windows NT 6.2.9200.0
Processor=Intel(R) Core(TM) i7-4712HQ CPU @ 2.30GHz, ProcessorCount=8
HostCLR=MS.NET 4.0.30319.42000, Arch=32-bit
Type=ThirtyTwentyCompetition Mode=Throughput Platform=HostPlatform Jit=HostJit .NET=HostFramework
```
Method | AvrTime | StdDev | op/s |
--------------------------- |----------- |---------- |------- |
PolynomialLazyThirtyTwenty | 1.2705 ms | 0.0114 ms | 787.08 |
PolynomialThirtyTwenty | 13.4990 ms | 0.2331 ms | 74.08 |
TranslatedThirtyTwenty | 9.0769 ms | 0.4108 ms | 110.17 |
Method | AvrTime | StdDev | op/s |
-------------------------- |-------------- |------------ |------ |
PolynomialLazyThirtySixty | 12.3559 ms | 0.2563 ms | 80.93 |
PolynomialThirtySixty | 2,610.2568 ms | 138.6060 ms | 0.38 |
TranslatedThirtySixty | 712.4170 m
Code Snippets
public struct Polynomial
{
private readonly IReadOnlyList<BigInteger> _coefficients;
public Polynomial(IEnumerable<BigInteger> coefficients)
{
if (coefficients == null)
{
throw new ArgumentNullException(nameof(coefficients));
}
_coefficients = coefficients.ToArray();
}
public static readonly Polynomial One = new Polynomial(new[] { BigInteger.One });
public static Polynomial operator *(Polynomial a, Polynomial b)
{
var degree = a.Degree + b.Degree;
var coefficients = new BigInteger[degree + 1];
for (var k = 0; k <= degree; k++)
{
coefficients[k] = Enumerable.Range(0, k + 1)
.Select(i => a[i] * b[k - i])
.Sum();
}
return new Polynomial(coefficients);
}
public BigInteger this[int n] => n >= 0 && n < _coefficients.Count ? _coefficients[n] : BigInteger.Zero;
private int Degree => _coefficients.Count - 1;
}internal static class EnumerableExtensions
{
public static BigInteger Sum(this IEnumerable<BigInteger> source)
{
return source.Aggregate(BigInteger.Zero, (x, y) => x + y);
}
public static Polynomial Product(this IEnumerable<Polynomial> source)
{
return source.Aggregate(Polynomial.One, (x, y) => x * y);
}
}public static BigInteger Solve(int n, IEnumerable<int> bags)
{
return bags.Select(bag => new Polynomial(Enumerable.Repeat(BigInteger.One, bag + 1)))
.Product()[n];
}public struct LazyPolynomial
{
private readonly IReadOnlyList<Lazy<BigInteger>> _coefficients;
public LazyPolynomial(IEnumerable<Lazy<BigInteger>> coefficients)
{
if (coefficients == null)
{
throw new ArgumentNullException(nameof(coefficients));
}
_coefficients = coefficients.ToArray();
}
public static readonly LazyPolynomial One = new LazyPolynomial(new[] { new Lazy<BigInteger>(() => BigInteger.One) });
public static LazyPolynomial operator *(LazyPolynomial a, LazyPolynomial b)
{
var degree = a.Degree + b.Degree;
var coefficients = new Lazy<BigInteger>[degree + 1];
for (var k = 0; k <= degree; k++)
{
var tmp = k;
coefficients[k] = new Lazy<BigInteger>(() => Enumerable.Range(0, tmp + 1)
.Select(i => a[i] * b[tmp - i])
.Sum());
}
return new LazyPolynomial(coefficients);
}
public BigInteger this[int n] => n >= 0 && n < _coefficients.Count ? _coefficients[n].Value : BigInteger.Zero;
private int Degree => _coefficients.Count - 1;
}```ini
BenchmarkDotNet=v0.7.8.0
OS=Microsoft Windows NT 6.2.9200.0
Processor=Intel(R) Core(TM) i7-4712HQ CPU @ 2.30GHz, ProcessorCount=8
HostCLR=MS.NET 4.0.30319.42000, Arch=32-bit
Type=ThirtyTwentyCompetition Mode=Throughput Platform=HostPlatform Jit=HostJit .NET=HostFramework
```
Method | AvrTime | StdDev | op/s |
--------------------------- |----------- |---------- |------- |
PolynomialLazyThirtyTwenty | 1.2705 ms | 0.0114 ms | 787.08 |
PolynomialThirtyTwenty | 13.4990 ms | 0.2331 ms | 74.08 |
TranslatedThirtyTwenty | 9.0769 ms | 0.4108 ms | 110.17 |
Method | AvrTime | StdDev | op/s |
-------------------------- |-------------- |------------ |------ |
PolynomialLazyThirtySixty | 12.3559 ms | 0.2563 ms | 80.93 |
PolynomialThirtySixty | 2,610.2568 ms | 138.6060 ms | 0.38 |
TranslatedThirtySixty | 712.4170 ms | 40.4284 ms | 1.40 |Context
StackExchange Code Review Q#108787, answer score: 16
Revisions (0)
No revisions yet.