diff --git a/src/tequila/simulators/simulator_api.py b/src/tequila/simulators/simulator_api.py index 4a93e592..8db872a9 100755 --- a/src/tequila/simulators/simulator_api.py +++ b/src/tequila/simulators/simulator_api.py @@ -9,6 +9,7 @@ from tequila.utils.exceptions import TequilaException, TequilaWarning from tequila.simulators.simulator_base import BackendCircuit, BackendExpectationValue from tequila.circuit.noise import NoiseModel +from tequila.wavefunction.qubit_wavefunction import QubitWaveFunction SUPPORTED_BACKENDS = ["qulacs", "qulacs_gpu", "qibo", "qiskit", "qiskit_gpu", "cirq", "pyquil", "symbolic", "qlm"] SUPPORTED_NOISE_BACKENDS = ["qiskit", "qiskit_gpu", "cirq", "pyquil"] # qulacs removed in v.1.9 @@ -22,7 +23,6 @@ from tequila.objective import Objective, Variable from tequila.circuit.gates import QCircuit import numbers.Real as RealNumber - from tequila.wavefunction.qubit_wavefunction import QubitWaveFunction """ Check which simulators are installed @@ -369,8 +369,9 @@ def simulate(objective: typing.Union['Objective', 'QCircuit', 'QTensor'], backend: str = None, noise: NoiseModel = None, device: str = None, + initial_state: Union[int, QubitWaveFunction] = 0, *args, - **kwargs) -> Union[RealNumber, 'QubitWaveFunction']: + **kwargs) -> Union[RealNumber, QubitWaveFunction]: """Simulate a tequila objective or circuit Parameters @@ -388,6 +389,8 @@ def simulate(objective: typing.Union['Objective', 'QCircuit', 'QTensor'], specify a noise model to apply to simulation/sampling device: a device upon which (or in emulation of which) to sample + initial_state: int or QubitWaveFunction: + the initial state of the circuit *args : **kwargs : @@ -409,7 +412,7 @@ def simulate(objective: typing.Union['Objective', 'QCircuit', 'QTensor'], compiled_objective = compile(objective=objective, samples=samples, variables=variables, backend=backend, noise=noise, device=device, *args, **kwargs) - return compiled_objective(variables=variables, samples=samples, *args, **kwargs) + return compiled_objective(variables=variables, samples=samples, initial_state=initial_state, *args, **kwargs) def draw(objective, variables=None, backend: str = None, name=None, *args, **kwargs): diff --git a/src/tequila/simulators/simulator_base.py b/src/tequila/simulators/simulator_base.py index dc0053bb..db67c5b3 100755 --- a/src/tequila/simulators/simulator_base.py +++ b/src/tequila/simulators/simulator_base.py @@ -7,6 +7,7 @@ from tequila import BitString from tequila.objective.objective import Variable, format_variable_dictionary from tequila.circuit import compiler +from typing import Union import numbers, typing, numpy, copy, warnings @@ -107,6 +108,11 @@ class BackendCircuit(): "cc_max": True } + # Can be overwritten by backends that allow basis state initialization when sampling + supports_sampling_initialization: bool = False + # Can be overwritten by backends that allow initializing arbitrary states + supports_generic_initialization: bool = False + @property def n_qubits(self) -> numbers.Integral: return len(self.qubit_map) @@ -328,7 +334,7 @@ def update_variables(self, variables): """ self.circuit = self.create_circuit(abstract_circuit=self.abstract_circuit, variables=variables) - def simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveFunction: + def simulate(self, variables, initial_state: Union[int, QubitWaveFunction] = 0, *args, **kwargs) -> QubitWaveFunction: """ simulate the circuit via the backend. @@ -348,13 +354,12 @@ def simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveFunc the wavefunction of the system produced by the action of the circuit on the initial state. """ + if isinstance(initial_state, QubitWaveFunction) and not self.supports_generic_initialization: + raise TequilaException("Backend does not support arbitrary initial states") + self.update_variables(variables) if isinstance(initial_state, BitString): initial_state = initial_state.integer - if isinstance(initial_state, QubitWaveFunction): - if initial_state.length() != 1: - raise TequilaException("only product states as initial states accepted") - initial_state = list(initial_state.keys())[0].integer all_qubits = list(range(self.abstract_circuit.n_qubits)) active_qubits = self.qubit_map.keys() @@ -362,11 +367,21 @@ def simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveFunc # Keymap is only necessary if not all qubits are active keymap_required = sorted(active_qubits) != all_qubits + # Combining keymap and general initial states is awkward, because it's not clear what should happen with + # different states on non-active qubits. For now, this is simply not allowed. + # A better solution might be to check if all components of the initial state differ only on the active qubits. + if keymap_required and isinstance(initial_state, QubitWaveFunction): + raise TequilaException("Can only set non-basis initial state if all qubits are used") + if keymap_required: # maps from reduced register to full register keymap = KeyMapSubregisterToRegister(subregister=active_qubits, register=all_qubits) - mapped_initial_state = keymap.inverted(initial_state).integer if keymap_required else int(initial_state) + if not isinstance(initial_state, QubitWaveFunction): + mapped_initial_state = keymap.inverted(initial_state).integer if keymap_required else int(initial_state) + else: + mapped_initial_state = initial_state + result = self.do_simulate(variables=variables, initial_state=mapped_initial_state, *args, **kwargs) @@ -375,7 +390,7 @@ def simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveFunc return result - def sample(self, variables, samples, read_out_qubits=None, circuit=None, *args, **kwargs): + def sample(self, variables, samples, read_out_qubits=None, circuit=None, initial_state=0, *args, **kwargs): """ Sample the circuit. If circuit natively equips paulistrings, sample therefrom. Parameters @@ -395,6 +410,12 @@ def sample(self, variables, samples, read_out_qubits=None, circuit=None, *args, The result of sampling, a recreated QubitWaveFunction in the sampled basis. """ + if initial_state != 0 and not self.supports_sampling_initialization: + raise TequilaException("Backend does not support initial states for sampling") + + if isinstance(initial_state, QubitWaveFunction) and not self.supports_generic_initialization: + raise TequilaException("Backend does not support arbitrary initial states") + self.update_variables(variables) if read_out_qubits is None: read_out_qubits = self.abstract_qubits @@ -406,7 +427,9 @@ def sample(self, variables, samples, read_out_qubits=None, circuit=None, *args, circuit = self.add_measurement(circuit=self.circuit, target_qubits=read_out_qubits) else: circuit = self.add_measurement(circuit=circuit, target_qubits=read_out_qubits) - return self.do_sample(samples=samples, circuit=circuit, read_out_qubits=read_out_qubits, *args, **kwargs) + + return self.do_sample(samples=samples, circuit=circuit, read_out_qubits=read_out_qubits, + initial_state=initial_state, *args, **kwargs) def sample_all_z_hamiltonian(self, samples: int, hamiltonian, variables, *args, **kwargs): """ @@ -511,7 +534,7 @@ def sample_paulistring(self, samples: int, paulistring, variables, *args, E = E / samples * paulistring.coeff return E - def do_sample(self, samples, circuit, noise, abstract_qubits=None, *args, **kwargs) -> QubitWaveFunction: + def do_sample(self, samples, circuit, noise, abstract_qubits=None, initial_state=0, *args, **kwargs) -> QubitWaveFunction: """ helper function for sampling. MUST be overwritten by inheritors. @@ -777,7 +800,7 @@ def __call__(self, variables, samples: int = None, *args, **kwargs): raise TequilaException( "BackendExpectationValue received not all variables. Circuit depends on variables {}, you gave {}".format( self._variables, variables)) - + if samples is None: data = self.simulate(variables=variables, *args, **kwargs) else: @@ -856,7 +879,7 @@ def sample(self, variables, samples, *args, **kwargs) -> numpy.array: samples = max(1, int(self.abstract_expectationvalue.samples * total_samples)) suggested = samples # samples are not necessarily set (either the user has to set it or some functions like optimize_measurements) - + if suggested is not None and suggested != samples: warnings.warn("simulating with samples={}, but expectationvalue carries suggested samples={}\nTry calling with samples='auto-total#ofsamples'".format(samples, suggested), TequilaWarning) diff --git a/src/tequila/simulators/simulator_qiskit.py b/src/tequila/simulators/simulator_qiskit.py index 6db20c2f..ee732518 100755 --- a/src/tequila/simulators/simulator_qiskit.py +++ b/src/tequila/simulators/simulator_qiskit.py @@ -5,9 +5,11 @@ from tequila import BitString, BitNumbering, BitStringLSB from tequila.utils.keymap import KeyMapRegisterToSubregister from tequila.utils import to_float +from typing import Union import warnings import numpy as np import qiskit, qiskit_aer, qiskit.providers.fake_provider +from qiskit import QuantumCircuit HAS_NOISE = True try: @@ -160,6 +162,9 @@ class BackendCircuitQiskit(BackendCircuit): numbering = BitNumbering.LSB + supports_sampling_initialization = True + supports_generic_initialization = True + def __init__(self, abstract_circuit: QCircuit, variables, qubit_map=None, noise=None, device=None, *args, **kwargs): """ @@ -262,6 +267,20 @@ def make_classical_map(self, qubit_map: dict): return classical_map + def add_state_init(self, circuit: QuantumCircuit, initial_state: Union[int, QubitWaveFunction]) -> QuantumCircuit: + if initial_state == 0: + return circuit + + if isinstance(initial_state, QubitWaveFunction): + statevector = initial_state.to_array(self.numbering) + else: + statevector = np.zeros(2 ** self.n_qubits) + statevector[reverse_int_bits(initial_state, self.n_qubits)] = 1.0 + + init_circuit = qiskit.QuantumCircuit(self.q, self.c) + init_circuit.set_statevector(statevector) + return init_circuit.compose(circuit) + def do_simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveFunction: """ Helper function for performing simulation. @@ -298,12 +317,7 @@ def do_simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveF circuit = self.circuit.assign_parameters(self.resolver) - if initial_state != 0: - state = np.zeros(2 ** self.n_qubits) - state[reverse_int_bits(initial_state, self.n_qubits)] = 1.0 - init_circuit = qiskit.QuantumCircuit(self.q, self.c) - init_circuit.set_statevector(state) - circuit = init_circuit.compose(circuit) + circuit = self.add_state_init(circuit, initial_state) circuit.save_statevector() @@ -311,7 +325,7 @@ def do_simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveF return QubitWaveFunction.from_array(array=backend_result.get_statevector(circuit).data, numbering=self.numbering) - def do_sample(self, circuit: qiskit.QuantumCircuit, samples: int, read_out_qubits, *args, + def do_sample(self, circuit: qiskit.QuantumCircuit, samples: int, read_out_qubits, initial_state=0, *args, **kwargs) -> QubitWaveFunction: """ Helper function for performing sampling. @@ -321,6 +335,8 @@ def do_sample(self, circuit: qiskit.QuantumCircuit, samples: int, read_out_qubit the circuit from which to sample. samples: the number of samples to take. + initial_state: + initial state of the circuit args kwargs @@ -371,14 +387,14 @@ def do_sample(self, circuit: qiskit.QuantumCircuit, samples: int, read_out_qubit else: use_basis = qiskit_backend.configuration().basis_gates circuit = circuit.assign_parameters(self.resolver) # this is necessary -- see qiskit-aer issue 1346 + circuit = self.add_state_init(circuit, initial_state) circuit = qiskit.transpile(circuit, backend=qiskit_backend, basis_gates=use_basis, optimization_level=optimization_level ) job = qiskit_backend.run(circuit, shots=samples) - return self.convert_measurements(job, - target_qubits=read_out_qubits) + return self.convert_measurements(job, target_qubits=read_out_qubits) def convert_measurements(self, backend_result, target_qubits=None) -> QubitWaveFunction: """ diff --git a/src/tequila/simulators/simulator_qulacs.py b/src/tequila/simulators/simulator_qulacs.py index 43d91ba1..71da737c 100755 --- a/src/tequila/simulators/simulator_qulacs.py +++ b/src/tequila/simulators/simulator_qulacs.py @@ -1,3 +1,5 @@ +from typing import Union + import qulacs import numbers, numpy import warnings @@ -63,6 +65,11 @@ class BackendCircuitQulacs(BackendCircuit): numbering = BitNumbering.LSB + quantum_state_class = qulacs.QuantumState + + supports_sampling_initialization = True + supports_generic_initialization = True + def __init__(self, abstract_circuit, noise=None, *args, **kwargs): """ @@ -110,10 +117,17 @@ def __init__(self, abstract_circuit, noise=None, *args, **kwargs): self.circuit=self.add_noise_to_circuit(noise) - def initialize_state(self, n_qubits:int=None) -> qulacs.QuantumState: + def initialize_state(self, n_qubits: int = None, initial_state: Union[int, QubitWaveFunction] = None) -> qulacs.QuantumState: if n_qubits is None: n_qubits = self.n_qubits - return qulacs.QuantumState(n_qubits) + + state = self.quantum_state_class(n_qubits) + if isinstance(initial_state, int): + state.set_computational_basis(reverse_int_bits(initial_state, self.n_qubits)) + elif isinstance(initial_state, QubitWaveFunction): + state.load(initial_state.to_array(self.numbering)) + + return state def update_variables(self, variables): """ @@ -130,7 +144,7 @@ def update_variables(self, variables): for k, angle in enumerate(self.variables): self.circuit.set_parameter(k, angle(variables)) - def do_simulate(self, variables, initial_state, *args, **kwargs): + def do_simulate(self, variables, initial_state: Union[int, QubitWaveFunction], *args, **kwargs): """ Helper function to perform simulation. @@ -148,8 +162,7 @@ def do_simulate(self, variables, initial_state, *args, **kwargs): QubitWaveFunction: QubitWaveFunction representing result of the simulation. """ - state = self.initialize_state(self.n_qubits) - state.set_computational_basis(reverse_int_bits(initial_state, self.n_qubits)) + state = self.initialize_state(self.n_qubits, initial_state) self.circuit.update_quantum_state(state) wfn = QubitWaveFunction.from_array(array=state.get_vector(), numbering=self.numbering) @@ -205,8 +218,7 @@ def do_sample(self, samples, circuit, noise_model=None, initial_state=0, *args, QubitWaveFunction: the results of sampling, as a Qubit Wave Function. """ - state = self.initialize_state(self.n_qubits) - state.set_computational_basis(reverse_int_bits(initial_state, self.n_qubits)) + state = self.initialize_state(self.n_qubits, initial_state) circuit.update_quantum_state(state) sampled = state.sampling(samples) return self.convert_measurements(backend_result=sampled, target_qubits=self.measurements) diff --git a/src/tequila/simulators/simulator_qulacs_gpu.py b/src/tequila/simulators/simulator_qulacs_gpu.py index 76dcb66e..41999b0b 100755 --- a/src/tequila/simulators/simulator_qulacs_gpu.py +++ b/src/tequila/simulators/simulator_qulacs_gpu.py @@ -1,16 +1,18 @@ import qulacs +from qulacs_core import QuantumStateGpu + from tequila import TequilaException from tequila.simulators.simulator_qulacs import BackendCircuitQulacs, BackendExpectationValueQulacs + class TequilaQulacsGpuException(TequilaException): def __str__(self): return "Error in qulacs gpu backend:" + self.message + class BackendCircuitQulacsGpu(BackendCircuitQulacs): - def initialize_state(self, n_qubits:int=None) -> qulacs.QuantumState: - if n_qubits is None: - n_qubits = self.n_qubits - return qulacs.QuantumStateGpu(n_qubits) + quantum_state_class = QuantumStateGpu + class BackendExpectationValueQulacsGpu(BackendExpectationValueQulacs): BackendCircuitType = BackendCircuitQulacsGpu diff --git a/src/tequila/wavefunction/qubit_wavefunction.py b/src/tequila/wavefunction/qubit_wavefunction.py index 6aa7199f..e7f48ad9 100644 --- a/src/tequila/wavefunction/qubit_wavefunction.py +++ b/src/tequila/wavefunction/qubit_wavefunction.py @@ -178,19 +178,24 @@ def dense(self) -> bool: """ return self._dense - def to_array(self, copy: bool = True) -> npt.NDArray[complex]: + def to_array(self, out_numbering: BitNumbering = BitNumbering.MSB, copy: bool = True) -> npt.NDArray[complex]: """ Returns array of amplitudes. + :param out_numbering: Whether the first qubit is the most or least significant in the output array indices. + For dense wavefunctions, this operation is significantly cheaper when this is the same as the numbering + of the wavefunction. :param copy: Whether to copy the array or use it directly for dense Wavefunctions. If False, changes to the array or wavefunction will affect each other. :return: Array of amplitudes. """ - if self._dense: + if self._dense and self._numbering == out_numbering: return self._state.copy() if copy else self._state else: result = np.zeros(2 ** self._n_qubits, dtype=complex) for k, v in self.raw_items(): + if self._numbering != out_numbering: + k = reverse_int_bits(k, self._n_qubits) result[k] = v return result @@ -246,6 +251,9 @@ def values(self) -> Generator[complex]: return (v for k, v in self.items()) def __eq__(self, other) -> bool: + if not isinstance(other, QubitWaveFunction): + return False + raise TequilaException("Wavefunction equality is not well-defined. Consider using isclose.") def isclose(self: QubitWaveFunction, @@ -261,8 +269,12 @@ def isclose(self: QubitWaveFunction, :return: Whether the wavefunctions are close. """ inner = self.inner(other) - cosine_similarity = inner / (self.norm() * other.norm()) - return np.isclose(abs(cosine_similarity), 1.0, rtol, atol) + self_norm = self.norm() + other_norm = other.norm() + cosine_similarity = inner / (self_norm * other_norm) + + return (np.isclose(abs(cosine_similarity), 1.0, rtol, atol) + and np.isclose(self_norm, other_norm, rtol, atol)) def __add__(self, other: QubitWaveFunction) -> QubitWaveFunction: if self._dense and other._dense and self._numbering == other._numbering: diff --git a/tests/test_simulator_backends.py b/tests/test_simulator_backends.py index 8d5cda8e..367aad40 100644 --- a/tests/test_simulator_backends.py +++ b/tests/test_simulator_backends.py @@ -4,12 +4,13 @@ import importlib import numpy +import numpy as np import pytest import random import numbers import tequila as tq import tequila.simulators.simulator_api -from tequila import BitString +from tequila import BitString, QubitWaveFunction """ Warn if Simulators are not installed @@ -356,6 +357,25 @@ def test_initial_state_from_integer(simulator, initial_state): assert numpy.isclose(wfn[BitString.from_int(initial_state, 6)], 1.0) +@pytest.mark.parametrize("simulator", tequila.simulators.simulator_api.INSTALLED_SIMULATORS.keys()) +def test_initial_state_from_wavefunction(simulator): + if not tequila.simulators.simulator_api.INSTALLED_SIMULATORS[simulator][0].supports_generic_initialization: + return + + U = tq.gates.H(target=0) + + state = QubitWaveFunction.from_array(np.array([1.0, 1.0])).normalize() + result = tq.simulate(U, initial_state=state, backend=simulator) + assert result.isclose(QubitWaveFunction.from_basis_state(n_qubits=1, basis_state=0)) + result = tq.simulate(U, initial_state=state, backend=simulator, samples=100) + assert result.isclose(QubitWaveFunction.from_array(np.array([100.0, 0.0]))) + + state = QubitWaveFunction.from_array(np.array([1.0, -1.0])).normalize() + result = tq.simulate(U, initial_state=state, backend=simulator, samples=100) + assert result.isclose(QubitWaveFunction.from_array(np.array([0.0, 100.0]))) + + + @pytest.mark.parametrize("backend", tequila.simulators.simulator_api.INSTALLED_SIMULATORS.keys()) def test_hamiltonian_reductions(backend): for q in [0, 1, 2, 3, 4]: