Visualizing Combinatorial Enumeration

We can not only count the number of cases resulting from the application of combinatorial functions — we can visualize them as well.

by Robert Dickau

Combinatorics is the study of counting things, typically the number of ways one can select or arrange objects taken from some finite collection, possibly given one or more constraints. One characteristic of many formulas derived in combinatorics is that the results quickly become quite large (hence the term “combinatorial explosion”); a common example that demonstrates one advantage of using Mathematica over a pocket calculator, abacus, heap of stones, etc., is computing an exact result for 100 factorial.

100!

It can be somewhat illuminating, as well as satisfying, to see the situations that (at least the smaller, down-to-earth siblings of) these unimaginably large numbers describe; it's one thing to be told there are thousands of stars visible to the naked eye, and quite another to go outside at night and look up.

The purpose of this article is to enumerate and visualize situations described by elementary combinatorial functions. Note that Steven Skiena's standard add-on package DiscreteMath`Combinatorica` defines hundreds of functions related to combinatorics and graph theory; our purpose here is not just to count cases, but also to visualize them.

Permutations and the factorial function

One of the first combinatorial functions that a student learns is the factorial function, which describes the number of ways to arrange a collection of distinct objects. For example, the letters A, B, and C can be arranged, or permuted, in six ways—ABC, ACB, BAC, BCA, CAB, and CBA. The derivation of a formula for computing the number of permutations of n objects typically goes something like this: assuming there are n items and n slots, there are n choices of item for the first slot; n - 1 choices of item for the second slot; n - 2 choices for the third slot, and so forth; until there is only one item left, which must be placed in the last ((n^th)) slot. Using the “multiplication rule”, we conclude that for a collection of n objects there are \(TraditionalForm
      \`n×\((n - 1)\)×\((n - 2)
          \)×⋯×1\), otherwise known as n factorial (and denoted n!), different permutations.

Texts usually go on to enumerate (or leave as an exercise to enumerate) the 4! = 24 permutations of a collection of four objects, after which the matter is rightly dropped.

Mathematica has a built-in function called Permutations which enumerates all the permutations of a list of items. Here is the answer to the traditional homework problem.

Permutations[{A, B, C, D}]

Of course, Mathematica allows us to list all the permutations of any number of objects. Here is an abbreviated list of the 5,040 permutations of the numbers 1 through 7.

Permutations[Range[7]] // Short

The trouble is, permutations of numbers all look pretty much alike, or at least are equally uninteresting to read. Our first goal is to find a better way to visualize permutations.

One possibility is to graphically connect elements that are transposed in a particular permutation. Suppose we decide to represent the n elements being permuted by a number of points along a line, and connect two transposed elements with a semicircular arc. Here is a function that accepts a list of two positions, and returns a graphics object representing an arc joining the points; the semicircle is centered between the two points, and has a radius of half the distance between the points.

Semicircle[{m_, n_}] := 
      Circle[{(1/2) (m + n), 0}, (1/2) Abs[m - n], {0, π}]

We know from Mathematica experience that circles and perpendicular lines look wrong using the default setting for the graphics option AspectRatio, so we'll use the setting Automatic.

SetOptions[Graphics, AspectRatio -> Automatic];

Here is our first arc.

Show[Graphics[Semicircle[{1, 2}]]];

We would now like to draw an arc connecting adjacent values in a given permutation. Let's first consider the identity permutation for the list {1, 2, 3, 4}. We can generate a list of adjacent pairs using Partition, giving as a third argument the overlap between successive partitions.

Partition[{1, 2, 3, 4}, 2, 1]

We can then map our function Semicircle over the partitioned list, and display the result.

Show[Graphics[Map[Semicircle, %]]];

Similarly, we can display a figure for an arbitrary permutation of a list.

Show[Graphics[Map[Semicircle, Partition[{1, 3, 4, 2}, 2, 1]]]];

One shortcoming of this way of visualizing permutations is that some different permutations have identical diagrams. Here is the diagram corresponding to {4, 3, 2, 1}, which is a permutation of {1, 2, 3, 4}.

Show[Graphics[Map[Semicircle, Partition[{4, 3, 2, 1}, 2, 1]]]];

The diagram for {4, 3, 2, 1} looks suspiciously similar to the diagram for {1, 2, 3, 4}. One way to distinguish between the two is to use different colors for the arcs in a diagram, in order to indicate which pair of adjacent elements is being displayed.

Here is a list of arcs.

arcs = Map[Semicircle, Partition[{1, 2, 3, 4}, 2, 1]]

MapIndexed lets us associate a gray level with each arc as a function of its position. Because GrayLevel expects a value between 0 and 1, we divide the position of each arc by the number of elements in the list being permuted. (Of course there are only n - 1 arcs involved in a list of n elements, but dividing by n keeps us from generating GrayLevel[1], or white, arcs, which don't show up against the default white background.)

MapIndexed[{GrayLevel[First[#2]/4], #}&, arcs]

Here is our new diagram corresponding to {1, 2, 3, 4}, where the darkest arc connects the elements 1 and 2, the next-darkest arc connects 2 and 3, and so forth.

Show[Graphics[%]];

And here is the diagram for {4, 3, 2, 1}.

\(arcs = Map[Semicircle, Partition[{4, 3, 2, 1}, 2, 1]]; 

    Show[Graphics[MapIndexed[{GrayLevel[First[#2]/4], #}&, arcs]]]; \)

(Note that this technique doesn't distinguish the permutations {1, 2} and {2, 1}, a complication that we'll just have to try to live with.)

A further improvement is to put all the diagrams corresponding to a list of permutations into an array. Since there are n! diagrams in each array, and n! is divisible by n, we can partition the complete list of permutations into rows of n elements. Actually, to keep the arrays from being too tall (as opposed to too wide), we'll partition each complete list into rows of 2  n elements, which is safe because n! is divisible by 2  n, provided n is greater than 2.

Here we put everything together into one function.

\(ShowPermutationDiagrams[n_Integer /; n > 2] := 
	
      With[{allpairs = Map[Partition[#, 2, 1]&, Permutations[Range[n]]]}, 
        \(Show@\(GraphicsArray@Partition[#, 2  n]\)&\)@
          Map[Graphics[#, 
                PlotRange -> {Automatic, {\(-0.1\), n/2 + 0.1}}]&, 
            Map[MapIndexed[{GrayLevel[First[#2]/n], #}&, #]&, 
              Map[Semicircle, allpairs, {2}]]]]\)

Here are the diagrams for permutations of lists containing 3, 4, and 5 elements.

ShowPermutationDiagrams[3];

ShowPermutationDiagrams[4];

ShowPermutationDiagrams[5];

Do-it-yourself permutations

A useful goal when programming combinatorial functions is to devise a method for generating results recursively, that is, computing results for case k using the results for case k - 1. While it is always possible to generate a (possibly huge) list of candidate sets and remove those not meeting our criteria, it is often inefficent and uninstructive. (However, the following section contains an example of this technique.)

For example, the relationship \(TraditionalForm
      \`\(n!\) = n×\(\((n - 1)\)!\)\) suggests that a procedure may exist for generating permutations of the numbers 1 through n using the permutations of the numbers 1 through n - 1.

There are two permutations of the list {1, 2}.

Permutations[{1, 2}]

And there are six permutations of the list {1, 2, 3}.

Permutations[{1, 2, 3}]

A bit of experimentation shows that inserting the number 3 into every possible location in the permutation {1, 2} leads to the lists {1, 2, 3}, {1, 3, 2}, and {3, 1, 2}, and inserting 3 into {2, 1} gives {2, 1, 3}, {2, 3, 1}, and {3, 2, 1}. As we can verify, this process generates a complete list of permutations of n elements, based on a complete list of permutations of n - 1 elements.

ReplaceList is a programming function new in Mathematica 3. ReplaceList takes as arguments an expression and a list of transformation rules, and returns all possible results of transforming the expression. Nontrivial results arise when the transformation rules include double- or triple-blank patterns, which respectively match one or more, or zero or more, elements.

The following example uses a list of expressions matched with a sequence of two triple-blank patterns and returns the matched sequences with a new element inserted between them; the result is a list containing the original list with the new element inserted in each possible location.

ReplaceList[{old, old, old}, {a___, b___} -> {a, NEW, b}]

In general, a list with n elements has n + 1 locations for a new element.

The following definitions generate a list of all possible permutations of the numbers 1 through n, using the list of permutations of the numbers 1 through n - 1; the elementary, or starting, case is the empty list, which represents a complete list of permutations of a list with zero elements.

PermutationRule[n_] = {a___, b___} -> {a, n, b};

PermutationN[0] = {{}};

Defining the recursive step using the form f[x_] := f[x] = ⋯ allows us to cache intermediate values, which speeds up subsequent calculations.

\(PermutationN[n_Integer?Positive] := 
      \(PermutationN[n] = 
		
        Flatten[Map[ReplaceList[#, PermutationRule[n]]&, 
            PermutationN[n - 1]], 1]\)\)

Here are some examples.

PermutationN[2]

PermutationN[4]

As an informal check of our results, we make certain the function returns the same results as the built-in Permutations function.

Sort[%] == Permutations[{1, 2, 3, 4}]

Derangements

One special type of permutation is the derangement, which is a permutation that leaves no element in its original position. Assuming that we'll permute only lists of consecutive integers 1 through n, the following predicate function detects derangements by subtracting the candidate permutation from the identity permutation, and checking to see that the result is free of any zeros.

DerangementQ[perm_List] := FreeQ[Range[Length[perm]] - perm, 0]

The following permutation is a derangement of {1, 2, 3}.

DerangementQ[{2, 3, 1}]

The permutation {2, 1, 3} leaves 3 in its original position, so it is not a derangement of {1, 2, 3}.

DerangementQ[{2, 1, 3}]

We can then enumerate the derangements of a list by enumerating the permutations and removing those that are not derangements.

Select[PermutationN[3], DerangementQ]

Select[PermutationN[4], DerangementQ]

At this point, we might as well define a function Derangements. (Note that this function does not generate the lists of derangements recursively; such a recursion is complicated enough to negate its educational benefits.)

Derangements[n_] := Select[PermutationN[n], DerangementQ]

We can then modify the body of PermutationDiagrams to visualize the list of derangements.

With[{allpairs = Map[Partition[#, 2, 1]&, Derangements[4]]}, 
      Show@GraphicsArray@Partition[#, 3]&@
        Map[Show[Graphics[#], DisplayFunction -> Identity, 
              PlotRange -> {Automatic, {-0.1, 2.1}}]&, 
          Map[MapIndexed[{GrayLevel[First[#2]/4], #}&, #]&, 
            Map[Semicircle, allpairs, {2}]]]];

The diagrams aren't particularly interesting, though, since they don't really show how derangements leave no element in its original position; and so it's time to think of another approach.

Another way to represent the elements in a list is as indexed points equally spaced around the perimeter of a circle. Since this is a representation that we'll use several times, we define a function Poly that generates such a list of coordinates.

Poly[n_] := 
      {Cos[#], Sin[#]}&/@Range[0, 2 π, ((2 π)/n)]

(Note that Poly creates a list whose first and last elements are the same, which comes in handy when we're drawing polygons using the Line graphics primitive; otherwise there would be an edge missing from the polygon, which gets to be a nuisance after you've had to correct for it a few times.)

Here we verify that Poly[5] represents a regular pentagon.

Show[Graphics[Line[Poly[5]]]];

We can then draw lines or arrows connecting the vertices corresponding to related elements in a list. For example, we might represent the relationship between the identity list {1, 2, 3} and its derangement {3, 1, 2} by drawing arrows connecting an element's position in the identity list to its position in the “deranged” list.

First we load the package that defines the Arrow graphics primitive, and make the arrowheads easier to see by giving new values to the HeadCenter and HeadWidth options.

Needs["Graphics`Arrow`"];

SetOptions[Arrow, HeadCenter -> 1, HeadWidth -> 1];

Here is our diagram for the derangement {3, 1, 2}, in which the arrows indicate that the element 3 moves to position 1, element 1 moves to position 2, and element 2 moves to position 3. (In this and subsequent diagrams, the leftmost vertex in a polygon represents element 1, and the other elements are arranged counterclockwise around the center of the polygon.)

tri = Poly[3]; 
    Show[Graphics[{
          Arrow[tri[[3]], 
            tri[[1]]], 
          Arrow[tri[[1]], 
            tri[[2]]], 
          Arrow[tri[[2]], 
            tri[[3]]]}]];

Here is a function that generates a graphics object for any derangement, additionally using a different color for each arrow.

DerangementGraphic[L_List] := 
      Module[{len = Length[L], thisPoly = Poly[Length[L]], endpoints}, 
        endpoints = Transpose[{L, Range[len]}]; 
        Graphics[
          MapIndexed[{Hue[First[#2]/len], Arrow@@Part[thisPoly, #]}&, 
            endpoints]]]

Here is our new picture for the list {3, 1, 2}.

Show[DerangementGraphic[{3, 1, 2}]];

Here is a function that shows all the derangements for the list {1, 2, …,  n}. (It turns out that the number of derangements of a list satisfies the recurrence relation \(TraditionalForm
      \`D\_n = \((n - 1)\) \((D\_\(n - 1\) + D\_\(n - 2\))\)\), so we can safely partition the list of derangements into an array with rows n - 1 elements wide.)

ShowDerangementDiagrams[n_Integer /; n > 2] := 
      Show[GraphicsArray[
          Partition[Map[DerangementGraphic, Derangements[n]], n - 1]]]

Here are the derangements of {1, 2, 3, 4}.

ShowDerangementDiagrams[4];

An improvement (arguably) is to redefine DerangementGraphic so that it places a border around the colored arrows, and to redefine ShowDerangementGraphics to transpose the graphics array, making it wide, rather than tall.

DerangementGraphic[L_List] := 
      Module[{len = Length[L], thisPoly = Poly[Length[L]], endpoints}, 
        endpoints = Transpose[{L, Range[len]}]; 
        Graphics[{GrayLevel[0.5], Line[1.25  Poly[len]], 
            MapIndexed[{Hue[First[#2]/len], Arrow@@Part[thisPoly, #]}&, 
              endpoints]}]]

ShowDerangementDiagrams[n_Integer /; n > 2] := 
      Show[GraphicsArray[
          Transpose[
            Partition[Map[DerangementGraphic, Derangements[n]], n - 1]]]]

ShowDerangementDiagrams[5];

Bell numbers

The Bell numbers (denoted (Subscript[B, n])) are the numbers of ways to partition a set into disjoint, non-empty subsets. For example, the set {1, 2} can be partitioned in the following two ways.

({{{{1}, {2}}}, {{{1, 2}}}})

Similarly, the set {1, 2, 3} can be partitioned in the following five ways.

({{{{1}, {2}, {3}}}, {{{1, 2}, {3}}}, {{{1, 3}, {2}}}, {{{1}, {2, 3}}}, {{{1, 2, 3}}}})

We can recursively compute the lists of Bell partitions by noting that the partitions of a set with n elements can be created by appending the singleton {n} to each list of partitions of n - 1 elements, or by appending the element n to each partition of n - 1 elements. An example: starting with the two partitions of {1, 2}, we can create the partitions {{1}, {2}, {3}} and {{1, 2}, {3}} using the first case, and the partitions {{1, 3}, {2}},  {{1}, {2, 3}}, and {{1, 2, 3}} using the second case.

First we translate our two cases into Mathematica rules.

\(BellRule1[n_] = {S__} -> {S, {n}}; 

    BellRule2[n_] = {b___, {S__}, a___} -> {b, {S, n}, a}; \)

Here we generate the list of partitions making up (Subscript[B, 3]) given the list of partitions making up (Subscript[B, 2]).

B2 = {{{1}, {2}}, {{1, 2}}};

These are the partitions generated from the first rule.

Map[ReplaceList[#, BellRule1[3]]&, B2]

And these are the partitions generated from the second rule.

Map[ReplaceList[#, BellRule2[3]]&, B2]

The following definitions generate the complete list of Bell partitions.

BellLists[1] = {{{1}}};

\(BellLists[n_Integer?Positive] := 
      \(BellLists[n] = 
		
        Flatten[
				
          Map[ReplaceList[#, {BellRule1[n], BellRule2[n]}]&, 
            BellLists[n - 1]], 1]\)\)

Here is the complete list of partitions making up (Subscript[B, 3]).

BellLists[3] // ColumnForm

And here is the complete list making up (Subscript[B, 4]).

BellLists[4] // ColumnForm

To visualize the partitions contributing to the (n^th) Bell number, we arrange n points (representing the n elements) along the circumference of a circle, and then connect points in the same subset with a line loop. (For a singleton {k} in a Bell list, we display a point.) Because the diagrams represent unordered lists, it is not necessary to use colors or arrows to distinguish different Bell lists.

Loop[ptL_List] := Line[Append[ptL, First[ptL]]]

BellDiagram[L_List] := 
      Module[{thisPoly = Poly[Max[L]], frm}, 
        frm = {GrayLevel[0.5], Line[1.5  thisPoly]}; 
        Graphics[{frm, 
            Map[Loop, 
              Map[thisPoly[[#]]&, L]]}]]

Here is the diagram corresponding to the Bell list {{1, 2}, {3}}.

Show[BellDiagram[{{1, 2}, {3}}]];

We next define a function to display the array of all Bell lists. In this case, it is more difficult to decide how to partition the graphics arrays because there are no simple divisors of the numbers of Bell lists. A technique that works for the arrays that we're likely to want to see is to factor the number of elements in the array, and then partition the array into rows whose length is equal to the largest prime divisor of the number of elements. (This technique fails for the 877 partitions in BellLists[7], but even we have our limits.)

PartitionGraphArray[
        L_List /; (PrimeQ[Length[L]] || Length[L] == 1)] := L

PartitionGraphArray[L_List] := 
      Partition[L, First[Last[FactorInteger[Length[L]]]]]

Here is our ShowBellDiagrams function.

ShowBellDiagrams[n_Integer /; n > 2] := 
      Show@GraphicsArray@
          PartitionGraphArray@Map[BellDiagram, BellLists[n]]

Here are the diagrams corresponding to the third, fourth, and fifth Bell numbers.

ShowBellDiagrams[3];

Bell4Plot = ShowBellDiagrams[4];

ShowBellDiagrams[5];

Stirling set numbers

The Stirling set numbers (often called Stirling numbers of the second kind, denoted {({{n}, {k}})}) describe the number of ways a set with n elements can be partitioned into k disjoint, non-empty subsets. For example, the set {1, 2, 3} can be partitioned into three subsets in one way — {{1}, {2}, {3}} — into two subsets in the following three ways — {{1, 2}, {3}},  {{1, 3}, {2}}, and {{1}, {2, 3}} — and into one subset in the following way — {{1, 2, 3}}. (The (n^th) Bell number is thus the sum of the (((n, k))^th) Stirling set numbers for k = 1, 2, …, n.)

To recursively generate the Stirling set numbers, we use a technique similar to that used for generating Bell numbers: for a given list of partitions of n - 1 elements, we add the singleton {n} to the partitions making up {({{n - 1}, {k - 1}})}, and add the element n to each list in the partitions making up {({{n - 1}, {k}})}. Here are the Mathematica translations of the two methods.

\(StirlingSetRule1[n_] = {S__} -> {S, {n}}; 

    StirlingSetRule2[n_] = {b___, {S__}, a___} -> {b, {S, n}, a}; \)

In this case, we give two “boundary conditions”.

StirlingSetLists[n_Integer?NonNegative, 1] := {{Range[n]}}

StirlingSetLists[n_Integer?NonNegative, n_] := {Partition[Range[n], 1]}

\(StirlingSetLists[n_, k_] := 
      \(StirlingSetLists[n, k] = 
		
        Flatten[Join[
			
            Map[ReplaceList[#, StirlingSetRule1[n]]&, 
              StirlingSetLists[n - 1, k - 1]], 
			
            Map[ReplaceList[#, StirlingSetRule2[n]]&, 
              StirlingSetLists[n - 1, k]]], 1]\)\)

Here are the partitions corresponding to {({{4}, {2}})}.

StirlingSetLists[4, 2] // ColumnForm

The functions for displaying a single Stirling number of the second kind diagram and for displaying the entire array are nearly identical to the code for displaying the Bell list diagrams.

StirlingSetDiagram[L_List] := 
      Module[{thisPoly = Poly[Max[L]], frm}, 
        frm = {GrayLevel[0.5], Line[1.5  thisPoly]}; 
        Graphics[{frm, 
            Map[Loop, 
              Map[thisPoly[[#]]&, L]]}]]

ShowStirlingSetDiagrams[n_, k_] := 
      Show[GraphicsArray[
          PartitionGraphArray[
            Map[StirlingSetDiagram, StirlingSetLists[n, k]]]]]

Here are our diagrams.

ShowStirlingSetDiagrams[4, 1];

ShowStirlingSetDiagrams[4, 2];

ShowStirlingSetDiagrams[4, 3];

ShowStirlingSetDiagrams[4, 4];

As we can see, (∑_n){({{4}, {n}})}=(Subscript[B, 4]).

S4kLists = Flatten[Map[StirlingSetLists[4, #]&, {1, 2, 3, 4}], 1];

S4kPlot = 
      GraphicsArray[
        Partition[Reverse[Map[StirlingSetDiagram, S4kLists]], 5]];

Show[GraphicsArray[{Bell4Plot, S4kPlot}]];

Bibliography and related reading

J. H. Conway and R. K. Guy, The Book of Numbers. Copernicus (Springer-Verlag), 1996.

R. L. Graham, D. E. Knuth, and O. Patashnik, Concrete Mathematics: A Foundation for Computer Science. Addison-Wesley, 1994.

F. S. Roberts, Applied Combinatorics. Prentice-Hall, 1984.

S. Skiena, Implementing Discrete Mathematics: Combinatorics and Graph Theory with Mathematica. Addison-Wesley, 1990.

E. W. Weisstein, CRC Concise Encyclopedia of Mathematics. CRC Press, 1998. Portions of the book are currently served from http://www.astro.virginia.edu/~eww6n/math/math0.html. This site is the reason the Web was invented.

About the author

Robert M. Dickau is a former member of technical staff at Wolfram Research, and is currently senior technical trainer for InstallShield Software Corporation in Schaumburg, Illinois. He holds degrees in mathematics and actuarial science from California Polytechnic University and University of Illinois at Urbana-Champaign.

E-Mail: pops@prairienet.org

ELECTRONIC SUBSCRIPTIONS

Included in the distribution for each electronic subscription is the file VisCombNum.nb containing Mathematica code for the material described in this article.