Merge pull request #372 from tequilahub/pr-1.9.8

merging to v.1.9.8

Author: kottmanj

Diff for: .github/workflows/ci_backends.yml

@@ -32,7 +32,7 @@ jobs:
             # myqlm does not work in python3.7
             pip install "cirq" "qiskit>=0.30" "qulacs" "qibo==0.1.1"
         elif [ $ver -eq 8 ]; then
-            pip install "cirq" "qiskit>=0.30" "qulacs" "qibo==0.1.1" "myqlm" "cirq-google"
+            pip install "cirq" "qiskit" "qulacs" "qibo==0.1.1" "myqlm" "cirq-google"
         fi
         pip install -e .
     - name: Lint with flake8

Diff for: README.md

@@ -342,12 +342,7 @@ You can avoid it by downgrading cirq and openfermion
 ```bash
 pip install --upgrade "openfermion<=1.0.0"
 pip install --upgrade "cirq<=0.9.1"
-```  
-
-
-## Qiskit backend
-Qiskit version 0.25 is not yet supported.
-`pip install --upgrade qiskit<0.25` fixes potential issues. If not: Please let us know.
+```
 
 ## Circuit drawing
 Standard graphical circuit representation within a Jupyter environment is often done using `tq.draw`.

Diff for: src/tequila/apps/unary_state_prep.py

@@ -89,7 +89,6 @@ def __init__(self, target_space: typing.List[BitString], max_repeat: int = 100,
         simulator.convert_to_numpy = False
         variables = None # {k:k.name.evalf() for k in self._abstract_circuit.extract_variables()}
         wfn = simulator.simulate(initial_state=BitString.from_int(0, nbits=self.n_qubits), variables=variables)
-        wfn.n_qubits = self._n_qubits
         equations = []
         for k in target_space:
             equations.append(wfn[k] - abstract_coefficients[k])
@@ -174,7 +173,7 @@ def __call__(self, wfn: QubitWaveFunction) -> QCircuit:
         :return:
         """
         try:
-            assert (len(wfn) == len(self._target_space))
+            assert wfn.length() == len(self._target_space)
             for key in wfn.keys():
                 try:
                     assert (key in self._target_space)

Diff for: src/tequila/hamiltonian/paulis.py

@@ -256,7 +256,7 @@ def Projector(wfn, threshold=0.0, n_qubits=None) -> QubitHamiltonian:
 
     """
 
-    wfn = QubitWaveFunction(state=wfn, n_qubits=n_qubits)
+    wfn = QubitWaveFunction.convert_from(n_qubits, wfn)
 
     H = QubitHamiltonian.zero()
     for k1, v1 in wfn.items():
@@ -304,8 +304,9 @@ def KetBra(ket: QubitWaveFunction, bra: QubitWaveFunction, hermitian: bool = Fal
 
     """
     H = QubitHamiltonian.zero()
-    ket = QubitWaveFunction(state=ket, n_qubits=n_qubits)
-    bra = QubitWaveFunction(state=bra, n_qubits=n_qubits)
+
+    ket = QubitWaveFunction.convert_from(n_qubits, ket)
+    bra = QubitWaveFunction.convert_from(n_qubits, bra)
 
     for k1, v1 in bra.items():
         for k2, v2 in ket.items():

Diff for: src/tequila/quantumchemistry/encodings.py

@@ -128,9 +128,9 @@ def map_state(self, state: list, *args, **kwargs) -> list:
         fop = openfermion.FermionOperator(string, 1.0)
         op = self(fop)
         from tequila.wavefunction.qubit_wavefunction import QubitWaveFunction
-        wfn = QubitWaveFunction.from_int(0, n_qubits=n_qubits)
+        wfn = QubitWaveFunction.from_basis_state(n_qubits, 0)
         wfn = wfn.apply_qubitoperator(operator=op)
-        assert (len(wfn.keys()) == 1)
+        assert wfn.length() == 1
         key = list(wfn.keys())[0].array
         return key
 

Diff for: src/tequila/quantumchemistry/orbital_optimizer.py

@@ -208,11 +208,16 @@ def kernel(self, h1, h2, *args, **kwargs):
             H = molecule.make_hardcore_boson_hamiltonian()
         else:
             H = molecule.make_hamiltonian()
+
+        rdm1 = None
+        rdm2 = None
         if self.vqe_solver is not None:
             vqe_solver_arguments = {}
             if self.vqe_solver_arguments is not None:
                 vqe_solver_arguments = self.vqe_solver_arguments
             result = self.vqe_solver(H=H, circuit=self.circuit, molecule=molecule, **vqe_solver_arguments)
+            if hasattr(self.vqe_solver, "compute_rdms"):
+                rdm1, rdm2 = self.vqe_solver.compute_rdms(U=self.circuit, variables=result.variables, molecule=molecule, use_hcb=restrict_to_hcb)
         elif self.circuit is None:
             raise Exception("Orbital Optimizer: Either provide a callable vqe_solver or a circuit")
         else:
@@ -233,8 +238,9 @@ def kernel(self, h1, h2, *args, **kwargs):
             # static ansatz
             U = self.circuit
 
-        rdm1, rdm2 = molecule.compute_rdms(U=U, variables=result.variables, spin_free=True, get_rdm1=True, get_rdm2=True, use_hcb=restrict_to_hcb)
-        rdm2 = self.reorder(rdm2, 'dirac', 'mulliken')
+        if rdm1 is None or rdm2 is None:
+            rdm1, rdm2 = molecule.compute_rdms(U=U, variables=result.variables, spin_free=True, get_rdm1=True, get_rdm2=True, use_hcb=restrict_to_hcb)
+            rdm2 = self.reorder(rdm2, 'dirac', 'mulliken')
         if not self.silent:
             print("{:20} : {}".format("energy", result.energy))
             if len(self.history) > 0:
@@ -259,3 +265,4 @@ def __str__(self):
             else:
                 result += "{:30} : {}\n".format(k, v)
         return result
+

Diff for: src/tequila/simulators/simulator_api.py

@@ -9,9 +9,10 @@
 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_gpu", "qulacs",'qibo', "qiskit", "cirq", "pyquil", "symbolic", "qlm"]
-SUPPORTED_NOISE_BACKENDS = ["qiskit", 'cirq', 'pyquil'] # qulacs removed in v.1.9
+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
 BackendTypes = namedtuple('BackendTypes', 'CircType ExpValueType')
 INSTALLED_SIMULATORS = {}
 INSTALLED_SAMPLERS = {}
@@ -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
@@ -43,6 +43,19 @@
     HAS_QISKIT = False
     HAS_QISKIT_NOISE = False
 
+try:
+    pkg_resources.require("qiskit-aer-gpu")
+    from tequila.simulators.simulator_qiskit_gpu import BackendCircuitQiskitGpu, BackendExpectationValueQiskitGpu
+    HAS_QISKIT_GPU = True
+    INSTALLED_SIMULATORS["qiskit_gpu"] = BackendTypes(BackendCircuitQiskitGpu, BackendExpectationValueQiskitGpu)
+    INSTALLED_SAMPLERS["qiskit_gpu"] = BackendTypes(BackendCircuitQiskitGpu, BackendExpectationValueQiskitGpu)
+    from tequila.simulators.simulator_qiskit import HAS_NOISE as HAS_QISKIT_GPU_NOISE
+    if HAS_QISKIT_GPU_NOISE:
+        INSTALLED_NOISE_SAMPLERS["qiskit_gpu"] = BackendTypes(BackendCircuitQiskitGpu, BackendExpectationValueQiskitGpu)
+except (ImportError, DistributionNotFound):
+    HAS_QISKIT_GPU = False
+    HAS_QISKIT_GPU_NOISE = False
+
 HAS_QIBO = True
 try:
     from tequila.simulators.simulator_qibo import BackendCircuitQibo, BackendExpectationValueQibo
@@ -82,8 +95,8 @@
     HAS_QULACS = False
 
 try:
-    pkg_resources.require("qulacs-gpu")
-    import qulacs
+    # pkg_resources.require("qulacs-gpu")
+    from qulacs import QuantumStateGpu
     from tequila.simulators.simulator_qulacs_gpu import BackendCircuitQulacsGpu, BackendExpectationValueQulacsGpu
 
     HAS_QULACS_GPU = True
@@ -350,14 +363,15 @@ def compile_circuit(abstract_circuit: 'QCircuit',
     return CircType(abstract_circuit=abstract_circuit, variables=variables, noise=noise, device=device, *args, **kwargs)
 
 
-def simulate(objective: typing.Union['Objective', 'QCircuit','QTensor'],
+def simulate(objective: typing.Union['Objective', 'QCircuit', 'QTensor'],
              variables: Dict[Union[Variable, Hashable], RealNumber] = None,
              samples: int = None,
              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
@@ -375,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 :
@@ -394,9 +410,9 @@ def simulate(objective: typing.Union['Objective', 'QCircuit','QTensor'],
                 objective.extract_variables()))
 
     compiled_objective = compile(objective=objective, samples=samples, variables=variables, backend=backend,
-                                 noise=noise,device=device, *args, **kwargs)
+                                 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 for: 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,34 +354,43 @@ 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 len(initial_state.keys()) != 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()
 
         # 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)
 
         if keymap_required:
-            result.apply_keymap(keymap=keymap, initial_state=initial_state)
+            result = QubitWaveFunction.from_wavefunction(result, keymap, n_qubits=len(all_qubits), initial_state=initial_state)
 
         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 for: src/tequila/simulators/simulator_cirq.py

@@ -173,7 +173,7 @@ def do_simulate(self, variables, initial_state=0, *args, **kwargs) -> QubitWaveF
         simulator = cirq.Simulator()
         backend_result = simulator.simulate(program=self.circuit, param_resolver=self.resolver,
                                             initial_state=initial_state)
-        return QubitWaveFunction.from_array(arr=backend_result.final_state_vector, numbering=self.numbering)
+        return QubitWaveFunction.from_array(array=backend_result.final_state_vector, numbering=self.numbering)
 
     def convert_measurements(self, backend_result: cirq.Result) -> QubitWaveFunction:
         """
@@ -186,18 +186,18 @@ def convert_measurements(self, backend_result: cirq.Result) -> QubitWaveFunction
         Returns
         -------
         QubitWaveFunction:
-            the result of sampling, as a tequila QubitWavefunction.
+            the result of sampling, as a tequila QubitWaveFunction.
 
         """
         assert (len(backend_result.measurements) == 1)
         for key, value in backend_result.measurements.items():
-            counter = QubitWaveFunction()
+            counter = QubitWaveFunction(self.n_qubits, self.numbering)
             for sample in value:
                 binary = BitString.from_array(array=sample.astype(int))
-                if binary in counter._state:
-                    counter._state[binary] += 1
+                if binary in counter.keys():
+                    counter[binary] += 1
                 else:
-                    counter._state[binary] = 1
+                    counter[binary] = 1
             return counter
 
     def do_sample(self, samples, circuit, *args, **kwargs) -> QubitWaveFunction:

Diff for: src/tequila/simulators/simulator_pyquil.py

@@ -439,7 +439,7 @@ def do_simulate(self, variables, initial_state, *args, **kwargs):
             if val > 0:
                 iprep += pyquil.gates.X(i)
         backend_result = simulator.wavefunction(iprep + self.circuit, memory_map=self.resolver)
-        return QubitWaveFunction.from_array(arr=backend_result.amplitudes, numbering=self.numbering)
+        return QubitWaveFunction.from_array(array=backend_result.amplitudes, numbering=self.numbering)
 
     def do_sample(self, samples, circuit, *args, **kwargs) -> QubitWaveFunction:
         """
@@ -495,7 +495,7 @@ def string_to_array(s):
                     listing.append(int(letter))
             return listing
 
-        result = QubitWaveFunction()
+        result = QubitWaveFunction(self.n_qubits, self.numbering)
         bit_dict = {}
         for b in backend_result:
             try:
@@ -505,7 +505,7 @@ def string_to_array(s):
 
         for k, v in bit_dict.items():
             arr = string_to_array(k)
-            result._state[BitString.from_array(arr)] = v
+            result[BitString.from_array(arr)] = v
         return result
 
     def no_translation(self, abstract_circuit):