diff --git a/README.md b/README.md index 9641bd0..b0cfd74 100644 --- a/README.md +++ b/README.md @@ -94,16 +94,6 @@ Documentation is stored in the [`docs/`](docs/) directory. It is organized acco The current version of this module has some important usage notes and limitations. -#### Encoding Limitations - -The [Magic Rounding scheme](https://github.com/qiskit-community/prototype-qrao/blob/main/qrao/magic_rounding.py) is currently only compatible with the (3,1,p) QRAC for encoding classical binary variables onto qubits. If using this rounding scheme, be sure to set `max_vars_per_qubit=3` when instantiating the encoding: -``` -encoding = QuantumRandomAccessEncoding(max_vars_per_qubit=3) -``` -Refer to the [tutorial on rounding schemes](https://github.com/qiskit-community/prototype-qrao/blob/main/docs/tutorials/02_advanced_usage.ipynb) for example usage. - -We plan to make Magic Rounding compatible with additional encoding options in the future. - #### Supported Problem Types As with `qiskit-optimization`, this module currently only supports input problems that are of type `QuadraticProgram` and that can be converted to a QUBO. Furthermore, the conversion to a QUBO must be made _before_ encoding the problem into a quantum Hamiltonian with `QuantumRandomAccessEncoding` (refer to the tutorial to [see an example of this](https://github.com/qiskit-community/prototype-qrao/blob/main/docs/tutorials/01_getting_started.ipynb)). diff --git a/qrao/magic_rounding.py b/qrao/magic_rounding.py index a45d461..da867da 100644 --- a/qrao/magic_rounding.py +++ b/qrao/magic_rounding.py @@ -24,7 +24,7 @@ from qiskit.opflow import PrimitiveOp from qiskit.utils import QuantumInstance -from .encoding import z_to_31p_qrac_basis_circuit +from .encoding import z_to_31p_qrac_basis_circuit, z_to_21p_qrac_basis_circuit from .rounding_common import ( RoundingSolutionSample, RoundingScheme, @@ -99,17 +99,22 @@ class MagicRounding(RoundingScheme): """ - # This information should really come from the encoding. - # Right now it's basically hard coded for one QRAC. - _ASSIGN = ( - {"0": (0, 0, 0), "1": (1, 1, 1)}, # IMI - {"0": (0, 1, 1), "1": (1, 0, 0)}, # XMX - {"0": (1, 0, 1), "1": (0, 1, 0)}, # YMY - {"0": (1, 1, 0), "1": (0, 0, 1)}, # ZMZ - ) - - # Pauli op string to label index in ops (assumes 3-QRAC) - _XYZ = {"X": 0, "Y": 1, "Z": 2} + _DECODING = { + 3: ( # Eq. (8) + {"0": [0, 0, 0], "1": [1, 1, 1]}, # I mu+ I, I mu- I + {"0": [0, 1, 1], "1": [1, 0, 0]}, # X mu+ X, X mu- X + {"0": [1, 0, 1], "1": [0, 1, 0]}, # Y mu+ Y, Y mu- Y + {"0": [1, 1, 0], "1": [0, 0, 1]}, # Z mu+ Z, Z mu- Z + ), + 2: ( # Sec. VII + {"0": [0, 0], "1": [1, 1]}, # I xi+ I, I xi- I + {"0": [0, 1], "1": [1, 0]}, # X xi+ X, X xi- X + ), + 1: ({"0": [0], "1": [1]},), + } + + # Pauli op string to label index in ops + _OP_INDICES = {1: {"Z": 0}, 2: {"X": 0, "Z": 1}, 3: {"X": 0, "Y": 1, "Z": 2}} def __init__( self, @@ -182,18 +187,19 @@ def _unpack_measurement_outcome( bits: str, basis: List[int], var2op: Dict[int, Tuple[int, PrimitiveOp]], + vars_per_qubit: int, ) -> List[int]: output_bits = [] # iterate in order over decision variables # (assumes variables are numbered consecutively beginning with 0) for var in range(len(var2op)): # pylint: disable=consider-using-enumerate q, op = var2op[var] - # get the index in [0,1,2] corresponding - # to each possible Pauli. - op_index = self._XYZ[str(op)] + # get the decoding outcome index for the variable + # corresponding to this Pauli op. + op_index = self._OP_INDICES[vars_per_qubit][str(op)] # get the bits associated to this magic basis' # measurement outcomes - bit_outcomes = self._ASSIGN[basis[q]] + bit_outcomes = self._DECODING[vars_per_qubit][basis[q]] # select which measurement outcome we observed # this gives up to 3 bits of information magic_bits = bit_outcomes[bits[q]] @@ -205,19 +211,26 @@ def _unpack_measurement_outcome( @staticmethod def _make_circuits( - circ: QuantumCircuit, bases: List[List[int]], measure: bool + circ: QuantumCircuit, bases: List[List[int]], measure: bool, vars_per_qubit: int ) -> List[QuantumCircuit]: circuits = [] for basis in bases: - qc = circ.compose( - z_to_31p_qrac_basis_circuit(basis).inverse(), inplace=False - ) + if vars_per_qubit == 3: + qc = circ.compose( + z_to_31p_qrac_basis_circuit(basis).inverse(), inplace=False + ) + elif vars_per_qubit == 2: + qc = circ.compose( + z_to_21p_qrac_basis_circuit(basis).inverse(), inplace=False + ) + elif vars_per_qubit == 1: + qc = circ.copy() if measure: qc.measure_all() circuits.append(qc) return circuits - def _evaluate_magic_bases(self, circuit, bases, basis_shots): + def _evaluate_magic_bases(self, circuit, bases, basis_shots, vars_per_qubit): """ Given a circuit you wish to measure, a list of magic bases to measure, and a list of the shots to use for each magic basis configuration. @@ -228,7 +241,7 @@ def _evaluate_magic_bases(self, circuit, bases, basis_shots): len(bases) == len(basis_shots) == len(basis_counts) """ measure = not _is_original_statevector_simulator(self.quantum_instance.backend) - circuits = self._make_circuits(circuit, bases, measure) + circuits = self._make_circuits(circuit, bases, measure, vars_per_qubit) # Execute each of the rotated circuits and collect the results @@ -281,7 +294,7 @@ def _evaluate_magic_bases(self, circuit, bases, basis_shots): return basis_counts - def _compute_dv_counts(self, basis_counts, bases, var2op): + def _compute_dv_counts(self, basis_counts, bases, var2op, vars_per_qubit): """ Given a list of bases, basis_shots, and basis_counts, convert each observed bitstrings to its corresponding decision variable @@ -294,7 +307,9 @@ def _compute_dv_counts(self, basis_counts, bases, var2op): for bitstr, count in counts.items(): # For each bit in the observed bitstring... - soln = self._unpack_measurement_outcome(bitstr, base, var2op) + soln = self._unpack_measurement_outcome( + bitstr, base, var2op, vars_per_qubit + ) soln = "".join([str(int(bit)) for bit in soln]) if soln in dv_counts: dv_counts[soln] += count @@ -302,17 +317,15 @@ def _compute_dv_counts(self, basis_counts, bases, var2op): dv_counts[soln] = count return dv_counts - def _sample_bases_uniform(self, q2vars): + def _sample_bases_uniform(self, q2vars, vars_per_qubit): bases = [ - self.rng.choice( - 4, size=len(q2vars), p=[1 / 4, 1 / 4, 1 / 4, 1 / 4] - ).tolist() + self.rng.choice(2 ** (vars_per_qubit - 1), size=len(q2vars)).tolist() for _ in range(self.shots) ] bases, basis_shots = np.unique(bases, axis=0, return_counts=True) return bases, basis_shots - def _sample_bases_weighted(self, q2vars, trace_values): + def _sample_bases_weighted(self, q2vars, trace_values, vars_per_qubit): """Perform weighted sampling from the expectation values. The goal is to make smarter choices about which bases to measure in @@ -327,23 +340,38 @@ def _sample_bases_weighted(self, q2vars, trace_values): # probability of picking the corresponding magic basis on that qubit. basis_probs = [] for dvars in q2vars: - x = 0.5 * (1 - tv[dvars[0]]) - y = 0.5 * (1 - tv[dvars[1]]) if (len(dvars) > 1) else 0 - z = 0.5 * (1 - tv[dvars[2]]) if (len(dvars) > 2) else 0 - # ppp: mu± = .5(I ± 1/sqrt(3)( X + Y + Z)) - # pmm: X mu± X = .5(I ± 1/sqrt(3)( X - Y - Z)) - # mpm: Y mu± Y = .5(I ± 1/sqrt(3)(-X + Y - Z)) - # mmp: Z mu± Z = .5(I ± 1/sqrt(3)(-X - Y + Z)) - # fmt: off - ppp_mmm = x * y * z + (1-x) * (1-y) * (1-z) - pmm_mpp = x * (1-y) * (1-z) + (1-x) * y * z - mpm_pmp = (1-x) * y * (1-z) + x * (1-y) * z - ppm_mmp = x * y * (1-z) + (1-x) * (1-y) * z - # fmt: on - basis_probs.append([ppp_mmm, pmm_mpp, mpm_pmp, ppm_mmp]) - + if vars_per_qubit == 3: + x = 0.5 * (1 - tv[dvars[0]]) + y = 0.5 * (1 - tv[dvars[1]]) if (len(dvars) > 1) else 0 + z = 0.5 * (1 - tv[dvars[2]]) if (len(dvars) > 2) else 0 + # ppp: mu± = .5(I ± 1/sqrt(3)( X + Y + Z)) + # pmm: X mu± X = .5(I ± 1/sqrt(3)( X - Y - Z)) + # mpm: Y mu± Y = .5(I ± 1/sqrt(3)(-X + Y - Z)) + # mmp: Z mu± Z = .5(I ± 1/sqrt(3)(-X - Y + Z)) + # fmt: off + ppp_mmm = x * y * z + (1-x) * (1-y) * (1-z) + pmm_mpp = x * (1-y) * (1-z) + (1-x) * y * z + mpm_pmp = (1-x) * y * (1-z) + x * (1-y) * z + ppm_mmp = x * y * (1-z) + (1-x) * (1-y) * z + # fmt: on + basis_probs.append([ppp_mmm, pmm_mpp, mpm_pmp, ppm_mmp]) + elif vars_per_qubit == 2: + x = 0.5 * (1 - tv[dvars[0]]) + z = 0.5 * (1 - tv[dvars[1]]) if (len(dvars) > 1) else 0 + # pp: xi± = .5(I ± 1/sqrt(2)( X + Z )) + # pm: X xi± X = .5(I ± 1/sqrt(2)( X - Z )) + # fmt: off + pp_mm = x * z + (1-x) * (1-z) + pm_mp = x * (1-z) + (1-x) * z + # fmt: on + basis_probs.append([pp_mm, pm_mp]) + elif vars_per_qubit == 1: + basis_probs.append([1.0]) bases = [ - [self.rng.choice(4, size=1, p=probs)[0] for probs in basis_probs] + [ + self.rng.choice(2 ** (vars_per_qubit - 1), p=probs) + for probs in basis_probs + ] for _ in range(self.shots) ] bases, basis_shots = np.unique(bases, axis=0, return_counts=True) @@ -352,13 +380,6 @@ def _sample_bases_weighted(self, q2vars, trace_values): def round(self, ctx: RoundingContext) -> MagicRoundingResult: """Perform magic rounding""" - if ctx._encoding is not None and ctx._encoding.max_vars_per_qubit != 3: - raise ValueError( - "Currently, MagicRounding only supports 3-QRACs, " - "but the passed encoding is a " - f"{ctx._encoding.max_vars_per_qubit}-QRAC." - ) - start_time = time.time() trace_values = ctx.trace_values circuit = ctx.circuit @@ -371,14 +392,18 @@ def round(self, ctx: RoundingContext) -> MagicRoundingResult: # We've already checked that it is one of these two in the constructor if self.basis_sampling == "uniform": - bases, basis_shots = self._sample_bases_uniform(ctx.q2vars) + bases, basis_shots = self._sample_bases_uniform( + ctx.q2vars, ctx._vars_per_qubit + ) elif self.basis_sampling == "weighted": if trace_values is None: raise NotImplementedError( "Magic rounding with weighted sampling requires the trace values " "to be available, but they are not." ) - bases, basis_shots = self._sample_bases_weighted(ctx.q2vars, trace_values) + bases, basis_shots = self._sample_bases_weighted( + ctx.q2vars, trace_values, ctx._vars_per_qubit + ) else: # pragma: no cover raise NotImplementedError( f'No such basis sampling method: "{self.basis_sampling}".' @@ -389,10 +414,14 @@ def round(self, ctx: RoundingContext) -> MagicRoundingResult: # the appropriate number of times (given by basis_shots) # and return the circuit results - basis_counts = self._evaluate_magic_bases(circuit, bases, basis_shots) + basis_counts = self._evaluate_magic_bases( + circuit, bases, basis_shots, ctx._vars_per_qubit + ) # keys will be configurations of decision variables # values will be total number of observations. - soln_counts = self._compute_dv_counts(basis_counts, bases, ctx.var2op) + soln_counts = self._compute_dv_counts( + basis_counts, bases, ctx.var2op, ctx._vars_per_qubit + ) soln_samples = [ RoundingSolutionSample( diff --git a/qrao/rounding_common.py b/qrao/rounding_common.py index c457cc1..15bbf35 100644 --- a/qrao/rounding_common.py +++ b/qrao/rounding_common.py @@ -44,26 +44,29 @@ def __init__( var2op: Optional[Dict[int, Tuple[int, PrimitiveOp]]] = None, q2vars: Optional[List[List[int]]] = None, trace_values=None, - circuit=None + circuit=None, + vars_per_qubit: Optional[int] = None, ): if encoding is not None: if var2op is not None or q2vars is not None: raise ValueError( "Neither var2op nor q2vars should be provided if encoding is" ) + if vars_per_qubit is not None: + raise ValueError("vars_per_qubit should not be provided if encoding is") self.var2op = encoding.var2op self.q2vars = encoding.q2vars - # We save the encoding here, as a temporary hack, so that the magic - # rounding backend can ensure that it is a 3-QRAC. Once the magic - # rounding backend supports all QRACs, we will remove this - # (protected) attribute. - self._encoding = encoding + self._vars_per_qubit = encoding.max_vars_per_qubit else: if var2op is None: raise ValueError("Either an encoding or var2ops must be provided") + if vars_per_qubit is None: + raise ValueError( + "vars_per_qubit must be provided if encoding is not provided" + ) self.var2op = var2op self.q2vars = q2vars_from_var2op(var2op) if q2vars is None else q2vars - self._encoding = None + self._vars_per_qubit = vars_per_qubit self.trace_values = trace_values # TODO: rename me self.circuit = circuit # TODO: rename me diff --git a/tests/test_magic_rounding.py b/tests/test_magic_rounding.py index ac0affd..84b2d92 100644 --- a/tests/test_magic_rounding.py +++ b/tests/test_magic_rounding.py @@ -97,7 +97,10 @@ def test_round_on_gate_and_sv_circs(self): magic_basis = 2 * (m[1] ^ m[2]) + (m[0] ^ m[2]) tv = self.deterministic_trace_vals[magic_basis] rounding_context = RoundingContext( - trace_values=tv, circuit=qrac_gate_circ, var2op=var2op + trace_values=tv, + circuit=qrac_gate_circ, + var2op=var2op, + vars_per_qubit=3, ) rounding_res = magic.round(rounding_context) self.assertEqual(rounding_res.samples[0].x.tolist(), list(m)) @@ -112,7 +115,10 @@ def test_round_on_gate_and_sv_circs(self): magic_basis = 2 * (m[1] ^ m[2]) + (m[0] ^ m[2]) tv = self.deterministic_trace_vals[magic_basis] rounding_context = RoundingContext( - trace_values=tv, circuit=qrac_sv_circ, var2op=var2op + trace_values=tv, + circuit=qrac_sv_circ, + var2op=var2op, + vars_per_qubit=3, ) rounding_res = magic.round(rounding_context) self.assertEqual(rounding_res.samples[0].x.tolist(), list(m)) @@ -128,7 +134,7 @@ def test_evaluate_magic_bases(self): qrac_state = qrac_state_prep_1q(*m).to_circuit() bases = [[2 * (m[1] ^ m[2]) + (m[0] ^ m[2])]] basis_counts = magic._evaluate_magic_bases( - qrac_state, bases=bases, basis_shots=[10] + qrac_state, bases=bases, basis_shots=[10], vars_per_qubit=3 ) self.assertEqual(len(basis_counts), 1) self.assertEqual(int(list(basis_counts[0].keys())[0]), m[0] ^ m[1] ^ m[2]) @@ -148,7 +154,7 @@ def test_dv_counts(self): for outcome in range(4): bases = [[b0, b1]] basis_counts = [{f"{outcome:02b}": 1}] - dv_counts = compute_dv_counts(basis_counts, bases, var2op) + dv_counts = compute_dv_counts(basis_counts, bases, var2op, 3) solns.append(list(dv_counts.keys())[0]) ref = [ "000000", @@ -242,16 +248,20 @@ def test_sample_bases_weighted(self): # pylint: disable=too-many-locals for b0, tv0 in enumerate(self.deterministic_trace_vals): for b1, tv1 in enumerate(self.deterministic_trace_vals): tv = tv0 + tv1 - bases, basis_shots = sample_bases_weighted(q2vars, tv) + bases, basis_shots = sample_bases_weighted(q2vars, tv, 3) self.assertTrue(np.all(np.array([b0, b1]) == bases)) self.assertEqual(basis_shots, (shots,)) self.assertEqual(bases.shape, (1, num_qubits)) # 1 == deterministic # Both trace values and a circuit must be provided with self.assertRaises(NotImplementedError): - magic.round(RoundingContext(trace_values=[1.0], var2op=var2op)) + magic.round( + RoundingContext(trace_values=[1.0], var2op=var2op, vars_per_qubit=3) + ) with self.assertRaises(NotImplementedError): - magic.round(RoundingContext(circuit=self.gate_circ, var2op=var2op)) + magic.round( + RoundingContext(circuit=self.gate_circ, var2op=var2op, vars_per_qubit=3) + ) def test_sample_bases_uniform(self): """ @@ -268,7 +278,7 @@ def test_sample_bases_uniform(self): quantum_instance=qi, basis_sampling="uniform", ) - bases, basis_shots = magic._sample_bases_uniform(q2vars) + bases, basis_shots = magic._sample_bases_uniform(q2vars, 3) self.assertEqual(basis_shots.shape, (4,)) self.assertEqual(np.sum(basis_shots), shots) self.assertEqual(bases.shape, (4, num_qubits)) @@ -276,18 +286,9 @@ def test_sample_bases_uniform(self): # A circuit must be provided, but trace values need not be circuit = QuantumCircuit(1) circuit.h(0) - magic.round(RoundingContext(circuit=circuit, var2op=var2op)) + magic.round(RoundingContext(circuit=circuit, var2op=var2op, vars_per_qubit=3)) with self.assertRaises(NotImplementedError): - magic.round(RoundingContext(var2op=var2op)) - - -def test_unsupported_qrac(): - qi = QuantumInstance(backend=Aer.get_backend("aer_simulator"), shots=1000) - encoding = QuantumRandomAccessEncoding(2) - rounding = MagicRounding(quantum_instance=qi) - ctx = RoundingContext(encoding=encoding) - with pytest.raises(ValueError): - rounding.round(ctx) + magic.round(RoundingContext(var2op=var2op, vars_per_qubit=3)) def test_unsupported_backend(): @@ -321,7 +322,9 @@ def test_magic_rounding_statevector_simulator(): circ.h(0) circ.h(1) circ.cx(0, 1) - ctx = RoundingContext(circuit=circ, trace_values=[1, 1, 1], var2op=var2op) + ctx = RoundingContext( + circuit=circ, trace_values=[1, 1, 1], var2op=var2op, vars_per_qubit=3 + ) res = magic.round(ctx) assert sum(s.probability for s in res.samples) == pytest.approx(1) @@ -343,7 +346,11 @@ def test_noisy_quantuminstance(): ops = [X, Y, Z] var2op = {i: (i // 3, ops[i % 3]) for i in range(3)} circuit = qrac_state_prep_1q(0, 1, 0).to_circuit() - magic.round(RoundingContext(trace_values=[1, 1, 1], var2op=var2op, circuit=circuit)) + magic.round( + RoundingContext( + trace_values=[1, 1, 1], var2op=var2op, circuit=circuit, vars_per_qubit=3 + ) + ) def test_magic_rounding_statistical(): @@ -366,18 +373,18 @@ def test_magic_rounding_statistical(): ((0, 1, 0), 2, 1), ((1, 1, 0), 3, 0), ((0, 0, 1), 3, 1), - # 2 variables encoded as a (3,1,p) QRAC (see issue #11) + # 2 variables encoded as a (2,1,p) QRAC ((0, 0), 0, 0), - ((1, 1), 3, 0), - ((0, 1), 2, 1), + ((1, 1), 0, 1), + ((0, 1), 1, 0), ((1, 0), 1, 1), - # 1 variable encoded as a (3,1,p) QRAC (see issue #11) + # 1 variable encoded as a (1,1,1) QRAC ((0,), 0, 0), - ((1,), 1, 1), + ((1,), 0, 1), ) for (m, basis, expected) in test_cases: - encoding = QuantumRandomAccessEncoding() + encoding = QuantumRandomAccessEncoding(len(m)) encoding._add_variables(list(range(len(m)))) circuit = encoding.state_prep(m).to_circuit() result = mr.round(RoundingContext(encoding=encoding, circuit=circuit))