Fix convention of PolynomialTensor basis change

src/openfermion/ops/representations/polynomial_tensor.py

@@ -31,8 +31,7 @@ def general_basis_change(general_tensor, rotation_matrix, key):
  r"""Change the basis of a general interaction tensor.
 
  M'^{p_1p_2...p_n} = R^{p_1}_{a_1} R^{p_2}_{a_2} ...
- R^{p_n}_{a_n} M^{a_1a_2...a_n} R^{p_n}_{a_n}^T ...
- R^{p_2}_{a_2}^T R_{p_1}_{a_1}^T
+ R^{p_n}_{a_n} M^{a_1a_2...a_n}
 
  where R is the rotation matrix, M is the general tensor, M' is the
  transformed general tensor, and a_k and p_k are indices. The formula uses
@@ -67,7 +66,7 @@ def general_basis_change(general_tensor, rotation_matrix, key):
 
  # Do the basis change through a single call of numpy.einsum. For example,
  # for the (1, 1, 0, 0) tensor, the call is:
- # numpy.einsum('abcd,aA,bB,cC,dD',
+ # numpy.einsum('abcd,Aa,Bb,Cc,Dd',
  # general_tensor,
  # rotation_matrix.conj(),
  # rotation_matrix.conj(),
@@ -79,7 +78,7 @@ def general_basis_change(general_tensor, rotation_matrix, key):
 
  # The 'Aa,Bb,Cc,Dd' part of the subscripts
  subscripts_rest = ','.join(
- chr(ord('a') + i) + chr(ord('A') + i) for i in range(order))
+ chr(ord('A') + i) + chr(ord('a') + i) for i in range(order))
 
  subscripts = subscripts_first + ',' + subscripts_rest
 

src/openfermion/ops/representations/polynomial_tensor_test.py

@@ -526,6 +526,47 @@ def test_rotate_basis_reverse(self):
  polynomial_tensor.rotate_basis(rotation_matrix_reverse)
  self.assertEqual(polynomial_tensor, want_polynomial_tensor)
 
+ def test_rotate_basis_90deg(self):
+ rotation_matrix_90deg = numpy.zeros((self.n_qubits, self.n_qubits))
+ rotation_matrix_90deg[0, 1] = -1
+ rotation_matrix_90deg[1, 0] = 1
+
+ one_body = numpy.zeros((self.n_qubits, self.n_qubits))
+ two_body = numpy.zeros(
+ (self.n_qubits, self.n_qubits, self.n_qubits, self.n_qubits))
+ one_body_90deg = numpy.zeros((self.n_qubits, self.n_qubits))
+ two_body_90deg = numpy.zeros(
+ (self.n_qubits, self.n_qubits, self.n_qubits, self.n_qubits))
+ i = 0
+ j = 0
+ i_90deg = pow(self.n_qubits, 2) - 1
+ j_90deg = pow(self.n_qubits, 4) - 1
+ for p in range(self.n_qubits):
+ for q in range(self.n_qubits):
+ one_body[p, q] = i
+ i = i + 1
+ one_body_90deg[p, q] = (-1)**([p, q].count(0))*i_90deg
+ i_90deg = i_90deg - 1
+ for r in range(self.n_qubits):
+ for s in range(self.n_qubits):
+ two_body[p, q, r, s] = j
+ j = j + 1
+ two_body_90deg[p, q, r, s] = \
+ (-1)**([p, q, r, s].count(0))*j_90deg
+ j_90deg = j_90deg - 1
+ polynomial_tensor = PolynomialTensor({
+ (): self.constant,
+ (1, 0): one_body,
+ (1, 1, 0, 0): two_body
+ })
+ want_polynomial_tensor = PolynomialTensor({
+ (): self.constant,
+ (1, 0): one_body_90deg,
+ (1, 1, 0, 0): two_body_90deg
+ })
+ polynomial_tensor.rotate_basis(rotation_matrix_90deg)
+ self.assertEqual(polynomial_tensor, want_polynomial_tensor)
+
  def test_rotate_basis_quadratic_hamiltonian_real(self):
  self.do_rotate_basis_quadratic_hamiltonian(True)
 
@@ -545,7 +586,7 @@ def do_rotate_basis_quadratic_hamiltonian(self, real):
  # Rotate a basis where the Hamiltonian is diagonal
  _, diagonalizing_unitary, _ = (
  quad_ham.diagonalizing_bogoliubov_transform())
- quad_ham.rotate_basis(diagonalizing_unitary.T)
+ quad_ham.rotate_basis(diagonalizing_unitary)
 
  # Check that the rotated Hamiltonian is diagonal with the correct
  # orbital energies
@@ -567,6 +608,8 @@ def test_rotate_basis_max_order(self):
  self.assertEqual(tensor, want_tensor)
  # I originally wanted to test 25 and 26, but it turns out that
  # numpy.einsum complains "too many subscripts in einsum" before 26.
+ # Currently, the sum of the number of output labels and combined labels
+ # can't not exceed NPY_MAXDIMS(=32).
 
  for order in [27, 28]:
  with self.assertRaises(ValueError):