## Computing a Stirling Number of the Second Kind from Scratch Using Python

A “Stirling number of the second kind” — S(n, k) — is the number of ways to partition n items into k subsets. For example, the number of ways to group n=4 items into k=2 subsets is 7. Suppose the 4 items are a, b, c, d. Then S(4,2) = 7. The seven partitions are: (a | bcd), (b | acd), (c | abd), (d | abc), (ab |cd), (ac | bd), and (ad | bc).

The value of S(n, k) gets astronomically large, very quickly. For example, S(100, 3) =

85,896,253,455,335,221,205,584,888,180,155,511,368,666,317,646

I’m sure there is certainly a Stirling function in one of the Python libraries (probably SciPy or NumPy) but I decided to implement a Stirling function from scratch using raw Python, without any libraries.

I got one of many math definitions for S(n,k) from Wikipedia:

Next, I had to worry about dealing with huge integer values. But one of the extremely nice things about the Python language is that is has built-in support for arbitrarily large integers. Dealing with huge integers in many languages is a huge pain point.

The math definition above uses the Factorial function (outside the summation), an integer Power function (used in two places), and a Choose function. The Choose(n,k) function returns the number of ways to choose k items from n items, where order does not matter. For example, Choose(4,2) = 6 because there are six ways to select 2 items from 4 items. Suppose the four items are q, r, s, t. Then the six ways are (q r), (q s), (q t), (r s), (r t), and (s t).

I implemented functions my_factorial(n), my_power(n, p), and my_choose(n, k) and then used them to implement a my_stirling_second(n, k) function.

The moral of the story is that you don’t have to be a software deity of some sort to implement functions. Except for a few rare situations, most functions can be implemented by mere mortals.

Three illustrations of ancient Egyptian deities by an artist named “Yliade”. From left to right: Seshat, Nephthys, Serqet. I like the modern style of art applied to an ancient theme.

```# stirling_second_demo.py
# note: // is integer division

def my_stirling_second(n, k):
sum = 0
for j in range(0, k+1):
a = my_power(-1, k-j)
b = my_choose(k, j)
c = my_power(j, n)
sum += a * b * c

return sum // my_factorial(k)

def my_factorial(n):
ans = 1
for i in range(1,n+1):
ans *= i
return ans

def my_power(n, p):
# n^p
ans = 1
for i in range(p):
ans *= n
return ans

def my_choose(n, k):
if n == k: return 1
if k "less-than" n-k:
delta = n - k
i_max = k
else:
delta = k
i_max = n - k

ans = delta + 1
for i in range(2, i_max+1):
ans = (ans * (delta + i)) // i
return ans

def main():
print("\nBegin Stirling numbers second kind demo ")

nf = my_factorial(10)
print("\nmy_factorial(10) = ")
print("{:,}".format(nf))

pwr = my_power(5, 3)
print("\nmy_power(5,3) = ")
print("{:,}".format(pwr))

c = my_choose(10, 3)
print("\nmy_choose(10,3) = ")
print("{:,}".format(c))

sn = my_stirling_second(100, 3)
print("\nmy_stirling_second(100,3) = ")
print("{:,}".format(sn))

print("\nEnd demo \n")

if __name__ == "__main__":
main()
```
This entry was posted in Miscellaneous. Bookmark the permalink.