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:

In [1]:
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
In [2]:
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

In [3]:
qubit_factors = []
cavity_factors = []
In [4]:
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)
In [5]:
qubit_factors
Out[5]:
[2, 2, 3]
In [6]:
cavity_factors
Out[6]:
[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

In [7]:
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
In [8]:
qubit_factors
Out[8]:
[2, 3]
In [9]:
cavity_factors
Out[9]:
[2, 2, 5]
In [10]:
nq = np.array(qubit_factors).prod()
nq
Out[10]:
6
In [11]:
nc = np.array(cavity_factors).prod()
nc
Out[11]:
20

The total algorithm looks like this (using prime_factors, and also allowing to specify the cavity dimension manually):

In [12]:
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
In [13]:
qubit_cavity_decomposition(5*5*10)
Out[13]:
(5, 10)
In [14]:
qubit_cavity_decomposition(5*5*25)
Out[14]:
(5, 25)
In [15]:
qubit_cavity_decomposition(10*10*20)
Out[15]:
(10, 20)
In [16]:
qubit_cavity_decomposition(8*8*25)
Out[16]:
(10, 16)
In [17]:
qubit_cavity_decomposition(8*8*25, nc=25)
Out[17]:
(8, 25)
In [18]:
qubit_cavity_decomposition(8*8*30)
Out[18]:
(8, 30)