This note was originally written as an IPython notebook.

We consider a Hilbert space that consists of two Transmon “qubits” (qudits, actually) with `nq`

levels, and a cavity with `nc`

levels. Given the dimension of the total Hilbert space, `n = nq*nq*nc`

, we want to find `nq`

and `nc`

, i.e. the tensor structure of the full Hilbert space. We demand that `nq <= nc`

.

The decomposition is based on the prime factors of `n`

:

```
def prime_factors(n):
"""Return the prime decomposition of n"""
i = 2
factors = []
while i * i <= n:
if n % i:
i += 1
else:
n //= i
factors.append(i)
if n > 1:
factors.append(n)
return factors
>>> n = 6*6*20
>>> factors = prime_factors(n)
>>> print factors
[2, 2, 2, 2, 3, 3, 5]
```

Now we put them in two list: the prime factors of the qubit Hilbert space, and the prime factors of the cavity Hilbert space

```
qubit_factors = []
cavity_factors = []
for factor in set(factors):
count = factors.count(factor)
while count >= 2:
qubit_factors.append(factor)
count += -2
if count > 0:
cavity_factors.append(factor)
>>> qubit_factors
[2, 2, 3]
>>> cavity_factors
[5]
```

We demand that the qubit dimension is not larger than the cavity dimension. In order to assure this, we start moving the smallest prime factors from the qubit to the cavity

```
import numpy as np
while np.array(qubit_factors).prod() > np.array(cavity_factors).prod():
factor = qubit_factors.pop(0)
cavity_factors = [factor, factor] + cavity_factors
>>> qubit_factors
[2, 3]
>>> cavity_factors
[2, 2, 5]
>>> nq = np.array(qubit_factors).prod()
>>> nq
6
>>> nc = np.array(cavity_factors).prod()
>> nc
20
```

The total algorithm looks like this (using `prime_factors`

, and also allowing to specify the cavity dimension manually):

```
def qubit_cavity_decomposition(n, nc=None):
"""
Return (nq, nc) of a decomposition n = nq*nq*nc, where nc >= nq
"""
if nc is not None:
n = n // nc
factors = prime_factors(n)
qubit_factors = []
cavity_factors = []
for factor in set(factors):
count = factors.count(factor)
while count >= 2:
qubit_factors.append(factor)
count += -2
if count > 0:
cavity_factors.append(factor)
# Ensure nc >= nq
if nc is None:
while np.array(qubit_factors).prod() > np.array(cavity_factors).prod():
factor = qubit_factors.pop(0)
cavity_factors = [factor, factor] + cavity_factors
nc = np.array(cavity_factors).prod()
nq = np.array(qubit_factors).prod()
return nq, nc
>>> qubit_cavity_decomposition(5*5*10)
(5, 10)
>>> qubit_cavity_decomposition(5*5*25)
(5, 25)
>>> qubit_cavity_decomposition(10*10*20)
(10, 20)
>>> qubit_cavity_decomposition(8*8*25)
(10, 16)
>>> qubit_cavity_decomposition(8*8*25, nc=25)
(8, 25)
>>> qubit_cavity_decomposition(8*8*30)
(8, 30)
```