Dimension mismatch with numpy

Question

I am trying to run a code for the following Hamiltoninan:

H = $\sum_{i<j} J_{ij}(\sigma_i^+ \sigma_j^- + \sigma_i^- \sigma_j^+) + \sum_j(B + B_j)\sigma_j^z$

The code is as follows:

import numpy as np
import qutip as qt
from scipy.constants import hbar

# Qubit-qubit coupling matrix:

alpha = 1.47
constant = 67.3765

def create_matrix(rows, columns):
    matrix = []
    for i in range(rows):
        row = []
        for j in range(columns):
            if i == j:
                entry = 0
                row.append(entry)
            else:
                entry = constant*np.power(abs(i - j), -alpha)
                row.append(entry)
        matrix.append(row)
    return matrix

num_qubits = 8

J_ij = create_matrix(num_qubits, num_qubits)
J_ij = np.matrix(J_ij)

# Magnetic field:

B_i = np.array([[2.85884931, 26.22601547, 48.5376065, 62.948092,
                 53.7557919, 37.31269595, 10.50862743, 0.1]])

B = 15707.96326795

# Hamiltonian:
# Define the spin operators for each ion
spin_operators = [qt.jmat(0.5, 'x'), qt.jmat(0.5, 'y'), qt.jmat(0.5, 'z')]

# Define the raising and lowering operators for each ion
s_plus = [qt.tensor([qt.qeye(2)] * i + [qt.create(2)] + [qt.qeye(2)] * (num_qubits - i - 1)) for i in range(num_qubits)]
s_minus = [qt.tensor([qt.qeye(2)] * i + [qt.destroy(2)] + [qt.qeye(2)] * (num_qubits - i - 1)) for i in range(num_qubits)]
s_z = [qt.tensor([qt.qeye(2)] * i + [spin_operators[0]] + [qt.qeye(2)] * (num_qubits - i - 1)) for i in range(num_qubits)]

# XY interaction term
interaction_term = sum([J_ij[i, j] * (s_plus[i] * s_minus[j] + s_minus[i] * s_plus[j]) for i in range(num_qubits-1) for j in range(i+1, num_qubits)])

# Transverse magnetic field term
#transverse_field_term = hbar*B * sum(s_z)
transverse_field_term = B*sum(s_z)

# Define the Hamiltonian
H = interaction_term + transverse_field_term
H = H.full()

I am providing the full code so it is reproducible. The problem is with the line transverse_field_term = B*sum(s_z) which I need to modify. I was working with H = $\sum_{i<j} J_{ij}(\sigma_i^+ \sigma_j^- + \sigma_i^- \sigma_j^+) + B\sum_j\sigma_j^z$ but the last term is now \sum_j(B + B_j)\sigma_j^z which I tried to write in the code as:

# Transverse magnetic field term
transverse_field_term = sum((B + B_i[i]) * s_z[i] for i in range(num_qubits))

But I got the following error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
 42 # Transverse magnetic field term
---> 43 transverse_field_term = sum((B + B_i[i]) * s_z[i] for i in range(num_qubits))
     45 # Define the Hamiltonian

<ipython-input-11-c550e8be0fe1> in <genexpr>(.0)
     41 
     42 # Transverse magnetic field term
---> 43 transverse_field_term = sum((B + B_i[i]) * s_z[i] for i in range(num_qubits))
     44 
     45 # Define the Hamiltonian

~/.virtualenvs/cgan_tf/lib/python3.10/site-packages/qutip/qobj.py in __rmul__(self, other)
     633                 return out
    634             else:
--> 635                 return other * self.data
    636 
    637         elif isinstance(other, list):

 ~/.virtualenvs/cgan_tf/lib/python3.10/site-packages/scipy/sparse/_base.py in      __rmul__(self, other)
    618 
    619     def __rmul__(self, other):  # other * self
--> 620         return self._rmul_dispatch(other)
    621 
    622     #######################

~/.virtualenvs/cgan_tf/lib/python3.10/site-packages/scipy/sparse/_base.py in _rmul_dispatch(self, other)
    612             except AttributeError:
    613                 tr = np.asarray(other).transpose()
--> 614             ret = self.transpose()._mul_dispatch(tr)
    615             if ret is NotImplemented:
    616                 return NotImplemented

~/.virtualenvs/cgan_tf/lib/python3.10/site-packages/scipy/sparse/_base.py in _mul_dispatch(self, other)
    557             # dense row or column vector
    558             if other.shape != (N,) and other.shape != (N, 1):
--> 559                 raise ValueError('dimension mismatch')
    560 
    561             result = self._mul_vector(np.ravel(other))

ValueError: dimension mismatch

I tried to correct it with

# Reshape B_i to match the shape of s_z
B_i_reshaped = np.reshape(B_i, (num_qubits, 1))

# Transverse magnetic field term
transverse_field_term = -B * sum((B_i_reshaped[i] + B) * s_z[i] for i in range(num_qubits))

But I got the same error.

How can I make it right?

Answer 1

From a quick read of qutib, I'm guessing it makes small square sparse matrices.

In [1]: from scipy import sparse
In [2]: M = sparse.random(4,4,.2, 'csr')

If I define another array, with a wrong shape, I get your error:

In [3]: B = numpy.array([1,2,3])
In [4]: B*M
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[4], line 1
----> 1 B*M

File ~/.local/lib/python3.10/site-packages/scipy/sparse/_matrix.py:51, in spmatrix.__rmul__(self, other)
     50 def __rmul__(self, other):
---> 51     return self._rmul_dispatch(other)

File ~/.local/lib/python3.10/site-packages/scipy/sparse/_base.py:608, in _spbase._rmul_dispatch(self, other)
    606 except AttributeError:
    607     tr = np.asarray(other).transpose()
--> 608 ret = self.transpose()._mul_dispatch(tr)
    609 if ret is NotImplemented:
    610     return NotImplemented

File ~/.local/lib/python3.10/site-packages/scipy/sparse/_base.py:553, in _spbase._mul_dispatch(self, other)
    550 if other.ndim == 1 or other.ndim == 2 and other.shape[1] == 1:
    551     # dense row or column vector
    552     if other.shape != (N,) and other.shape != (N, 1):
--> 553         raise ValueError('dimension mismatch')
    555     result = self._mul_vector(np.ravel(other))
    557     if isinstance

(other, np.matrix):

ValueError: dimension mismatch

Matching shape,

In [5]: B = numpy.array([1,2,3,4])
In [6]: B*M
Out[6]: array([0.90670377, 2.04396223, 0.        , 1.71075771])

s_z is a list of qt terms. sum(s_z) is presumably a simple sparse matrix; or atleast one that multiplies as such.

B is a scalar, so B*sum(s_z) works

transverse_field_term = sum((B + B_i[i]) * s_z[i] for i in range(num_qubits))

(B + B_i[i]) * s_z[i]

B_i is (1,8), so B_i[0] is (8,). So unless s_z[0] is also (8,8) you'll get the shape mismatch.

pseudo math expressions (latex?) are not displayed nicely with SO:

H = $\sum_{i<j} J_{ij}(\sigma_i^+ \sigma_j^- + \sigma_i^- \sigma_j^+) + B\sum_j\sigma_j^z$