Skip to content

Instantly share code, notes, and snippets.

@fredRos fredRos/partitions.cxx
Last active Mar 23, 2017

Embed
What would you like to do?
Create integer partitions in multiplicity representation. Either all partitions of n or into exactly k parts. The basis is algorithm Z from A. Zoghbi: Algorithms for generating integer partitions, Ottawa (1993) http://dx.doi.org/10.20381/ruor-11312
// Copyright (c) 2017 Frederik Beaujean (Frederik.Beaujean@lmu.de)
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to deal
// in the Software without restriction, including without limitation the rights
// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
// The above copyright notice and this permission notice shall be included in all
// copies or substantial portions of the Software.
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
// SOFTWARE.
#include <cassert>
#include <algorithm>
#include <iostream>
#include <vector>
using namespace std;
using Int_t = int;
using vec = std::vector<Int_t>;
/*
* The output has to be filtered because Algorithm Z keeps higher
* coefficients so the naive \sum_i c_i*y_i may exceed `n`.
*/
void filter_partition(vec& y, vec& c, const Int_t n)
{
assert(y.size() == c.size());
assert(y.size() > 0);
Int_t sum = 0;
for(size_t i=1; i < y.size(); ++i) {
// check if n already reached, then zero the rest
sum += c[i]*y[i];
if (sum > n) {
fill(y.begin() + i, y.end(), 0.0);
fill(c.begin() + i, c.end(), 0.0);
break;
}
}
}
/*
* Print a partition of `n` in multiplicity representation.
*
* The output has to be filtered because Algorithm Z keeps higher
* coefficients so the naive sum may exceed `n`.
*/
void print_partition(vec& y, vec& c, bool welcome=false)
{
assert(y.size() == c.size());
assert(y.size() > 0);
if (welcome)
cout << "sum_i c_i*y_i" << endl;
bool init = false;
for(size_t i=1; i < y.size(); ++i) {
// only print valid entries
if (c[i] > 0 && y[i] > 0) {
std::cout << (init? " + " : "") << c[i] << "*" << y[i];
init = true;
}
}
std::cout << std::endl;
}
/*!
* Code up the update step inside the while loop of algorithm Z from
*
* A. Zoghbi: Algorithms for generating integer partitions, Ottawa (1993)
*
* to generate the next integer partition of `n` in multiplicity
* representation `\sum_i c_i * y_i` such that `c_i` is the
* multiplicity of part `y_i`.
*
*/
void algorithmZ_update(const Int_t n, vec& y, vec& c, Int_t& h)
{
assert(y.size() == c.size());
assert(y.size() > 0);
// running index
Int_t i = h-1;
// calculated part
Int_t k = c[h];
// calculated remainder
Int_t r = c[h] * y[h];
// calculate the remainder
while (y[h] - y[i] < 2) {
k += c[i];
r += c[i] * y[i];
--i;
}
// update current part when it equals 1
if (c[i] == 1) {
if (i != 0) {
r += c[i] * y[i];
++y[i];
} else {
i = 1;
y[i] = 1;
}
}
// update current part != 1
else {
--c[i];
r += y[i];
++i;
y[i] = y[i-1] + 1;
}
// calculate next parts based on remainder left from previous update
c[i] = k;
r -= c[i] * y[i];
h = i+1;
// update last modified part if it's the remainder
if (r == y[i]) {
++c[i];
h = i;
}
// add new part with multiplicity 1
else {
y[h] = r;
c[h] = 1;
}
filter_partition(y, c, n);
}
/*!
* Print all partitions of `n` in multiplicity representation where
* `c_i` is the multiplicity and `y_i` is the part.
*
* Based on the pseudocode of Algorithm Z in A. Zoghbi: Algorithms for
* generating integer partitions, Ottawa (1993), which generates all
* partitions of `n`.
*
* Example: partition(5)
sum_i c_i*y_i
1*5
1*1 + 1*4
1*2 + 1*3
2*1 + 1*3
1*1 + 2*2
3*1 + 1*2
5*1
*/
void partition(const Int_t n)
{
assert(n > 0);
vec y(n+1, 0), c(n+1, 0);
// initialize
y[0] = -1;
c.at(0) = 1;
y.at(1) = n;
c[1] = 1;
// start with only one part
Int_t h = 1;
print_partition(y, c, true);
// quit if all parts are 1
while (c[1] != n) {
algorithmZ_update(n, y, c, h);
print_partition(y, c);
}
}
/*!
* Print all partitions of `n` into exactly `s` parts in multiplicity
* representation where `c_i` is the multiplicity and `y_i` is the
* part.
*
* Based on the pseudocode of Algorithm Z in A. Zoghbi: Algorithms for
* generating integer partitions, Ottawa (1993), which generates all
* partitions of `n`. I made the necessary modifications to partition
* into `s` parts.
*
* Example: partition(8, 3)
sum_i c_i*y_i
2*1 + 1*6
1*1 + 1*2 + 1*5
1*1 + 1*3 + 1*4
2*2 + 1*4
1*2 + 2*3
*/
void partition(const Int_t n, const Int_t s)
{
assert(n > 0);
assert(s > 0);
assert(n >= s);
// have at most s distinct parts + one buffer value at index 0
vec y(s+1, 0), c(s+1, 0);
y[0] = -1;
c[0] = 1;
// h is the number of distinct parts
Int_t h = 0;
// the largest part of any partition of n into s parts
const auto maxPart = n-s+1;
// define initial partition
// If there is only one distinct part then maxPart divides n. But
// part(6,3) = 5+1 = 4+2 = 3+3 shows that it doesn't have to be
// the first partition, rather it is always the last in our order
// because it cannot be mutated further. So if first and last
// agree, there is only one partition. From a quick glance at the
// partition table or the recurrence relation defining it T(n, 1)
// = T(n, n) = 1, T(n, k) = 0 (k>n), T(n, k) = T(n-1, k-1) +
// T(n-k, k), this reduces the choices to s=1,n-1,n. But for
// s=n-1, the partition is n = (n-2)*1 + 1*2 so has two distinct
// parts if n > 2. For n=2, s=n-1=1 is no special case, and for
// n=1 all options collapse.
if (s == 1 || s == n) {
y[1] = maxPart;
c[1] = n / maxPart;
h = 1;
} else {
// initialize to the partition where the largest part has the
// maximum value and everything else is 1. There are two
// distinct parts
y[1] = 1;
c[1] = s-1;
y[2] = maxPart;
c[2] = 1;
h = 2;
}
print_partition(y, c, true);
// quit if all parts differ by at most one
while (y[h] - y[1] > 1) {
algorithmZ_update(n, y, c, h);
print_partition(y, c);
}
}
int main()
{
partition(6);
cout << "-------------\n";
for (Int_t i=1; i <= 6; ++i)
partition(6, i);
return 0;
}
// Local Variables:
// compile-command:"g++ -std=c++11 -g -O2 partitions.cxx -o partitions && ./partitions"
// End:
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.