diff --git a/jupyter_notebooks/Examples/MCFE-Example.ipynb b/jupyter_notebooks/Examples/MCFE-Example.ipynb new file mode 100644 index 000000000..dc2fc9681 --- /dev/null +++ b/jupyter_notebooks/Examples/MCFE-Example.ipynb @@ -0,0 +1,417 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Running MCFE on target circuits" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "import pygsti\n", + "from collections import defaultdict\n", + "import numpy as np\n", + "import time" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "# Create pyGSTi circs\n", + "unmapped_circs = [pygsti.circuits.Circuit([[\"Gxpi2\", \"Q0\"], [\"Gypi2\", \"Q1\"]]),pygsti.circuits.Circuit([[\"Gypi2\", \"Q0\"], [\"Gxpi2\", \"Q1\"]])]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Map circuits to device connectivity and U3-CX gate set" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This step will be different depending on what architecture you are using. For this example, we are using an IBM device. You need to end up with pyGSTi circuits in a U3-CX gate set so that circuit mirroring can be performed." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/ndsieki/pygsti_development/src/pygsti/pygsti/circuits/circuit.py:4154: UserWarning: pyGSTi circuit mapping discards classical registers.\n", + " _warnings.warn('pyGSTi circuit mapping discards classical registers.')\n", + "/home/ndsieki/pygsti_development/src/pygsti/pygsti/circuits/circuit.py:4210: UserWarning: skipping measure\n", + " _warnings.warn('skipping measure')\n" + ] + } + ], + "source": [ + "mapped_circs = defaultdict(list)\n", + "\n", + "import qiskit\n", + "from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager as _pass_manager\n", + "from qiskit_ibm_runtime.fake_provider import FakeSherbrooke, FakeAthensV2\n", + "\n", + "fake_backend = FakeAthensV2()\n", + "\n", + "pm = _pass_manager(coupling_map=fake_backend.coupling_map, basis_gates=['u3', 'cx'], optimization_level=0)\n", + "\n", + "\n", + "for i, circ in enumerate(unmapped_circs):\n", + " # Convert from pyGSTi to Qiskit\n", + " # Comment these lines out and do qiskit_circ = circ if passing in Qiskit\n", + " pygsti_openqasm_circ = circ.convert_to_openqasm(block_between_layers=True, include_delay_on_idle=False)\n", + " # print(pygsti_openqasm_circ)\n", + " qiskit_circ = qiskit.QuantumCircuit.from_qasm_str(pygsti_openqasm_circ)\n", + "\n", + " # print(qiskit_circ.draw())\n", + "\n", + " mapped_qiskit_circ = pm.run(qiskit_circ)\n", + "\n", + " # print(mapped_qiskit_circ.draw())\n", + " pygsti_circ, _ = pygsti.circuits.Circuit.from_qiskit(mapped_qiskit_circ)\n", + " # print(pygsti_circ)\n", + "\n", + " mapped_circ = pygsti_circ\n", + "\n", + " metadata = {'width': len(mapped_circ.line_labels), 'depth': mapped_circ.depth, 'dropped_gates': 0, 'id': i}\n", + " mapped_circs[mapped_circ] += [metadata]\n", + "\n", + "\n", + "unmirrored_design = pygsti.protocols.FreeformDesign(mapped_circs)\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Mirror circuit generation" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We use Pauli random compiling (`pauli_rc`) here. Central Pauli (`central_pauli`) is also an option." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "using provided edesign for both reference and test compilations\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Sampling mirror circuits: 100%|##########| 2/2 [00:00<00:00, 3.54it/s]\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Sampling reference circuits for width 5 with 1 subsets\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Sampling reference circuits for subset ('Q0', 'Q1', 'Q2', 'Q3', 'Q4'): 100%|##########| 100/100 [00:00<00:00, 4501.20it/s]" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Mirroring time: 0.6060056686401367\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\n" + ] + } + ], + "source": [ + "# Highly recommended to seed all RNG\n", + "mcfe_rand_state = np.random.RandomState(20240718)\n", + "\n", + "start = time.time()\n", + "mirror_design = pygsti.protocols.mirror_edesign.make_mirror_edesign(\n", + " unmirrored_design,\n", + " account_for_routing=False,\n", + " num_mcs_per_circ=100,\n", + " num_ref_per_qubit_subset=100,\n", + " mirroring_strategy='pauli_rc',\n", + " rand_state=mcfe_rand_state)\n", + "print(f'Mirroring time:', time.time() - start)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We have created the MCFE experiment design." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Run the Edesign" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This example will run the edesign on a fake IBM backend, but this is not strictly required. This step needs to generate a `ProtocolData(edesign=mirror_edesign, dataset=circuit_counts_data)` where `mirror_edesign` is the variable defined earlier and `circuit_counts_data` is a `DataSet` that contains the outcomes for each circuit." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "0.1233065128326416\n" + ] + } + ], + "source": [ + "from pygsti.extras.devices import ExperimentalDevice\n", + "from pygsti.extras import devices, ibmq\n", + "\n", + "device = ExperimentalDevice.from_qiskit_backend(fake_backend)\n", + "pspec = device.create_processor_spec(['Gc{}'.format(i) for i in range(24)] + ['Gcnot'])\n", + "\n", + "start = time.time()\n", + "exp = ibmq.IBMQExperiment(mirror_design, pspec, circuits_per_batch=300, num_shots=1024, seed=20240718, checkpoint_override=True)\n", + "print(time.time() - start)" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "transpiling to basis gates ['measure', 'id', 'rz', 'reset', 'delay', 'sx', 'x', 'cx']\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "100%|██████████| 2/2 [00:13<00:00, 6.84s/it]" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Total transpilation time: 13.707915306091309\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\n" + ] + } + ], + "source": [ + "from qiskit_aer import AerSimulator\n", + "\n", + "sim_backend = AerSimulator.from_backend(fake_backend)\n", + "\n", + "qiskit_convert_kwargs={}\n", + "\n", + "start = time.time()\n", + "exp.transpile(sim_backend, direct_to_qiskit=True, qiskit_convert_kwargs=qiskit_convert_kwargs)\n", + "end = time.time()\n", + "print(f'Total transpilation time: {end - start}')" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Submitting batch 1\n", + " - Job ID is 7a324892-5a58-4fb1-a545-3183ba678622\n", + "Submitting batch 2\n", + " - Job ID is e6888670-fe9d-4cf1-bea8-6b047927b9ec\n" + ] + } + ], + "source": [ + "exp.submit(sim_backend)" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Querying IBMQ for results objects for batch 1...\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Querying IBMQ for results objects for batch 2...\n", + "11.363933324813843\n" + ] + } + ], + "source": [ + "start = time.time()\n", + "exp.batch_results = []\n", + "exp.retrieve_results()\n", + "end = time.time()\n", + "print(end - start)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [], + "source": [ + "data = exp.data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Compute process fidelity for each circuit" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Calculating effective polarizations: 67%|######7 | 336/500 [00:00<00:00, 3355.23it/s]" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Calculating effective polarizations: 100%|##########| 500/500 [00:00<00:00, 3217.16it/s]\n" + ] + } + ], + "source": [ + "from pygsti.protocols.vbdataframe import VBDataFrame\n", + "\n", + "df = VBDataFrame.from_mirror_experiment(unmirrored_design, data)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "If you used Central Pauli instead, you can swap `'RC Process Fidelity'` for `'CP Process Fidelity'` in the cell below." + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "process_fidelities = df.dataframe['RC Process Fidelity']" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "[0.9973274623920335, 0.9969363722899848]" + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "process_fidelities.tolist()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "qiskit2_venv", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.11" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/pygsti/baseobjs/label.py b/pygsti/baseobjs/label.py index db1d92745..5d3b93ee1 100644 --- a/pygsti/baseobjs/label.py +++ b/pygsti/baseobjs/label.py @@ -1893,7 +1893,7 @@ def to_native(self): """ return tuple(self) - def replacename(self, oldname, newname): + def replace_name(self, oldname, newname): """ Returns a label with `oldname` replaced by `newname`. diff --git a/pygsti/baseobjs/qubitgraph.py b/pygsti/baseobjs/qubitgraph.py index 46eba6231..397a0e8aa 100644 --- a/pygsti/baseobjs/qubitgraph.py +++ b/pygsti/baseobjs/qubitgraph.py @@ -853,7 +853,7 @@ def shortest_path_predecessor_matrix(self): self._refresh_dists_and_predecessors() return self._predecessors.copy() - def subgraph(self, nodes_to_keep, reset_nodes=False): + def subgraph(self, nodes_to_keep, reset_nodes=False, include_directions=True): """ Return a graph that includes only `nodes_to_keep` and the edges between them. @@ -867,6 +867,10 @@ def subgraph(self, nodes_to_keep, reset_nodes=False): be the integers starting at 0 (in 1-1 correspondence with the ordering in `nodes_to_keep`). + include_directions : bool, optional + Whether the subgraph should retain the direction information + of the parent graph. Defaults to True. + Returns ------- QubitGraph @@ -878,7 +882,7 @@ def subgraph(self, nodes_to_keep, reset_nodes=False): qubit_labels = nodes_to_keep edges = [] - for edge in self.edges(): + for edge in self.edges(include_directions=include_directions): if edge[0] in nodes_to_keep and edge[1] in nodes_to_keep: if reset_nodes: edges.append((labelmap[edge[0]], labelmap[edge[1]])) diff --git a/pygsti/circuits/__init__.py b/pygsti/circuits/__init__.py index 31b2217d5..c4adb7a7c 100644 --- a/pygsti/circuits/__init__.py +++ b/pygsti/circuits/__init__.py @@ -19,3 +19,5 @@ from .cloudcircuitconstruction import * from .gstcircuits import * # Unused: from rpecircuits import * + +from .subcircuit_selection import * diff --git a/pygsti/circuits/circuit.py b/pygsti/circuits/circuit.py index 2432cb692..52ae8e11f 100644 --- a/pygsti/circuits/circuit.py +++ b/pygsti/circuits/circuit.py @@ -9,7 +9,16 @@ # in compliance with the License. You may obtain a copy of the License at # http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. #*************************************************************************************************** + from __future__ import annotations +from typing import Dict, Tuple, Union, Optional, List, TYPE_CHECKING +if TYPE_CHECKING: + try: + import qiskit + import stim + except: + pass + import collections as _collections import itertools as _itertools import warnings as _warnings @@ -22,9 +31,6 @@ from pygsti.tools import slicetools as _slct from pygsti.tools.legacytools import deprecate as _deprecate_fn -from typing import Union, Optional, TYPE_CHECKING -if TYPE_CHECKING: - import stim #Externally, we'd like to do thinks like: # c = Circuit( LabelList ) @@ -4183,6 +4189,196 @@ def from_cirq(cls, circuit, qubit_conversion=None, cirq_gate_conversion= None, else: return cls(circuit_layers) + @classmethod + def from_qiskit(cls, + circuit: qiskit.QuantumCircuit, + qubit_conversion: Optional[Dict[qiskit.circuit.Qubit, str]] = None, + qiskit_gate_conversion: Optional[Dict[str, str]] = None, + use_standard_gate_conversion_as_backup: bool = True, + allow_different_gates_in_same_layer: bool = True, + verbose: bool = False + ) -> Tuple[Circuit, Dict[int, str]]: + """ + Converts and instantiates a pyGSTi Circuit object from a Qiskit QuantumCircuit object. + + Parameters + ---------- + circuit : Qiskit QuantumCircuit + The Qiskit QuantumCircuit object to parse into a pyGSTi circuit. + + qubit_conversion : dict, optional (default None) + A dictionary specifying a mapping between Qiskit qubit indices and + pyGSTi qubit labels (either integers or strings). + If None, then a default mapping is created. + + qiskit_gate_conversion : dict, optional (default None) + If specified a dictionary with keys given by Qiskit gate objects, + and values given by pyGSTi gate names which overrides the built-in + conversion dictionary used by default. + + use_standard_gate_conversion_as_backup : boolean (default True) + Determines how the circuit conversion will be handled when the custom + Qiskit gate conversion dict does not have an entry for the encountered + gate. If True, this method will fall back on the standard conversions + found in qiskit_gatenames_standard_conversions(). If False, the method + will fail. + + allow_different_gates_in_same_layer: boolean (default True) + determines if gates with different names can be in the same layer. + For instance, if a CZ gate and CX gate can both fit in the same layer, + they will either be placed in the same layer (if True) or split into + separate layers (if False). + + Returns + ------- + Tuple: + pygsti_circuit + A pyGSTi Circuit instance equivalent to the specified Qiskit one. + + Dict {qiskit_qubit_idx, pyGSTi_qubit} + Dictionary that contains the mapping from the Qiskit qubit index + to the corresponding pyGSTi qubit. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("Circuit class method `from_qiskit()` is designed for qiskit version 2.1.1 and may not \ + function properly for your qiskit version, which is " + qiskit.__version__) + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + # mapping between qiskit gates and pygsti gate names: + if qiskit_gate_conversion is not None: + if use_standard_gate_conversion_as_backup == True: + qiskit_to_gate_name_mapping = _itgs.qiskit_gatenames_standard_conversions() + for key, value in qiskit_gate_conversion.items(): + qiskit_to_gate_name_mapping[key] = value + else: + qiskit_to_gate_name_mapping = qiskit_gate_conversion + else: + qiskit_to_gate_name_mapping = _itgs.qiskit_gatenames_standard_conversions() + + # get all of the qubits in the Qiskit circuit + if len(circuit.qregs) > 1: + _warnings.warn('pyGSTi circuit mapping does not preserve Qiskit qreg structure.') + + if len(circuit.cregs): + _warnings.warn('pyGSTi circuit mapping discards classical registers.') + + qubits = circuit.qubits + + #ensure all of these have a conversion available. + if qubit_conversion is not None: + unmapped_qubits = set(qubits).difference(set(qubit_conversion.keys())) + assert len(unmapped_qubits) == 0, f'Missing Qiskit to pygsti conversions for some qubits: {unmapped_qubits}' + + qubit_idx_conversion = {i: qubit_conversion[circuit._qbit_argument_conversion(i)[0]] for i in range(circuit.num_qubits)} + + #if it is None, build a default mapping. + else: + # default mapping is the identity mapping: qubit i in the Qiskit circuit maps to qubit i in the pyGSTi circuit + qubit_conversion = {circuit._qbit_argument_conversion(i)[0]: f'Q{i}' for i in range(circuit.num_qubits)} # in Qiskit 1.1.1, the method is called qbit_argument_conversion. In Qiskit >=1.2 (as far as Noah can tell), the method is called _qbit_argument_conversion. + + qubit_idx_conversion = {i: f'Q{i}' for i in range(circuit.num_qubits)} + + line_labels = tuple(sorted([qubit_conversion[qubit] for qubit in qubits])) + + if verbose: + print(f'qubit_conversion: {qubit_conversion}') + print(f'line_labels: {line_labels}') + + layer_indices = {} + for line_label in line_labels: + layer_indices[line_label] = 0 + + instructions = circuit.data + + pygsti_circ_layers = [] + + if allow_different_gates_in_same_layer == False: + layer_names = [] + + for instruction in instructions: + if verbose: + print(f'instruction: {instruction}') + + if allow_different_gates_in_same_layer == False: + assert len(pygsti_circ_layers) == len(layer_names), "there must be one layer name for each layer!" + + name = instruction.operation.name + num_qubits = instruction.operation.num_qubits + params = instruction.operation.params + if len(params) == 0: + params = None + + if verbose: + print(f'qiskit instruction name: {name}') + print(f'number of qubits in instruction: {num_qubits}') + print(f'instruction params: {params}') + + pygsti_gate_qubits = [qubit_conversion[qubit] for qubit in instruction.qubits] + + if name == 'measure': + _warnings.warn('skipping measure') + continue + + if name == 'barrier': + next_index = max(layer_indices[qubit] for qubit in pygsti_gate_qubits) + for qubit in pygsti_gate_qubits: + layer_indices[qubit] = next_index + continue + + pygsti_gate_name = qiskit_to_gate_name_mapping[name][0] + + label = _Label(pygsti_gate_name, pygsti_gate_qubits, args=params) + + next_index = max(layer_indices[qubit] for qubit in pygsti_gate_qubits) + if verbose: + print(f'computed next index for {pygsti_gate_name} gate on lines {pygsti_gate_qubits}: {next_index}') + + if allow_different_gates_in_same_layer == True: # layers are not separated by the type of gate they contain + + if next_index < len(pygsti_circ_layers): # there is an existing layer in the circuit where the gate can be inserted + pygsti_circ_layers[next_index].append(label) + if verbose: + print(f"inserting {pygsti_gate_name} gate in layer {next_index} on lines {pygsti_gate_qubits}") + + for pygsti_qubit in pygsti_gate_qubits: # update where a gate on these qubits can be placed + layer_indices[pygsti_qubit] = next_index + 1 + else: + if verbose: + print(f"inserting {pygsti_gate_name} gate at end of circuit, which is layer {len(pygsti_circ_layers)} on lines {pygsti_gate_qubits}") + + pygsti_circ_layers.append([label]) # need to append gate at the end of the circuit, thus creating a new layer + for pygsti_qubit in pygsti_gate_qubits: + layer_indices[pygsti_qubit] = len(pygsti_circ_layers) + + else: + for i in range(next_index, len(pygsti_circ_layers)): # searching for a layer where the gate can be inserted. Layer name needs to match the name of the gate being inserted + if name == layer_names[i]: + if verbose: + print(f'inserting gate {pygsti_gate_name} on qubits {pygsti_gate_qubits} in layer {i}') + + pygsti_circ_layers[i].append(label) + for pygsti_qubit in pygsti_gate_qubits: + layer_indices[pygsti_qubit] = i + 1 + break + + else: # no place to put the gate in the existing layers. New layer is added with corresponding name + if verbose: + print(f'inserting gate {pygsti_gate_name} on qubits {pygsti_gate_qubits} in layer {len(pygsti_circ_layers)}') + + pygsti_circ_layers.append([label]) + layer_names.append(name) + for pygsti_qubit in pygsti_gate_qubits: + layer_indices[pygsti_qubit] = len(pygsti_circ_layers) + + circuit = cls(pygsti_circ_layers, line_labels=line_labels) + + return (circuit, qubit_idx_conversion) + + def convert_to_quil(self, num_qubits=None, gatename_conversion=None, @@ -4350,6 +4546,115 @@ def convert_to_quil(self, return quil + + def convert_to_qiskit(self, + num_qubits: int = None, + qubit_conversion: Optional[Union[str, Dict[str, Union[int, qiskit.circuit.Qubit]]]] = None, + gatename_conversion: Optional[Dict[str, qiskit.circuit.Instruction]] = None, + block_between_layers: bool = True, + qubits_to_measure: Optional[Union[str, List[str]]] = None, + ) -> qiskit.QuantumCircuit: + + """ + Convert a pyGSTi circuit to a Qiskit QuantumCircuit. + + Parameters + ---------- + num_qubits : int, optional + size of Qiskit QuantumCircuit to create from the pyGSTi circuit. If None, + the number of line labels will set the size of the QuantumCircuit. + It is often useful to provide this field if the pyGSTi circuit is to be + executed on a backend that has more qubits than the pyGSTi circuit. + + qubit_conversion : dict, optional + mapping from pyGSTi line labels to Qiskit qubits, either as indices + or Qiskit Qubit objects. If none, a literal mapping is used. If 'remove-Q' is set, + then the 'Q' at the beginning of the line label is removed: e.g., 'Q53' becomes 53 (integer). + + gatename_conversion : dict, optional + A dictionary mapping gate names contained in this circuit to the corresponding + gate names used in the rendered Qiskit QuantumCircuit. If None, a standard set of conversions + is used (see :func:`standard_gatenames_qiskit_conversions`). + + block_between_layers : bool, optional + Set whether or not pyGSTi layer structure is maintained in the Qiskit circuit. Default is True. + If set to False, Qiskit will move all gates to their earliest possible execution point. + + qubits_to_measure : str or list, optional + Set which pyGSTi qubits should be measured in the rendered Qiskit QuantumCircuit. + If None, no qubits are measured. + If 'all', all qubits in the length-`num_qubits` Qiskit QuantumRegister are measured. + If 'active', only the qubits for which a qubit conversion is specified are measured. + If a list of pyGSTi line labels is provided, then only the corresponding Qiskit qubits are measured. + + + Returns + --------- + qiskit.QuantumCircuit + a Qiskit QuantumCircuit corresponding to the pyGSTi circuits. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("Circuit class method `convert_to_qiskit()` is designed for qiskit version 2.1.1 and may not \ + function properly for your qiskit version, which is " + qiskit.__version__) + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + depth = self.depth + + if num_qubits is None: + num_qubits = self.width + + if qubit_conversion is None: + qubit_conversion = {label: label for label in self.line_labels} + elif qubit_conversion == 'remove-Q': + qubit_conversion = {label: (int(label[1:]) if (isinstance(label, str) and label[0]=='Q' and label[1:].isnumeric()) else label) for label in self.line_labels} + + + qiskit_qc = qiskit.QuantumCircuit(num_qubits) + + qiskit_gate_conversion = _itgs.standard_gatenames_qiskit_conversions() + + for i in range(depth): + layer = self.layer_label(i).components + for gate in layer: + qiskit_gate, qiskit_gate_name, is_standard_gate = qiskit_gate_conversion[gate.name] + qiskit_qubits = [qubit_conversion[qubit] for qubit in gate.qubits] + qiskit_qc.append(qiskit_gate(*(gate.args)), qiskit_qubits, copy=False) + + if block_between_layers: + qiskit_qc.barrier() + + if qubits_to_measure is not None: + if isinstance(qubits_to_measure, str): + if qubits_to_measure == 'all': + qiskit_qc.measure_all() + + elif qubits_to_measure == 'active': + qiskit_qubits_to_measure = [v for v in qubit_conversion.values()] + new_creg = qiskit_qc._create_creg(len(qiskit_qubits_to_measure), "cr") + qiskit_qc.add_register(new_creg) + qiskit_qc.barrier() + qiskit_qc.measure(qiskit_qubits_to_measure, new_creg) + + else: + raise ValueError(f"unknown string option for 'qubits_to_measure': {qubits_to_measure}") + + elif isinstance(qubits_to_measure, list): + qiskit_qubits_to_measure = [qubit_conversion[qubit] for qubit in qubits_to_measure] + new_creg = qiskit_qc._create_creg(len(qiskit_qubits_to_measure), "cr") + qiskit_qc.add_register(new_creg) + qiskit_qc.barrier() + qiskit_qc.measure(qiskit_qubits_to_measure, new_creg) + + else: + raise ValueError(f"could not parse argument for 'qubits_to_measure': {qubits_to_measure}") + + return qiskit_qc + + def convert_to_openqasm(self, num_qubits=None, standard_gates_version='u3', gatename_conversion=None, qubit_conversion=None, diff --git a/pygsti/circuits/subcircuit_selection.py b/pygsti/circuits/subcircuit_selection.py new file mode 100644 index 000000000..dfebc7888 --- /dev/null +++ b/pygsti/circuits/subcircuit_selection.py @@ -0,0 +1,811 @@ +""" +Utility functions for subcircuit selection +""" +#*************************************************************************************************** +# Copyright 2015, 2019, 2025 National Technology & Engineering Solutions of Sandia, LLC (NTESS). +# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights +# in this software. +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +# in compliance with the License. You may obtain a copy of the License at +# http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. +#*************************************************************************************************** + +from __future__ import annotations +from typing import Dict, List, Tuple, Callable, Union, Optional, Any, Protocol, Literal, Set, TYPE_CHECKING + +if TYPE_CHECKING: + try: + import qiskit + except: + pass + +import warnings as _warnings +from collections import Counter as _Counter, defaultdict as _defaultdict +import itertools as _itertools +import networkx as _nx + +import tqdm as _tqdm + +import numpy as _np +import json as _json + +from pygsti.circuits.circuit import Circuit as _Circuit +from pygsti.baseobjs.label import Label as _Label +from pygsti.protocols.protocol import FreeformDesign as _FreeformDesign +from pygsti.tools import internalgates as _itgs + +MAX_STARTING_LAYER_ATTEMPTS = 1000 + + +class HasGetMethod(Protocol): + def get(gate_name: str, qubits: List[int]) -> float: + pass + +def sample_subcircuits(full_circs: Union[_Circuit, List[_Circuit]], + width_depths: Dict[int, List[int]], + coupling_map: Union[str, List[int], List[str], qiskit.transpiler.CouplingMap], + instruction_durations: Union[qiskit.transpiler.InstructionDurations, HasGetMethod], + use_qiskit_for_instruction_durations: Optional[bool] = None, + num_samples_per_width_depth: int = 10, + strategy: Union[Literal['simple', 'greedy'], Callable[..., Any]] = 'simple', + strategy_args: Optional[dict] = None, + depth_metric: Literal['layer_count', 'falcon_depth'] = 'layer_count', + num_test_samples: Optional[int] = None, + rand_state: Optional[_np.random.RandomState] = None, + ) -> _FreeformDesign: + """ + Samples subcircuits from a full circuit based on specified width and depth pairs. + + Parameters + ------------- + full_circs : Union[pygsti.circuits.Circuit, List[pygsti.circuits.Circuit]] + The full circuit(s) from which to sample subcircuits. + + width_depths : Dict[int, List[int]] + A dictionary whose keys are subcircuit widths and whose values are lists of depths + to sample for that subcircuit depth. Defines a list of (width, depth) pairs for + subcircuits. + + coupling_map : Union[str, qiskit.transpiler.CouplingMap] + The coupling map defining the connectivity of qubits. Can be 'all-to-all', 'linear', + or a qiskit CouplingMap object. + + instruction_durations : qiskit.transpiler.InstructionDurations + A qiskit InstructionDurations object used to determine delay times for + idle subcircuit layers. + + use_qiskit_for_instruction_durations : bool, optional + Whether to use qiskit gatenames instead of pyGSTi gatenames when looking + up instruction durations in `instruction_durations`. If this argument + is not provided, then the lookup gatename is based on the type of + `instruction_durations`. If `instruction_durations` is an + InstructionDurations object, the qiskit gatenames are used, else + the pyGSTi gatenames are used. + + num_samples_per_width_depth : int, optional + The number of subcircuits to sample for each width-depth combination. Default is 10. + + strategy : Union[str, Callable[..., Any]], optional + The subcircuit sampling strategy to use ('simple', 'greedy', or a custom function). + Default is 'simple'. + + strategy_args : Dict, optional + Additional arguments if a custom sampling strategy function is used. Default is None. + + depth_metric : str, optional + The metric to use for measuring depth ('layer_count' or 'falcon_depth'). If `layer_count`, + the depth is calculated as the number of layers in the subcircuit. If `falcon_depth`, + the gate set must be U3-CX; U3 gates contribute 2 to the depth and CX gates contribute 1 + to the depth. The `falcon_depth` metric is based on the Falcon generation of IBMQ devices. + Default is 'layer_count'. + + num_test_samples : int, optional + The number of test samples to use if the 'greedy' subcircuit sampling strategy is employed. + If `strategy` is not `greedy`, this argument is ignored. + Default is None. + + rand_state : _np.random.RandomState, optional + A random state for reproducibility. Default is None. + + Returns + -------- + pygsti.protocols.FreeformDesign + A FreeformDesign object containing the sampled subcircuits and auxiliary + information, including a circuit ID and depth. + """ + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("Subcircuit selection with a qiskit CouplingMap and/or InstructionDurations object " + "is designed for qiskit version 2.1.1 " + "and may not function properly for your qiskit version, which is " + qiskit.__version__) + except: + _warnings.warn('Qiskit is required if using a CouplingMap or InstructionDurations object ' + 'for subcircuit selection, and does not appear to be installed.') + + if rand_state is None: + rand_state = _np.random.RandomState() + + # Build unmirrored circuit + subcircuits = _defaultdict(list) + counter = 0 + + if not isinstance(full_circs, list): # package pygsti circuit into list if list was not provided. + full_circs = [full_circs] + + for full_circ in _tqdm.tqdm(full_circs, ascii=True, desc='sampling subcircuits'): + for w, ds in width_depths.items(): + print(f'Width: {w}, Depth: ', end='') + for d in ds: + print(d, end=' ') + if strategy == 'simple': + + subcircs, drops = simple_weighted_subcirc_selection(full_circ, + w, d, + num_subcircs=num_samples_per_width_depth, + depth_metric=depth_metric, + coupling_map=coupling_map, + instruction_durations=instruction_durations, + use_qiskit_for_instruction_durations=use_qiskit_for_instruction_durations, + rand_state=rand_state, + verbosity=0) + + elif strategy == 'greedy': + num_test_samples = 50 * num_samples_per_width_depth + subcircs, drops = greedy_growth_subcirc_selection(full_circ, w, d, num_subcircs=num_samples_per_width_depth, + num_test_subcircs=num_test_samples, rand_state=rand_state, verbosity=0) + elif callable(strategy): + subcircs, drops = strategy(full_circ, w, d, num_subcircs=num_samples_per_width_depth, **strategy_args) + else: + raise ValueError("'strategy' is not a function or known string") + + for subcirc, drop in zip(subcircs, drops): + subcircuits[subcirc].append({'width': w, 'depth': d, 'dropped_gates': drop, 'id': counter}) + counter += 1 + + print() + + return _FreeformDesign(subcircuits) + +def simple_weighted_subcirc_selection(full_circ: _Circuit, + width: int, + depth: int, + num_subcircs: int, + coupling_map: Union[str, List[int], List[str], qiskit.transpiler.CouplingMap], + instruction_durations: Union[qiskit.transpiler.InstructionDurations, HasGetMethod], + use_qiskit_for_instruction_durations: Optional[bool] = None, + depth_metric: Literal['layer_count', 'falcon_depth'] = 'layer_count', + rand_state: Optional[_np.random.RandomState] = None, + return_depth_info: bool = False, + stochastic_2q_drops: bool = False, + verbosity: Literal[0, 1] = 1, + ) -> List[List[_Circuit], + List[int], + Optional[List[int]], + Optional[List[Tuple[int, int]]], + Optional[List[int]], + Optional[List[List[int]]]]: + """ + Samples subcircuits from a full circuit using a simple approach. The simple approach + is to identify a starting layer, along with a connected subset of active qubits, and + snip out a subcircuit with the desired width and depth. + + Parameters + ------------- + full_circs : Union[pygsti.circuits.Circuit, List[pygsti.circuits.Circuit]] + The full circuit(s) from which to sample subcircuits. + + width : int + width of subcircuit to snip out. + + depth: int + depth of subcircuit to snip out. + + num_subcircs : int + The number of subcircuits to snip out for the given width and depth. + + coupling_map : Union[str, qiskit.transpiler.CouplingMap] + The coupling map defining the connectivity of qubits. Can be 'all-to-all', 'linear', + or a qiskit CouplingMap object. + + instruction_durations : qiskit.transpiler.InstructionDurations + A qiskit InstructionDurations object used to determine delay times for + idle subcircuit layers. + + use_qiskit_for_instruction_durations : bool, optional + Whether to use qiskit gatenames instead of pyGSTi gatenames when looking + up instruction durations in `instruction_durations`. If this argument + is not provided, then the lookup gatename is based on the type of + `instruction_durations`. If `instruction_durations` is an + InstructionDurations object, the qiskit gatenames are used, else + the pyGSTi gatenames are used. + + rand_state : _np.random.RandomState, optional + A random state for reproducibility. Default is None. + + return_depth_info : bool, optional + Whether to include compiled depths and the start and end layer for + each subcircuit. Default is False. + + stochastic_2q_drops : bool, optional + Whether to apply stochastic dropping of 2-qubit gates. Default is False, + in which case all dangling 2-qubit gates are dropped. A gate is considered + dangling if it has support on at least one qubit that is not in the + selected subset of qubits snipped out for the subcircuit. + + verbosity : int, optional + Level of verbosity for logging. Default is 1. + + Returns + -------- + Tuple[List[pygsti.circuits.Circuit], List[int], + Optional[List[int], Optional[List[Tuple[int,int]]]], + Optional[List[int]], Optional[List[int]]] + + A tuple containing the selected subcircuits and the counts of dropped gates. + If `return_depth_info` is set to True, then the returns are extended to include + the compiled depth of each subcircuit and the start and end layers of each subcircuit. + If `stochastic_2q_drops` is set to True, then the returns are extended to include + the number of dangling gates in each subcircuit and the indices of added layers. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("Simple subcircuit selection with a qiskit CouplingMap and/or InstructionDurations object " + "is designed for qiskit version 2.1.1 " + "and may not function properly for your qiskit version, which is " + qiskit.__version__) + except: + _warnings.warn('Qiskit is required if using a CouplingMap or InstructionDurations object ' + 'for simple subcircuit selection, and does not appear to be installed.') + + full_width = len(full_circ.line_labels) + full_depth = len(full_circ) + assert width > 1 and depth > 1, "Target width and depth must be greater than 1" + assert width <= full_width, f"Target width has to be less than full circuit width ({full_width})" + assert depth <= full_depth, f"Target depth has to be less than full circuit depth ({full_depth})" + + if rand_state is None: + rand_state = _np.random.RandomState() + + subcircs = [] + failures = 0 + + possible_starts = list(range(full_depth - depth + 1)) + + while (len(subcircs) < num_subcircs) and (failures < MAX_STARTING_LAYER_ATTEMPTS): + # Sample depth with cumulative layer weights + start = rand_state.choice(possible_starts) + + # Identify layer subset to snip out. + + # Calculate physical depth + compiled_depth = 0 + end = start - 1 + while compiled_depth < depth: + end += 1 + if depth_metric == 'layer_count': + layer_depth = 1 + + # qubits = set(full_circ.line_labels) + elif depth_metric == 'falcon_depth': + layer_depth = 1 + for comp in full_circ._layer_components(end): + if comp.name == 'Gu3': + layer_depth = 2 + + else: + raise RuntimeError('Invalid layer type!') + + else: + raise ValueError(f"Unknown value for parameter `depth_metric`: {depth_metric}") + + compiled_depth += layer_depth + + if compiled_depth > depth: + # We overshot (added a Gu3 at the end with only one space) + failures += 1 + # Skip this and try again + continue + + # calculate layer durations to be used for any idle layers + layer_durations = [] + + if use_qiskit_for_instruction_durations is None: + try: + if isinstance(instruction_durations, qiskit.transpiler.InstructionDurations): + use_qiskit_for_instruction_durations = True + gate_name_conversions = _itgs.standard_gatenames_qiskit_conversions() + else: + use_qiskit_for_instruction_durations = False + except: + use_qiskit_for_instruction_durations = False + + + for layer_idx in range(start, end+1): + max_comp_duration = 0 + for comp in full_circ._layer_components(layer_idx): + if use_qiskit_for_instruction_durations: + backend_comp_name = gate_name_conversions[comp.name][1] + backend_comp_qubits = [int(q[1:]) for q in comp.qubits] + else: + backend_comp_name = comp.name + backend_comp_qubits = comp.qubits + + comp_duration = instruction_durations.get(backend_comp_name, backend_comp_qubits) + if comp_duration > max_comp_duration: + max_comp_duration = comp_duration + + layer_durations.append(max_comp_duration) + + # Determine a width subset to snip out. + if coupling_map == 'all-to-all': + qubit_subset = set(rand_state.choice(full_circ.line_labels, size=width, replace=False)) + + elif coupling_map == 'linear': + q_start = rand_state.choice(full_width - width) + q_end = q_start + width + qubit_subset = set(full_circ.line_labels[q_start:q_end]) + + elif isinstance(coupling_map, list): # this list needs to already include the 'Q' prefix on the qubits + G = _nx.Graph() + G.add_edges_from(coupling_map) + qubit_subset = random_connected_subgraph(G, width, rand_state) + + else: # likely the coupling_map is a CouplingMap instance + try: + if not isinstance(coupling_map, qiskit.transpiler.CouplingMap): + raise RuntimeError("If `coupling_map` is not 'all-to-all, 'linear', or a list," \ + "it must be a qiskit CouplingMap object.") + except: + raise ImportError('Qiskit is required when providing a CouplingMap for subcircuit selection,' \ + 'and does not appear to be installed.') + + edges = [] + for cs in coupling_map: + qubits = [f'Q{c}' for c in cs] + if all([q in full_circ.line_labels for q in qubits]): + edges.append(qubits) + + G = _nx.Graph() + G.add_edges_from(edges) + qubit_subset = random_connected_subgraph(G, width, rand_state) + + qubit_subset = set(str(qubit) for qubit in qubit_subset) + + # We have identified a width/depth to snip out, so do so now + subcirc_layers = [] + possible_dropped_gates = [] + for layer_idx in range(start, end+1): + new_layer = [] + for op in full_circ._layer_components(layer_idx): + if all([q in qubit_subset for q in op.qubits]): + # This is inside our width/depth, so add it + new_layer.append(op) + elif any([q in qubit_subset for q in op.qubits]): + # This is dangling, account for it and deal with it later + possible_dropped_gates.append((op, len(subcirc_layers))) + # else we are totally outside, so not dropped + + subcirc_layers.append(new_layer) + + # Handle dropped gates + total_dropped_gates = 0 + total_dangling_gates = 0 + added_layer_indices = [] + if stochastic_2q_drops: + # Split gates into two random sets of almost equal size + op_idxs_to_drop = rand_state.choice(list(range(len(possible_dropped_gates))), len(possible_dropped_gates)//2, replace=False) + total_dropped_gates = len(op_idxs_to_drop) + + # Add some dangling gates to make up for drops + ops_to_add = [op for i,op in enumerate(possible_dropped_gates) if i not in op_idxs_to_drop] + total_dangling_gates = 2*len(ops_to_add) + + offset = 0 + last_added_layer = -1 + new_layer = [] + for op, layer_idx in ops_to_add: + # If we're adding to a new layer, add the previously built additional layer + if layer_idx != last_added_layer: + if len(new_layer): + print(f'Adding new layer at layer {last_added_layer} with offset {offset}') + subcirc_layers.insert(last_added_layer + offset + 1, new_layer) + added_layer_indices.append(last_added_layer + offset + 1) + offset += 1 + new_layer = [] + + last_added_layer = layer_idx + + # Add the dangling gate back to the current subcirc layer + print(f'Adding op to layer {layer_idx} with offset {offset}') + subcirc_layers[layer_idx + offset].append(op) + # Add the dangling gate to an additional layer as well + new_layer.append(op) + + # Ensure we add the last additional layer + if len(new_layer): + print(f'Adding new layer at layer {last_added_layer} with offset {offset}') + subcirc_layers.insert(last_added_layer + offset + 1, new_layer) + added_layer_indices.append(last_added_layer + offset + 1) + else: + total_dropped_gates = len(possible_dropped_gates) + + + for i in range(len(subcirc_layers)): + scl = subcirc_layers[i] + if len(scl) == 0: # layer is empty. add delay of appropriate duration + delay_layer = [_Label('Gdelay', qubit, args=[layer_durations[i]]) for qubit in qubit_subset] + subcirc_layers[i] = delay_layer + + # Build subcircuit + subcirc = _Circuit(subcirc_layers, line_labels=sorted(qubit_subset)) + + subcircs.append((subcirc, total_dropped_gates, compiled_depth, (start, end), total_dangling_gates, added_layer_indices)) + + if verbosity > 0: + print(f'Found subcircuit with {total_dropped_gates} dropped gates, {compiled_depth} depth, and {total_dangling_gates} dangling gates') + + if (failures == MAX_STARTING_LAYER_ATTEMPTS): + raise RuntimeError(f"Failed to find a valid starting layer {MAX_STARTING_LAYER_ATTEMPTS} times!") + + # Unpacking to match greedy growth function + selected_subcircs, dropped_counts, compiled_depths, start_ends, dangling_counts, added_layers = list(zip(*subcircs)) + + if verbosity > 0: + print(f'Dropped gate counts for selected circuits: {dropped_counts}') + print(f'Compiled depths for selected circuits: {compiled_depths}') + print(f'Dangling gate counts for selected circuits: {dangling_counts}') + + returns = [list(selected_subcircs), dropped_counts] + + if return_depth_info: + returns.extend([compiled_depths, start_ends]) + if stochastic_2q_drops: + returns.extend([dangling_counts, added_layers]) + + return returns + +def greedy_growth_subcirc_selection(full_circ: _Circuit, + width: int, + depth: int, + num_subcircs: int = 1, + num_test_subcircs: int = 10, + rand_state: Optional[_np.random.RandomState] = None, + verbosity: Literal[0, 1] = 1, + return_depth_info: bool = False + ) -> Tuple[List[_Circuit], List[int], + Optional[List[int]], Optional[List[Tuple[int,int]]]]: + """ + Selects subcircuits using a greedy growth strategy, starting with one gate in one layer + and adding gates with overlapping support in subsequent layers. + + Parameters + ------------- + full_circ : pygsti.circuits.Circuit + The full circuit from which to select subcircuits. + + width : int + The target width of the subcircuits to be selected. + + depth : int + The target depth of the subcircuits to be selected. + + num_subcircs : int, optional + The number of subcircuits to select. Default is 1. + + num_test_subcircs : int, optional + The number of test subcircuits to generate. Default is 10. + + rand_state : Optional[_np.random.RandomState], optional + A random state for reproducibility. Default is None. + + verbosity : int, optional + Level of verbosity for logging. Default is 1. + + return_depth_info : bool, optional + Whether to include compiled depths and the start and end layer for + each subcircuit. Default is False. + + Returns + -------- + Tuple[List[_Circuit], List[int], + Optional[List[int]], Optional[List[Tuple[int,int]]]] + A tuple containing the selected subcircuits and the counts of dropped gates. + If `return_depth_info` is set to True, then the returns are extended to include + the compiled depth of each subcircuit and the start and end layers of each subcircuit. + """ + + # TODO: sample until success wuth some upper bound on the number of samples that can be taken? + + full_width = len(full_circ.line_labels) + full_depth = len(full_circ) + assert width > 1 and depth > 1, "Target width and depth must be greater than 1" + assert width <= full_width, f"Target width has to be less than full circuit width ({full_width})" + assert depth <= full_depth, f"Target depth has to be less than full circuit depth ({full_depth})" + assert num_subcircs <= num_test_subcircs, f"Must try at least {num_subcircs} test subcircuits" + + if rand_state is None: + rand_state = _np.random.RandomState() + + # Get test subcircuits + # TODO: Could be parallelized + test_subcircs = [_greedy_growth_subcirc(full_circ, width, depth, rand_state, verbosity-1) + for _ in range(num_test_subcircs)] + + # Drop any with lower physical depth (possible error mode) + pruned_test_subcircs = [] + seen_circuits = set() + for sc in test_subcircs: + # If not unique or not exact depth (possible error or trailing Gu3), skip + if sc[0] in seen_circuits or sc[2] != depth: + continue + + seen_circuits.add(sc[0]) + pruned_test_subcircs.append(sc) + + # Sort on lowest number of dropped gate counts, and then physical depth + sorted_test_subcircs = sorted(pruned_test_subcircs, key=lambda x: (x[1], x[2])) + + print(f'greedy num_subcircs arg: {num_subcircs}') + + print(f"number of subcircuits in 'sorted_test_subcircs': {len(sorted_test_subcircs)}") + + # Select + unzip subcircs & dropped counts + if len(sorted_test_subcircs) < num_subcircs: + raise ValueError(f"Not enough subcircuits, only found {len(sorted_test_subcircs)}. Try increasing 'num_test_subcircs'") + + selected_subcircs, dropped_counts, physical_depths, start_ends = list(zip(*sorted_test_subcircs[:num_subcircs])) + + + print(f'number of subcircuits (greedy): {len(selected_subcircs)}') + + if verbosity > 0: + print(f'Dropped gate counts for selected circuits: {dropped_counts}') + print(f'Physical depths for selected circuits: {physical_depths}') + + if return_depth_info: + return list(selected_subcircs), dropped_counts, physical_depths, start_ends + + return list(selected_subcircs), dropped_counts + +# Workhorse function for greedy growth subcircuit selection +def _greedy_growth_subcirc(circ: _Circuit, + width: int, + depth: int, + rand_state: _np.random.RandomState, + verbosity: Literal[0, 1] = 0 + ) -> Tuple[_Circuit, int, int, Tuple[int, int]]: + """ + Workhorse function for greedy generation of candidate subcircuits. + + Parameters + ------------- + circ : pygsti.circuits.Circuit + The full circuit from which to grow the subcircuit. + + width : int + The target width of the subcircuit to be grown. + + depth : int + The target depth of the subcircuit to be grown. + + rand_state : Optional[_np.random.RandomState] + A random state for reproducibility. Default is None. + + verbosity : int, optional + Level of verbosity for logging. Default is 0. + + Returns + -------- + Tuple[pygsti.circuits.Circuit, int, int, Tuple[int, int]] + A tuple containing the grown subcircuit, the count of dropped gates, + the physical depth, and the start-end indices of the subcircuit. + """ + + # Pick an initial layer + start = end = rand_state.randint(circ.depth) + + # Pick an operation in layer + ops = circ._layer_components(start) + op_idx = rand_state.randint(len(ops)) + + # Start tracking our subcircuit + qubit_subset = set(ops[op_idx].qubits) + physical_depth = 2 if ops[op_idx].name == 'Gu3' else 1 + + # Helper function for adding new layers to our subcircuit + def add_new_layer(layer_idx: int) -> Tuple[Set[_Label], int, int]: + """ + Helper function for adding new layers to the subcircuit. + + Parameters + ------------ + layer_idx : int + Index of layer to be added to the subcircuit. + + Returns + ----------- + Tuple[Set[pygsti.baseobjs.Label], int, int] + Tuple containing labels to be added, the number of dropped gates, + and the Falcon depth of the layer. The Falcon depth of a U3 layer + is 2; the Falcon depth of a CX layer is 1. + """ + + labels_to_add = set() + drops = 0 + new_depth = 0 + for op in circ._layer_components(layer_idx): + if any([q in qubit_subset for q in op.qubits]): + # This operation overlaps with our current subset + # We need to see if it extends our width or must be dropped + new_qubits = set(op.qubits) - qubit_subset + + if new_qubits and len(new_qubits) + len(labels_to_add) + len(qubit_subset) <= width: + # This operation expands our qubit subset, but we still have room the new qubit + labels_to_add.update(new_qubits) + new_depth = 2 if op.name == 'Gu3' else 1 + elif new_qubits: + # This operation expands our qubit subset but we don't have room, so must drop it + drops += 1 + else: + # This op is in our subset, get the physical depth + new_depth = 2 if op.name == 'Gu3' else 1 + + return labels_to_add, drops, new_depth + + # Now try to grow the circuit + # Technically this means physical_depth = depth or depth + 1 when we terminate + while physical_depth < depth: + # Look at previous layer + prev_failed = False + try: + prev_labels, prev_drops, prev_depth = add_new_layer(start - 1) + except IndexError: + prev_failed = True + + # Look at next layer + next_failed = False + try: + next_labels, next_drops, next_depth = add_new_layer(end + 1) + except IndexError: + next_failed = True + + if prev_failed and next_failed: + # Both failed, exit early + break + elif prev_failed: + labels = next_labels + new_depth = next_depth + end += 1 + elif next_failed: + labels = prev_labels + new_depth = prev_depth + start -= 1 + else: + # Prefer direction with fewer drops, otherwise choose randomly + if prev_drops < next_drops or (prev_drops == next_drops and rand_state.randint(2)): + labels = prev_labels + new_depth = prev_depth + start -= 1 + else: + labels = next_labels + new_depth = next_depth + end += 1 + + # Add new width/depth + qubit_subset.update(labels) + physical_depth += new_depth + + # We have the right depth, but not necessarily the right width + if len(qubit_subset) < width: + # We can look at same depth/complement width circuit to find qubits to add + remaining_width_circ = _Circuit(circ[start:end+1], editable=True) + remaining_width_circ.delete_lines(list(qubit_subset), delete_straddlers=True) + + # We can iterate through and check for 2-qubit gate counts as a heuristic for what qubits + # are likely to not add too many dangling gates + # Can probably do something even better here, but this may be good enough + remaining_qubits = set(circ.line_labels) - qubit_subset + remaining_2q_gate_counts = {q: 0 for q in remaining_qubits} + for layer_idx in range(remaining_width_circ.depth): + for op in remaining_width_circ._layer_components(layer_idx): + for q in op.qubits: + remaining_2q_gate_counts[q] += 1 + + # Sort and add qubits with lowest 2Q gate counts + sorted_gate_counts = sorted(list(remaining_2q_gate_counts.items()), key=lambda x: x[1]) + qubits_to_add = [gc[0] for gc in sorted_gate_counts[:(width - len(qubit_subset))]] + + qubit_subset.update(qubits_to_add) + + # We have identified a width/depth to snip out, so do so now + # Technically this can have more dropped ops since we didn't know our full width for all layers above + # But under the assumption that we fill out width pretty quickly, it should be close + subcirc_layers = [] + total_dropped_gates = 0 + for layer_idx in range(start, end+1): + new_layer = [] + + for op in circ._layer_components(layer_idx): + if all([q in qubit_subset for q in op.qubits]): + # This is inside our width/depth, so add it + new_layer.append(op) + elif any([q in qubit_subset for q in op.qubits]): + # This has some overlap but dangles, so drop it + total_dropped_gates += 1 + # else we are totally outside, so not dropped + + if len(new_layer): + subcirc_layers.append(new_layer) + + if verbosity > 0: + print(f'Found subcircuit with {total_dropped_gates} dropped gates and {physical_depth} depth') + + # Return circuit + dropped gates/physical depth for external selection + return _Circuit(subcirc_layers, line_labels=qubit_subset), total_dropped_gates, physical_depth, (start, end) + + +def random_connected_subgraph(G: _nx.Graph, + width: int, + rand_state: Optional[_np.random.RandomState] = None + ) -> Union[Set[int], Set[str]]: + """ + Generates a random set of nodes that form a + connected subgraph of a specified width from a given graph. + + Parameters + ------------- + G : networkx.Graph + The graph from which to generate the connected subgraph. + width : int + The target width of the subgraph to be generated. + rand_state : Optional[_np.random.RandomState] + A random state for reproducibility. Default is None. + + Returns + -------- + set[int or str] + A set of nodes representing the generated connected subgraph. + """ + + if rand_state is None: + rand_state = _np.random.RandomState() + + # pick a random starting node + starting_node = rand_state.choice(G.nodes) + + used_nodes = set([starting_node]) + + plausible_growth_nodes = set([starting_node]) + + for i in range(width - 1): + valid_node_found = False + while valid_node_found == False and len(plausible_growth_nodes): + growth_node = rand_state.choice(list(plausible_growth_nodes)) + neighbors = set(G.neighbors(growth_node)) + new_neighbors = neighbors.difference(used_nodes) + + if len(new_neighbors): + new_node = rand_state.choice(list(new_neighbors)) + used_nodes.add(new_node) + plausible_growth_nodes.add(new_node) + valid_node_found = True + + else: + plausible_growth_nodes.remove(growth_node) + + if valid_node_found == False: # we exited because all possible nodes where the graph could be extended have been exhausted + raise RuntimeError(f'Could not generate a subgraph with {width} nodes') + # failure should only occur if the initial node is on a connected component of the graph with less nodes than 'width'. + + + assert len(used_nodes) == width, f'set of selected nodes has length {used_nodes} but should have length {width}' + assert _nx.is_connected(G.subgraph(used_nodes)), f'subgraph on nodes {used_nodes} should be connected but is not' + + return used_nodes + + + + + diff --git a/pygsti/extras/ibmq/ibmqexperiment.py b/pygsti/extras/ibmq/ibmqexperiment.py index 38a370c74..ceb714b72 100644 --- a/pygsti/extras/ibmq/ibmqexperiment.py +++ b/pygsti/extras/ibmq/ibmqexperiment.py @@ -17,21 +17,32 @@ import time as _time import tqdm as _tqdm import warnings as _warnings +from collections import defaultdict as _defaultdict + # Try to load Qiskit try: import qiskit as _qiskit from qiskit.providers import JobStatus as _JobStatus + from qiskit.providers.fake_provider import GenericBackendV2 as _GenericBackendV2 except: _qiskit = None + _GenericBackendV2 = None # Try to load IBM Runtime try: from qiskit_ibm_runtime import SamplerV2 as _Sampler from qiskit_ibm_runtime import Session as _Session from qiskit_ibm_runtime import RuntimeJobV2 as _RuntimeJobV2 + from qiskit_ibm_runtime import IBMBackend as _IBMBackend from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager as _pass_manager -except: _Sampler = None +except: + _Sampler = None + +try: + from qiskit_aer import AerSimulator as _AerSimulator +except: + _AerSimulator = None # Tim updated to Qiskit 2.1.0 # OLD COMMENT FROM STEFAN? @@ -55,18 +66,23 @@ # Needs to be defined first for multiprocessing reasons -def _transpile_batch(circs, pass_manager, qasm_convert_kwargs): +def _transpile_batch(circs, pass_manager, direct_to_qiskit, qasm_convert_kwargs=None, qiskit_convert_kwargs=None): batch = [] for circ in circs: # TODO: Replace this with direct to qiskit - pygsti_openqasm_circ = circ.convert_to_openqasm(**qasm_convert_kwargs) - qiskit_qc = _qiskit.QuantumCircuit.from_qasm_str(pygsti_openqasm_circ) + if direct_to_qiskit: + qiskit_qc = circ.convert_to_qiskit(**qiskit_convert_kwargs) + else: + pygsti_openqasm_circ = circ.convert_to_openqasm(**qasm_convert_kwargs) + qiskit_qc = _qiskit.QuantumCircuit.from_qasm_str(pygsti_openqasm_circ) + batch.append(qiskit_qc) # Run pass manager on batch return pass_manager.run(batch) + class IBMQExperiment(_TreeNode, _HasPSpec): """ A object that converts pyGSTi ExperimentDesigns into jobs to be submitted to IBM Q, submits these @@ -276,15 +292,17 @@ def reverse_dict_key_bits(counts_dict): # NOTE: This is probably duplicative of some other code in pyGSTi def partial_trace(ordered_target_indices, input_dict): - output_dict = {} + output_dict = _defaultdict(int) for bitstring in input_dict.keys(): new_string = '' for index in ordered_target_indices: new_string += bitstring[index] - try: - output_dict[new_string] += input_dict[bitstring] - except: - output_dict[new_string] = input_dict[bitstring] + + output_dict[new_string] += input_dict[bitstring] + # try: + # output_dict[new_string] += input_dict[bitstring] + # except: + # output_dict[new_string] = input_dict[bitstring] return output_dict #if len(self.batch_results): @@ -297,16 +315,23 @@ def partial_trace(ordered_target_indices, input_dict): qjob = self.qjobs[exp_idx] print(f"Querying IBMQ for results objects for batch {exp_idx + 1}...") batch_result = qjob.result() + + for i, circ in enumerate(self.pygsti_circuit_batches[exp_idx]): + #ordered_target_indices = _np.argsort(_np.argsort([int(label[1:]) for label in circ.line_labels])) + # ordered_target_indices = [self.processor_spec.qubit_labels.index(q) for q in circ.line_labels] + ordered_target_indices = list(range(len(circ.line_labels))) + # pygsti labels are sorted lexicographically, qiskit ordering is numeric. #TODO: this assumes that we only measure the pyGSTi qubits on the device, e.g., if the pyGSTi circuit only uses Q100 and Q101, then we only measure device qubits 100 and 101. It may be useful to make this more flexible and associate some kind of 'measured_qubits' field with a metadata dict on each pyGSTi circuit (or transpiled qiskit circ). + + counts_data = partial_trace(ordered_target_indices, reverse_dict_key_bits(batch_result[i].data.cr.get_counts()))#TODO: make the name of the measurement register a kwarg (often it is named 'meas' and not 'cr') + # print(counts_data) + + ds.add_count_dict(circ, counts_data) + self.batch_results.append(batch_result) if not self.disable_checkpointing: self._write_checkpoint() - for i, circ in enumerate(self.pygsti_circuit_batches[exp_idx]): - ordered_target_indices = [self.processor_spec.qubit_labels.index(q) for q in circ.line_labels] - counts_data = partial_trace(ordered_target_indices, reverse_dict_key_bits(batch_result[i].data.cr.get_counts())) - ds.add_count_dict(circ, counts_data) - self.data = _ProtocolData(self.edesign, ds) if not self.disable_checkpointing: @@ -319,7 +344,7 @@ def submit(self, ibmq_backend, start=None, stop=None, ignore_job_limit=True, wai Parameters ---------- - ibmq_backend: qiskit.providers.ibmq.ibmqbackend.IBMQBackend + ibmq_backend: qiskit.providers.ibmq.ibmqbackend.IBMQBackend or qiskit_aer.AerSimulator The IBM Q backend to submit the jobs to. Should be the backend corresponding to the processor that this experiment has been designed for. @@ -360,9 +385,20 @@ def submit(self, ibmq_backend, start=None, stop=None, ignore_job_limit=True, wai "Transpilation missing! Either run .transpile() first, or if loading from file, " + \ "use the regen_qiskit_circs=True option in from_dir()." + #Get the backend version backend_version = ibmq_backend.version assert backend_version >= 2, "IBMQExperiment no longer supports v1 backends due to their deprecation by IBM" + + + if isinstance(ibmq_backend, _IBMBackend): + is_simulator = ibmq_backend.simulator + elif _AerSimulator is not None and isinstance(ibmq_backend, _AerSimulator): # without lazy evaluation, this expression would throw an error if _AerSimulator is None + is_simulator = True + elif _GenericBackendV2 is not None and isinstance(ibmq_backend, _GenericBackendV2): # without lazy evaluation, this expression would throw an error if _AerSimulator is None + is_simulator = True + else: + raise ValueError(f'ibmq_backend has an unsupported type, {type(ibmq_backend)}. Support is currently only available for arguments which are instances of `IBMBackend`, `AerSimulator`, or `GenericBackendV2`.') total_waits = 0 self.qjobs = [] if self.qjobs is None else self.qjobs @@ -393,6 +429,7 @@ def submit(self, ibmq_backend, start=None, stop=None, ignore_job_limit=True, wai print(f"Submitting batch {batch_idx + 1}") submit_status = False batch_waits = 0 + while not submit_status and batch_waits < max_attempts: try: #If submitting to a real device, get calibration data @@ -406,59 +443,66 @@ def submit(self, ibmq_backend, start=None, stop=None, ignore_job_limit=True, wai # Submit job self.qjobs.append(sampler.run(batch, shots = self.num_shots)) + + if not is_simulator: + status = self.qjobs[-1].status() + initializing = True + initializing_steps = 0 + while initializing and initializing_steps < max_attempts: + if status in [_JobStatus.INITIALIZING, "INITIALIZING", _JobStatus.VALIDATING, "VALIDATING"]: + status = self.qjobs[-1].status() + print(f' - {status} (query {initializing_steps})') + _time.sleep(wait_time) + initializing_steps += 1 + else: + initializing = False + + try: + job_id = self.qjobs[-1].job_id() + print(f' - Job ID is {job_id}') + self.job_ids.append(job_id) + except Exception: + print(' - Failed to get job_id.') + self.job_ids.append(None) - status = self.qjobs[-1].status() - initializing = True - initializing_steps = 0 - while initializing and initializing_steps < max_attempts: - if status in [_JobStatus.INITIALIZING, "INITIALIZING", _JobStatus.VALIDATING, "VALIDATING"]: - status = self.qjobs[-1].status() - print(f' - {status} (query {initializing_steps})') - _time.sleep(wait_time) - initializing_steps += 1 - else: - initializing = False - - try: + try: + print(f' - Queue position is {self.qjobs[-1].queue_position()}') + except Exception: + print(f' - Failed to get queue position for batch {batch_idx + 1}') + if isinstance(self.qjobs[-1], _RuntimeJobV2): + print(' (because queue position not available in RuntimeJobV2)') + try: + metrics = self.qjobs[-1].metrics() + start_time = _datetime.fromisoformat(metrics["estimated_start_time"]) + print(f' - Estimated start time: {start_time.astimezone()} (local timezone)') + except Exception: + print(f' - Unable to retrieve estimated start time') + + else: job_id = self.qjobs[-1].job_id() print(f' - Job ID is {job_id}') self.job_ids.append(job_id) - except Exception: - - print(' - Failed to get job_id.') - self.job_ids.append(None) - - try: - print(f' - Queue position is {self.qjobs[-1].queue_position()}') - except Exception: - print(f' - Failed to get queue position for batch {batch_idx + 1}') - if isinstance(self.qjobs[-1], _RuntimeJobV2): - print(' (because queue position not available in RuntimeJobV2)') - try: - metrics = self.qjobs[-1].metrics() - start_time = _datetime.fromisoformat(metrics["estimated_start_time"]) - print(f' - Estimated start time: {start_time.astimezone()} (local timezone)') - except Exception: - print(f' - Unable to retrieve estimated start time') - + submit_status = True + except Exception as ex: template = " An exception of type {0} occurred. Arguments:\n{1!r}" message = template.format(type(ex).__name__, ex.args) print(message) - try: - print(' Machine status is {}.'.format(ibmq_backend.status().status_msg)) - except Exception as ex1: - print(' Failed to get machine status!') - template = " An exception of type {0} occurred. Arguments:\n{1!r}" - message = template.format(type(ex).__name__, ex1.args) - print(message) - total_waits += 1 - batch_waits += 1 - print(f"This batch has failed {batch_waits} times and there have been {total_waits} total failures") - print('Waiting', end='') - _time.sleep(wait_time) + if not is_simulator: + try: + print(' Machine status is {}.'.format(ibmq_backend.status().status_msg)) + except Exception as ex1: + print(' Failed to get machine status!') + template = " An exception of type {0} occurred. Arguments:\n{1!r}" + message = template.format(type(ex).__name__, ex1.args) + print(message) + total_waits += 1 + batch_waits += 1 + print(f"This batch has failed {batch_waits} times and there have been {total_waits} total failures") + print('Waiting', end='') + _time.sleep(wait_time) finally: # Checkpoint calibration and job id data if not self.disable_checkpointing: @@ -483,7 +527,14 @@ def submit(self, ibmq_backend, start=None, stop=None, ignore_job_limit=True, wai if submit_status is False: raise RuntimeError("Ran out of max attempts and job was still not submitted successfully") - def transpile(self, ibmq_backend, qiskit_pass_kwargs=None, qasm_convert_kwargs=None, num_workers=1): + def transpile(self, + ibmq_backend, + direct_to_qiskit=False, + qiskit_pass_kwargs=None, + qasm_convert_kwargs=None, + qiskit_convert_kwargs=None, + num_workers=1, + ): """Transpile pyGSTi circuits into Qiskit circuits for submission to IBMQ. Parameters @@ -491,36 +542,61 @@ def transpile(self, ibmq_backend, qiskit_pass_kwargs=None, qasm_convert_kwargs=N ibmq_backend: IBM backend to use during Qiskit transpilation - opt_level: int, optional - Optimization level for Qiskit `generate_preset_pass_manager`. + direct_to_qiskit: bool, optional + Whether to use OpenQASM as an intermediary in the conversion from pyGSTi circut + to Qiskit circuit. If True, `qiskit_convert_kwargs` is passed to + `Circuit.convert_to_qiskit`. If False, the pyGSTi circuits are first converted to + OpenQASM and then converted into Qiskit circuits. qiskit_pass_kwargs: dict, optional + Only used if direct_to_qiskit is False. Additional kwargs to pass in to `generate_preset_pass_manager`. If not defined, the default is {'seed_transpiler': self.seed, 'optimization_level': 0, 'basis_gates': ibmq_backend.operation_names} Note that "optimization_level" is a required argument to the pass manager. qasm_convert_kwargs: dict, optional + Only used if direct_to_qiskit is False. Additional kwargs to pass in to `Circuit.convert_to_openqasm`. If not defined, the default is {'num_qubits': self.processor_spec.num_qubits, 'standard_gates_version': 'x-sx-rz'} - + + qiskit_convert_kwargs: dict, optional + Only used if direct_to_qiskit is True. + Additional kwargs to pass to `Circuit.convert_to_qiskit`. + If not defined, the default is {'num_qubits': self.processor_spec.num_qubits, + 'qubit_conversion': 'remove-Q', 'block_between_layers': True, + 'qubits_to_measure': 'active'}. + num_workers: int, optional Number of workers to use for parallel (by batch) transpilation """ circuits = self.edesign.all_circuits_needing_data.copy() num_batches = int(_np.ceil(len(circuits) / self.circuits_per_batch)) + if direct_to_qiskit: + if qiskit_convert_kwargs is None: + qiskit_convert_kwargs = {} + qiskit_convert_kwargs['num_qubits'] = qiskit_convert_kwargs.get('num_qubits', self.processor_spec.num_qubits) + qiskit_convert_kwargs['qubit_conversion'] = qiskit_convert_kwargs.get('qubit_conversion', 'remove-Q') + qiskit_convert_kwargs['block_between_layers'] = qiskit_convert_kwargs.get('block_between_layers', True) + qiskit_convert_kwargs['qubits_to_measure'] = qiskit_convert_kwargs.get('qubits_to_measure', 'active') + + else: + if qasm_convert_kwargs is None: + qasm_convert_kwargs = {} + qasm_convert_kwargs['num_qubits'] = qasm_convert_kwargs.get('num_qubits', self.processor_spec.num_qubits) + qasm_convert_kwargs['standard_gates_version'] = qasm_convert_kwargs.get('standard_gates_version', 'x-sx-rz') + if qiskit_pass_kwargs is None: qiskit_pass_kwargs = {} qiskit_pass_kwargs['seed_transpiler'] = qiskit_pass_kwargs.get('seed_transpiler', self.seed) - qiskit_pass_kwargs['optimization_level'] = qiskit_pass_kwargs.get('optimization_level', 0) + qiskit_pass_kwargs['layout_method'] = qiskit_pass_kwargs.get('layout_method', 'trivial') + qiskit_pass_kwargs['routing_method'] = qiskit_pass_kwargs.get('routing_method', 'none') + qiskit_pass_kwargs['optimization_level'] = qiskit_pass_kwargs.get('optimization_level', 1) qiskit_pass_kwargs['basis_gates'] = qiskit_pass_kwargs.get('basis_gates', ibmq_backend.operation_names) - - if qasm_convert_kwargs is None: - qasm_convert_kwargs = {} - qasm_convert_kwargs['num_qubits'] = qasm_convert_kwargs.get('num_qubits', self.processor_spec.num_qubits) - qasm_convert_kwargs['standard_gates_version'] = qasm_convert_kwargs.get('standard_gates_version', 'x-sx-rz') + + print(f"transpiling to basis gates {qiskit_pass_kwargs['basis_gates']}") if not len(self.pygsti_circuit_batches): rand_state = _np.random.RandomState(self.seed) # TODO: Should this be a different seed as transpiler? @@ -553,30 +629,38 @@ def transpile(self, ibmq_backend, qiskit_pass_kwargs=None, qasm_convert_kwargs=N print(f'Already completed transpilation of {len(self.qiskit_isa_circuit_batches)}/{num_batches} circuit batches') if len(self.qiskit_isa_circuit_batches) == num_batches: return - - pm = _pass_manager(backend=ibmq_backend, **qiskit_pass_kwargs) - + # Set up parallel tasks tasks = [self.pygsti_circuit_batches[i] for i in range(len(self.qiskit_isa_circuit_batches), num_batches)] - # We want to use transpile_batch and it's the same pm/convert kwargs, so create a new function with partially applied kwargs + # we are running transpile_batch where the only arg that changes + # across ranks is the circuits. Thus, create a new function with partially applied kwargs # This function now only takes circs as an argument (which are our task elements above) - task_fn = _partial(_transpile_batch, pass_manager=pm, qasm_convert_kwargs=qasm_convert_kwargs) + + pm = _pass_manager(**qiskit_pass_kwargs) + + task_fn = _partial(_transpile_batch, + pass_manager=pm, + qasm_convert_kwargs=qasm_convert_kwargs, + direct_to_qiskit=direct_to_qiskit, + qiskit_convert_kwargs=qiskit_convert_kwargs) + for task in _tqdm.tqdm(tasks): self.qiskit_isa_circuit_batches.append(task_fn(task)) - # Save single batch - chkpt_path = _pathlib.Path(self.checkpoint_path) / "ibmqexperiment" - with open(chkpt_path / 'meta.json', 'r') as f: - metadata = _json.load(f) + # Save single batch + if not self.disable_checkpointing: + chkpt_path = _pathlib.Path(self.checkpoint_path) / "ibmqexperiment" + with open(chkpt_path / 'meta.json', 'r') as f: + metadata = _json.load(f) - filenm = f"qiskit_isa_circuit_batches{len(self.qiskit_isa_circuit_batches)-1}" - _metadir._write_auxfile_member(chkpt_path, filenm, 'qpy', self.qiskit_isa_circuit_batches[-1]) - if 'qiskit_isa_circuit_batches' in metadata: - metadata['qiskit_isa_circuit_batches'].append(None) - - with open(chkpt_path / 'meta.json', 'w') as f: - _json.dump(metadata, f) + filenm = f"qiskit_isa_circuit_batches{len(self.qiskit_isa_circuit_batches)-1}" + _metadir._write_auxfile_member(chkpt_path, filenm, 'qpy', self.qiskit_isa_circuit_batches[-1]) + if 'qiskit_isa_circuit_batches' in metadata: + metadata['qiskit_isa_circuit_batches'].append(None) + + with open(chkpt_path / 'meta.json', 'w') as f: + _json.dump(metadata, f) def write(self, dirname=None): """ diff --git a/pygsti/processors/__init__.py b/pygsti/processors/__init__.py index 413dad846..547b12f8a 100644 --- a/pygsti/processors/__init__.py +++ b/pygsti/processors/__init__.py @@ -10,5 +10,6 @@ # http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. #*************************************************************************************************** -from .processorspec import ProcessorSpec, QubitProcessorSpec, QuditProcessorSpec -from .compilationrules import CompilationRules, CliffordCompilationRules, CompilationError +from .processorspec import * +from .compilationrules import * +from .random_compilation import * \ No newline at end of file diff --git a/pygsti/processors/random_compilation.py b/pygsti/processors/random_compilation.py new file mode 100644 index 000000000..223ade0a9 --- /dev/null +++ b/pygsti/processors/random_compilation.py @@ -0,0 +1,689 @@ +""" +Randomized circuit compilation via random compiling and central Pauli propagation. +""" +#*************************************************************************************************** +# Copyright 2015, 2019, 2025 National Technology & Engineering Solutions of Sandia, LLC (NTESS). +# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights +# in this software. +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +# in compliance with the License. You may obtain a copy of the License at +# http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. +#*************************************************************************************************** + +from __future__ import annotations +from typing import Literal, Optional, Union, List, Iterable, Dict, Tuple + +import numpy as _np + +from pygsti.circuits.circuit import Circuit as _Circuit +from pygsti.baseobjs.label import Label as _Label + +class RandomCompilation(object): + """ + A class for performing randomized circuit compilation. + + Attributes + ---------- + rc_strategy : str + The strategy used for randomized compilation. Currently, + 'pauli_rc' (pauli randomized compiling on a U3-CX-CZ gate set) and + 'central_pauli' (central Pauli propagation for a U3-CX-CZ gate set) + are supported. + + return_bs : bool + If True, the `compile` method will return the target bitstring for + the randomly compiled circuit. + + testing : bool + Flag for unit testing. If True, the user can provide test Pauli layers for + random compilation instead of the layers being randomly generated. + + rand_state : np.random.RandomState + A random state for reproducibility of random operations. + """ + + def __init__(self, + rc_strategy: Optional[Literal['rc', 'cp']] = None, + return_bs: Optional[bool] = False, + testing: Optional[bool] = False, + rand_state: Optional[_np.random.RandomState] = None): + """ + Initialize the RandomCompilation object. + + Parameters + ---------- + rc_strategy : str + The strategy used for randomized compilation. Currently, + 'pauli_rc' (pauli randomized compiling on a U3-CX-CZ gate set, see + https://arxiv.org/abs/2204.07568) and 'central_pauli' + (central Pauli propagation for a U3-CX-CZ gate set, see + https://www.nature.com/articles/s41567-021-01409-7) + are supported. Defaults to 'pauli_rc'. + + return_bs : bool + If True, the `compile` method will return the target bitstring for + the randomly compiled circuit. Default is False. + + testing : bool + Flag for unit testing. If True, the user can provide test Pauli layers for + random compilation instead of the layers being randomly generated. Default is False. + + rand_state : np.random.RandomState + A random state for reproducibility of random operations. Default is None. + """ + + self.rc_strategy = rc_strategy if rc_strategy is not None else "pauli_rc" + self.return_bs = return_bs + self.testing = testing + + + if isinstance(rand_state, _np.random.RandomState): + self.rand_state = rand_state + elif isinstance(rand_state, int): + self.rand_state = _np.random.RandomState(seed=rand_state) + else: + self.rand_state = _np.random.RandomState() + + + def compile(self, + circ: _Circuit, + test_layers: Optional[Union[List[_np.ndarray], _np.ndarray]] = None + ) -> _Circuit: + """ + Compiles the given circuit using the specified randomized compilation strategy. + + Parameters + ------------- + circ : pygsti.circuits.Circuit + The n-qubit circuit to be compiled. + + test_layers : list[np.ndarray[int]], optional + A list of test layers to be used in the random compilation + if `self.testing` is True. Layers are specified by a length-2*n array + whose entries are either 0 or 2. Indices 0:n correspond to Pauli Z errors: + a 2 indicates the presence an error. Likewise indices n:2*n indicate a + Pauli Z error. If using central Pauli, only one layer must be provided. + If using random compilation, a number of layers equal to the number of + layers of single-qubit gates must be provided. Default is None. + + Returns + -------- + list[pygsti.circuits.Circuit, str (optional), np.ndarray (optional)] + A list containing the randomized circuit, and optionally the bitstring and target Pauli vector. + """ + + if self.rc_strategy == 'pauli_rc': + return_bs = False + return_target_pauli = False + insert_test_layers = False + if self.return_bs: + return_bs = True + if self.testing: + insert_test_layers = True + return_target_pauli = True + return_bs = True + + return pauli_randomize_circuit(circ=circ, + rand_state=self.rand_state, + return_bs=return_bs, + return_target_pauli=return_target_pauli, + insert_test_layers=insert_test_layers, + test_layers=test_layers + ) + + elif self.rc_strategy == 'central_pauli': + return_bs = False + return_target_pauli = False + insert_test_layer = False + if self.return_bs: + return_bs = True + if self.testing: + insert_test_layer = True + return_target_pauli = True + return_bs = True + + return randomize_central_pauli(circ=circ, + rand_state=self.rand_state, + return_bs=return_bs, + return_target_pauli=return_target_pauli, + insert_test_layer=insert_test_layer, + test_layer=test_layers + ) + else: + raise ValueError(f"unknown compilation strategy '{self.rc_strategy}'!") + + +def pauli_randomize_circuit(circ: _Circuit, + rand_state: Optional[_np.random.RandomState] = None, + return_bs: bool = False, + return_target_pauli: bool = False, + insert_test_layers: bool = False, + test_layers: Optional[List[_np.ndarray]] = None + ) -> _Circuit: + """ + Performs random compilation on a given circuit by inserting Pauli gates between layers. + + Parameters + ------------- + circ : pygsti.circuits.Circuit + The circuit to be randomized. + + rand_state : np.random.RandomState, optional + A random state for reproducibility. Default is None, which initializes a new random state. + + return_bs : bool, optional + If True, returns the target bitstring for the randomly compiled circuit. Default is False. + + return_target_pauli : bool, optional + If True, returns the target Pauli vector for the circuit. Default is False. + + insert_test_layers : bool, optional + If True, uses `test_layers` as the Pauli layers to randomly compile instead of + randomly generating Pauli layers. + + test_layers : list[np.ndarray[int]], optional + A list of length-2*n arrays representing the test layers to be inserted + if `insert_test_layers `is True. The number of test layers must equal + the number of U3 layers in the circuit. Default is None. + + Returns + -------- + list[pygsti.circuits.Circuit, str (optional), np.ndarray (optional)] + A list containing the randomized circuit, and optionally the bitstring and target Pauli vector if specified. + """ + + if rand_state is None: + rand_state = _np.random.RandomState() + + d = circ.depth + n = circ.width + p = _np.zeros(2*n, int) + + + if insert_test_layers: + num_u3_layers = 0 + for i in range(d): + layer = circ.layer_label(i).components + if layer[0].name == 'Gu3': + num_u3_layers += 1 + assert len(test_layers) == num_u3_layers, f'expected {num_u3_layers} Pauli vectors but got {len(test_layers)} instead' + + if return_bs: + layer_0 = circ.layer_label(0).components + if layer_0[0].name == 'Gcnot': + return_0_bs = True + + qubit_map = {j:i for i,j in enumerate(circ.line_labels)} + + layers = [] + + for i in range(d): + layer = circ.layer_label(i).components + + if layer[0].name in ['Gi', 'Gdelay']: # making this explicit for the sake of clarity + layers.append(layer) + + elif len(layer) == 0 or layer[0].name == 'Gu3': + if insert_test_layers: + q = test_layers.pop(0) + if len(q) != 2*n: + raise ValueError(f"test layer {q} should have length {2*n} but has length {len(q)} instead") + + else: + q = 2 * rand_state.randint(0, 2, 2*n) + # update_u3_parameters now twirls/adds label for empty qubits, so don't prepad for speed + rc_layer = update_u3_parameters(layer, p, q, qubit_map) + layers.append(rc_layer) + p = q + + else: # we must have a layer of 2Q gates. this implementation allows for Gcnot and Gcphase gates in the same layer. + layers.append(layer) + for g in layer: + if g.name == 'Gcnot': + (control, target) = g.qubits + p[qubit_map[control]] = (p[qubit_map[control]] + p[qubit_map[target]]) % 4 + p[n + qubit_map[target]] = (p[n + qubit_map[control]] + p[n + qubit_map[target]]) % 4 + elif g.name == 'Gcphase': + (control, target) = g.qubits + p[qubit_map[control]] = (p[qubit_map[control]] + p[n + qubit_map[target]]) % 4 + p[qubit_map[target]] = (p[n + qubit_map[control]] + p[qubit_map[target]]) % 4 + else: + raise ValueError("Circuit can only contain Gcnot, Gcphase, Gu3, and Gi gates in separate layers!") + + bs = ''.join([str(b // 2) for b in p[n:]]) + + # Avoid checks for speed + rc_circ = _Circuit(layers, line_labels=circ.line_labels, check=False, expand_subcircuits=False) + + out = [rc_circ] + + if return_bs: + out.append(bs) + if return_target_pauli: + out.append(p) + + return out + + +def randomize_central_pauli(circ: _Circuit, + rand_state: Optional[_np.random.RandomState] = None, + return_bs: bool = False, + return_target_pauli: bool = False, + insert_test_layer: bool = False, + test_layer: Optional[_np.ndarray] = None + ) -> _Circuit: + """ + Perform circuit randomization by propagating a central Pauli layer through the circuit. + This function is designed to handle the "back half" of the mirror circuit: i.e., given a circuit C + whose fidelity is to be estimated using central Pauli mirroring, this function should be passed + C_inv + L_inv, where L_inv is Haar-random U3 layer. Refer to `make_mirror_edesign` in `protocols/mirror_design.py` + for more information. + + Parameters + ------------- + circ : Circuit + The circuit through which the central Pauli layer is to be propagated. + + rand_state : np.random.RandomState, optional + A random state for reproducibility. Default is None, which initializes a new random state. + + return_bs : bool, optional + If True, returns the target bitstring for the full mirror central Pauli circuit. Default is False. + + return_target_pauli : bool, optional + If True, returns the target Pauli vector that has been propagated + through the circuit. Default is False. + + insert_test_layer : bool, optional + If True, uses `test_layer` as the central Pauli layer + instead of randomly generating a central Pauli layer. + + test_layer : np.ndarray[int], optional + A length-2*n array representing the test layer to be inserted + if `insert_test_layer `is True. Default is None. + + Returns + -------- + list[pygsti.circuits.Circuit, str (optional), np.ndarray (optional)] + A list containing the randomized circuit, and optionally the bitstring and target Pauli vector if specified. + """ + + if rand_state is None: + rand_state = _np.random.RandomState() + + qubits = circ.line_labels + + qubit_map = {j:i for i,j in enumerate(circ.line_labels)} + + n = circ.width + d = circ.depth + + if insert_test_layer: + assert len(test_layer) == 2*n, f"Central Pauli vector must be length {2*n} but test_layer has length {len(test_layer)}" + central_pauli = test_layer + else: + central_pauli = 2 * rand_state.randint(0, 2, 2*n) + + central_pauli_layer = pauli_vector_to_u3_layer(central_pauli, qubits) + p = central_pauli.copy() + + layers = [central_pauli_layer] + + for i in range(d): + + layer = circ.layer_label(i).components + + if layer[0].name in ['Gi', 'Gdelay']: #making this explicit for the sake of clarity + layers.append(layer) + + elif len(layer) == 0 or layer[0].name == 'Gu3': + rc_layer = update_u3_parameters(layer, p, p, qubit_map) + layers.append(rc_layer) + + else: + layers.append(layer) + for g in layer: + if g.name == 'Gcnot': + (control, target) = g.qubits + p[qubit_map[control]] = (p[qubit_map[control]] + p[qubit_map[target]]) % 4 + p[n + qubit_map[target]] = (p[n + qubit_map[control]] + p[n + qubit_map[target]]) % 4 + elif g.name == 'Gcphase': + (control, target) = g.qubits + p[qubit_map[control]] = (p[qubit_map[control]] + p[n + qubit_map[target]]) % 4 + p[qubit_map[target]] = (p[n + qubit_map[control]] + p[qubit_map[target]]) % 4 + else: + raise ValueError("Circuit can only contain Gcnot, Gcphase, Gu3, and Gi gates in separate layers!") + + bs = ''.join([str(b // 2) for b in p[n:]]) + + # Avoid checks for speed + cp_circ = _Circuit(layers, line_labels=circ.line_labels, check=False, expand_subcircuits=False) + + out = [cp_circ] + + if return_bs: + out.append(bs) + if return_target_pauli: + out.append(p) + + return out + + +def update_u3_parameters(layer: Iterable[_Label], + p: _np.ndarray, + q: _np.ndarray, + qubit_map: Union[Dict[str, int], Dict[int, int]] + ) -> List[_Label]: + """ + Updates the parameters of U3 gates in a given layer based on the provided Pauli random compiling vectors. + + Parameters + ------------- + layer : iterable[pygsti.baseobjs.Label] + A list of gate labels representing the layer containing U3 gates. + + p : np.ndarray[int] + A vector describing the Pauli gates preceding the layer. + For an n-qubit layer, p is a length-2n array. + p[0:n] indicates Pauli-Z (2 is yes Z, 0 is no Z), p[n:2*n] is Pauli-X (2 yes, 0 no). + E.g., if n = 5, p[3] = 2, and p[8] = 2, then there is a Y gate on qubit 3. + + p : np.ndarray[int] + A vector describing the Pauli gates foloowing the layer. + For an n-qubit layer, q is a length-2n array. + q[0:n] indicates Pauli-Z (2 is yes Z, 0 is no Z), q[n:2*n] is Pauli-X (2 yes, 0 no). + E.g., if n = 5, p[1] = 0, and p[6] = 2, then there is an X gate on qubit 3. + + qubit_map : dict[str, int] + A mapping of qubit labels to their corresponding indices. + + Returns + -------- + list + A new layer containing updated U3 gates based on the applied Pauli gates. + """ + used_qubits = set() + + new_layer = [] + n = len(qubit_map) + + for g in layer: + assert(g.name == 'Gu3') + (theta, phi, lamb) = (float(g.args[0]), float(g.args[1]), float(g.args[2])) + qubit = g.qubits[0] + qubit_index = qubit_map[qubit] + if p[qubit_index] == 2: # Z gate preceeding the layer + lamb = lamb + _np.pi + if q[qubit_index] == 2: # Z gate following the layer + phi = phi + _np.pi + if p[n + qubit_index] == 2: # X gate preceeding the layer + theta = theta - _np.pi + phi = phi + lamb = -lamb - _np.pi + if q[n + qubit_index] == 2: # X gate following the layer + theta = theta - _np.pi + phi = -phi - _np.pi + lamb = lamb + + new_args = (mod_2pi(theta), mod_2pi(phi), mod_2pi(lamb)) + new_label = _Label('Gu3', qubit, args=new_args) + + new_layer.append(new_label) + used_qubits.add(qubit) + + for qubit, qubit_index in qubit_map.items(): + if qubit in used_qubits: + continue + + # Insert twirled idle on unpadded qubit + (theta, phi, lamb) = (0.0, 0.0, 0.0) + if p[qubit_index] == 2: # Z gate preceeding the layer + lamb = lamb + _np.pi + if q[qubit_index] == 2: # Z gate following the layer + phi = phi + _np.pi + if p[n + qubit_index] == 2: # X gate preceeding the layer + theta = theta - _np.pi + phi = phi + lamb = -lamb - _np.pi + if q[n + qubit_index] == 2: # X gate following the layer + theta = theta - _np.pi + phi = -phi - _np.pi + lamb = lamb + + if _np.allclose((theta, phi, lamb), (0.0, 0.0, 0.0)): + new_label = _Label('Gi', qubit, args=None) + else: + new_args = (mod_2pi(theta), mod_2pi(phi), mod_2pi(lamb)) + new_label = _Label('Gu3', qubit, args=new_args) + new_layer.append(new_label) + used_qubits.add(qubit) + + assert(set(used_qubits) == set(qubit_map.keys())) + + return new_layer + +def mod_2pi(theta: float) -> float: + """ + Modifies an angle to be within the range of -π to π. + + Parameters + ------------- + theta : float + The angle in radians to be modified. + + Returns + -------- + float + The modified angle within the range of -π to π. + """ + + while (theta > _np.pi or theta <= -1 * _np.pi): + if theta > _np.pi: + theta = theta - 2 * _np.pi + elif theta <= -1 * _np.pi: + theta = theta + 2 * _np.pi + return theta + + +def pauli_vector_to_u3_layer(p: _np.ndarray, + qubits: Union[List[str], List[int]] + ) -> _Label: + """ + Converts a Pauli vector into a corresponding layer of U3 gates. + + Parameters + ------------- + p : np.ndarray[int] + A vector representing the Pauli gates to be converted. + For an n-qubit layer, p is a length-2n array. + p[0:n] indicates Pauli-Z (2 is yes Z, 0 is no Z), p[n:2*n] is Pauli-X (2 yes, 0 no). + E.g., if n = 5, p[3] = 2, and p[8] = 2, then there is a Y gate on qubit 3. + + qubits : list[str] + A list of qubit labels corresponding to the Pauli vector. + + Returns + -------- + pygsti.baseobjs.Label + Label containing the layer of U3 gates derived from the Pauli vector. + """ + + n = len(qubits) + layer = [] + for i, q in enumerate(qubits): + + if p[i] == 0 and p[i+n] == 0: # I + theta = 0.0 + phi = 0.0 + lamb = 0.0 + if p[i] == 2 and p[i+n] == 0: # Z + theta = 0.0 + phi = _np.pi / 2 + lamb = _np.pi / 2 + if p[i] == 0 and p[i+n] == 2: # X + theta = _np.pi + phi = 0.0 + lamb = _np.pi + if p[i] == 2 and p[i+n] == 2: # Y + theta = _np.pi + phi = _np.pi / 2 + lamb = _np.pi / 2 + + layer.append(_Label('Gu3', q, args=(theta, phi, lamb))) + + return _Label(layer) + +def haar_random_u3_layer(qubits: Union[List[str], List[int]], + rand_state: Optional[_np.random.RandomState] = None + ) -> _Label: + """ + Generates a layer of Haar-random U3 gates. + + Parameters + ------------- + qubits : list[str] + A list of qubit labels for which to generate U3 gates. + rand_state : np.random.RandomState, optional + A random state for reproducibility. Default is None, which initializes a new random state. + + Returns + -------- + pygsti.baseobjs.Label + A label containing the layer of randomly generated U3 gates. + """ + + return _Label([haar_random_u3(q, rand_state) for q in qubits]) + +def haar_random_u3(q: Union[str, int], + rand_state: Optional[_np.random.RandomState] = None + ) -> _Label: + """ + Generates a Haar-random U3 gate. + + Parameters + ------------- + q : str + The qubit label for which to generate the U3 gate. + rand_state : np.random.RandomState, optional + A random state for reproducibility. Default is None, which initializes a new random state. + + Returns + -------- + pygsti.baseobjs.Label + A label representing the randomly generated U3 gate for the specified qubit. + """ + + if rand_state is None: + rand_state = _np.random.RandomState() + + a, b = 2 * _np.pi * rand_state.rand(2) + theta = mod_2pi(2 * _np.arcsin(_np.sqrt(rand_state.rand(1)))[0]) + phi = mod_2pi(a - b + _np.pi) + lamb = mod_2pi(-1 * (a + b + _np.pi)) + return _Label('Gu3', q, args=(theta, phi, lamb)) + + +def u3_cx_cz_inv(circ: _Circuit) -> _Circuit: + """ + Computes the inverse of a circuit composed of U3, CX and CZ gates. + + Parameters + ------------- + circ : pygsti.circuits.Circuit + The circuit for which to compute the inverse. + + Returns + -------- + pygsti.circuits.Circuit + A new circuit representing the inverse of the input circuit. + """ + + inverse_layers = [] + d = circ.depth + + for j in range(d): + layer = circ.layer(j) + inverse_layer = [gate_inverse(gate_label) for gate_label in layer] + inverse_layers.insert(0, inverse_layer) + + inverse_circ = _Circuit(inverse_layers, line_labels = circ.line_labels, check=False, expand_subcircuits=False) + + return inverse_circ + +def gate_inverse(label: _Label) -> _Label: + """ + Computes the inverse of a given gate label. + + Parameters + ------------- + label : pygsti.baseobjs.Label + The gate label for which to compute the inverse. + + Returns + -------- + pygsti.baseobjs.Label + A new label representing the inverse of the input gate. + """ + + if label.name == 'Gcnot': + return label + elif label.name == 'Gcphase': + return label + elif label.name == 'Gu3': + return _Label('Gu3', label.qubits, args=inverse_u3(label.args)) + elif label.name in ['Gi', 'Gdelay']: + return label + else: + raise RuntimeError(f'cannot compute gate inverse for {label}') + +def inverse_u3(args: Tuple[float, float, float]) -> Tuple[float, float, float]: + """ + Computes the inverse parameters for a U3 gate given its parameters. + + Parameters + ------------- + args : tuple[float] + A tuple containing the parameters (theta, phi, lambda) of the U3 gate. + + Returns + -------- + tuple[float] + A tuple containing the parameters of the inverse U3 gate. + """ + + theta_inv = mod_2pi(-float(args[0])) + phi_inv = mod_2pi(-float(args[2])) + lambda_inv = mod_2pi(-float(args[1])) + return (theta_inv, phi_inv, lambda_inv) + + +def pad_layer(layer: Iterable[_Label], + qubits: Union[List[str], List[int]] + ) -> List[_Label]: + """ + Pads a layer of gates with idle gates for any unused qubits. + + Parameters + ------------- + layer : list[pygsti.baseobjs.Label] + A list of gate labels representing the layer to be padded. + qubits : list[str] + A list of qubit labels to ensure all qubits are represented in the padded layer. + + Returns + -------- + list[pygsti.baseobjs.Label] + A new layer containing the original gates and idle gates for unused qubits. + """ + + padded_layer = list(layer) + used_qubits = [] + for g in layer: + for q in g.qubits: + used_qubits.append(q) + + for q in qubits: + if q not in used_qubits: + padded_layer.append(_Label('Gu3', (q,), args=(0.0, 0.0, 0.0))) + + return padded_layer \ No newline at end of file diff --git a/pygsti/protocols/__init__.py b/pygsti/protocols/__init__.py index 815daba32..4d0bd0d1b 100644 --- a/pygsti/protocols/__init__.py +++ b/pygsti/protocols/__init__.py @@ -19,3 +19,4 @@ from .stability import * from .vb import * from .vbdataframe import * +from .mirror_edesign import * \ No newline at end of file diff --git a/pygsti/protocols/mirror_edesign.py b/pygsti/protocols/mirror_edesign.py new file mode 100644 index 000000000..b50c77911 --- /dev/null +++ b/pygsti/protocols/mirror_edesign.py @@ -0,0 +1,984 @@ +""" +Functions for creating mirror experiment designs for mirror circuit fidelity estimation (MCFE) +""" +#*************************************************************************************************** +# Copyright 2015, 2019, 2025 National Technology & Engineering Solutions of Sandia, LLC (NTESS). +# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights +# in this software. +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +# in compliance with the License. You may obtain a copy of the License at +# http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. +#*************************************************************************************************** + +from __future__ import annotations +from typing import Callable, Union, Literal, Optional, Tuple, List, Dict, Any, TYPE_CHECKING +if TYPE_CHECKING: + try: + import qiskit + except: + pass + + +import numpy as _np +import tqdm as _tqdm +import warnings as _warnings +import time + +from collections import defaultdict + +from pygsti.circuits.circuit import Circuit as _Circuit +from pygsti.circuits import subcircuit_selection as _subcircsel +from pygsti.protocols.protocol import FreeformDesign as _FreeformDesign +from pygsti.protocols.protocol import CombinedExperimentDesign as _CombinedExperimentDesign +from pygsti.processors import random_compilation as _rc + + +def qiskit_circuits_to_mirror_edesign(qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + mirroring_kwargs_dict: Dict[str, Any] = {} + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a noise benchmark from transpiled Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a noise benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring. If an arg + is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + Returns + --------- + Tuple + pygsti.protocols.FreeformDesign + Experiment design containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + + pygsti.protocols.CombinedExperimentDesign + Experiment design containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. + """ + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'qiskit_circuits_to_mirror_edesign' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + + from qiskit import transpile + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + test_circ_dict = defaultdict(list) + ref_circ_dict = defaultdict(list) + + ref_circ_id_lookup_dict = {} + + if isinstance(qk_circs, list): + qk_circs = {i: qk_circ for i, qk_circ in enumerate(qk_circs)} + + for k, qk_test_circ in qk_circs.items(): + + qk_ref_circ = transpile(qk_test_circ, + basis_gates=['u3', 'cz'], + layout_method='trivial', + routing_method='none', + optimization_level=1, + approximation_degree=1.0, + seed_transpiler=0) + + ps_test_circ, _ = _Circuit.from_qiskit(qk_test_circ, allow_different_gates_in_same_layer=True) + ps_ref_circ, _ = _Circuit.from_qiskit(qk_ref_circ, allow_different_gates_in_same_layer=False) + + ps_test_circ = ps_test_circ.delete_idling_lines() + ps_ref_circ = ps_ref_circ.delete_idling_lines() + + test_circ_metadata = { + 'id': k, + 'width': ps_test_circ.width, + 'depth': ps_test_circ.depth, + # other information as it becomes necessary + } + + test_circ_dict[ps_test_circ] += [test_circ_metadata] + + ref_circ_metadata = { + 'id': k, + } + + ref_circ_dict[ps_ref_circ] += [ref_circ_metadata] + + ref_circ_id_lookup_dict[k] = ps_ref_circ + + test_edesign = _FreeformDesign(test_circ_dict) + ref_edesign = _FreeformDesign(ref_circ_dict) + + start = time.time() + + mirror_edesign = make_mirror_edesign(test_edesign=test_edesign, + ref_edesign=ref_edesign, + ref_id_lookup_dict=ref_circ_id_lookup_dict, + account_for_routing=False, + **mirroring_kwargs_dict + ) + + elapsed = time.time() - start + print(f'mirroring time: {elapsed}') + + return test_edesign, mirror_edesign + + +def qiskit_circuits_to_fullstack_mirror_edesign( + qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + qk_backend: Optional[qiskit.providers.BackendV2] = None, + coupling_map: Optional[qiskit.transpiler.CouplingMap] = None, + basis_gates: Optional[List[str]] = None, + transpiler_kwargs_dict: Dict[str, Any] = {}, + mirroring_kwargs_dict: Dict[str, Any] = {}, + num_transpilation_attempts: int = 100, + return_qiskit_time: bool = False + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a full-stack benchmark from high-level Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a full-stack benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + qk_backend : qiskit-ibm-runtime.IBMBackend, optional + IBM backend whose native gate set, connectivity, and error rates + are to be targeted when doing a full-stack transpilation. Fake backends + are also acceptable. If not provided, both `coupling_map` and `basis_gates` + must be provided, and certain transpiler optimizations that depend on + backend error rates may not be accessible. + + coupling_map : qiskit.transpiler.CouplingMap, optional + couplinp map for a target backend that must be provided if `qk_backend` is None. + This argument is ignored if `qk_backend` is not None. + + basis_gates : list[str], optional + list of native gates for a target backend that must be provided + if `qk_backend` is None. + This argument is ignored if `qk_backend` is not None. + + transpiler_kwargs_dict : dict, optional + keyword arguments that are passed to the Qiskit transpiler for full-stack + transpilation. Please see the Qiskit transpiler documentation + for a comprehensive list of options. If any or all of the following keys + are not provided, the listed defaults are used: + + 'optimization_level': default 3 + 'approximation_degree': default 1.0 + 'seed_transpiler': default None + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring. If an arg + is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + num_transpilation_attempts : int, optional + number of times to attempt full-stack circuit transpilation. Circuit mirroring requires + that the transpilation not use ancilla qubits, which is difficult to enforce with the + other transpiler options. Instead, we adopt a try-until-success strategy, which may fail + if an insufficient number of attempts are allowed. Increase this number if you are + having issues with the transpilation. If not provided, the default is 100 attempts. + + return_qiskit_time : bool, optional + Debug flag that sets whether or not to report the time spent in the Qiskit transpiler. + Qiskit transpilation is often the most costly part of benchmark creation and it can be + helpful to know how much time it is consuming. + + Returns + --------- + Tuple + pygsti.protocols.FreeformDesign + Experiment design containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + + pygsti.protocols.CombinedExperimentDesign + Experiment design containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. + + float, optional + amount of time spent in Qiskit transpiler. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'qiskit_circuits_to_fullstack_mirror_edesign' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + from qiskit import transpile + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + if qk_backend is None: + assert coupling_map is not None and basis_gates is not None, "'coupling_map' and 'basis_gates' must be provided if 'qk_backend is not provided." + + + def get_active_qubits_no_dag(circ: qiskit.QuantumCircuit, + ignore: set = set() + ) -> List[qiskit.circuit.Qubit]: + """ + Utility function for determining how many qubits are + active. A qubit is defined as active if there exists at + least one gate, excluding those in `ignore`, + that has support on that qubit. + + Parameters + ------------ + circ: pygsti.circuits.Circuit + Qiskit circuit whose active qubits are to be determined. + + ignore: set, optional + set of Qiskit instruction names that are ignored for the + purposes of determining activity. For instance, it may be + reasonable to exclude barriers and measures. + + Returns + ---------- + List + """ + + active_qubits = set() + for instruction, qargs, _ in circ.data: + if instruction.name in ignore: + continue + for qubit in qargs: + active_qubits.add(qubit) + + return active_qubits + + test_circ_dict = defaultdict(list) + ref_circ_dict = defaultdict(list) + + ref_circ_id_lookup_dict = {} + + if isinstance(qk_circs, list): + qk_circs = {i: qk_circ for i, qk_circ in enumerate(qk_circs)} + + # set default transpilation kwargs + transpiler_kwargs_dict['optimization_level'] = transpiler_kwargs_dict.get('optimization_level', 3) + transpiler_kwargs_dict['approximation_degree'] = transpiler_kwargs_dict.get('approximation_degree', 1.0) + transpiler_kwargs_dict['seed_transpiler'] = transpiler_kwargs_dict.get('seed_transpiler', None) + + if return_qiskit_time: + qiskit_time = 0.0 + + for k, qk_circ in qk_circs.items(): + num_virtual_qubits = qk_circ.num_qubits + for _ in range(num_transpilation_attempts): + + if return_qiskit_time: start = time.time() + + if qk_backend is not None: + if coupling_map is not None: + _warnings.warn("'coupling_map' is ignored when 'qk_backend' is provided.") + if basis_gates is not None: + _warnings.warn("'basis_gates' is ignored when 'qk_backend' is provided.") + + qk_test_circ = transpile(qk_circ, + backend=qk_backend, + coupling_map=coupling_map, + **transpiler_kwargs_dict) + else: + qk_test_circ = transpile(qk_circ, + coupling_map=coupling_map, + basis_gates=basis_gates, + **transpiler_kwargs_dict + ) + + if return_qiskit_time: qiskit_time += time.time() - start + + active_qubits = get_active_qubits_no_dag(qk_test_circ) + + uses_ancilla = False + + for qubit in active_qubits: + idx = qubit._index + init_register = qk_test_circ.layout.initial_layout[idx]._register + if init_register.name == 'ancilla': + uses_ancilla = True + print(qubit) + break + + if not uses_ancilla: + break + else: + raise RuntimeError('Could not generate transpilation that does not use ancilla qubits. Maybe try increasing `num_transpilation_attempts?`') + + qk_init_opt_layout = qk_test_circ.layout.initial_index_layout(filter_ancillas=False)[:num_virtual_qubits] + qk_final_opt_layout = qk_test_circ.layout.final_index_layout(filter_ancillas=False)[:num_virtual_qubits] + + assert set(qk_init_opt_layout) == set(qk_final_opt_layout) + + ps_test_circ, qubit_mapping = _Circuit.from_qiskit(qk_test_circ, allow_different_gates_in_same_layer=False) + ps_test_circ = ps_test_circ.delete_idling_lines() + + assert len(ps_test_circ.line_labels) <= num_virtual_qubits + + ps_init_opt_layout = [qubit_mapping[i] for i in qk_init_opt_layout] + ps_final_opt_layout = [qubit_mapping[i] for i in qk_final_opt_layout] + ps_routing_perm = {ps_init_opt_layout[i]: ps_final_opt_layout[i] for i in range(num_virtual_qubits)} + + assert set(ps_init_opt_layout) == set(ps_final_opt_layout) + + test_metadata = { + 'id': k, + 'width': ps_test_circ.width, + 'depth': ps_test_circ.depth, + 'initial_layout': ps_init_opt_layout, + 'final_layout': ps_final_opt_layout, + 'routing_permutation': ps_routing_perm, + } + + test_circ_dict[ps_test_circ] += [test_metadata] + + if qk_backend is not None: + reduced_coupling_map = qk_backend.coupling_map.reduce(qk_final_opt_layout) + else: # there must be a coupling map that was separately provided + reduced_coupling_map = coupling_map.reduce(qk_final_opt_layout) + + # generate exact reference circuit, ensuring that the line labels will align with the test circuit when inverted. + ref_seed = 0 + for _ in range(num_transpilation_attempts): + + qk_ref_inv_circ = transpile(qk_circ.inverse(), + coupling_map=reduced_coupling_map, + basis_gates=['u3', 'cz'], + initial_layout=list(range(num_virtual_qubits)), + optimization_level=1, + approximation_degree=1.0, + seed_transpiler=ref_seed + ) + + active_qubits = get_active_qubits_no_dag(qk_ref_inv_circ) + + uses_ancilla = False + + for qubit in active_qubits: + idx = qubit._index + init_register = qk_ref_inv_circ.layout.initial_layout[idx]._register + if init_register.name == 'ancilla': + uses_ancilla = True + print(qubit) + break + + if not uses_ancilla: + break + + ref_seed += 1 + + else: + raise RuntimeError('Could not generate transpilation that does not use ancilla qubits.') + + qk_init_ref_inv_reduced_layout = qk_ref_inv_circ.layout.initial_index_layout(filter_ancillas=False)[:num_virtual_qubits] + + qk_final_ref_inv_reduced_layout = qk_ref_inv_circ.layout.final_index_layout(filter_ancillas=False)[:num_virtual_qubits] + + qk_ref_circ = qk_ref_inv_circ.inverse() + + qubit_map = {qk_ref_inv_circ._qbit_argument_conversion(i)[0]: f'Q{qk_final_opt_layout[i]}' + for i in qk_init_ref_inv_reduced_layout} + + ps_ref_circ, qubit_idx_map = _Circuit.from_qiskit(qk_ref_circ, qubit_conversion=qubit_map, allow_different_gates_in_same_layer=False) + ps_ref_circ = ps_ref_circ.delete_idling_lines() + + assert len(ps_ref_circ.line_labels) <= num_virtual_qubits + + ps_init_ref_inv_layout = [qubit_idx_map[i] for i in qk_init_ref_inv_reduced_layout] + ps_final_ref_inv_layout = [qubit_idx_map[i] for i in qk_final_ref_inv_reduced_layout] + routing_perm_ref_inv_ps = {ps_init_ref_inv_layout[i]: ps_final_ref_inv_layout[i] for i in range(num_virtual_qubits)} + + ref_metadata = { + 'id': k, + 'width': ps_ref_circ.width, + 'depth': ps_ref_circ.depth, + 'initial_layout_inv': ps_init_ref_inv_layout, + 'final_layout_inv': ps_final_ref_inv_layout, + 'routing_permutation_inv': routing_perm_ref_inv_ps + } + + ref_circ_dict[ps_ref_circ] += [ref_metadata] + ref_circ_id_lookup_dict[k] = ps_ref_circ + + test_edesign = _FreeformDesign(test_circ_dict) + ref_edesign = _FreeformDesign(ref_circ_dict) + + mirror_edesign = make_mirror_edesign(test_edesign=test_edesign, + ref_edesign=ref_edesign, + ref_id_lookup_dict=ref_circ_id_lookup_dict, + account_for_routing=True, + **mirroring_kwargs_dict + ) + if return_qiskit_time: + return test_edesign, mirror_edesign, qiskit_time + else: + return test_edesign, mirror_edesign + + +def qiskit_circuits_to_subcircuit_mirror_edesign( + qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + aggregate_subcircs: bool, + width_depth_dict: Dict[int, List[int]], + coupling_map: qiskit.transpiler.CouplingMap, + instruction_durations: qiskit.transpiler.InstructionDurations, + subcirc_kwargs_dict: Dict[str, Any] = {}, + mirroring_kwargs_dict: Dict[str, Any] = {} + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a subcircuit benchmark from transpiled Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a subcircuit benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + aggregate_subcircs : bool + Whether or not the provided Qiskit circuits should be used to create + one combined subcircuit experiment design or kept separate. + Circuit aggregration can be useful if the provided circuits are all + instances of the same 'family' of the circuit, e.g., all + Bernstein-Vazirani circuits with different secret keys. + + width_depth_dict : dict[int, list[int]] + dictionary whose keys are subcircuit widths and whose values are + lists of depths to snip out for that width. + + coupling_map : str or qiskit.transpiler.CouplingMap + coupling map for the device the Qiskit circuits were transpiled to. + If 'all-to-all', an all-to-all coupling map is used. + If 'linear', a linear topology is used. + + instruction_durations : qiskit.transpiler.InstructionDurations + instruction durations for each gate in the target device. These + durations are needed to calculate the appropriate delay time + when only idling qubits are sampled out from a full circuit layer. + + subcirc_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in subcircuit selection. + If an arg is not provided, a default value is used. + The args are: + 'num_samples_per_width_depth': default 10. number of subcircuits to sample + from each full circuit. if `aggregate_subcircuits` is set to + True, `num_samplse_per_width_depth` subcircuits are drawn from + each full circuit and combined into one experiment design. + + 'rand_state': default None. np.random.RandomState to be used for subcircuit selection. + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring + If an arg is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + Returns + --------- + Tuple + dict[hashable, pygsti.protocols.FreeformDesign] or pygsti.protocols.FreeformDesign + Experiment design(s) containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + A dictionary is returned if `aggregate_subcircs` is False, otherwise a FreeformDesign + is returned. + + dict[hashable, pygsti.protocols.CombinedExperimentDesign] or pygsti.protocols.CombinedExperimentDesign + Experiment design(s) containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. A dictionary is returned + if `aggregate_subcircs` is False, otherwise a FreeformDesign is returned. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'qiskit_circuits_to_subcircuit_mirror_edesign' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + + from qiskit import transpile + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + if isinstance(qk_circs, list): + qk_circs = {i: qk_circ for i, qk_circ in enumerate(qk_circs)} + + + ps_circs = {} + for k, qk_circ in qk_circs.items(): + ps_circ, _ = _Circuit.from_qiskit(qk_circ, allow_different_gates_in_same_layer=True) + ps_circ = ps_circ.delete_idling_lines() + + ps_circs[k] = ps_circ + + # subcirc sampling in device native gate set + subcirc_kwargs_dict['num_samples_per_width_depth'] = subcirc_kwargs_dict.get('num_samples_per_width_depth', 10) + + test_edesigns = {} + if aggregate_subcircs: + test_edesign = _subcircsel.sample_subcircuits(full_circs=list(ps_circs.values()), + width_depths=width_depth_dict, + instruction_durations=instruction_durations, + coupling_map=coupling_map, + **subcirc_kwargs_dict) + test_edesigns[0] = test_edesign + + else: + for k, ps_circ in ps_circs.items(): + test_edesign = _subcircsel.sample_subcircuits(full_circs=ps_circ, + width_depths=width_depth_dict, + instruction_durations=instruction_durations, + coupling_map=coupling_map, + **subcirc_kwargs_dict) + + test_edesigns[k] = test_edesign + + mirror_edesigns = {} + + for k, test_edesign in test_edesigns.items(): + ref_circ_dict = defaultdict(list) + ref_circ_id_lookup_dict = {} + + for ps_test_circ, auxlist in test_edesign.aux_info.items(): + qubit_mapping_dict = {sslbl: i for i, sslbl in enumerate(ps_test_circ.line_labels)} # avoid mapping small circuits to large backends since connectivity is already ensured + + # convert back to qiskit to perform reference transpilations in u3-cz gate set (no layer blocking required, just need logical equivalence) + qk_test_circ = ps_test_circ.convert_to_qiskit(qubit_conversion=qubit_mapping_dict, + ) + qk_ref_circ = transpile(qk_test_circ, + basis_gates=['u3', 'cz'], + layout_method='trivial', + routing_method='none', + optimization_level=1, + approximation_degree=1.0, + seed_transpiler=0 + ) + + # convert those reference circuits back to pyGSTi (layer blocking required to separate u3 and cz) + qubit_mapping_dict_2 = {qk_test_circ._qbit_argument_conversion(i)[0]: sslbl for sslbl, i in qubit_mapping_dict.items()} # now map back + # in Qiskit 1.1.1, the method is called qbit_argument_conversion. In Qiskit >=1.2 (as far as Noah can tell), the method is called _qbit_argument_conversion. + + ps_ref_circ, _ = _Circuit.from_qiskit(qk_ref_circ, + qubit_conversion=qubit_mapping_dict_2, allow_different_gates_in_same_layer=False) + ps_ref_circ = ps_ref_circ.delete_idling_lines() + + for aux in auxlist: + ref_circ_dict[ps_ref_circ] += [aux] + ref_circ_id_lookup_dict[aux['id']] = ps_ref_circ + + ref_edesign = _FreeformDesign(ref_circ_dict) + + mirror_edesign = make_mirror_edesign(test_edesign=test_edesign, + ref_edesign=ref_edesign, + ref_id_lookup_dict=ref_circ_id_lookup_dict, + account_for_routing=False, + **mirroring_kwargs_dict + ) + + mirror_edesigns[k] = mirror_edesign + + if aggregate_subcircs: + return test_edesigns[0], mirror_edesigns[0] + else: + return test_edesigns, mirror_edesigns + + +def make_mirror_edesign(test_edesign: _FreeformDesign, + account_for_routing: bool, # handles permutation compiler optimization + ref_edesign: Optional[_FreeformDesign] = None, + ref_id_lookup_dict: Optional[dict] = None, + num_mcs_per_circ: int = 10, + num_ref_per_qubit_subset: int = 10, + mirroring_strategy: Literal['pauli_rc', 'central_pauli'] = 'pauli_rc', + gate_set: str = 'u3_cx_cz', + inverse: Optional[Callable[[_Circuit], _Circuit]] = None, + inv_kwargs: Optional[dict] = None, + rc_function: Optional[Callable[[_Circuit], Tuple[_Circuit, str]]] = None, + rc_kwargs: Optional[dict] = None, + state_initialization: Optional[Union[str, Callable[..., _Circuit]]] = None, + state_init_kwargs: Optional[dict] = None, + rand_state: Optional[_np.random.RandomState] = None) -> _CombinedExperimentDesign: + """ + Creates an experiment design containing the mirror circuits needed for mirror circuit fidelity estimation (MCFE). + + Parameters + ---------- + test_edesign : pygsti.protocols.FreeformDesign + The experiment design containing the test circuits. + + account_for_routing : bool + Indicates whether to account for routing in the design. + + ref_edesign : Optional[pygsti.protocols.FreeformDesign], optional + The experiment design containing the reference circuits. Default is None. + + ref_id_lookup_dict : Optional[dict], optional + A lookup dictionary for matching test circuits with reference circuits. Default is None. + + num_mcs_per_circ : int, optional + The number of mirror circuits to generate for each test circuit. Default is 10. + + num_ref_per_qubit_subset : int, optional + The number of SPAM reference circuits to use for each qubit subset. Default is 10. + + mirroring_strategy : Literal['pauli_rc', 'central_pauli'], optional + The strategy to use for mirroring ('pauli_rc' or 'central_pauli'). Default is 'pauli_rc'. + + gate_set : str, optional + The set of gates to be used in the design. Default is 'u3_cx_cz'. + + inverse : Optional[Callable[[_Circuit], _Circuit]], optional + A custom function to compute the inverse of a circuit. Default is None. + Signature: inverse(circ: pygsti.circuits.Circuit, ...) -> pygsti.circuits.Circuit + If providing a custom inverse function, 'circ' must be the circuit parameter name. + + inv_kwargs : Optional[dict], optional + Additional keyword arguments for a custom inverse function. Default is None. + + rc_function : Optional[Callable[[_Circuit], Tuple[_Circuit, str]]], optional + A custom function for random compilation. Default is None. + + Signature: rc_function(circ: pygsti.circuits.Circuit, rand_state: Optional[_np.random.RandomState] = None, ...) + -> Tuple[pygsti.circuits.Circuit, str] + The user-defined function must return the randomized circuit, + along with the expected bitstring measurement given the randomization. + + This function is called twice: + 1) On the reverse half of the circuit, when creating the init-test-ref_inv-init_inv circuit. ref_inv and init_inv are randomized. The bitstring should be the expected measurement *for the full circuit*, *not* ref_inv-init_inv in isolation from the forward half of the circuit. Pass the forward half of the circuit as a kwarg in rc_kwargs if necessary. + + 2) On the entire circuit, when creating the init-ref-ref_inv-init_inv circuit. All four pieces are randomized. The bitstring should again be the expected measurement for the full circuit, but it is more clear in this case than in the previous. + + rc_kwargs : Optional[dict], optional + Additional keyword arguments for the random compilation function. Default is None. + + state_initialization : Optional[Union[str, Callable[..., _Circuit]]], optional + A function or string for state initialization. Default is None. + Signature: state_initialization(qubits, rand_state: Optional[_np.random.RandomState] = None, ...) + -> pygsti.circuits.Circuit + If providing a custom state initialization, + the parameter names must be 'qubits' (list of the qubits being used) and 'rand_state'. + + state_init_kwargs : Optional[dict], optional + Additional keyword arguments for a custom state initialization function. Default is None. + + rand_state : Optional[_np.random.RandomState], optional + A random state for reproducibility. Default is None. + + Returns + -------- + pygsti.protocols.CombinedExperimentDesign + A combined experiment design containing the MCFE circuits. + """ + + # Power user could supply their own functions. + + if rand_state is None: + rand_state = _np.random.RandomState() + + qubit_subsets = defaultdict(list) + + test_ref_invs = defaultdict(list) + ref_ref_invs = defaultdict(list) + spam_refs = defaultdict(list) + + + central_pauli_allowed = True + + if ref_edesign is not None: + assert ref_id_lookup_dict is not None, "when providing separate test and reference compilations, you must provide a lookup dictionary for the reference circuits so they can be matched with the correct test circuits." + central_pauli_allowed = False # forward and reverse compilations must match for central pauli to be a valid fidelity estimation method. + + else: + print("using provided edesign for both reference and test compilations") + + for c, auxlist in _tqdm.tqdm(test_edesign.aux_info.items(), ascii=True, desc='Sampling mirror circuits'): + test_aux = auxlist[0] + # add a qubit subset to the pool if there's a circuit that uses it. if multiple identical circuits use it, do not double count; if distinct circuits use the same qubit subset, then do add it again. + qubits = c.line_labels + + qubit_subsets[test_aux['width']].append(qubits) + + if ref_edesign is not None: + # find the corresponding exact circuit for the inexact circuit 'c' using the look-up table generated + circ_id = test_aux['id'] + + exact_circ = ref_id_lookup_dict[circ_id] + + valid_test_ids = set(aux['id'] for aux in ref_edesign.aux_info[exact_circ]) + + assert circ_id in valid_test_ids, f"Invalid test ID {circ_id} for ref circuit corresponding to test IDs {valid_test_ids}" + + else: + exact_circ = c + + # R for "Reference" circuit, T for "Test" circuit + R = exact_circ + T = c + + R_inv = compute_inverse(circ=R, gate_set=gate_set, inverse=inverse, inv_kwargs=inv_kwargs) + + for j in range(num_mcs_per_circ): + + random_compiler = _rc.RandomCompilation(rc_strategy=mirroring_strategy, return_bs=True, + rand_state=rand_state) + + L_bareref = init_layer(qubits=qubits, gate_set=gate_set, state_initialization=state_initialization, rand_state=rand_state, state_init_kwargs=state_init_kwargs) + + L_refref = init_layer(qubits=qubits, gate_set=gate_set, state_initialization=state_initialization, rand_state=rand_state, state_init_kwargs=state_init_kwargs) + + L_bareref_inv = compute_inverse(circ=L_bareref, gate_set=gate_set, inverse=inverse, inv_kwargs=inv_kwargs) + + L_refref_inv = compute_inverse(circ=L_refref, gate_set=gate_set, inverse=inverse, inv_kwargs=inv_kwargs) + + if mirroring_strategy == 'pauli_rc': + if account_for_routing: # only the bare-ref circuit has non-trivial routing + assert ref_edesign is not None, "'account_for_routing' set to True but no ref_edesign has been provided to match routing. If you are not providing a 'ref_edesign' that needs matched with a 'test_edesign', please set 'account_for_routing' to False." + + T_routing = test_aux['routing_permutation'] + + ref_aux = ref_edesign.aux_info[R][0] + + R_routing = ref_aux['routing_permutation_inv'] + + mc_routing_perm = {k: R_routing[v] for k,v in T_routing.items()} # compute overall routing permutation for the composition of the test circuit and the inverse of the reference circuit + + L_bareref_inv = L_bareref_inv.map_state_space_labels(mc_routing_perm) + + + if rc_function is not None: + try: + # bare-ref + Rinv_Linv, L_T_Rinv_Linv_bs = rc_function(circ=R_inv + L_bareref_inv, rand_state=rand_state, **rc_kwargs) + L_T_Rinv_Linv = L_bareref + T + Rinv_Linv + + # ref-ref + L_R_Rinv_Linv, L_R_Rinv_Linv_bs = rc_function(circ=L_refref + R + R_inv + L_refref_inv, rand_state=rand_state, **rc_kwargs) + + # the assertions below check if the circuit addition in the function calls above have caused a line label reordering. If so, the MCFE code will compare bitstrings incorrectly, which is bad. This issue is fixable with enough Circuit.reorder_lines() calls, but the better approach is simply to ensure that all circuits in test_edesign and ref_edesign obey a lexicographical ordering. This is easily done by using something like 'c = c.reorder_lines(sorted(c.line_labels)) on the circuits *prior* to creating the mirror edesign. + + assert L_T_Rinv_Linv.line_labels == qubits, f'line labels have been modified/permuted: should be {qubits} but is {L_T_Rinv_Linv.line_labels} instead.' + + assert L_R_Rinv_Linv.line_labels == qubits, f'line labels have been modified/permuted: should be {qubits} but is {L_R_Rinv_Linv.line_labels} instead.' + + except Exception as e: + raise RuntimeError(f"User-provided RC function for gate set '{gate_set}' returned an error: {e}") + elif gate_set == 'u3_cx_cz': + + # bare-ref + Rinv_Linv, L_T_Rinv_Linv_bs = random_compiler.compile(R_inv + L_bareref_inv) + L_T_Rinv_Linv = L_bareref + T + Rinv_Linv + + # ref-ref + L_R_Rinv_Linv, L_R_Rinv_Linv_bs = random_compiler.compile(L_refref + R + R_inv + L_refref_inv) + + # the assertions below check if the circuit addition in the function calls above have caused a line label reordering. If so, the MCFE code will compare bitstrings incorrectly, which is bad. This issue is fixable with enough Circuit.reorder_lines() calls, but the better approach is simply to ensure that all circuits in test_edesign and ref_edesign obey a lexicographical ordering. This is easily done by using something like 'c = c.reorder_lines(sorted(c.line_labels)) on the circuits *prior* to creating the mirror edesign. + + assert L_T_Rinv_Linv.line_labels == qubits, f'line labels have been permuted: should be {qubits} but is {L_T_Rinv_Linv.line_labels} instead.' + + assert L_R_Rinv_Linv.line_labels == qubits, f'line labels have been permuted: should be {qubits} but is {L_R_Rinv_Linv.line_labels} instead.' + + elif gate_set == 'clifford': + #TODO: add clifford RC function + raise NotImplementedError("Clifford RC is not yet supported!") + elif gate_set == 'clifford_rz': + #TODO: add clifford_rz RC function + raise NotImplementedError("Clifford_rz RC is not yet supported!") + else: #we have no support for this gate set at this point, the user must provide a custom RC function + raise RuntimeError(f"No default RC function for gate set '{gate_set}' exists, you must provide your own!") + + elif mirroring_strategy == 'central_pauli': + assert central_pauli_allowed, "Central Pauli is not allowed when 'ref_edesign' is provided." + + if gate_set == 'u3_cx_cz': + + CP_Rinv_Linv, L_T_Rinv_Linv_bs = random_compiler.compile(circ=R_inv+L_refref_inv) + L_T_Rinv_Linv = L_refref + T + CP_Rinv_Linv + + # the assertion belows check if the circuit addition in the function call above has caused a line label reordering. If so, the MCFE code will compare bitstrings incorrectly, which is bad. This issue is fixable with enough Circuit.reorder_lines() calls, but the better approach is simply to ensure that all circuits in test_edesign and ref_edesign obey a lexicographical ordering. This is easily done by using something like 'c = c.reorder_lines(sorted(c.line_labels)) on the circuits *prior* to creating the mirror edesign. + + assert L_T_Rinv_Linv.line_labels == qubits, f'line labels have been permuted: should be {qubits} but is {L_T_Rinv_Linv.line_labels} instead.' + + elif gate_set == 'clifford': + #TODO: add clifford CP function + raise NotImplementedError("Clifford CP is not yet supported!") + elif gate_set == 'clifford_rz': + #TODO: add clifford_rz CP function + raise NotImplementedError("Clifford_rz CP is not yet supported!") + else: #we have no support for this gate set at this point, the user must provide a custom CP function + raise RuntimeError(f"No default CP function for gate set '{gate_set}' exists, you must provide your own!") + + else: + raise RuntimeError("'mirroring_strategy' must be either 'pauli_rc' or 'central_pauli'") + + + L_T_Rinv_Linv_aux = [{'base_aux': a, + 'idealout': L_T_Rinv_Linv_bs, + 'qs_to_measure': L_T_Rinv_Linv.line_labels, + 'id': j} for a in auxlist] + + test_ref_invs[L_T_Rinv_Linv] = test_ref_invs[L_T_Rinv_Linv] + L_T_Rinv_Linv_aux + + if mirroring_strategy == 'pauli_rc': + L_R_Rinv_Linv_aux = [{'base_aux': a, + 'idealout': L_R_Rinv_Linv_bs, + 'qs_to_measure': L_R_Rinv_Linv.line_labels, + 'id': j} for a in auxlist] + + ref_ref_invs[L_R_Rinv_Linv] = ref_ref_invs[L_R_Rinv_Linv] + L_R_Rinv_Linv_aux + + + for w, width_subsets in qubit_subsets.items(): + + unique_subsets = set(width_subsets) + print(f'Sampling reference circuits for width {w} with {len(unique_subsets)} subsets') + for unique_subset in unique_subsets: + for j in _tqdm.tqdm(range(num_ref_per_qubit_subset), ascii=True, desc=f'Sampling reference circuits for subset {unique_subset}'): + L = init_layer(qubits=unique_subset, gate_set=gate_set, state_initialization=state_initialization, rand_state=rand_state, state_init_kwargs=state_init_kwargs) + spam_ref = L + compute_inverse(circ=L, gate_set=gate_set, inverse=inverse, inv_kwargs=inv_kwargs) + spam_refs[spam_ref].append({'idealout': '0'*w, 'id': j, 'qs_to_measure': spam_ref.line_labels, 'width': w}) + + + edesigns = {} + if mirroring_strategy == 'pauli_rc': + edesigns['br'] = _FreeformDesign(test_ref_invs) + edesigns['rr'] = _FreeformDesign(ref_ref_invs) + edesigns['ref'] = _FreeformDesign(spam_refs) + elif mirroring_strategy == 'central_pauli': + edesigns['cp'] = _FreeformDesign(test_ref_invs) + edesigns['cpref'] = _FreeformDesign(spam_refs) + else: + raise RuntimeError("unknown mirroring strategy!") + + return _CombinedExperimentDesign(edesigns) + + +def compute_inverse(circ: _Circuit, + gate_set: str, + inverse: Optional[Callable[[_Circuit], _Circuit]] = None, + inv_kwargs: Optional[dict] = None) -> _Circuit: + """ + Computes the inverse of a given circuit based on the specified gate set. + + Parameters + ---------- + circ : pygsti.circuits.Circuit + The circuit for which to compute the inverse. + + gate_set : str + The set of gates used in the circuit (e.g., 'u3_cx_cz'). + + inverse : Optional[Callable[[_Circuit], _Circuit]], optional + A custom function to compute the inverse of the circuit. Default is None. + + inv_kwargs : Optional[dict], optional + Additional keyword arguments for the custom inverse function. Default is None. + + Returns + -------- + pygsti.circuits.Circuit + A new circuit that is the inverse of the input circuit. + """ + + if inverse is not None: + try: + circ_inv = inverse(circ=circ, **inv_kwargs) + except Exception as e: + raise RuntimeError(f"User-provided inverse function for gate set '{gate_set}' returned an error: {e}") + if gate_set == 'u3_cx_cz': + circ_inv = _rc.u3_cx_cz_inv(circ) + elif gate_set == 'clifford': + #TODO: add clifford inverse function to circuit_inverse.py + raise NotImplementedError("Clifford inversion is not yet supported!") + elif gate_set == 'clifford_rz': + #TODO: add clifford_rz inverse function to circuit_inverse.py + raise NotImplementedError("Clifford_rz inversion is not yet supported!") + else: #we have no support for this gate set at this point, the user must provide a custom inverse function + raise RuntimeError(f"No default inverse function for gate set '{gate_set}' exists, you must provide your own!") + + return circ_inv + +def init_layer(qubits: List[str], + gate_set: str, + state_initialization: Optional[Union[str, Callable[..., _Circuit]]] = None, + state_init_kwargs: Optional[dict] = None, + rand_state: Optional[_np.random.RandomState] = None, + ) -> _Circuit: + """ + Create initial layer for mirror circuit. + + Parameters + ---------- + qubits : list + A list of qubit labels for the layer. + + gate_set : str + The gate set for the mirror circuit (e.g., 'u3_cx_cz'). + + state_initialization : Optional[Union[str, Callable[..., _Circuit]]], optional + Custom function for creating initial layer, or the string 'none'. If + 'none' is provided, then an empty initial layer is created, i.e., there + is no state preparation layer. Default is None, which prepares an initial layer + according to the gate set. + + state_init_kwargs : Optional[dict], optional + Additional keyword arguments for the state initialization function. + Default is None. + + rand_state : Optional[_np.random.RandomState], optional + A random state for reproducibility. Default is None. + + + Returns + -------- + pygsti.circuits.Circuit + A new circuit representing the initial layer of gates. + """ + + if state_initialization == 'none': + L = _Circuit([], qubits) + elif state_initialization is not None: + try: + L = state_initialization(qubits=qubits, rand_state=rand_state, **state_init_kwargs) + except Exception as e: + raise RuntimeError(f"User-provided state_initialization function for gate set '{gate_set}' returned an error: {e}") + elif gate_set == 'u3_cx_cz': + prep_layer = _rc.haar_random_u3_layer(qubits, rand_state) + L = _Circuit([prep_layer], qubits) + elif gate_set == 'clifford': + #TODO: add clifford default initialization + raise NotImplementedError("Clifford state initialization is not yet supported!") + elif gate_set == 'clifford_rz': + #TODO: add clifford_rz default initialization + raise NotImplementedError("Clifford_rz state initialization is not yet supported!") + else: #we have no support for this gate set at this point, the user must provide a custom state initialization + raise RuntimeError(f"No default state_initialization function for gate set '{gate_set}' exists, you must provide your own!") + + return L \ No newline at end of file diff --git a/pygsti/protocols/vbdataframe.py b/pygsti/protocols/vbdataframe.py index ef7998f9f..21b64f288 100644 --- a/pygsti/protocols/vbdataframe.py +++ b/pygsti/protocols/vbdataframe.py @@ -9,10 +9,34 @@ # in compliance with the License. You may obtain a copy of the License at # http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. #*************************************************************************************************** +from __future__ import annotations +from typing import Callable, Literal, Optional, TYPE_CHECKING +if TYPE_CHECKING: + try: + import matplotlib as _mp + except: + pass + + import pandas as _pandas import numpy as _np +import tqdm as _tqdm from scipy.stats import chi2 as _chi2 +from pygsti.protocols import ( + ExperimentDesign as _ExperimentDesign, + FreeformDesign as _FreeformDesign, + ProtocolData as _ProtocolData +) + +from pygsti.tools.mcfetools import (hamming_distance_counts, + effective_polarization, + predicted_process_fidelity_for_central_pauli_mcs, + fidelity_to_polarization, + polarization_to_success_probability, + rc_predicted_process_fidelity, + rc_bootstrap_predicted_pfid + ) def _calculate_summary_statistic(x, statistic, lower_cutoff=None): """ @@ -33,21 +57,6 @@ def _calculate_summary_statistic(x, statistic, lower_cutoff=None): return _np.max([v, lower_cutoff]) -def polarization_to_success_probability(p, n): - """ - Inverse of success_probability_to_polarization. - """ - return p * (1 - 1 / 2**n) + 1 / 2**n - - -def success_probability_to_polarization(s, n): - """ - Maps a success probablity s for an n-qubit circuit to - the polarization, defined by `p = (s - 1/2^n)/(1 - 1/2^n)` - """ - return (s - 1 / 2**n) / (1 - 1 / 2**n) - - def classify_circuit_shape(success_probabilities, total_counts, threshold, significance=0.05): """ @@ -154,8 +163,13 @@ class VBDataFrame(object): A class for storing a DataFrame that contains volumetric benchmarking data, and that has methods for manipulating that data and creating volumetric-benchmarking-like plots. """ - def __init__(self, df, x_axis='Depth', y_axis='Width', x_values=None, y_values=None, - edesign=None): + def __init__(self, + df: _pandas.DataFrame, + x_axis: str = 'Depth', + y_axis: str = 'Width', + x_values: Optional[str] = None, + y_values: Optional[str] = None, + edesign: Optional[_ExperimentDesign] = None): """ Initialize a VBDataFrame object. @@ -200,6 +214,177 @@ def __init__(self, df, x_axis='Depth', y_axis='Width', x_values=None, y_values=N else: self.y_values = y_values + @classmethod + def from_mirror_experiment(cls, + unmirrored_design: _FreeformDesign, + mirrored_data: _ProtocolData, + include_dropped_gates: bool = False, + bootstrap: bool = True, + num_bootstraps: int = 50, + rand_state: Optional[_np.random.RandomState] = None, + verbose: bool = False, + ) -> VBDataFrame: + """ + Create a dataframe from MCFE data and edesigns. + + Parameters + ---------- + unmirrored_design: pygsti.protocols.protocol.FreeformDesign + Edesign containing the circuits whose process fidelities are to be estimated. + + mirrored_data: pygsti.protocols.protocol.ProtocolData + Data object containing the full mirror edesign and the outcome counts for each + circuit in the full mirror edesign. + + include_dropped_gates: bool + Whether to include the number of gates dropped from each subcircuit during + subcircuit creation. This flag should be set to False for noise benchmark and + fullstack benchmark analysis, but can be set to True for subcircuit benchmark + analysis. Default is False. + + bootstrap: bool + Toggle the calculation of error bars from bootstrapped process fidelity calculations. If True, + error bars are calculated. If False, error bars are not calculated. + + num_bootstraps: int + Number of samples to draw from the bootstrapped process fidelity calculations. This argument + is ignored if 'bootstrap' is False. + + rand_state: np.random.RandomState + random state used to seed bootstrapping. If 'bootstrap' is set to False, this argument is ignored. + + verbose: bool + Toggle print statements with debug information. If True, print statements are + turned on. If False, print statements are omitted. + + Returns + --------- + VBDataFrame + A VBDataFrame whose dataframe contains calculated MCFE values and circuit statistics. + + """ + + eff_pols = {k: {} for k in mirrored_data.edesign.keys()} + + # Stats about the base circuit + dropped_gates = {} + occurrences = {} + + # Get a dict going from id to circuit for easy lookup + reverse_circ_ids = {} + for k,v in unmirrored_design.aux_info.items(): + if isinstance(v, dict): + reverse_circ_ids[v['id']] = k + else: + for entry in v: + reverse_circ_ids[entry['id']] = k + seen_keys = set() + + for c in _tqdm.tqdm(mirrored_data.dataset, ascii=True, desc='Calculating effective polarizations'): + for edkey, ed in mirrored_data.edesign.items(): + auxlist = ed.aux_info.get(c, None) + if auxlist is None: + continue + elif isinstance(auxlist, dict): + auxlist = [auxlist] + + for aux in auxlist: + if edkey.endswith('ref'): + # For reference circuits, only width matters, so aggregate on that now + key = (aux['width'], c.line_labels) + else: + key = (aux['base_aux']['width'], aux['base_aux']['depth'], aux['base_aux']['id'], c.line_labels) + + + if edkey == 'br': + base_line_labels = reverse_circ_ids[aux['base_aux']['id']].line_labels + if c.line_labels != base_line_labels: + raise RuntimeError('line labels permuted') + + # Calculate effective polarization + hdc = hamming_distance_counts(mirrored_data.dataset[c], c, (aux['idealout'],)) + ep = effective_polarization(hdc) + + # Append to other mirror circuit samples + eps = eff_pols[edkey].get(key, []) + eps.append(ep) + eff_pols[edkey][key] = eps + + if edkey.endswith('ref') or key in seen_keys: + # Skip statistic gathering for reference circuits or base circuits we've seen already + continue + + # TODO: add 2Q gate density metric + + if include_dropped_gates: + dropped_gates[key] = aux['base_aux']['dropped_gates'] + occurrences[key] = len(auxlist) + + seen_keys.add(key) + + # Calculate process fidelities + df_data = {} + for i, key in enumerate(sorted(list(seen_keys), key=lambda x: x[2])): + cp_pfid = cp_pol = cp_success_prob = None + if 'cp' in mirrored_data.edesign and 'cpref' in mirrored_data.edesign: + if verbose and i == 0: print('Central pauli data detected, computing CP process fidelity') + cp_pfid = predicted_process_fidelity_for_central_pauli_mcs(eff_pols['cp'][key], + eff_pols['cpref'][(key[0], key[3])], + key[0]) + + cp_pol = fidelity_to_polarization(cp_pfid, key[0]) + cp_success_prob = polarization_to_success_probability(cp_pol, key[0]) + + rc_pfid = rc_pol = rc_success_prob = rc_pfid_stdev = None + if 'rr' in mirrored_data.edesign and 'br' in mirrored_data.edesign and 'ref' in mirrored_data.edesign: + if verbose and i == 0: print('Random compilation data detected, computing RC process fidelity') + + + if verbose: + spam_ref_key = (key[0], key[3]) + print(f'key for SPAM reference circuits: {spam_ref_key}') + print(f"number of SPAM reference circuits: {len(eff_pols['ref'][spam_ref_key])}") + + rc_pfid = rc_predicted_process_fidelity(eff_pols['br'][key], + eff_pols['rr'][key], + eff_pols['ref'][(key[0], key[3])], + key[0]) + + rc_pol = fidelity_to_polarization(rc_pfid, key[0]) # should this be fidelity_to_polarization? + rc_success_prob = polarization_to_success_probability(rc_pol, key[0]) + + if bootstrap: + # do bootstrapping to obtain confidence intervals + rc_pfid_stdev = rc_bootstrap_predicted_pfid(brs=eff_pols['br'][key], + rrs=eff_pols['rr'][key], + refs=eff_pols['ref'][(key[0], key[3])], + n=key[0], + num_bootstraps=num_bootstraps, + rand_state=rand_state + ) + + data_dict = {'Width': key[0], 'Depth': key[1], 'Circuit Id': key[2], + 'CP Process Fidelity': cp_pfid, + 'CP Polarization': cp_pol, 'CP Success Probability': cp_success_prob, + 'RC Process Fidelity': rc_pfid, 'RC Process Fidelity stdev': rc_pfid_stdev, + 'RC Polarization': rc_pol, 'RC Success Probability': rc_success_prob, + 'Occurrences': occurrences[key], + } + + if include_dropped_gates: + data_dict['Dropped Gates'] = dropped_gates[key] + + # Depth is doubled for conventions (same happens for RMC/PMC) + df_data[i] = data_dict + + + df = _pandas.DataFrame.from_dict(df_data, orient='index') + df = df.sort_values(by='Circuit Id') + df = df.reset_index(drop=True) + + return cls(df, x_axis='Depth', y_axis='Width') + + def select_column_value(self, column_label, column_value): """ Filters the dataframe, by discarding all rows of the @@ -480,3 +665,136 @@ def capability_regions(self, metric='polarization', threshold=1 / _np.e, signifi capreg = trimmed_capreg return capreg + + + + def create_vb_plot(self, + title: str, + accumulator: Callable = _np.mean, + cp_or_rc: Literal['cp', 'rc'] = 'rc', + show_dropped_gates: bool = False, + dg_accumulator: Callable = _np.mean, + cmap: _mp.colors.Colormap = None, + margin: float = 0.15, + save_fig: bool = False, + fig_path: str = None, + fig_format: str = None): + + """ + Generate process fidelity volumetric benchmarking (VB) plot from dataframe. + This function is designed with subcircuit volumetric benchmarking in mind, + where the x-axis is the subicrcuit depth and the y-axis is the subcircuit + width. + + Parameters + ------------ + title : str + The title of the plot. + + accumulator : callable, optional + Function used to accumulate process fidelities for a (width, depth) pair + on the VB plot. Default is np.mean. + + cp_or_rc : str, optional + Whether the process fidelities were computed via randomly compiled + circuits ('rc') or central Pauli mirroring ('cp'). Default is 'rc'. + + show_dropped_gates : bool, optional + Whether the plot should visualize the average (see dg_accumulator) + number of dropped gates for each subcircuit width-depth pair. + Subcircuit sampling can drop gates when a gate has only partial support + on the qubits in the selected width subset. + + dg_accumulator : callable, optional + Function used to accumulate the dropped gate counts + for a (width, depth) pair on the VB plot. Default is np.mean. + + cmap : matplotlib.colors.Colormap, optional + Colormap to use for plotting process fidelities. Default is spectral. + + margin : float, optional + Margin between adjacent width-depth pairs in the VB plot. + Default is 0.15. + + save_fig : bool, optional + Whether to write the generated VB plot to file. Default is False. + + fig_path : str, optional + If `save_fig` is set to True, this argument is used as the path + the figure is saved to. If `save_fig` is False, this argument + is ignored. Default is None + + fig_format : str, optional + If `save_fig` is set to True, this argument is the file format + for the generated VB plot. If `save_fig` is False, this argument + is ignored. Acceptable values are any file format recognized by matplotlib. + """ + try: + import matplotlib.pyplot as _plt + import matplotlib as _mp + except: + raise RuntimeError('matplotlib is required for dataframe plotting, and does not appear to be installed.') + + if cmap is None: + cmap = _mp.colormaps['Spectral'] + + + fig = _plt.figure(figsize=(5.5,8)) + ax = _plt.gca() + ax.set_aspect('equal') + + x_axis = 'Depth' + y_axis = 'Width' + + if cp_or_rc == 'rc': + c_axis = 'RC Process Fidelity' + elif cp_or_rc == 'cp': + c_axis = 'CP Process Fidelity' + else: + raise ValueError(f"invalid argument passed to 'cp_or_rc': {cp_or_rc}") + + xticks = self.dataframe[x_axis].unique() + yticks = self.dataframe[y_axis].unique() + + xmap = {j: i for i, j in enumerate(xticks)} + ymap = {j: i for i, j in enumerate(yticks)} + + acc_values = self.dataframe.groupby([x_axis, y_axis])[c_axis].apply(accumulator) + if show_dropped_gates: + dg_values = self.dataframe.groupby([x_axis, y_axis])["Dropped Gates"].apply(dg_accumulator) + + for xlbl, ylbl in acc_values.keys(): + x = xmap[xlbl] + y = ymap[ylbl] + + side = 0.5-margin + + cval = acc_values[(xlbl, ylbl)] + if show_dropped_gates: + dgval = dg_values[(xlbl, ylbl)] + # Inner square: Smaller by ratio of dropped gates + inside = side * (xlbl-dgval)/xlbl + + else: + inside = side + + inpoints = [(x-side, y-side), (x-side, y+inside), (x+inside , y+inside), (x+inside, y-side)] + tin = _plt.Polygon(inpoints, color=cmap(cval)) + ax.add_patch(tin) + + # Outer square + outer_points = [(x-side, y-side), (x-side, y+side), (x+side, y+side), (x+side, y-side)] + tout = _plt.Polygon(outer_points, edgecolor='k', fill=None) + ax.add_patch(tout) + + _plt.xlabel(x_axis, {'size': 20}) + _plt.ylabel(y_axis, {'size': 20}) + _plt.xticks(list(range(len(xticks))), labels=xticks) + _plt.yticks(list(range(len(yticks))), labels=yticks) + _plt.xlim([-0.5, len(xticks)-0.5]) + _plt.ylim([-0.5, len(yticks)-0.5]) + ax.tick_params(axis='both', which='major', labelsize=14) + _plt.title(title, {"size": 24}) + + if save_fig: + _plt.savefig(fig_path, format=fig_format) \ No newline at end of file diff --git a/pygsti/tools/bsqale.py b/pygsti/tools/bsqale.py new file mode 100644 index 000000000..812a9ce3b --- /dev/null +++ b/pygsti/tools/bsqale.py @@ -0,0 +1,343 @@ +""" +B-Sqale is a software tool for creating scalable, robust benchmarks from any quantum circuit. +Please see for more information. +""" +#*************************************************************************************************** +# Copyright 2015, 2019, 2025 National Technology & Engineering Solutions of Sandia, LLC (NTESS). +# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights +# in this software. +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +# in compliance with the License. You may obtain a copy of the License at +# http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. +#*************************************************************************************************** + +from __future__ import annotations +from typing import Tuple, Optional, Dict, List, Union, Any, TYPE_CHECKING +if TYPE_CHECKING: + try: + import qiskit + except: + pass + + +import warnings as _warnings + +from pygsti.protocols import ( + ProtocolData as _ProtocolData, + FreeformDesign as _FreeformDesign, + CombinedExperimentDesign as _CombinedExperimentDesign, + VBDataFrame as _VBDataFrame, + mirror_edesign as _mirror +) + +import numpy as _np + + +def noise_mirror_benchmark(qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + mirroring_kwargs_dict: Dict[str, Any] = {} + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a noise benchmark from transpiled Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a noise benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring. If an arg + is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + Returns + --------- + Tuple + pygsti.protocols.FreeformDesign + Experiment design containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + + pygsti.protocols.CombinedExperimentDesign + Experiment design containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'noise_mirror_benchmark' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + return _mirror.qiskit_circuits_to_mirror_edesign(qk_circs, mirroring_kwargs_dict) + + +def fullstack_mirror_benchmark(qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + qk_backend: Optional[qiskit.providers.BackendV2] = None, + coupling_map: Optional[qiskit.transpiler.CouplingMap] = None, + basis_gates: Optional[List[str]] = None, + transpiler_kwargs_dict: Dict[str, Any] = {}, + mirroring_kwargs_dict: Dict[str, Any] = {}, + num_transpilation_attempts: int = 100, + return_qiskit_time: bool = False + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a full-stack benchmark from high-level Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a full-stack benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + qk_backend : qiskit-ibm-runtime.IBMBackend, optional + IBM backend whose native gate set, connectivity, and error rates + are to be targeted when doing a full-stack transpilation. Fake backends + are also acceptable. If not provided, both `coupling_map` and `basis_gates` + must be provided, and certain transpiler optimizations that depend on + backend error rates may not be accessible. + + coupling_map : qiskit.transpiler.CouplingMap, optional + couplinp map for a target backend that must be provided if `qk_backend` is None. + This argument is ignored if `qk_backend` is not None. + + basis_gates : list[str], optional + list of native gates for a target backend that must be provided + if `qk_backend` is None. + This argument is ignored if `qk_backend` is not None. + + transpiler_kwargs_dict : dict, optional + keyword arguments that are passed to the Qiskit transpiler for full-stack + transpilation. Please see the Qiskit transpiler documentation + for a comprehensive list of options. If any or all of the following keys + are not provided, the listed defaults are used: + + 'optimization_level': default 3 + 'approximation_degree': default 1.0 + 'seed_transpiler': default None + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring. If an arg + is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + num_transpilation_attempts : int, optional + number of times to attempt full-stack circuit transpilation. Circuit mirroring requires + that the transpilation not use ancilla qubits, which is difficult to enforce with the + other transpiler options. Instead, we adopt a try-until-success strategy, which may fail + if an insufficient number of attempts are allowed. Increase this number if you are + having issues with the transpilation. If not provided, the default is 100 attempts. + + return_qiskit_time : bool, optional + Debug flag that sets whether or not to report the time spent in the Qiskit transpiler. + Qiskit transpilation is often the most costly part of benchmark creation and it can be + helpful to know how much time it is consuming. + + Returns + --------- + Tuple + pygsti.protocols.FreeformDesign + Experiment design containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + + pygsti.protocols.CombinedExperimentDesign + Experiment design containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. + + float, optional + amount of time spent in Qiskit transpiler. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'fullstack_mirror_benchmark' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + return _mirror.qiskit_circuits_to_fullstack_mirror_edesign(qk_circs, #not transpiled + qk_backend, + coupling_map, + basis_gates, + transpiler_kwargs_dict, + mirroring_kwargs_dict, + num_transpilation_attempts, + return_qiskit_time + ) + + +def subcircuit_mirror_benchmark(qk_circs: Union[Dict[Any, qiskit.QuantumCircuit], List[qiskit.QuantumCircuit]], + aggregate_subcircs: bool, + width_depth_dict: Dict[int, List[int]], + coupling_map: qiskit.transpiler.CouplingMap, + instruction_durations: qiskit.transpiler.InstructionDurations, + subcirc_kwargs_dict: Dict[str, Any] = {}, + mirroring_kwargs_dict: Dict[str, Any] = {} + ) -> Tuple[_FreeformDesign, _CombinedExperimentDesign]: + """ + Create a subcircuit benchmark from transpiled Qiskit circuits. + + Parameters + ---------- + qk_circs : List[qiskit.QuantumCircuit] | Dict[qiskit.QuantumCircuit] + Qiskit QuantumCircuits from which a subcircuit benchmark is to be created. + If a dictionary is provided, those keys are used. + If a list is provided, default integer keys are used. + + aggregate_subcircs : bool + Whether or not the provided Qiskit circuits should be used to create + one combined subcircuit experiment design or kept separate. + Circuit aggregration can be useful if the provided circuits are all + instances of the same 'kind' of the circuit, e.g., all + Bernstein-Vazirani circuits with different secret keys. + + width_depth_dict : dict[int, list[int]] + dictionary whose keys are subcircuit widths and whose values are + lists of depths to snip out for that width. + + coupling_map : str or qiskit.transpiler.CouplingMap + coupling map for the device the Qiskit circuits were transpiled to. + If 'all-to-all', an all-to-all coupling map is used. + If 'linear', a linear topology is used. + + instruction_durations : qiskit.transpiler.InstructionDurations + instruction durations for each gate in the target device. These + durations are needed to calculate the appropriate delay time + when only idling qubits are sampled out from a full circuit layer. + + + subcirc_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in subcircuit selection. + If an arg is not provided, a default value is used. + The args are: + 'num_samples_per_width_depth': default 10. number of subcircuits to sample + from each full circuit. if `aggregate_subcircuits` is set to + True, `num_samplse_per_width_depth` subcircuits are drawn from + each full circuit and combined into one experiment design. + + 'rand_state': default None. np.random.RandomState to be used for subcircuit selection. + + mirroring_kwargs_dict : dict, optional + dictionary of keyword arguments to be used in circuit mirroring + If an arg is not provided, a default value is used. + The args are: + 'mirror_circuits_per_circ': default 10. The number of mirror circuits of the + test-exact and exact-exact varieties to be used for the process fidelity estimation + of each provided Qiskit circuit. + + 'num_ref_per_qubit_subset': default 10. The number of SPAM reference circuits to use + for each qubit subset that is represented among the provided Qiskit circuits. + + 'rand_state': default None. np.random.RandomState to be used for circuit mirroring. + + Returns + --------- + Tuple + dict[hashable, pygsti.protocols.FreeformDesign] or pygsti.protocols.FreeformDesign + Experiment design(s) containing the pyGSTi conversion of all Qiskit circuits that + were passed in. Does not need executed, but is needed for fidelity calculations. + A dictionary is returned if `aggregate_subcircs` is False, otherwise a FreeformDesign + is returned. + + dict[hashable, pygsti.protocols.CombinedExperimentDesign] or pygsti.protocols.CombinedExperimentDesign + Experiment design(s) containing all mirror circuits that must be executed + in order to perform mirror circuit fidelity estimation. A dictionary is returned + if `aggregate_subcircs` is False, otherwise a FreeformDesign is returned. + """ + + try: + import qiskit + if qiskit.__version__ != '2.1.1': + _warnings.warn("The function 'noise_mirror_benchmark' is designed for qiskit 2.1.1." \ + "Your version is " + qiskit.__version__) + + except: + raise RuntimeError('Qiskit is required for this operation, and does not appear to be installed.') + + return _mirror.qiskit_circuits_to_subcircuit_mirror_edesign(qk_circs, + aggregate_subcircs, + width_depth_dict, + coupling_map, + instruction_durations, + subcirc_kwargs_dict, + mirroring_kwargs_dict + ) # qk_circs must already be transpiled to the device + + +def calculate_mirror_benchmark_results(unmirrored_design: _FreeformDesign, + mirrored_data: _ProtocolData, + include_dropped_gates: bool = False, + bootstrap: bool = True, + num_bootstraps: int = 50, + rand_state: _np.random.RandomState = None, + verbose: bool = False, + ) -> _VBDataFrame: + """ + Create a dataframe from MCFE data and edesigns. + + Parameters + ---------- + unmirrored_design: pygsti.protocols.protocol.FreeformDesign + Edesign containing the circuits whose process fidelities are to be estimated. + + mirrored_data: pygsti.protocols.protocol.ProtocolData + Data object containing the full mirror edesign and the outcome counts for each + circuit in the full mirror edesign. + + include_dropped_gates: bool + Whether to include the number of gates dropped from each subcircuit during + subcircuit creation. This flag should be set to False for noise benchmark and + fullstack benchmark analysis, but can be set to True for subcircuit benchmark + analysis. Default is False. + + bootstrap: bool + Toggle the calculation of error bars from bootstrapped process fidelity calculations. If True, + error bars are calculated. If False, error bars are not calculated. + + num_bootstraps: int + Number of samples to draw from the bootstrapped process fidelity calculations. This argument + is ignored if 'bootstrap' is False. + + rand_state: np.random.RandomState + random state used to seed bootstrapping. If 'bootstrap' is set to False, this argument is ignored. + + verbose: bool + Toggle print statements with debug information. If True, print statements are + turned on. If False, print statements are omitted. + + Returns + --------- + VBDataFrame + A VBDataFrame whose dataframe contains calculated MCFE values and circuit statistics. + + """ + + return _VBDataFrame.from_mirror_experiment(unmirrored_design, mirrored_data, + include_dropped_gates, + bootstrap, + num_bootstraps, + rand_state, + verbose + ) \ No newline at end of file diff --git a/pygsti/tools/internalgates.py b/pygsti/tools/internalgates.py index b9dab7420..1b088a2d1 100644 --- a/pygsti/tools/internalgates.py +++ b/pygsti/tools/internalgates.py @@ -10,6 +10,16 @@ # http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. #*************************************************************************************************** +from __future__ import annotations +from typing import Dict, List, Tuple, TYPE_CHECKING +if TYPE_CHECKING: + try: + import qiskit + except: + pass + +import warnings as _warnings + import numpy as _np import scipy.linalg as _spl @@ -18,12 +28,11 @@ from pygsti.baseobjs.unitarygatefunction import UnitaryGateFunction as _UnitaryGateFunction from pygsti.tools.gatetools import sigmax, sigmay, sigmaz, sigmaxz - class Gzr(_UnitaryGateFunction): shape = (2, 2) def __call__(self, theta): - return _np.array([[1., 0.], [0., _np.exp(-1j * float(theta[0]))]]) + return _np.array([[1., 0.], [0., _np.exp(1j * float(theta[0]))]]) @classmethod def _from_nice_serialization(cls, state): @@ -35,13 +44,26 @@ class Gczr(_UnitaryGateFunction): def __call__(self, theta): return _np.array([[1., 0., 0., 0.], [0., 1., 0., 0.], [0., 0., 1., 0.], - [0., 0., 0., _np.exp(-1j * float(theta[0]))]], complex) + [0., 0., 0., _np.exp(1j * float(theta[0]))]], complex) @classmethod def _from_nice_serialization(cls, state): return super(Gczr, cls)._from_nice_serialization(state) + +class Gu3(_UnitaryGateFunction): + shape = (2,2) + def __call__(self, arg): + theta, phi, lamb = (float(arg[0]), float(arg[1]), float(arg[2])) + return _np.array([[_np.cos(theta/2), -_np.exp(1j*lamb)*_np.sin(theta/2)], + [_np.exp(1j*phi)*_np.sin(theta/2), _np.exp(1j*(phi + lamb))*_np.cos(theta/2)]]) + + @classmethod + def _from_nice_serialization(cls, state): + return super(Gu3, cls)._from_nice_serialization(state) + + def internal_gate_unitaries(): """ The unitaries for the *internally* defined gates. @@ -203,6 +225,7 @@ def standard_gatename_unitaries(): * 'Gt', 'Gtdag' : the T and inverse T gates (T is a Z rotation by pi/4). * 'Gzr' : a parameterized gate that is a Z rotation by an angle, where when the angle = pi then it equals Z. * 'Gn' : N gate, pi/2 rotation about the (np.sqrt(3)/2, 0, -1/2) axis of the Bloch sphere, native gate in some spin qubit systems. + * 'Gu3' : parameterized gate that can encode any single-qubit rotation with three parameters (theta, phi, lambda). Mostly, pyGSTi does not assume that a gate with one of these names is indeed the unitary specified here. Instead, these names are intended as short-hand @@ -308,6 +331,8 @@ def u_op(exp): std_unitaries['Gzr'] = Gzr() std_unitaries['Gczr'] = Gczr() + std_unitaries['Gu3'] = Gu3() + #Add these at the end, since we don't want unitary_to_standard_gatenames to return these "shorthand" names std_unitaries['Gx'] = std_unitaries['Gxpi2'] std_unitaries['Gy'] = std_unitaries['Gypi2'] @@ -365,8 +390,8 @@ def unitary_to_standard_gatename(unitary, up_to_phase = False, return_phase = Fa return std_name, phase else: return std_name - return None + def standard_gatenames_stim_conversions(): """ A dictionary converting the gates with standard names to stim tableus for these gates. Currently is only capable of converting @@ -566,6 +591,7 @@ def cirq_gatenames_standard_conversions(): return cirq_to_standard_mapping + def standard_gatenames_quil_conversions(): """ A dictionary converting the gates with standard names to the QUIL names for these gates. @@ -867,6 +893,104 @@ def standard_gatenames_openqasm_conversions(version='u3'): return std_gatenames_to_qasm, std_gatenames_to_argmap +def qiskit_gatenames_standard_conversions() -> Dict[str, List[str, bool]]: + """ + A dictionary converting Qiskit gates (based on Instruction.name) + to built-in pyGSTi names for these gates. + Additional flag in dict value indicates if the gate has params. + + Does not currently support conversion of all qiskit gate types. + + Returns + --------- + Dictionary mapping qiskit gate names to pyGSTi gate names and a flag + for whether or not the gate has parameters. + """ + + qiskit_to_standard_mapping = {} + + qiskit_to_standard_mapping['id'] = ['Gi', False] + qiskit_to_standard_mapping['x'] = ['Gxpi', False] + qiskit_to_standard_mapping['y'] = ['Gypi', False] + qiskit_to_standard_mapping['z'] = ['Gzpi', False] + qiskit_to_standard_mapping['sx'] = ['Gxpi2', False] + qiskit_to_standard_mapping['sxdg'] = ['Gxmpi2', False] + qiskit_to_standard_mapping['t'] = ['Gt', False] + qiskit_to_standard_mapping['h'] = ['Gh', False] + + qiskit_to_standard_mapping['rz'] = ['Gzr', True] + qiskit_to_standard_mapping['ry'] = ['Gyr', True] + qiskit_to_standard_mapping['rx'] = ['Gxr', True] + qiskit_to_standard_mapping['u'] = ['Gu3', True] + qiskit_to_standard_mapping['u3'] = ['Gu3', True] + + qiskit_to_standard_mapping['cx'] = ['Gcnot', False] + qiskit_to_standard_mapping['cz'] = ['Gcphase', False] + qiskit_to_standard_mapping['ecr'] = ['Gecres', False] + qiskit_to_standard_mapping['swap'] = ['Gswap', False] + + qiskit_to_standard_mapping['delay'] = ['Gdelay', True] + + return qiskit_to_standard_mapping + + +def standard_gatenames_qiskit_conversions() -> Dict[str, Tuple[qiskit.circuit.Instruction, str, bool]]: + """ + A dictionary converting the gates with standard names to the Qiskit gates/names. + + See :func:`standard_gatename_unitaries`. + + Note that throughout pyGSTi the standard gatenames (e.g., 'Gh' for Hadamard) + are not enforced to correspond to the expected unitaries. So, if the user + as, say, defined 'Gh' to be something other than the Hadamard gate this + conversion dictionary will be incorrect. + + Later versions of qiskit may be able to utilize an accelerated circuit append. + + Returns + ------- + dict[str, tuple[qiskit.circuit.Instruction, str, bool]] + Maps pyGSTi gatenames to tuples of the qiskit gate class, gatename, and a boolean + that indicates whether the append of the gate to the circuit can be accelerated in + a later version of qiskit. + """ + + try: + import qiskit + from qiskit.circuit import Delay + from qiskit.circuit.library import standard_gates + + if qiskit.__version__ != '2.1.1': + _warnings.warn("function 'standard_gatenames_qiskit_conversions()' is designed for qiskit version 2.1.1 \ + and may not function properly for your qiskit version, which is " + qiskit.__version__) + + except ImportError: + _warnings.warn("This operation requires qiskit, which does not appear to be installed.") + + std_gatenames_to_qiskit = {} + + # third entry of tuple indicates whether the gate is a accelerateable standard Qiskit gate. + + std_gatenames_to_qiskit['Gcphase'] = (standard_gates.CZGate, 'cz', False) + std_gatenames_to_qiskit['Gcnot'] = (standard_gates.CXGate, 'cx', False) + std_gatenames_to_qiskit['Gecres'] = (standard_gates.ECRGate, 'ecr', False) + std_gatenames_to_qiskit['Gxpi2'] = (standard_gates.SXGate, 'sx', False) + std_gatenames_to_qiskit['Gzr'] = (standard_gates.RZGate, 'rz', False) + std_gatenames_to_qiskit['Gyr'] = (standard_gates.RYGate, 'ry', False) + std_gatenames_to_qiskit['Gxr'] = (standard_gates.RXGate, 'rx', False) + std_gatenames_to_qiskit['Gxpi'] = (standard_gates.XGate, 'x', False) + std_gatenames_to_qiskit['Gzpi'] = (standard_gates.ZGate, 'z', False) + std_gatenames_to_qiskit['Gypi'] = (standard_gates.YGate, 'y', False) + std_gatenames_to_qiskit['Gi'] = (standard_gates.IGate, 'id', False) + std_gatenames_to_qiskit['Gu3'] = (standard_gates.UGate, 'u', False) + std_gatenames_to_qiskit['Gh'] = (standard_gates.HGate, 'h', False) + + std_gatenames_to_qiskit['Gdelay'] = (Delay, 'delay', False) + + + return std_gatenames_to_qiskit + + def qasm_u3(theta, phi, lamb, output='unitary'): """ The u3 1-qubit gate of QASM, returned as a unitary. diff --git a/pygsti/tools/mcfetools.py b/pygsti/tools/mcfetools.py new file mode 100644 index 000000000..584bfb6ab --- /dev/null +++ b/pygsti/tools/mcfetools.py @@ -0,0 +1,435 @@ +""" +Tools for analyzing MCFE data +""" +#*************************************************************************************************** +# Copyright 2015, 2019, 2025 National Technology & Engineering Solutions of Sandia, LLC (NTESS). +# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights +# in this software. +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +# in compliance with the License. You may obtain a copy of the License at +# http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory. +#*************************************************************************************************** + +from __future__ import annotations +from typing import Optional + +import numpy as _np + +from pygsti.tools.rbtools import hamming_distance + +from pygsti.data.dataset import _DataSetRow +from pygsti.circuits import Circuit as _Circuit + +def success_probability_to_polarization(s: float, n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Maps a success probablity `s` for an n-qubit circuit to + the polarization `s`, defined by `p = (s - 1/2^n)/(1 - 1/2^n)`. + For large n, the difference between `p` and `s` is negligible + and the calculation of 2**n is prohibitive, so we impose + a cutoff above which we assert `p = s`. + + Parameters + ------------ + s : float + success probability for the circuit. + + n : int + number of qubits in the circuit. + + + Returns + ----------- + float + circuit polarization. + """ + + if n < 20: + return (s - 1 / 2**n) / (1 - 1 / 2**n) + else: + return s + + +def polarization_to_success_probability(p: float, n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Maps a polarization `p` for an n-qubit circuit to + the success probability `s`, defined by `s = p * (1 - 1/2**n) + 1/2**n`. + For large n, the difference between `p` and `s` is negligible + and the calculation of 2**n is prohibitive, so we impose + a cutoff above which we assert `s = p`. + + Parameters + ------------ + p : float + circuit polarization + + n : int + number of qubits in the circuit. + + + Returns + ----------- + float + circuit success probability. + """ + + if n < 20: + return p * (1 - 1 / 2**n) + 1 / 2**n + else: + return p + + +def polarization_to_fidelity(p: float, n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Maps a polarization `p` for an n-qubit circuit to + the process fidelity `f`, defined by `f = 1 - (4**n - 1)*(1 - p)/4**n`. + For large n, the difference between `p` and `f` is negligible + and the calculation of 4**n is prohibitive, so we impose + a cutoff above which we assert `f = p`. + + Parameters + ------------ + p : float + circuit polarization + + n : int + number of qubits in the circuit. + + + Returns + ----------- + float + circuit process fidelity. + """ + + if n < 20: + return 1 - (4**n - 1)*(1 - p)/4**n + else: + return p + + +def fidelity_to_polarization(f: float, n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Maps a process fidelity `f` for an n-qubit circuit to + the polarization `p`, defined by `p = 1 - (4**n)*(1 - f)/(4**n - 1)`. + For large n, the difference between `p` and `f` is negligible + and the calculation of 4**n is prohibitive, so we impose + a cutoff above which we assert `p = f`. + + Parameters + ------------ + p : float + circuit polarization + + n : int + number of qubits in the circuit. + + + Returns + ----------- + float + circuit process fidelity. + """ + + if n < 20: + return 1 - (4**n)*(1 - f)/(4**n - 1) + else: + return f + + +def hamming_distance_counts(dsrow: _DataSetRow, circ: _Circuit, idealout: str) -> _np.ndarray: + """ + Utility function for MCFE VBDataFrame creation. + + Compute the hamming distance counts for the outcomes + of a circuit on `n` qubits. + + Hamming distance is defined as the number of bits + that are different between two bitstrings of the same + length. E.g., the Hamming distance between 0010 and 1110 + is 2. + + Parameters + ------------ + dsrow : pygsti.data.dataset.DataSetRow + Row in dataset containing outcome information + for all mirror circuits. + + circ : pygsti.circuits.Circuit + pyGSTi circuit to which the outcome counts in + `dsrow` pertain. + + idealout : string + length-`n` string that defines the target + bitstring for the circuit. All Hamming distances + are calculated relative to this string. + + + Returns + -------------- + np.ndarray[float] + Array whose indices correspond to Hamming distances + from `idealout`, and whose values at each index + correspond to the number of counts in `dsrow` that + have that Hamming distance. + """ + + nQ = len(circ.line_labels) # number of qubits + assert nQ == len(idealout[-1]), f'{nQ} != {len(idealout[-1])}' + hamming_distance_counts = _np.zeros(nQ + 1, float) + if dsrow.total > 0: + for outcome_lbl, counts in dsrow.counts.items(): + outbitstring = outcome_lbl[-1] + hamming_distance_counts[hamming_distance(outbitstring, idealout[-1])] += counts + return hamming_distance_counts + + +def adjusted_success_probability(hamming_distance_counts: _np.ndarray) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Compute the adjusted success probability `adjSP` of a circuit + from its Hamming distance counts, according to the formula + + adjSP = \sum_{k=0}^n (-1/2)^k * h_k, + + where h_k is the probability of Hamming distance k. + + Parameters + ------------ + hamming_distance_counts : np.ndarray[float] + Array whose indices correspond to Hamming distances, + and whose values at each index correspond to the number of + counts of that Hamming distance in the circuit outcome + data from which `hamming_distance_counts` was derived. + + Returns + ----------- + float + adjusted success probability. + """ + + if _np.sum(hamming_distance_counts) == 0.: + return 0. + else: + hamming_distance_pdf = _np.array(hamming_distance_counts) / _np.sum(hamming_distance_counts) + adjSP = _np.sum([(-1 / 2)**n * hamming_distance_pdf[n] for n in range(len(hamming_distance_pdf))]) + return adjSP + + +def effective_polarization(hamming_distance_counts: _np.ndarray) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Compute the effective polarization `p` of a circuit + from its Hamming distance counts, according to the formula + + `p = (4**n * asp - 1)/(4**n - 1)` + + where `asp` is the adjusted success probability of + the circuit. + + Parameters + ------------ + hamming_distance_counts : np.ndarray[float] + Array whose indices correspond to Hamming distances, + and whose values at each index correspond to the number of + counts of that Hamming distance in the circuit outcome + data from which `hamming_distance_counts` was derived. + + Returns + ----------- + float + polarization. + """ + + n = len(hamming_distance_counts) - 1 + asp = adjusted_success_probability(hamming_distance_counts) + + if n < 20: + return (4**n * asp - 1)/(4**n - 1) + else: + return asp + + +def rc_predicted_process_fidelity(bare_rc_effective_pols: _np.ndarray, + rc_rc_effective_pols: _np.ndarray, + reference_effective_pols: _np.ndarray, + n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Compute the process fidelity `f` for a circuit on `n` qubits according to the formula + + `f = 1 - (4**n - 1)/4**n * (1 - p)`, + + Where `p` is the effective polarization of the circuit, given by + + `p = E[p(M_1)] / sqrt( E[p(M_2)] * E[p(M_3)] )`. + + Here, M_i refers to the i-th family of mirror circuits used in mirror circuit fidelity + estimation with randomized compiling. See https://arxiv.org/pdf/2204.07568 for more + information. + + The process fidelity estimate is clamped to the [0.0, 1.0] range. + + Parameters + ------------ + bare_rc_effective_pols : np.ndarray[float] + list of effective polarizations for bare_rc (M_1) mirror circuits. + + rc_rc_effective_pols : np.ndarray[float] + list of effective polarizations for rc_rc (M_2) mirror circuits. + + reference_effective_pols : np.ndarray[float] + list of effective polarizations for SPAM reference (M_3) mirror circuits. + + n : int + number of qubits in the quantum circuit. + + + Returns: + ---------- + float + process fidelity estimate. + """ + + a = _np.mean(bare_rc_effective_pols) + b = _np.mean(rc_rc_effective_pols) + c = _np.mean(reference_effective_pols) + + # print(a) + + # print(b) + + # print(c) + + if c <= 0.: + return _np.nan # raise ValueError("Reference effective polarization zero or smaller! Cannot estimate the process fidelity") + elif b <= 0: + return 0. + else: + pfid = polarization_to_fidelity(a / _np.sqrt(b * c), n) + if pfid < 0.0: + return 0.0 + elif pfid > 1.0: + return 1.0 + else: + return pfid + + +def predicted_process_fidelity_for_central_pauli_mcs(central_pauli_effective_pols: _np.ndarray, + reference_effective_pols: _np.ndarray, + n: int) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Compute the process fidelity `f` for a circuit on `n` qubits according to the formula + + `f = 1 - (4**n - 1)/4**n * (1 - p)`, + + Where `p` is the effective polarization of the circuit, given by + + `p = sqrt( E[p(M_1)] / E[p(M_2)]) + + Here, M_1 refers to central Pauli quasi-mirror circuits and + M_2 refers to SPAM reference circuits. + + Parameters + ------------ + central_effective_pols : np.ndarray[float] + list of effective polarizations for central Pauli (M_1) mirror circuits. + + + reference_effective_pols : np.ndarray[float] + list of effective polarizations for SPAM reference (M_2) mirror circuits. + + n : int + number of qubits in the quantum circuit. + + + Returns: + ---------- + float + process fidelity estimate. + """ + + a = _np.mean(central_pauli_effective_pols) + c = _np.mean(reference_effective_pols) + if c <= 0.: + return _np.nan # raise ValueError("Reference effective polarization zero or smaller! Cannot estimate the process fidelity") + elif a <= 0: + return 0. + else: + return polarization_to_fidelity(_np.sqrt(a / c), n) + + +def rc_bootstrap_predicted_pfid(brs: _np.ndarray, + rrs: _np.ndarray, + refs: _np.ndarray, + n: int, + num_bootstraps: int = 50, + rand_state: Optional[_np.random.RandomState] = None + ) -> float: + """ + Utility function for MCFE VBDataFrame creation. + + Compute bootstrapped error bars for the MCFE random compilation process fidelity + estimate. + + Parameters + ------------ + brs : np.ndarray[float] + list of effective polarizations for bare_rc (M_1) mirror circuits. + + rrs : np.ndarray[float] + list of effective polarizations for rc_rc (M_2) mirror circuits. + + refs : np.ndarray[float] + list of effective polarizations for SPAM reference (M_3) mirror circuits. + + n : int + number of qubits in the quantum circuit. + + num_bootstraps : int + number of samples to take from the bootstrapped distribution. + + rand_state : optional, np.random.RandomState + random state to be used for bootstrapepd distribution sampling. + + + Returns: + ---------- + float + estimate for standard deviation of process fidelity estimate. + """ + + if rand_state is None: + rand_state = _np.random.RandomState() + + # print(num_bootstraps) + + pfid_samples = [] + for _ in range(num_bootstraps): + br_sample = rand_state.choice(brs, len(brs), replace=True) + rr_sample = rand_state.choice(rrs, len(rrs), replace=True) + ref_sample = rand_state.choice(refs, len(refs), replace=True) + + pfid = rc_predicted_process_fidelity( + br_sample, + rr_sample, + ref_sample, + n) + + pfid_samples.append(pfid) + + bootstrapped_stdev = _np.std(pfid_samples) + + return bootstrapped_stdev \ No newline at end of file diff --git a/test/performance/mpi_2D_scaling/mpi_timings.py b/test/performance/mpi_2D_scaling/mpi_timings.py deleted file mode 100644 index b4b966c52..000000000 --- a/test/performance/mpi_2D_scaling/mpi_timings.py +++ /dev/null @@ -1,67 +0,0 @@ -#!/usr/bin/env python - -import cProfile -import os -from pathlib import Path -import pickle -import time - -from mpi4py import MPI - -import pygsti -from pygsti.modelpacks import smq2Q_XYICNOT as std - -comm = MPI.COMM_WORLD -resource_alloc = pygsti.baseobjs.ResourceAllocation(comm) -mdl = std.target_model() - -exp_design = std.get_gst_experiment_design(64) - -mdl_datagen = mdl.depolarize(op_noise=0.01, spam_noise=0.01) - -# First time running through, generate reference dataset -if not Path('reference_ds.pkl').exists(): - if comm.rank == 0: - ds = pygsti.data.simulate_data(mdl_datagen, exp_design, 1000, seed=1234, comm=resource_alloc.comm) - pickle.dump(ds, open('reference_ds.pkl','wb')) - else: - time.sleep(2) - -ds_ref = pickle.load(open('reference_ds.pkl','rb')) -ds = ds_ref - -MINCLIP = 1e-4 -chi2_builder = pygsti.objectivefns.Chi2Function.builder( - 'chi2', regularization={'min_prob_clip_for_weighting': MINCLIP}, penalties={'cptp_penalty_factor': 0.0}) -mle_builder = pygsti.objectivefns.PoissonPicDeltaLogLFunction.builder( - 'logl', regularization={'min_prob_clip': MINCLIP, 'radius': MINCLIP}) -iteration_builders = [chi2_builder]; final_builders = [mle_builder] -builders = pygsti.protocols.GSTObjFnBuilders(iteration_builders, final_builders) - -tol = 1e-6 - -opt = None # default - -#GST TEST -data = pygsti.protocols.ProtocolData(exp_design, ds) -#mdl.sim = pygsti.baseobjs.MatrixForwardSimulator(num_atoms=1) -mdl.sim = pygsti.forwardsims.MapForwardSimulator(num_atoms=1, max_cache_size=0) -gst = pygsti.protocols.GateSetTomography(mdl, - objfn_builders=builders, optimizer=opt, verbosity=4) - -profiler = cProfile.Profile() -profiler.enable() - -results = gst.run(data, comm=comm) - -profiler.disable() - -num_procs = comm.Get_size() -num_procs_host = os.environ.get('PYGSTI_MAX_HOST_PROCS', num_procs) -os.makedirs(f'{num_procs}_{num_procs_host}.profile', exist_ok=True) -profiler.dump_stats(f'{num_procs}_{num_procs_host}.profile/{comm.rank}.prof') - -results = None # Needed to cause shared mem to be freed by garbage collection *before* python shuts down shared mem "system" - -comm.barrier() -if comm.rank == 0: print("DONE") diff --git a/test/unit/extras/ibmq/test_ibmqexperiment.py b/test/unit/extras/ibmq/test_ibmqexperiment.py index 4e4c25ed4..bff25d2ef 100644 --- a/test/unit/extras/ibmq/test_ibmqexperiment.py +++ b/test/unit/extras/ibmq/test_ibmqexperiment.py @@ -1,9 +1,9 @@ +from ...util import BaseCase + import pygsti from pygsti.extras.devices.experimentaldevice import ExperimentalDevice from pygsti.extras import ibmq from pygsti.processors import CliffordCompilationRules as CCR -from qiskit.providers.fake_provider import GenericBackendV2 -from qiskit_ibm_runtime import QiskitRuntimeService from pygsti.protocols import MirrorRBDesign as RMCDesign from pygsti.protocols import PeriodicMirrorCircuitDesign as PMCDesign from pygsti.protocols import ByDepthSummaryStatistics @@ -11,9 +11,26 @@ from pygsti.protocols import StandardGSTDesign import numpy as np +try: + from qiskit.providers.fake_provider import GenericBackendV2 +except: + GenericBackendV2 = None + +try: + from qiskit_ibm_runtime import QiskitRuntimeService +except: + QiskitRuntimeService = None + +import pytest + class IBMQExperimentTester(): @classmethod def setup_class(cls): + if GenericBackendV2 is None: + pytest.skip('Qiskit is required for this operation, and does not appear to be installed.') + elif QiskitRuntimeService is None: + pytest.skip('Qiskit Runtime is required for this operation, and does not appear to be installed.') + cls.backend = GenericBackendV2(num_qubits=4) cls.device = ExperimentalDevice.from_qiskit_backend(cls.backend) cls.pspec = cls.device.create_processor_spec(['Gc{}'.format(i) for i in range(24)] + ['Gcphase']) diff --git a/test/unit/objects/test_circuit.py b/test/unit/objects/test_circuit.py index 57143152e..cfd7dc3af 100644 --- a/test/unit/objects/test_circuit.py +++ b/test/unit/objects/test_circuit.py @@ -602,6 +602,64 @@ def test_from_cirq(self): self.assertEqual(ckt_global_idle_custom, converted_pygsti_circuit_global_idle_custom) self.assertEqual(ckt_global_idle_custom, converted_pygsti_circuit_global_idle_custom_1) self.assertEqual(ckt_global_idle_none, converted_pygsti_circuit_global_idle_none) + + + def test_from_qiskit(self): + try: + import qiskit + except ImportError: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + qc = qiskit.QuantumCircuit(4) + + qc.x(0) + qc.cx(0,1) + qc.h(3) + qc.cz(2,3) + qc.rz(0.5, 2) + qc.rz(0.5, 3) + qc.rz(0.5, 0) + qc.barrier(0,1) + qc.z(3) + qc.cz(0,2) + + expected_ps_circ = circuit.Circuit([Label([Label('Gxpi', 'Q0'), Label('Gh', 'Q3')]), + Label([Label('Gcnot', ('Q0', 'Q1')), Label('Gcphase', ('Q2','Q3'))]), + Label([Label('Gzr', 'Q2', args=[0.5]), Label('Gzr', 'Q3', args=[0.5]), Label('Gzr', 'Q0', args=[0.5])]), + Label([Label('Gzpi', 'Q3'), Label('Gcphase', ('Q0', 'Q2'))]) + ], line_labels=('Q0', 'Q1', 'Q2', 'Q3')) + + expected_mapping = {i: f'Q{i}' for i in range(4)} + + observed_ps_circ, observed_mapping = circuit.Circuit.from_qiskit(qc) + + self.assertEqual(expected_ps_circ, observed_ps_circ) + + self.assertEqual(expected_mapping, observed_mapping) + + + def test_to_from_qiskit_round_trip(self): + try: + import qiskit + except ImportError: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + + expected_ps_circ = circuit.Circuit([Label([Label('Gxpi', 'Q0'), Label('Gh', 'Q3')]), + Label([Label('Gcnot', ('Q0', 'Q1')), Label('Gcphase', ('Q2','Q3'))]), + Label([Label('Gzr', 'Q2', args=[0.5]), Label('Gzr', 'Q3', args=[0.5]), Label('Gzr', 'Q0', args=[0.5])]), + Label([Label('Gzpi', 'Q3'), Label('Gcphase', ('Q0', 'Q2'))]) + ], line_labels=('Q0', 'Q1', 'Q2', 'Q3')) + + + qk_circ = expected_ps_circ.convert_to_qiskit(qubit_conversion='remove-Q') + observed_ps_circ, observed_mapping = circuit.Circuit.from_qiskit(qk_circ) + + expected_mapping = {i: f'Q{i}' for i in range(4)} + + self.assertEqual(expected_ps_circ, observed_ps_circ) + self.assertEqual(expected_mapping, observed_mapping) + def test_convert_to_stim_tableau(self): #TODO: Add correctness checks for generated Tableaus. diff --git a/test/unit/objects/test_subcirc_selection.py b/test/unit/objects/test_subcirc_selection.py new file mode 100644 index 000000000..158f0107a --- /dev/null +++ b/test/unit/objects/test_subcirc_selection.py @@ -0,0 +1,304 @@ +from ..util import BaseCase + +from pygsti.circuits import subcircuit_selection as _subcircsel, Circuit as C +from pygsti.baseobjs import Label as L + +import numpy as _np +import networkx as _nx + +from collections import defaultdict + +class TestSubcircuitSelection(BaseCase): + def test_random_subgraph_creation(self): + + rand_state = _np.random.RandomState(0) + + n = 100 + g = _nx.gnp_random_graph(n, 2 * _np.log(n)/n, seed=0) + + widths = list(range(1,21)) + + subgraphs = [_subcircsel.random_connected_subgraph(g, width, rand_state) for width in widths] + + + self.assertTrue(all(len(subgraphs[i]) == widths[i] for i in range(len(widths)))) + self.assertTrue(all(_nx.is_connected(g.subgraph(subgraph)) for subgraph in subgraphs)) + + + def test_random_subgraph_reproducibility(self): + + rand_state = _np.random.RandomState(0) + + n = 100 + g = _nx.gnp_random_graph(n, 2 * _np.log(n)/n, seed=0) + + widths = list(range(1,21)) + + subgraphs_1 = [_subcircsel.random_connected_subgraph(g, width, rand_state) for width in widths] + + rand_state = _np.random.RandomState(0) + subgraphs_2 = [_subcircsel.random_connected_subgraph(g, width, rand_state) for width in widths] + + self.assertTrue(len(widths) == len(subgraphs_1)) + self.assertTrue(len(widths) == len(subgraphs_2)) + + self.assertTrue(all(subgraphs_1[i] == subgraphs_2[i] for i in range(len(widths)))) + + + def test_simple_subcirc_selection_linear_connectivity(self): + + rand_state = _np.random.RandomState(0) + + class NoDelayInstructions(object): + def get(self, *args): + return 0.0 + + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [_np.pi,0,_np.pi] + Y_args = [_np.pi,_np.pi/2,_np.pi/2] + Z_args = [0,0,_np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + + width_depths = {2: [2,4,6], + 3: [3,6,9]} + + num_subcircs_per_width_depth = 5 + + + design = _subcircsel.sample_subcircuits(circ, width_depths=width_depths, + instruction_durations=NoDelayInstructions(), + coupling_map='linear', + num_samples_per_width_depth=num_subcircs_per_width_depth, + rand_state=rand_state + ) + + # check that the width and depth metadata agree with the circuit width and depth for each circuit in the design. + for subcirc, auxlist in design.aux_info.items(): + self.assertTrue(all(subcirc.width == aux['width'] for aux in auxlist)) + self.assertTrue(all(subcirc.depth == aux['depth'] for aux in auxlist)) + + # check that the correct number of circuits of each width and depth exist by looking at aux_info + reshaped_width_depths = [] + for width, depths in width_depths.items(): + for depth in depths: + reshaped_width_depths.append((width, depth)) + + created_width_depths = defaultdict(int) + for subcirc, auxlist in design.aux_info.items(): + created_width_depths[(subcirc.width, subcirc.depth)] += len(auxlist) + + self.assertListEqual(sorted(reshaped_width_depths), sorted(list(created_width_depths.keys()))) + self.assertTrue(all(created_width_depths[k] == num_subcircs_per_width_depth for k in reshaped_width_depths)) + + + + def test_simple_subcirc_selection_full_connectivity(self): + + rand_state = _np.random.RandomState(0) + + class NoDelayInstructions(object): + def get(self, *args): + return 0.0 + + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [_np.pi,0,_np.pi] + Y_args = [_np.pi,_np.pi/2,_np.pi/2] + Z_args = [0,0,_np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + + width_depths = {2: [2,4,6], + 3: [3,6,9]} + + num_subcircs_per_width_depth = 20 + + + design = _subcircsel.sample_subcircuits(circ, width_depths=width_depths, + instruction_durations=NoDelayInstructions(), + coupling_map='all-to-all', + num_samples_per_width_depth=num_subcircs_per_width_depth, + rand_state=rand_state + ) + + + # check that the width and depth metadata agree with the circuit width and depth for each circuit in the design. + for subcirc, auxlist in design.aux_info.items(): + self.assertTrue(all(subcirc.width == aux['width'] for aux in auxlist)) + self.assertTrue(all(subcirc.depth == aux['depth'] for aux in auxlist)) + + # check that the correct number of circuits of each width and depth exist by looking at aux_info + reshaped_width_depths = [] + for width, depths in width_depths.items(): + for depth in depths: + reshaped_width_depths.append((width, depth)) + + created_width_depths = defaultdict(int) + for subcirc, auxlist in design.aux_info.items(): + created_width_depths[(subcirc.width, subcirc.depth)] += len(auxlist) + + self.assertListEqual(sorted(reshaped_width_depths), sorted(list(created_width_depths.keys()))) + self.assertTrue(all(created_width_depths[k] == num_subcircs_per_width_depth for k in reshaped_width_depths)) + + + def test_simple_subcirc_selection_full_connectivity(self): + + rand_state = _np.random.RandomState(0) + + class CustomInstructions(object): + def get(self, *args): + if args[0] == 'Gu3': + return 5.0 + elif args[0] == 'Gcnot': + return 12.0 + elif args[0] == 'Gcphase': + return 13.0 + else: + raise RuntimeError(f'No duration known for instruction {args[0]}') + + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [_np.pi,0,_np.pi] + Y_args = [_np.pi,_np.pi/2,_np.pi/2] + Z_args = [0,0,_np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + + width_depths = {2: [2,4,6], + 4: [3,6,9]} + + num_subcircs_per_width_depth = 31 + + graph = _nx.barbell_graph(10, 3) + + int_coupling_list = list(graph.edges()) + coupling_list = [(f'Q{i}', f'Q{j}') for i,j in int_coupling_list] + + + design = _subcircsel.sample_subcircuits(circ, width_depths=width_depths, + instruction_durations=CustomInstructions(), + coupling_map=coupling_list, + num_samples_per_width_depth=num_subcircs_per_width_depth, + rand_state=rand_state + ) + + + # check that the width and depth metadata agree with the circuit width and depth for each circuit in the design. + for subcirc, auxlist in design.aux_info.items(): + self.assertTrue(all(subcirc.width == aux['width'] for aux in auxlist)) + self.assertTrue(all(subcirc.depth == aux['depth'] for aux in auxlist)) + + # check that the correct number of circuits of each width and depth exist by looking at aux_info + reshaped_width_depths = [] + for width, depths in width_depths.items(): + for depth in depths: + reshaped_width_depths.append((width, depth)) + + created_width_depths = defaultdict(int) + for subcirc, auxlist in design.aux_info.items(): + created_width_depths[(subcirc.width, subcirc.depth)] += len(auxlist) + + self.assertListEqual(sorted(reshaped_width_depths), sorted(list(created_width_depths.keys()))) + self.assertTrue(all(created_width_depths[k] == num_subcircs_per_width_depth for k in reshaped_width_depths)) + + + def test_simple_subcirc_selection_qiskit_coupling_map(self): + try: + import qiskit + except: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + try: + import qiskit_ibm_runtime + except: + self.skipTest('Qiskit Runtime is required for this operation, and does not appear to be installed.') + + + backend = qiskit_ibm_runtime.fake_provider.FakeFez() + + qk_circ = qiskit.circuit.library.QFT(4) + qk_circ = qiskit.transpile(qk_circ, backend=backend) + + ps_circ, _ = C.from_qiskit(qk_circ) + ps_circ = ps_circ.delete_idling_lines() + + rand_state = _np.random.RandomState(0) + + width_depths = {2: [2,4,6], + 3: [3,6,9]} + + num_subcircs_per_width_depth = 15 + + design = _subcircsel.sample_subcircuits(ps_circ, width_depths=width_depths, + instruction_durations=backend.instruction_durations, + coupling_map=backend.coupling_map, + num_samples_per_width_depth=num_subcircs_per_width_depth, + rand_state=rand_state + ) + + + # check that the width and depth metadata agree with the circuit width and depth for each circuit in the design. + for subcirc, auxlist in design.aux_info.items(): + self.assertTrue(all(subcirc.width == aux['width'] for aux in auxlist)) + self.assertTrue(all(subcirc.depth == aux['depth'] for aux in auxlist)) + + # check that the correct number of circuits of each width and depth exist by looking at aux_info + reshaped_width_depths = [] + for width, depths in width_depths.items(): + for depth in depths: + reshaped_width_depths.append((width, depth)) + + created_width_depths = defaultdict(int) + for subcirc, auxlist in design.aux_info.items(): + created_width_depths[(subcirc.width, subcirc.depth)] += len(auxlist) + + + self.assertListEqual(sorted(reshaped_width_depths), sorted(list(created_width_depths.keys()))) + self.assertTrue(all(created_width_depths[k] == num_subcircs_per_width_depth for k in reshaped_width_depths)) + + + diff --git a/test/unit/processors/__init__.py b/test/unit/processors/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/test/unit/processors/test_cp.py b/test/unit/processors/test_cp.py new file mode 100644 index 000000000..6c3cb5458 --- /dev/null +++ b/test/unit/processors/test_cp.py @@ -0,0 +1,630 @@ +from ..util import BaseCase + +import numpy as np + +from pygsti.baseobjs.label import Label as L +from pygsti.circuits.circuit import Circuit as C +from pygsti.processors.random_compilation import RandomCompilation +from pygsti.tools.internalgates import Gu3, standard_gatename_unitaries as _standard_gatename_unitaries + + +def get_clifford_from_unitary(U): + clifford_unitaries = {k: v for k, v in _standard_gatename_unitaries().items() + if 'Gc' in k and v.shape == (2, 2)} + for k,v in clifford_unitaries.items(): + for phase in [1, -1, 1j, -1j]: + if np.allclose(U, phase*v): + return k + + raise RuntimeError(f'Failed to look up Clifford for unitary:\n{U}') + +class TestCentralPauli(BaseCase): + + def test_u3_x_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([0,0,0,0,2,2,2,2,]) # all X errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # This is passing but with (global) phase differences + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc3','Gc3','Gc3','Gc3'],['Gc0', 'Gc3', 'Gc6', 'Gc9']]) + + def test_u3_z_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,2,2,0,0,0,0]) # all Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc9','Gc9','Gc9','Gc9'],['Gc0', 'Gc3', 'Gc6', 'Gc9']]) + + def test_u3_x_and_z_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,2,2,2,2,2,2]) # all X and Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc6','Gc6','Gc6','Gc6'],['Gc0', 'Gc3', 'Gc6', 'Gc9']]) + + def test_u3_x_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([0,2]) # all X errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc3'], ['Gc0'], ['Gc3'], ['Gc6'], ['Gc9']]) + + def test_u3_z_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,0]) # all Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc9'], ['Gc0'], ['Gc3'], ['Gc6'], ['Gc9']]) + + def test_u3_x_and_z_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2]) # all X and Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layer) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc6'], ['Gc0'], ['Gc3'], ['Gc6'], ['Gc9']]) + + # def test_u3_assorted(self): + # #TODO: implement + + #CNOT + def test_cnot_x_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([0,0,2,2]) # all X errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc3', 'Gc3'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,0,2,0]))) # X I is propagated through + + def test_cnot_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,0,0]) # all Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc9', 'Gc9'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,2,0,0]))) # Z I is propagated through + + + def test_cnot_x_and_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,2,2]) # all X and Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc6', 'Gc6'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,2,2,0]))) # X Z is propagated through + + + def test_cphase_x_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([0,0,2,2]) # all X errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc3', 'Gc3'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,2,2,2]))) # XZ ZX is propagated through + + def test_cphase_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,0,0]) # all Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc9', 'Gc9'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,2,0,0]))) # Z Z is propagated through + + + def test_cphase_x_and_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,2,2,2]) # all X and Z errors + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc6', 'Gc6'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,0,2,2]))) # X Z is propagated through + + + def test_big_circuit(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [np.pi,0,np.pi] + Y_args = [np.pi,np.pi/2,np.pi/2] + Z_args = [0,0,np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + + test_layer = np.array([2,0,0,2,2,2,2,2,2,0]) #YXXYZ + + rc = RandomCompilation(rc_strategy='central_pauli', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layer) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase' or gate.name == 'Gcnot': + pauli_layer.append(gate.name) + else: + print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + correct_output_circ = [['Gc6', 'Gc3', 'Gc3', 'Gc6', 'Gc9'], + ['Gc9', 'Gc6', 'Gc3', 'Gc0', 'Gc0'], + ['Gcnot', 'Gcnot'], + ['Gc6', 'Gc9', 'Gc3', 'Gc0', 'Gc3'], + ['Gcnot', 'Gcnot'], + ['Gcphase', 'Gcphase'], + ['Gcnot', 'Gcnot'], + ['Gc3', 'Gc6', 'Gc0', 'Gc0', 'Gc6'], + ['Gc3', 'Gc6', 'Gc0', 'Gc6', 'Gc6'], + ['Gc6', 'Gc3', 'Gc9', 'Gc6', 'Gc0'], + ['Gcphase', 'Gcphase']] + + self.assertTrue(paulis == correct_output_circ) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,0,0,2,0,2,0,0,2,0]))) \ No newline at end of file diff --git a/test/unit/processors/test_rc.py b/test/unit/processors/test_rc.py new file mode 100644 index 000000000..8658e2b3c --- /dev/null +++ b/test/unit/processors/test_rc.py @@ -0,0 +1,644 @@ +from ..util import BaseCase + +import numpy as np + +from pygsti.baseobjs.label import Label as L +from pygsti.circuits.circuit import Circuit as C +from pygsti.processors.random_compilation import RandomCompilation +from pygsti.tools.internalgates import Gu3, standard_gatename_unitaries as _standard_gatename_unitaries + + + + +def get_clifford_from_unitary(U): + clifford_unitaries = {k: v for k, v in _standard_gatename_unitaries().items() + if 'Gc' in k and v.shape == (2, 2)} + for k,v in clifford_unitaries.items(): + for phase in [1, -1, 1j, -1j]: + if np.allclose(U, phase*v): + return k + + raise RuntimeError(f'Failed to look up Clifford for unitary:\n{U}') + +class TestRandomCompilation(BaseCase): + + def test_u3_x_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([0,0,0,0,2,2,2,2,])] # all X errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + # output should be X, I, Z, Y + + # print(rc_circ) + + # This is passing but with (global) phase differences + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc3', 'Gc0', 'Gc9', 'Gc6']]) + + def test_u3_z_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,2,2,0,0,0,0])] # all Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc9', 'Gc6', 'Gc3', 'Gc0']]) + + def test_u3_x_and_z_errors_parallel(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), #I + L('Gu3', ['Q2'], args=[np.pi,0,np.pi]), #X + L('Gu3', ['Q3'], args=[np.pi,np.pi/2,np.pi/2]), #Y + L('Gu3', ['Q4'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,2,2,2,2,2,2])] # all X and Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc6', 'Gc9', 'Gc0', 'Gc3']]) + + def test_u3_x_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([0,2]), np.array([0,2]), np.array([0,2]), np.array([0,2])] # all X errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc3'], ['Gc3'], ['Gc6'], ['Gc9']]) + + def test_u3_z_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,0]), np.array([2,0]), np.array([2,0]), np.array([2,0])] # all Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc9'], ['Gc3'], ['Gc6'], ['Gc9']]) + + def test_u3_x_and_z_errors_serial(self): + line_labels = ['Q1'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0])], #I + [L('Gu3', ['Q1'], args=[np.pi,0,np.pi])], #X + [L('Gu3', ['Q1'], args=[np.pi,np.pi/2,np.pi/2])], #Y + [L('Gu3', ['Q1'], args=[0,0,np.pi])]] #Z + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2]), np.array([2,2]), np.array([2,2]), np.array([2,2])] # all X and Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, _] = rc.compile(circ, test_layers) + + # print(rc_circ) + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + # print(paulis) + + self.assertTrue(paulis == [['Gc6'], ['Gc3'], ['Gc6'], ['Gc9']]) + + # def test_u3_assorted(self): + # #TODO: implement + + #CNOT + def test_cnot_x_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([0,0,2,2])] # all X errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc3', 'Gc3'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,0,2,0]))) # X I is propagated through + + def test_cnot_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,0,0])] # all Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc9', 'Gc9'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,2,0,0]))) # Z I is propagated through + + + def test_cnot_x_and_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcnot', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,2,2])] # all X and Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcnot': + pauli_layer.append('Gcnot') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc6', 'Gc6'], ['Gcnot']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,2,2,0]))) # X Z is propagated through + + + def test_cphase_x_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([0,0,2,2])] # all X errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc3', 'Gc3'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,2,2,2]))) # XZ ZX is propagated through + + def test_cphase_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,0,0])] # all Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + # print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + # print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc9', 'Gc9'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,2,0,0]))) # Z Z is propagated through + + + def test_cphase_x_and_z_errors(self): + + line_labels = ['Q1', 'Q2'] + layers = [[L('Gu3', ['Q1'], args=[0,0,0]), L('Gu3', ['Q2'], args=[0,0,0])], + [L('Gcphase', ['Q1', 'Q2'], args=None)]] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,2,2,2])] # all X and Z errors + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase': + pauli_layer.append('Gcphase') + continue + print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + self.assertTrue(paulis == [['Gc6', 'Gc6'], ['Gcphase']]) + self.assertTrue(np.array_equal(pauli_frame, np.array([0,0,2,2]))) # X Z is propagated through + + + def test_big_circuit(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [np.pi,0,np.pi] + Y_args = [np.pi,np.pi/2,np.pi/2] + Z_args = [0,0,np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + + test_layers = [np.array([2,0,0,2,2,2,2,2,2,0]), #YXXYZ + np.array([0,0,2,0,2,2,0,0,2,0]), #XIZXZ + np.array([0,0,2,2,2,2,2,2,2,2]), #XXYYY + np.array([2,0,2,2,2,2,2,2,0,2]), #YXYZY + np.array([2,2,0,2,2,0,0,2,0,0]) #ZZXZZ + ] + + rc = RandomCompilation(rc_strategy='pauli_rc', testing=True) + + [rc_circ, _, pauli_frame] = rc.compile(circ, test_layers) + + # print(rc_circ) + + # Output should be + + paulis = [] + + U = Gu3() + + # Gc0 is I + # Gc3 is X + # Gc6 is Y + # Gc9 is Z + + for i in range(len(rc_circ)): + layer = rc_circ.layer_label(i).components + print(layer) + pauli_layer = [] + for gate in layer: + if gate.name == 'Gcphase' or gate.name == 'Gcnot': + pauli_layer.append(gate.name) + else: + print(gate) + unitary = U(gate.args) + clifford = get_clifford_from_unitary(unitary) + pauli_layer.append(clifford) + paulis.append(pauli_layer) + + print(paulis) + + print(pauli_frame) + + correct_output_circ = [['Gc3', 'Gc9', 'Gc0', 'Gc6', 'Gc9'], + ['Gcnot', 'Gcnot'], + ['Gc3', 'Gc9', 'Gc6', 'Gc9', 'Gc3'], + ['Gcnot', 'Gcnot'], + ['Gcphase', 'Gcphase'], + ['Gcnot', 'Gcnot'], + ['Gc3', 'Gc9', 'Gc6', 'Gc9', 'Gc9'], + ['Gc6', 'Gc6', 'Gc0', 'Gc9', 'Gc6'], + ['Gc9', 'Gc9', 'Gc0', 'Gc6', 'Gc3'], + ['Gcphase', 'Gcphase']] + + self.assertTrue(paulis == correct_output_circ) + self.assertTrue(np.array_equal(pauli_frame, np.array([2,2,0,2,2,0,0,2,0,0]))) \ No newline at end of file diff --git a/test/unit/protocols/test_mirror_edesign.py b/test/unit/protocols/test_mirror_edesign.py new file mode 100644 index 000000000..7d0eedfce --- /dev/null +++ b/test/unit/protocols/test_mirror_edesign.py @@ -0,0 +1,334 @@ +from ..util import BaseCase + +import numpy as np + +from pygsti.circuits import Circuit as C +from pygsti.baseobjs import Label as L +from pygsti.protocols import FreeformDesign +from pygsti.protocols.mirror_edesign import * +from pygsti.models import modelconstruction as mc +from pygsti.data import simulate_data +from pygsti.processors import QubitProcessorSpec + +class TestMirrorEDesign(BaseCase): + def test_circuit_count_and_target_outcome_rc(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [np.pi,0,np.pi] + Y_args = [np.pi,np.pi/2,np.pi/2] + Z_args = [0,0,np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + circ_dict = {circ: [{'id': 0, 'width': circ.width, 'depth': circ.depth}]} + test_edesign = FreeformDesign(circ_dict) + + num_mcs_per_circ = 2 + num_ref_per_qubit_subset = 2 + + mirror_edesign = make_mirror_edesign(test_edesign=test_edesign, account_for_routing=False, + ref_edesign=None, + num_mcs_per_circ=num_mcs_per_circ, + num_ref_per_qubit_subset=num_ref_per_qubit_subset, + mirroring_strategy='pauli_rc', + gate_set='u3_cx_cz' + ) + + self.assertEqual(2*num_mcs_per_circ + num_ref_per_qubit_subset, len(mirror_edesign.all_circuits_needing_data)) + + flat_test_auxlist = [aux for auxlist in test_edesign.aux_info.values() for aux in auxlist] + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['rr'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['br'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + pspec = QubitProcessorSpec(5, ['Gu3', 'Gcphase', 'Gcnot', 'Gi'], + availability={'Gcnot':'all-permutations', 'Gcphase': 'all-permutations'}, + qubit_labels=['Q1', 'Q2', 'Q3', 'Q4', 'Q5']) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + for edkey, ed in mirror_edesign.items(): + for circ, auxlist in ed.aux_info.items(): + target_bs = auxlist[0]['idealout'] + result = ideal_model.probabilities(circ) + self.assertAlmostEqual(result[(target_bs,)], 1.0) + + + def test_circuit_count_and_target_outcome_cp(self): + line_labels = ['Q1', 'Q2', 'Q3', 'Q4', 'Q5'] + + I_args = [0,0,0] + X_args = [np.pi,0,np.pi] + Y_args = [np.pi,np.pi/2,np.pi/2] + Z_args = [0,0,np.pi] + + layers = [ + [L('Gu3', ['Q1'], args=Z_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcnot', ['Q1', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q3'], args=None)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=Z_args), L('Gu3', ['Q3'], args=X_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=X_args)], + [L('Gcnot', ['Q2', 'Q3'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q3', 'Q4'], args=None)], + [L('Gcnot', ['Q3', 'Q2'], args=None), L('Gcnot', ['Q4', 'Q5'], args=None)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=I_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=X_args), L('Gu3', ['Q2'], args=Y_args), L('Gu3', ['Q3'], args=I_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=Y_args)], + [L('Gu3', ['Q1'], args=Y_args), L('Gu3', ['Q2'], args=X_args), L('Gu3', ['Q3'], args=Z_args), L('Gu3', ['Q4'], args=Y_args), L('Gu3', ['Q5'], args=I_args)], + [L('Gcphase', ['Q1', 'Q2'], args=None), L('Gcphase', ['Q4', 'Q5'], args=None)] + ] + + circ = C(layers, line_labels=line_labels) + circ_dict = {circ: [{'id': 0, 'width': circ.width, 'depth': circ.depth}]} + test_edesign = FreeformDesign(circ_dict) + + num_mcs_per_circ = 4 + num_ref_per_qubit_subset = 10 + + mirror_edesign = make_mirror_edesign(test_edesign=test_edesign, account_for_routing=False, + ref_edesign=None, + num_mcs_per_circ=num_mcs_per_circ, + num_ref_per_qubit_subset=num_ref_per_qubit_subset, + mirroring_strategy='central_pauli', + gate_set='u3_cx_cz' + ) + + self.assertEqual(num_mcs_per_circ + num_ref_per_qubit_subset, len(mirror_edesign.all_circuits_needing_data)) + + flat_test_auxlist = [aux for auxlist in test_edesign.aux_info.values() for aux in auxlist] + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['cp'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + pspec = QubitProcessorSpec(5, ['Gu3', 'Gcphase', 'Gcnot', 'Gi'], + availability={'Gcnot':'all-permutations', 'Gcphase': 'all-permutations'}, + qubit_labels=['Q1', 'Q2', 'Q3', 'Q4', 'Q5']) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + + for edkey, ed in mirror_edesign.items(): + for circ, auxlist in ed.aux_info.items(): + target_bs = auxlist[0]['idealout'] + result = ideal_model.probabilities(circ) + self.assertAlmostEqual(result[(target_bs,)], 1.0) + + +class TestScalableBenchmarks(BaseCase): + def test_noise_mirror_benchmark(self): + try: + import qiskit + except: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + try: + import qiskit_ibm_runtime + except: + self.skipTest('Qiskit Runtime is required for this operation, and does not appear to be installed.') + + backend = qiskit_ibm_runtime.fake_provider.FakeFez() + + num_mcs_per_circ = 5 + num_ref_per_qubit_subset = 7 + mirroring_kwargs_dict = {'num_mcs_per_circ': num_mcs_per_circ, + 'num_ref_per_qubit_subset': num_ref_per_qubit_subset} + + qk_circ = qiskit.circuit.library.QFT(6) + qk_circ = qiskit.transpile(qk_circ, backend=backend) + + test_edesign, mirror_edesign = qiskit_circuits_to_mirror_edesign([qk_circ], + mirroring_kwargs_dict=mirroring_kwargs_dict) + + self.assertEqual(2*num_mcs_per_circ + num_ref_per_qubit_subset, len(mirror_edesign.all_circuits_needing_data)) + + flat_test_auxlist = [aux for auxlist in test_edesign.aux_info.values() for aux in auxlist] + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['br'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['rr'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + test_circ = test_edesign.all_circuits_needing_data[0] + + pspec = QubitProcessorSpec(test_circ.width, ['Gu3','Gxpi2','Gxpi','Gzr','Gi','Gcphase'], + availability={'Gcnot':'all-permutations','Gcphase':'all-permutations'}, + qubit_labels=test_circ.line_labels) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + for edkey, ed in mirror_edesign.items(): + for circ, auxlist in ed.aux_info.items(): + target_bs = auxlist[0]['idealout'] + result = ideal_model.probabilities(circ) + + ###### + # qk_circ = circ.convert_to_qiskit(qubit_conversion={q:i for i,q in enumerate(circ.line_labels)}, + # block_between_layers=True) + + # statevec = qiskit.quantum_info.Statevector.from_instruction(qk_circ) + # probs = statevec.probabilities_dict() + # probs = {k[::-1]: v for k,v in probs.items()} # Qiskit endianness + + self.assertAlmostEqual(result[(target_bs,)], 1.0) + + + + def test_fullstack_mirror_benchmark(self): + try: + import qiskit + except: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + try: + import qiskit_ibm_runtime + except: + self.skipTest('Qiskit Runtime is required for this operation, and does not appear to be installed.') + + backend = qiskit_ibm_runtime.fake_provider.FakeFez() + + num_mcs_per_circ = 3 + num_ref_per_qubit_subset = 4 + + mirroring_kwargs_dict = {'num_mcs_per_circ': num_mcs_per_circ, + 'num_ref_per_qubit_subset': num_ref_per_qubit_subset} + transpiler_kwargs_dict = {'optimization_level': 2} + + qk_circ = qiskit.circuit.library.QFT(6) + + test_edesign, mirror_edesign = qiskit_circuits_to_fullstack_mirror_edesign([qk_circ], + qk_backend=backend, + transpiler_kwargs_dict=transpiler_kwargs_dict, + mirroring_kwargs_dict=mirroring_kwargs_dict) + + self.assertEqual(2*num_mcs_per_circ + num_ref_per_qubit_subset, len(mirror_edesign.all_circuits_needing_data)) + + flat_test_auxlist = [aux for auxlist in test_edesign.aux_info.values() for aux in auxlist] + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['br'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['rr'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + test_circ = test_edesign.all_circuits_needing_data[0] + + pspec = QubitProcessorSpec(test_circ.width, ['Gu3','Gxpi2','Gxpi','Gzr','Gi','Gcphase'], + availability={'Gcnot':'all-permutations','Gcphase':'all-permutations'}, + qubit_labels=test_circ.line_labels) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + for edkey, ed in mirror_edesign.items(): + for circ, auxlist in ed.aux_info.items(): + target_bs = auxlist[0]['idealout'] + result = ideal_model.probabilities(circ) + + ###### + # qk_circ = circ.convert_to_qiskit(qubit_conversion={q:i for i,q in enumerate(circ.line_labels)}, + # block_between_layers=True) + + # statevec = qiskit.quantum_info.Statevector.from_instruction(qk_circ) + # probs = statevec.probabilities_dict() + # probs = {k[::-1]: v for k,v in probs.items()} # Qiskit endianness + + self.assertAlmostEqual(result[(target_bs,)], 1.0) + + def test_subcircuit_mirror_benchmark(self): + try: + import qiskit + except: + self.skipTest('Qiskit is required for this operation, and does not appear to be installed.') + + try: + import qiskit_ibm_runtime + except: + self.skipTest('Qiskit Runtime is required for this operation, and does not appear to be installed.') + + from pygsti.baseobjs.unitarygatefunction import UnitaryGateFunction + + class Gdelay(UnitaryGateFunction): + shape = (2, 2) + def __call__(self, dt): + return np.eye(2) + + + backend = qiskit_ibm_runtime.fake_provider.FakeFez() + + width_depths = {2: [2,4,6], + 4: [4,8,10]} + + num_samples_per_width_depth = 3 + + subcirc_kwargs_dict = {'num_samples_per_width_depth': num_samples_per_width_depth} + + num_mcs_per_circ = 3 + num_ref_per_qubit_subset = 4 + + mirroring_kwargs_dict = {'num_mcs_per_circ': num_mcs_per_circ, + 'num_ref_per_qubit_subset': num_ref_per_qubit_subset} + + + + + qk_circ = qiskit.circuit.library.QFT(6) + qk_circ = qiskit.transpile(qk_circ, backend=backend) + + test_edesign, mirror_edesign = qiskit_circuits_to_subcircuit_mirror_edesign( + [qk_circ], + aggregate_subcircs=True, + width_depth_dict=width_depths, + coupling_map=backend.coupling_map, + instruction_durations=backend.instruction_durations, + subcirc_kwargs_dict=subcirc_kwargs_dict, + mirroring_kwargs_dict=mirroring_kwargs_dict + ) + + + flat_test_auxlist = [aux for auxlist in test_edesign.aux_info.values() for aux in auxlist] + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['br'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + flat_mirror_auxlist = [aux['base_aux'] for auxlist in mirror_edesign['rr'].aux_info.values() for aux in auxlist] + self.assertTrue(all(flat_mirror_auxlist.count(test_aux) == num_mcs_per_circ for test_aux in flat_test_auxlist)) + + test_circ = test_edesign.all_circuits_needing_data[0] + + pspec = QubitProcessorSpec(test_circ.width, ['Gu3','Gxpi2','Gxpi','Gzr','Gi','Gcphase','Gdelay'], + nonstd_gate_unitaries={'Gdelay': Gdelay()}, + availability={'Gcnot':'all-permutations','Gcphase':'all-permutations'}, + qubit_labels=test_circ.line_labels) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + for edkey, ed in mirror_edesign.items(): + for circ, auxlist in ed.aux_info.items(): + + pspec = QubitProcessorSpec(circ.width, ['Gu3','Gxpi2','Gxpi','Gzr','Gi','Gcphase','Gdelay'], + nonstd_gate_unitaries={'Gdelay': Gdelay()}, + availability={'Gcnot':'all-permutations','Gcphase':'all-permutations'}, + qubit_labels=circ.line_labels) + + ideal_model = mc.create_crosstalk_free_model(pspec) + + target_bs = auxlist[0]['idealout'] + result = ideal_model.probabilities(circ) + + ###### + # qk_circ = circ.convert_to_qiskit(qubit_conversion={q:i for i,q in enumerate(circ.line_labels)}, + # block_between_layers=True) + + # statevec = qiskit.quantum_info.Statevector.from_instruction(qk_circ) + # probs = statevec.probabilities_dict() + # probs = {k[::-1]: v for k,v in probs.items()} # Qiskit endianness + + self.assertAlmostEqual(result[(target_bs,)], 1.0)