# Efficient algorithm for getting number of partitions of integer with distinct parts (Partition function Q)

**Tested two algorithms**

Simple recurrence relation

WolframMathword algorithm (based upon Georgiadis, Kediaya, Sloane)

Both implemented with Memoization using LRUCache.

Results: WolframeMathword approach orders of magnitude faster.

**1. Simple recurrence relation (with Memoization)**

Reference

*Code*

```
@lru_cache(maxsize=None)
def p(n, d=0):
if n:
return sum(p(n-k, n-2*k+1) for k in range(1, n-d+1))
else:
return 1
```

*Performance*

```
n Time (sec)
10 time elapsed: 0.0020
50 time elapsed: 0.5530
100 time elapsed: 8.7430
200 time elapsed: 168.5830
```

**2. WolframMathword algorithm**

(based upon Georgiadis, Kediaya, Sloane)

Reference

**Code**

```
# Implementation of q recurrence
# https://mathworld.wolfram.com/PartitionFunctionQ.html
class PartitionQ():
def __init__(self, MAXN):
self.MAXN = MAXN
self.j_seq = self.calc_j_seq(MAXN)
@lru_cache
def q(self, n):
" Q strict partition function "
assert n < self.MAXN
if n == 0:
return 1
sqrt_n = int(sqrt(n)) + 1
temp = sum(((-1)**(k+1))*self.q(n-k*k) for k in range(1, sqrt_n))
return 2*temp + self.s(n)
def s(self, n):
if n in self.j_seq:
return (-1)**self.j_seq[n]
else:
return 0
def calc_j_seq(self, MAX_N):
""" Used to determine if n of form j*(3*j (+/-) 1) / 2
by creating a dictionary of n, j value pairs "
result = {}
j = 0
valn = -1
while valn <= MAX_N:
jj = 3*j*j
valp, valn = (jj - j)//2, (jj+j)//2
result[valp] = j
result[valn] = j
j += 1
return result
```

**Performance**

```
n Time (sec)
10 time elapsed: 0.00087
50 time elapsed: 0.00059
100 time elapsed: 0.00125
200 time elapsed: 0.10933
```

Conclusion: This algorithm is orders of magnitude faster than the simple recurrence relationship

**Algorithm**

Reference

I think a straightforward and efficient way to solve this is to explicitly compute the coefficient of the generating function from the Wolfram PartitionsQ link in the original post.

This is a pretty illustrative example of how to construct generating functions and how they can be used to count solutions. To start, we recognize that the problem may be posed as follows:

```
Let m_1 + m_2 + ... + m_{n-1} = n where m_j = 0 or m_j = j for all j.
Q(n) is the number of solutions of the equation.
```

We can find `Q(n)`

by constructing the following polynomial (i.e. the generating function)

```
(1 + x)(1 + x^2)(1 + x^3)...(1 + x^(n-1))
```

The number of solutions is the number of ways the terms combine to make `x^n`

, i.e. the coefficient of `x^n`

after expanding the polynomial. Therefore, we can solve the problem by simply performing the polynomial multiplication.

```
def Q(n):
# Represent polynomial as a list of coefficients from x^0 to x^n.
# G_0 = 1
G = [int(g_pow == 0) for g_pow in range(n + 1)]
for k in range(1, n):
# G_k = G_{k-1} * (1 + x^k)
# This is equivalent to adding G shifted to the right by k to G
# Ignore powers greater than n since we don't need them.
G = [G[g_pow] if g_pow - k < 0 else G[g_pow] + G[g_pow - k] for g_pow in range(n + 1)]
return G[n]
```

Timing (average of 1000 iterations)

```
import time
print("n Time (sec)")
for n in [10, 50, 100, 200, 300, 500, 1000]:
t0 = time.time()
for i in range(1000):
Q(n)
elapsed = time.time() - t0
print('%-5d%.08f'%(n, elapsed / 1000))
```

```
n Time (sec)
10 0.00001000
50 0.00017500
100 0.00062900
200 0.00231200
300 0.00561900
500 0.01681900
1000 0.06701700
```