1# Licensed under the Apache License, Version 2.0 (the "License"); 2# you may not use this file except in compliance with the License. 3# You may obtain a copy of the License at 4# 5# http://www.apache.org/licenses/LICENSE-2.0 6# 7# Unless required by applicable law or agreed to in writing, software 8# distributed under the License is distributed on an "AS IS" BASIS, 9# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 10# See the License for the specific language governing permissions and 11# limitations under the License. 12"""Tests _bravyi_kitaev_fast_test.py.""" 13 14import os 15import unittest 16 17import numpy 18 19from openfermion.config import DATA_DIRECTORY 20from openfermion.chem import MolecularData 21from openfermion.ops.operators import FermionOperator, QubitOperator 22from openfermion.ops.representations import InteractionOperator 23from openfermion.transforms.opconversions import (bksf, get_fermion_operator, 24 normal_ordered) 25from openfermion.transforms.opconversions.jordan_wigner import ( 26 jordan_wigner, jordan_wigner_one_body) 27from openfermion.linalg import get_sparse_operator, eigenspectrum 28from openfermion.utils import count_qubits 29 30 31class bravyi_kitaev_fastTransformTest(unittest.TestCase): 32 33 def setUp(self): 34 geometry = [('H', (0., 0., 0.)), ('H', (0., 0., 0.7414))] 35 basis = 'sto-3g' 36 multiplicity = 1 37 filename = os.path.join(DATA_DIRECTORY, 'H2_sto-3g_singlet_0.7414') 38 self.molecule = MolecularData(geometry, 39 basis, 40 multiplicity, 41 filename=filename) 42 self.molecule.load() 43 44 # Get molecular Hamiltonian. 45 self.molecular_hamiltonian = self.molecule.get_molecular_hamiltonian() 46 47 # Get FCI RDM. 48 self.fci_rdm = self.molecule.get_molecular_rdm(use_fci=1) 49 # Get explicit coefficients. 50 self.nuclear_repulsion = self.molecular_hamiltonian.constant 51 self.one_body = self.molecular_hamiltonian.one_body_tensor 52 self.two_body = self.molecular_hamiltonian.two_body_tensor 53 54 # Get fermion Hamiltonian. 55 self.fermion_hamiltonian = normal_ordered( 56 get_fermion_operator(self.molecular_hamiltonian)) 57 58 # Get qubit Hamiltonian. 59 self.qubit_hamiltonian = jordan_wigner(self.fermion_hamiltonian) 60 61 # Get the sparse matrix. 62 self.hamiltonian_matrix = get_sparse_operator( 63 self.molecular_hamiltonian) 64 65 def test_bad_input(self): 66 with self.assertRaises(TypeError): 67 bksf.bravyi_kitaev_fast(FermionOperator()) 68 69 def test_bravyi_kitaev_fast_edgeoperator_Bi(self): 70 # checking the edge operators 71 edge_matrix = numpy.triu(numpy.ones((4, 4))) 72 edge_matrix_indices = numpy.array( 73 numpy.nonzero( 74 numpy.triu(edge_matrix) - numpy.diag(numpy.diag(edge_matrix)))) 75 76 correct_operators_b0 = ((0, 'Z'), (1, 'Z'), (2, 'Z')) 77 correct_operators_b1 = ((0, 'Z'), (3, 'Z'), (4, 'Z')) 78 correct_operators_b2 = ((1, 'Z'), (3, 'Z'), (5, 'Z')) 79 correct_operators_b3 = ((2, 'Z'), (4, 'Z'), (5, 'Z')) 80 81 qterm_b0 = QubitOperator(correct_operators_b0, 1) 82 qterm_b1 = QubitOperator(correct_operators_b1, 1) 83 qterm_b2 = QubitOperator(correct_operators_b2, 1) 84 qterm_b3 = QubitOperator(correct_operators_b3, 1) 85 self.assertTrue( 86 qterm_b0 == bksf.edge_operator_b(edge_matrix_indices, 0)) 87 self.assertTrue( 88 qterm_b1 == bksf.edge_operator_b(edge_matrix_indices, 1)) 89 self.assertTrue( 90 qterm_b2 == bksf.edge_operator_b(edge_matrix_indices, 2)) 91 self.assertTrue( 92 qterm_b3 == bksf.edge_operator_b(edge_matrix_indices, 3)) 93 94 def test_bravyi_kitaev_fast_edgeoperator_Aij(self): 95 # checking the edge operators 96 edge_matrix = numpy.triu(numpy.ones((4, 4))) 97 edge_matrix_indices = numpy.array( 98 numpy.nonzero( 99 numpy.triu(edge_matrix) - numpy.diag(numpy.diag(edge_matrix)))) 100 correct_operators_a01 = ((0, 'X'),) 101 correct_operators_a02 = ((0, 'Z'), (1, 'X')) 102 correct_operators_a03 = ((0, 'Z'), (1, 'Z'), (2, 'X')) 103 correct_operators_a12 = ((0, 'Z'), (1, 'Z'), (3, 'X')) 104 correct_operators_a13 = ((0, 'Z'), (2, 'Z'), (3, 'Z'), (4, 'X')) 105 correct_operators_a23 = ((1, 'Z'), (2, 'Z'), (3, 'Z'), (4, 'Z'), (5, 106 'X')) 107 108 qterm_a01 = QubitOperator(correct_operators_a01, 1) 109 qterm_a02 = QubitOperator(correct_operators_a02, 1) 110 qterm_a03 = QubitOperator(correct_operators_a03, 1) 111 qterm_a12 = QubitOperator(correct_operators_a12, 1) 112 qterm_a13 = QubitOperator(correct_operators_a13, 1) 113 qterm_a23 = QubitOperator(correct_operators_a23, 1) 114 115 self.assertTrue( 116 qterm_a01 == bksf.edge_operator_aij(edge_matrix_indices, 0, 1)) 117 self.assertTrue( 118 qterm_a02 == bksf.edge_operator_aij(edge_matrix_indices, 0, 2)) 119 self.assertTrue( 120 qterm_a03 == bksf.edge_operator_aij(edge_matrix_indices, 0, 3)) 121 self.assertTrue( 122 qterm_a12 == bksf.edge_operator_aij(edge_matrix_indices, 1, 2)) 123 self.assertTrue( 124 qterm_a13 == bksf.edge_operator_aij(edge_matrix_indices, 1, 3)) 125 self.assertTrue( 126 qterm_a23 == bksf.edge_operator_aij(edge_matrix_indices, 2, 3)) 127 128 def test_bravyi_kitaev_fast_jw_number_operator(self): 129 # bksf algorithm allows for even number of particles. So, compare the 130 # spectrum of number operator from jordan-wigner and bksf algorithm 131 # to make sure half of the jordan-wigner number operator spectrum 132 # can be found in bksf number operator spectrum. 133 bravyi_kitaev_fast_n = bksf.number_operator(self.molecular_hamiltonian) 134 jw_n = QubitOperator() 135 n_qubits = count_qubits(self.molecular_hamiltonian) 136 for i in range(n_qubits): 137 jw_n += jordan_wigner_one_body(i, i) 138 jw_eig_spec = eigenspectrum(jw_n) 139 bravyi_kitaev_fast_eig_spec = eigenspectrum(bravyi_kitaev_fast_n) 140 evensector = 0 141 for i in range(numpy.size(jw_eig_spec)): 142 if bool( 143 numpy.size( 144 numpy.where( 145 jw_eig_spec[i] == bravyi_kitaev_fast_eig_spec))): 146 evensector += 1 147 self.assertEqual(evensector, 2**(n_qubits - 1)) 148 149 def test_bravyi_kitaev_fast_jw_hamiltonian(self): 150 # make sure half of the jordan-wigner Hamiltonian eigenspectrum can 151 # be found in bksf Hamiltonian eigenspectrum. 152 n_qubits = count_qubits(self.molecular_hamiltonian) 153 bravyi_kitaev_fast_H = bksf.bravyi_kitaev_fast( 154 self.molecular_hamiltonian) 155 jw_H = jordan_wigner(self.molecular_hamiltonian) 156 bravyi_kitaev_fast_H_eig = eigenspectrum(bravyi_kitaev_fast_H) 157 jw_H_eig = eigenspectrum(jw_H) 158 bravyi_kitaev_fast_H_eig = bravyi_kitaev_fast_H_eig.round(5) 159 jw_H_eig = jw_H_eig.round(5) 160 evensector = 0 161 for i in range(numpy.size(jw_H_eig)): 162 if bool( 163 numpy.size( 164 numpy.where(jw_H_eig[i] == bravyi_kitaev_fast_H_eig))): 165 evensector += 1 166 self.assertEqual(evensector, 2**(n_qubits - 1)) 167 168 def test_bravyi_kitaev_fast_generate_fermions(self): 169 n_qubits = count_qubits(self.molecular_hamiltonian) 170 # test for generating two fermions 171 edge_matrix = bksf.bravyi_kitaev_fast_edge_matrix( 172 self.molecular_hamiltonian) 173 edge_matrix_indices = numpy.array( 174 numpy.nonzero( 175 numpy.triu(edge_matrix) - numpy.diag(numpy.diag(edge_matrix)))) 176 fermion_generation_operator = bksf.generate_fermions( 177 edge_matrix_indices, 2, 3) 178 fermion_generation_sp_matrix = get_sparse_operator( 179 fermion_generation_operator) 180 fermion_generation_matrix = fermion_generation_sp_matrix.toarray() 181 bksf_vacuum_state_operator = bksf.vacuum_operator(edge_matrix_indices) 182 bksf_vacuum_state_sp_matrix = get_sparse_operator( 183 bksf_vacuum_state_operator) 184 bksf_vacuum_state_matrix = bksf_vacuum_state_sp_matrix.toarray() 185 vacuum_state = numpy.zeros((2**(n_qubits), 1)) 186 vacuum_state[0] = 1. 187 bksf_vacuum_state = numpy.dot(bksf_vacuum_state_matrix, vacuum_state) 188 two_fermion_state = numpy.dot(fermion_generation_matrix, 189 bksf_vacuum_state) 190 # using total number operator to check the number of fermions generated 191 tot_number_operator = bksf.number_operator(self.molecular_hamiltonian) 192 number_operator_sp_matrix = get_sparse_operator(tot_number_operator) 193 number_operator_matrix = number_operator_sp_matrix.toarray() 194 tot_fermions = numpy.dot( 195 two_fermion_state.conjugate().T, 196 numpy.dot(number_operator_matrix, two_fermion_state)) 197 # checking the occupation number of site 2 and 3 198 number_operator_2 = bksf.number_operator(self.molecular_hamiltonian, 2) 199 number_operator_3 = bksf.number_operator(self.molecular_hamiltonian, 3) 200 number_operator_23 = number_operator_2 + number_operator_3 201 number_operator_23_sp_matrix = get_sparse_operator(number_operator_23) 202 number_operator_23_matrix = number_operator_23_sp_matrix.toarray() 203 tot_23_fermions = numpy.dot( 204 two_fermion_state.conjugate().T, 205 numpy.dot(number_operator_23_matrix, two_fermion_state)) 206 self.assertTrue(2.0 - float(tot_fermions.real) < 1e-13) 207 self.assertTrue(2.0 - float(tot_23_fermions.real) < 1e-13) 208 209 def test_bravyi_kitaev_fast_excitation_terms(self): 210 # Testing on-site and excitation terms in Hamiltonian 211 constant = 0 212 one_body = numpy.array([[1, 2, 0, 3], [2, 1, 2, 0], [0, 2, 1, 2.5], 213 [3, 0, 2.5, 1]]) 214 # No Coloumb interaction 215 two_body = numpy.zeros((4, 4, 4, 4)) 216 molecular_hamiltonian = InteractionOperator(constant, one_body, 217 two_body) 218 n_qubits = count_qubits(molecular_hamiltonian) 219 # comparing the eigenspectrum of Hamiltonian 220 bravyi_kitaev_fast_H = bksf.bravyi_kitaev_fast(molecular_hamiltonian) 221 jw_H = jordan_wigner(molecular_hamiltonian) 222 bravyi_kitaev_fast_H_eig = eigenspectrum(bravyi_kitaev_fast_H) 223 jw_H_eig = eigenspectrum(jw_H) 224 bravyi_kitaev_fast_H_eig = bravyi_kitaev_fast_H_eig.round(5) 225 jw_H_eig = jw_H_eig.round(5) 226 evensector_H = 0 227 for i in range(numpy.size(jw_H_eig)): 228 if bool( 229 numpy.size( 230 numpy.where(jw_H_eig[i] == bravyi_kitaev_fast_H_eig))): 231 evensector_H += 1 232 233 # comparing eigenspectrum of number operator 234 bravyi_kitaev_fast_n = bksf.number_operator(molecular_hamiltonian) 235 jw_n = QubitOperator() 236 n_qubits = count_qubits(molecular_hamiltonian) 237 for i in range(n_qubits): 238 jw_n += jordan_wigner_one_body(i, i) 239 jw_eig_spec = eigenspectrum(jw_n) 240 bravyi_kitaev_fast_eig_spec = eigenspectrum(bravyi_kitaev_fast_n) 241 evensector_n = 0 242 for i in range(numpy.size(jw_eig_spec)): 243 if bool( 244 numpy.size( 245 numpy.where( 246 jw_eig_spec[i] == bravyi_kitaev_fast_eig_spec))): 247 evensector_n += 1 248 self.assertEqual(evensector_H, 2**(n_qubits - 1)) 249 self.assertEqual(evensector_n, 2**(n_qubits - 1)) 250 251 def test_bravyi_kitaev_fast_number_excitation_operator(self): 252 # using hydrogen Hamiltonian and introducing some number operator terms 253 constant = 0 254 one_body = numpy.zeros((4, 4)) 255 one_body[(0, 0)] = .4 256 one_body[(1, 1)] = .5 257 one_body[(2, 2)] = .6 258 one_body[(3, 3)] = .7 259 two_body = self.molecular_hamiltonian.two_body_tensor 260 # initiating number operator terms for all the possible cases 261 two_body[(1, 2, 3, 1)] = 0.1 262 two_body[(1, 3, 2, 1)] = 0.1 263 two_body[(1, 2, 1, 3)] = 0.15 264 two_body[(3, 1, 2, 1)] = 0.15 265 two_body[(0, 2, 2, 1)] = 0.09 266 two_body[(1, 2, 2, 0)] = 0.09 267 two_body[(1, 2, 3, 2)] = 0.11 268 two_body[(2, 3, 2, 1)] = 0.11 269 two_body[(2, 2, 2, 2)] = 0.1 270 molecular_hamiltonian = InteractionOperator(constant, one_body, 271 two_body) 272 # comparing the eigenspectrum of Hamiltonian 273 n_qubits = count_qubits(molecular_hamiltonian) 274 bravyi_kitaev_fast_H = bksf.bravyi_kitaev_fast(molecular_hamiltonian) 275 jw_H = jordan_wigner(molecular_hamiltonian) 276 bravyi_kitaev_fast_H_eig = eigenspectrum(bravyi_kitaev_fast_H) 277 jw_H_eig = eigenspectrum(jw_H) 278 bravyi_kitaev_fast_H_eig = bravyi_kitaev_fast_H_eig.round(5) 279 jw_H_eig = jw_H_eig.round(5) 280 evensector_H = 0 281 for i in range(numpy.size(jw_H_eig)): 282 if bool( 283 numpy.size( 284 numpy.where(jw_H_eig[i] == bravyi_kitaev_fast_H_eig))): 285 evensector_H += 1 286 287 # comparing eigenspectrum of number operator 288 bravyi_kitaev_fast_n = bksf.number_operator(molecular_hamiltonian) 289 jw_n = QubitOperator() 290 n_qubits = count_qubits(molecular_hamiltonian) 291 for i in range(n_qubits): 292 jw_n += jordan_wigner_one_body(i, i) 293 jw_eig_spec = eigenspectrum(jw_n) 294 bravyi_kitaev_fast_eig_spec = eigenspectrum(bravyi_kitaev_fast_n) 295 evensector_n = 0 296 for i in range(numpy.size(jw_eig_spec)): 297 if bool( 298 numpy.size( 299 numpy.where( 300 jw_eig_spec[i] == bravyi_kitaev_fast_eig_spec))): 301 evensector_n += 1 302 self.assertEqual(evensector_H, 2**(n_qubits - 1)) 303 self.assertEqual(evensector_n, 2**(n_qubits - 1)) 304