Tag Archives: Combinatorial Algorithms

Recursive integer partitions in C

Here is a recursive algorithm to generate integer partitions in antilexicographic order. The function partitions() takes an integer to partition, and a callback function to call for each partition found.

#include <stdlib.h>

unsigned int min(int a, int b)
{
    return a < b ? a : b;
}

typedef void (*partitionfn)(const unsigned int *, size_t);

static void partitions_recursive(unsigned int n, unsigned int maximum, unsigned int *partition,
        size_t length, partitionfn fun)
{
    unsigned int i;
    if (n == 0) {
        fun(partition, length);
    }
    for (i = min(maximum, n); i >= 1; i--) {
        partition[length] = i;
        partitions_recursive(n - i, i, partition, length + 1, fun);
    }
}

void partitions(unsigned int n, partitionfn fun)
{
    unsigned int *partition = malloc(n * sizeof(unsigned int));
    if (partition) {
        partitions_recursive(n, n, partition, 0, fun);
        free(partition);
    }
}

An example program:

#include <stdio.h>

void print(const unsigned int *partition, size_t length)
{
    unsigned int i;
    for (i = 0; i < length; i++) {
        printf("%u ", partition[i]);
    }
    putchar('\n');
}

int main(void)
{
    partitions(6, print);
    return 0;
}

The output:

6
5 1
4 2
4 1 1
3 3
3 2 1
3 1 1 1
2 2 2
2 2 1 1
2 1 1 1 1
1 1 1 1 1 1

Reference: Partition.java

Permutation cycles in C

A cycle of a permutation is a subset of the elements that replace one another in sequence, until the last element is replaced by the first. For example, consider the permutation below:

\(\sigma=\begin{pmatrix}
0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
1 & 3 & 7 & 0 & 4 & 8 & 2 & 6 & 5\end{pmatrix}\)

We can find the cycles:
\(0 \rightarrow 1, 1 \rightarrow 3, 3 \rightarrow 0\)
\(2 \rightarrow 7, 7 \rightarrow 6, 6 \rightarrow 2\)
\(4 \rightarrow 4\)
\(5 \rightarrow 8, 8 \rightarrow 5\)

These can be written as:
\((0, 1, 3)(2, 7, 6)(4)(5, 8)\)

It’s customary to omit cycles of length 1, so this would usually be written as
\((0, 1, 3)(2, 7, 6)(5, 8).\)

To find the cycle decomposition of a permutation, we can use an algorithm that is very similar to depth-first search (DFS) on a graph. We begin a new search for each unvisited vertex (number), and visit its neighbour (image in the permutation) until we get back to the first vertex again.

Below is an implementation in C. The function permutation_cycles() takes a permutation in the form of integers starting from 0, its length, and a callback function that is called for each cycle found.

#include <stdlib.h>

typedef void(*cyclefn)(const unsigned int *, size_t);

void permutation_cycles_recursive(const unsigned int *permutation, unsigned int *visited, 
        unsigned int start, unsigned int current, unsigned int *path,
        size_t length, cyclefn fun)
{
    visited[current] = 1;
    path[length] = current;
    if (start == current && length > 0) {
        fun(path, length);
    }
    else {
        permutation_cycles_recursive(permutation, visited, start, permutation[current],
                path, length + 1, fun);
    }
}

void permutation_cycles(const unsigned int *permutation, size_t n, cyclefn fun)
{
    unsigned int i;
    unsigned int *visited = calloc(n, sizeof(unsigned int));
    unsigned int *path = malloc(n * sizeof(unsigned int));
    if (!(visited && path)) {
        free(visited);
        free(path);
        return;
    }
    for (i = 0; i < n; i++) {
        if (!visited[i]) {
            permutation_cycles_recursive(permutation, visited, i, i, 
                    path, 0, fun);
        }
    }
    free(visited);
    free(path);
}

Here is an example program that finds the cycle decomposition of the permutation shown above:

#include <stdio.h>

void print(const unsigned int *cycle, size_t length)
{
    if (length > 1) {
        unsigned int i;
        putchar('(');
        for (i = 0; i < length; i++) {
            printf("%u", cycle[i]);
            if (i < length - 1) {
                printf(", ");
            }
        }
        putchar(')');
    }
}

int main(void)
{
    unsigned int permutation[] = {1, 3, 7, 0, 4, 8, 2, 6, 5};
    const size_t n = sizeof(permutation) / sizeof(permutation[0]);
    permutation_cycles(permutation, n, print);
    putchar('\n');
    return 0;
}

The output:

(0, 1, 3)(2, 7, 6)(5, 8)

Permutations in lexicographic order in C

This recursive algorithm produces the permutations in the most natural order, and is also the easiest to understand. It uses two buffers, one containing the permutation being built, and another for the remaining, unused letters.

The algorithm is:

if (there are no more characters to arrange) {
    print (the current permutation);
}
else {
    choose a remaining character;
    remove it from the remaining characters;
    add it to the permutation;
    use recursion to add the remaining characters to the permutation in the same way;
}

Here it is in C:

#include <stdlib.h>
#include <string.h>

typedef void(*permutationfn)(const char *);

static void permutations_recursive(char *perm, size_t used, char *remaining, size_t len,
        permutationfn fun)
{
    if (used == len) {
        perm[used] = '\0';
        fun(perm);
    } 
    else {
        const size_t unused = len - used;
        unsigned int i;
        for (i = 0; i < unused; i++) {
            /* Remove character at i from remaining */
            char c = remaining[i];
            memmove(remaining + i, remaining + i + 1, unused - i);
            /* Append it to perm */
            perm[used] = c;
            permutations_recursive(perm, used + 1, remaining, len, fun);
            /* Put character back */
            memmove(remaining + i + 1, remaining + i, unused - i);
            remaining[i] = c;
        }
    }
}

void permutations(const char *str, permutationfn fun)
{
    const size_t len = strlen(str);
    char *perm = malloc(len + 1);
    char *remaining = malloc(len + 1);
    if (perm && remaining) {
        strcpy(remaining, str);
        permutations_recursive(perm, 0, remaining, len, fun);
    }
    free(perm);
    free(remaining);

Example program:

#include <stdio.h>

int main(void)
{
    char str[] = "ABCD";
    permutations(str, (permutationfn)puts);
    return 0;
}

The output:

ABCD
ABDC
ACBD
ACDB
ADBC
ADCB
BACD
BADC
BCAD
BCDA
BDAC
BDCA
CABD
CADB
CBAD
CBDA
CDAB
CDBA
DABC
DACB
DBAC
DBCA
DCAB
DCBA

Reference: Exhaustive recursion and backtracking

Permutations in transposition order in C

This algorithm is similar to the direct insertion order permutation algorithm, except that instead of inserting a new element into the spaces between existing elements, the new element is inserted into a space already occupied by an existing element or at the end. If the space is already occupied, the existing element is bumped to the end of the permutation.

#include <string.h>
#include <stdlib.h>

typedef void(*permutationfn)(const char *);

static void permutations_recursive(const char *str, char *perm, unsigned int level, size_t len, 
        permutationfn fun)
{
    if (level == len) {
        perm[level] = '\0';
        fun(perm);
    }
    else {
        /* Insert character str[level] in every position */
        unsigned int i;
        for (i = 0; i <= level; i++) {
            if (level > 0 && i < level) {
                /* Bump the existing element to the end */
                perm[level] = perm[i];
            }
            perm[i] = str[level];
            permutations_recursive(str, perm, level + 1, len, fun);
            if (level > 0 && i < level) {
                /* Move the existing element back */
                perm[i] = perm[level];
            }
        } 
    }
}

void permutations(const char *str, permutationfn fun)
{
    char *perm = malloc(strlen(str) + 1);
    permutations_recursive(str, perm, 0, strlen(str), fun);
    free(perm);
}

Example program:

#include <stdio.h>

int main(void)
{
    char str[] = "ABCD";
    permutations(str, (permutationfn)puts);
    return 0;
}

The output:

DABC
CDBA
CADB
CABD
DCAB
BDAC
BCDA
BCAD
DACB
BDCA
BADC
BACD
DBAC
CDAB
CBDA
CBAD
DCBA
ADBC
ACDB
ACBD
DBCA
ADCB
ABDC
ABCD

Reference: Decision Trees

Permutations in direct insertion order in C

To form a permutation in direct insertion order, we take an existing partial permutation containing \(k\) elements, and then insert the next element into each of the \(k+1\) possible places immediately before an element or at the end. So for example if we start with the permutation “A”, we can insert “B” before the “A” or at the end, producing “BA” and “AB”. Taking “BA”, we can then insert “C” at the beginning, between the “B” and the “A”, or at the end, producing “CBA”, “BCA”, “BAC”. The process continues until every element has been added.

Below is an implementation in C. I used a recursive function that takes a callback to process each permutation as it is found. Notice the use of memmove() to move existing elements out of the way, and then to move them back after the recursive call. We can’t use memcpy() for this, because the source and destination ranges overlap.

typedef void(*permutationfn)(const char *);

static void permutations_recursive(const char *str, char *perm, unsigned int level, size_t len, 
        permutationfn fun)
{
    if (level == len) {
        perm[level] = '\0';
        fun(perm);
    }
    else {
        /* Insert character str[level] in every position */
        unsigned int i;
        for (i = 0; i <= level; i++) {
            if (level > 0 && i < level) {
                /* Move the existing elements out of the way */
                memmove(perm + i + 1, perm + i, level - i);
            }
            perm[i] = str[level];
            permutations_recursive(str, perm, level + 1, len, fun);
            if (level > 0 && i < level) {
                /* Move the existing elements back */
                memmove(perm + i, perm + i + 1, len - i);
            }
        } 
    }
}

void permutations(const char *str, permutationfn fun)
{
    char *perm = malloc(strlen(str) + 1);
    permutations_recursive(str, perm, 0, strlen(str), fun);
    free(perm);
}

An example program:

int main(void)
{
    char str[] = "ABCD";
    permutations(str, (permutationfn)puts);
    return 0;
}

The output:

DCBA
CDBA
CBDA
CBAD
DBCA
BDCA
BCDA
BCAD
DBAC
BDAC
BADC
BACD
DCAB
CDAB
CADB
CABD
DACB
ADCB
ACDB
ACBD
DABC
ADBC
ABDC
ABCD

Reference: Decision Trees

Recursive permutations in C

Recursion tree for a permutation algorithm

Another permutation algorithm in C, this time using recursion. I got this algorithm from Eitan Gurari’s CIS 680 lecture notes, which sadly are no longer online, although they are available on the Wayback Machine here: CIS 680: DATA STRUCTURES. I’ve stolen the image above, which shows a partial recursion tree, from him.

Here is the implementation:

#include <string.h>

typedef void (*permutationfn)(const char *);

static void swap(char *a, char *b)
{
    char temp = *a;
    *a = *b;
    *b = temp;
}

static void permutations_recursive(char *str, unsigned int left, unsigned int right,
        permutationfn fun)
{
    if (left == right) {
        fun(str);
    }
    else {
        unsigned int i;
        for (i = left; i <= right; i++) {
          swap(str + left, str + i);
          permutations_recursive(str, left + 1, right, fun);
          swap(str + left, str + i);
       }
   }
}

void permutations(char *str, permutationfn fun)
{
    permutations_recursive(str, 0, strlen(str) - 1, fun);
}

Example program:

#include <stdio.h>

int main(void)
{
    char str[] = "ABCD";
    permutations(str, (permutationfn)puts);
    return 0;
}

The output:

ABCD
ABDC
ACBD
ACDB
ADCB
ADBC
BACD
BADC
BCAD
BCDA
BDCA
BDAC
CBAD
CBDA
CABD
CADB
CDAB
CDBA
DBCA
DBAC
DCBA
DCAB
DACB
DABC

xkcd 287

General solutions get you a 50% tip

I’m not sure what this problem is called – I’m going to call it "multicombination sum" – but I don’t doubt that it is NP-complete, as it’s a variety of knapsack problem in which the values of the items are the same as their weight.

Below are three methods of solving it: a brute force method, using backtracking, and using dynamic programming.

Brute force method

The brute force method is just to construct all of the possible orders that might total $15.05.
The combinatorial algorithm we want is combinations with duplicates, or multicombinations.
Since 3 of the most expensive appetizer – SAMPLER PLATE – exceeds the target, and 7 of the cheapest appetizer – MIXED FRUIT – equals the target (so that’s one of the solutions), we want to iterate over all multicombinations with k ranging from 3 to 7.

unsigned int multiset_sum(const unsigned int *multiset, const unsigned int *values, unsigned int k)
{
    unsigned int i;
    unsigned int sum = 0;
    for (i = 0; i < k; i++) {
        sum += values[multiset[i]];
    }
    return sum;
}

typedef void (*multiset1fn)(const unsigned int *, unsigned int);

void multicombination_sum(const unsigned int *items, size_t len, unsigned int target,
        multiset1fn fun)
{
    unsigned int first = target / array_max(items, len);
    unsigned int last = target / array_min(items, len);
    unsigned int *multiset = calloc(last + 1, sizeof(unsigned int));
    unsigned int k;
    if (!multiset) {
        return;
    }
    for (k = first; k <= last; k++) {
        do {
            if (multiset_sum(multiset, items, k) == target) {
                fun(multiset, k);
            }
        } while (next_multicombination(multiset, len, k));
    }
    free(multiset);
}

void order_print(const unsigned int *numbers, unsigned int k)
{
	const char *appetizers[] = {"MIXED FRUIT", "FRENCH FRIES", "SIDE SALAD", 
        "HOT WINGS", "MOZARELLA STICKS", "SAMPLER PLATE"};
    unsigned int i, item, count;
    for (i = 0; i < k; i++) {
        if (i == 0 || numbers[i] != item) {
            if (i > 0) {
                printf("%s (x%d) ", appetizers[item], count);
            }
            count = 1;
            item = numbers[i];
        }
        else {
            count++;
        }
    }
    printf("%s (x%d)\n", appetizers[item], count);
}

int main(void)
{
    unsigned int prices[] = {215, 275, 335, 355, 420, 580};
	const size_t n = sizeof(prices) / sizeof(prices[0]);
    const unsigned int target = 1505;
    multicombination_sum(prices, n, target, (multiset1fn)order_print);

	return 0;
}

Output:

MIXED FRUIT (x1) HOT WINGS (x2) SAMPLER PLATE (x1)
MIXED FRUIT (x7)

This took 1709 iterations to come up with the answer.

Backtracking

We can drastically reduce the search space by building the orders up one item at a time, and backtracking if the target is exceeded.


typedef void(*multiset2fn)(unsigned int *, const unsigned int *, size_t); 

static void multicombination_sum_recursive(int i, unsigned int *combination, 
       const unsigned int *items, size_t len, unsigned int target, multiset2fn fun,
       unsigned int sum)
{
    if (i == (int)len - 1) {
        if (sum == target) {
            fun(combination, items, i);
        }
    }
    else  {
        unsigned int quantity;
        unsigned int max_quantity = (target - sum) / items[i + 1];
        for (quantity = 0; quantity  <= max_quantity; quantity++) {
            combination[i + 1] = quantity;
            multicombination_sum_recursive(i + 1, combination, items, len, target, fun, 
                    sum + quantity * items[i + 1]);
        }
    }
}

void multicombination_sum(const unsigned int *items, size_t len, unsigned int target,
        multiset2fn fun)
{
    unsigned int *combination = malloc(len * sizeof(unsigned int));
    multicombination_sum_recursive(-1, combination, items, len, target, fun, 0);
    free(combination);
}

void order_print(const unsigned int *combination, const unsigned int *items, size_t len)
{
    const char *appetizers[] = {"MIXED FRUIT", "FRENCH FRIES", "SIDE SALAD", 
        "HOT WINGS", "MOZARELLA STICKS", "SAMPLER PLATE"};
    unsigned int i;
    for (i = 0; i <= len; i++) {
        if (combination[i] > 0) {
            printf("%s (x%d) ", appetizers[i], combination[i]);
        }
    }
    putchar('\n');
}

int main(void)
{
    unsigned int prices[] = {215, 275, 335, 355, 420, 580};
    const unsigned int len = sizeof(prices) / sizeof(unsigned int);
    const unsigned int target = 1505;
    multicombination_sum(prices, len, target, (multiset2fn)order_print);
    return 0;
}

Output:

MIXED FRUIT (x1) HOT WINGS (x2) SAMPLER PLATE (x1)
MIXED FRUIT (x7)

This took just 211 iterations.

Dynamic Programming

As I said at the beginning, this problem is a special case of the knapsack problem in which the values of the items are the same as their weight. This means that we can use the classic dynamic programming algorithm for knapsack to solve it. The algorithm works by calculating the most profitable knapsack for each capacity up to the target capacity.

Here is an implementation of the algorithm:

struct knapsack {
    unsigned int profit;
    struct knapsack *prev;
};
typedef struct knapsack knapsack;

/* Find the minimum weight item with this profit */
int min_weight_item(unsigned int profit, const unsigned int *weights, const unsigned int *profits,
        size_t len)
{
    int item = -1;
    unsigned int i;
    for (i = 0; i < len; i++) {
        if (profits[i] == profit) {
            if (item == -1 || weights[i] < weights[item]) {
                item = i;
            }
        }
    }
    return item;
}

unsigned int unbounded_knapsack(unsigned int capacity, unsigned int *weights, 
       unsigned int *profits, unsigned int *counts, size_t len)
{
    knapsack *z = malloc((capacity + 1) * sizeof(knapsack));
    unsigned int c, i;
    unsigned int solution, profit;
    z[0].profit = 0;
    z[0].prev = NULL;
    knapsack *current;
    /* Fill in the array */
    for (c = 1; c <= capacity; c++) {
        z.profit = z.profit;
        z.prev = &(z);
        for (i = 0; i < len; i++) {
            if (weights[i] <= c) {
                knapsack *prev = z + (c - weights[i]);
                if (prev->profit + profits[i] > z.profit) {
                    z.profit = prev->profit + profits[i];
                    z.prev = prev;
                }
            }
        }
    }
    /* Read back the best solution */
    for (profit = z[capacity].profit, current = z[capacity].prev;
            current != NULL;
            profit = current->profit, current = current->prev) {
        counts[min_weight_item(profit - current->profit, weights, profits, len)]++;
        
    }
    solution = z[capacity].profit;
    free(z);
    return solution;
}

We just need to use this algorithm and pass the prices of the menu items for both weights and profits:

int main(void)
{
    unsigned int values[] = {215, 275, 335, 355, 420, 580};
    const size_t len = sizeof(values) / sizeof(values[0]);
    unsigned int counts[6] = {0};
    const unsigned int target = 1505;
    unbounded_knapsack(target, values, values, counts, len);
    order_print(counts, len);
    return 0;
}

This produces one of the solutions:

MIXED FRUIT (x7)

We could modify the algorithm to produce all of them.

K-Permutations

The k-permutations of a set are the permutations of the combinations of size k. They are also known as sequences without repetitions.

There are 24 3-permutations of the 4-set {0, 1, 2, 3}:

[0, 1, 2]
[0, 2, 1]
[1, 0, 2]
[1, 2, 0]
[2, 0, 1]
[2, 1, 0]
[0, 1, 3]
[0, 3, 1]
[1, 0, 3]
[1, 3, 0]
[3, 0, 1]
[3, 1, 0]
[0, 2, 3]
[0, 3, 2]
[2, 0, 3]
[2, 3, 0]
[3, 0, 2]
[3, 2, 0]
[1, 2, 3]
[1, 3, 2]
[2, 1, 3]
[2, 3, 1]
[3, 1, 2]
[3, 2, 1]

The algorithm simply finds the next permutation of the array. If the array is at the last permutation, the next combination from the set is constructed.

Note that this algorithm does not produce the k-permutations in lexicographic order.

unsigned int next_k_permutation(unsigned int *ar, size_t n, unsigned int k)
{
	unsigned int result = next_permutation(ar, k);
	if (result == 0) {
		result = next_combination(ar, n, k);
	}
	return result;
}

Subsets of a multiset

This is the set of all subsets of a multiset, which are themselves multisets, i.e., the power set of a multiset.

For example, the power set of the multiset [a, b, b, c] consists of the sets:

()
(c)
(b)
(b, c)
(b, b)
(b, b, c)
(a)
(a, c)
(a, b)
(a, b, c)
(a, b, b)
(a, b, b, c)

The set of all subsets of a particular size are combinations of a multiset.

The multiset and its subsets are represented as a vector containing, for each element, the count of its occurrences in the multiset. For example, the multiset [a, b, b, c] is represented as [1, 2, 1]. This is similar to the characteristic vector used for power set, but with counts rather than boolean values.

The correspondences for the subsets are then as follows:

[0, 0, 0] = ()
[0, 0, 1] = (c)
[0, 1, 0] = (b)
[0, 1, 1] = (b, c)
[0, 2, 0] = (b, b)
[0, 2, 1] = (b, b, c)
[1, 0, 0] = (a)
[1, 0, 1] = (a, c)
[1, 1, 0] = (a, b)
[1, 1, 1] = (a, b, c)
[1, 2, 0] = (a, b, b)
[1, 2, 1] = (a, b, b, c)

The algorithm can then simply count from [0, 0, 0] to [1, 2, 1], using the values in the source multiset as the upper limit for the value of an element.

Note that this algorithm does not produce the subsets in lexicographic order.

unsigned int next_multiset_subset(const unsigned int *multiset, unsigned int *ar, size_t n)
{
	unsigned int changed = 0;
	int i;

	for (i = n - 1; i >= 0 && !changed; i--) {
		if (ar[i] < multiset[i]) {
			/* Increment */
			ar[i]++;
			changed = 1;
		}
		else {
			/* Roll over */
			ar[i] = 0;
		}
	}
	if (!changed) {
		/* Reset to first combination */
		for (i = 0; i < n; i++) {
			ar[i] = 0;
		}
	}
	return changed;
}

Here is an example program:

#include <stdio.h>

#include <multiset-subset.h>

static void print_list(const unsigned int *ar, size_t len, FILE *fptr)
{
	unsigned int i;
	fputc('(', fptr);
	for (i = 0; i < len; i++) {
		fprintf(fptr, "%d", ar[i]);
		if (i < len - 1) {
			fputs(", ", fptr);
		}
	}
	fputc(')', fptr);
}

int main(void)
{
	unsigned int multiset[] = {1, 2, 1};
	const size_t n = sizeof(multiset) / sizeof(unsigned int);
	unsigned int numbers[] = {0, 0, 0};

	do {
		print_list(numbers, n, stdout);
		putchar('\n');
	} while (next_multiset_subset(multiset, numbers, n));

	return 0;
}

Here is a second example that prints the elements of the multisets:

#include <stdio.h>

#include <multiset-subset.h>

static void print_array(const unsigned int *ar, size_t len, FILE *fptr)
{
	unsigned int i;
	fputc('[', fptr);
	for (i = 0; i < len; i++) {
		fprintf(fptr, "%d", ar[i]);
		if (i < len - 1) {
			fputs(", ", fptr);
		}
	}
	fputc(']', fptr);
}


static void print_multiset_subset(const unsigned int *ar, size_t len, 
		const void **elements, printfn print, FILE *fptr)
{
	unsigned int i, started = 0;
	fputc('(', fptr);
	for (i = 0; i < len; i++) {
		unsigned int j;
		for (j = 0; j < ar[i]; j++) {
			if (started) {
				fputs(", ", fptr);
			}
			print(elements[i], fptr);
			started = 1;
		}
	}
	fputc(')', fptr); 
}


int main(void)
{
	unsigned int multiset[] = {1, 2, 1};
	const size_t n = sizeof(multiset) / sizeof(unsigned int);
	char *elements[] = {"a", "b", "c"};
	unsigned int numbers[] = {0, 0, 0};

	do {
        print_array(numbers, n, stdout);
        printf(" = "); 
		print_multiset_subset(numbers, n, (void*)elements, (printfn)fputs, stdout);
		putchar('\n');
	} while (next_multiset_subset(multiset, numbers, n));

	return 0;
}