In the last post I ended by coming up with an algorithm that would randomly select one out of N qubits. Given this we could actually implement a shuffling algorithm which would be my first real quantum algorithm.
The principle of the algorithm is this:
It should be possible to implement this as an entire quantum program that does not measure any internal state until the very end.
This works with 8-bit values, so I need a way to initialize them. What I can do is allow the “program” to occupy the high qubits and just copy the bit representation of the ints over to the low bits.
It’s easy enough to do this using the initialize method which can set multiple qubits at once. Let’s try that with a fixed binary string:
from qiskit import QuantumCircuit, Aer, execute
qc = QuantumCircuit(8, 2)
qc.initialize("01010101")
backend = Aer.get_backend('qasm_simulator')
job = execute(qc, backend, memory=True)
job.result().get_counts(qc)Simulation failed and returned the following error message:
ERROR: Failed to load qobj: Unable to cast Python instance to C++ type (compile in debug mode for details)QiskitError: 'Data for experiment "circuit-0" could not be found.'Boo, boo. This looks like a bug as the underlying issue is the failure to load the qobj.
from qiskit import QuantumCircuit, Aer, execute
qc = QuantumCircuit(4, 2)
qc.initialize([float(bit) for bit in "0101"])
backend = Aer.get_backend('qasm_simulator')
job = execute(qc, backend, memory=True)
job.result().get_counts(qc)QiskitError: 'Sum of amplitudes-squared does not equal one.'The problem here is that it is tying together these separate qubits into a single initialization. It’s easier to just flip the bits that are set:
from qiskit import QuantumCircuit
from qiskit.result.result import Result
from qiskit.visualization import plot_histogram
def run(qc: QuantumCircuit, shots: int = 1024) -> Result:
# execute is a super generic term
from qiskit import Aer, execute
backend = Aer.get_backend('qasm_simulator')
job = execute(qc, backend, shots=shots, memory=True)
return job.result()
def show_results(qc: QuantumCircuit, shots: int = 1024) -> Result:
output = run(qc, shots=shots)
display(qc.draw(output="mpl"))
display(
plot_histogram(output.get_counts(qc))
)
return outputfrom qiskit import QuantumCircuit
value = 42 # 00101010
# get the value, which looks like 0b101010
# then reverse it so that enumerate matches bit index
binary_value = bin(value)[2:][::-1]
qc = QuantumCircuit(8, 8)
for index, value in enumerate(binary_value):
if value == "1":
qc.x(index)
qc.barrier() # makes it render nicely
qc.measure(range(8), range(8))
display(qc.draw(output="mpl"))
binary_result = run(qc, shots=1).get_memory()[0]
print(f"output: {binary_result}")
print(f"as int: {int(binary_result, base=2)}")
output: 00101010
as int: 42That was more painful than I expected. Anyway with this done I can initialize the qubits and read the results. Let’s try unconditionally swapping two values.
There is an operation available on the Quantum Circuit called swap which is associated with the SwapGate. It does exactly what it sounds like - it swaps the two states of the provided qubits without involving an intermediate qubit. The mechanics of this is interesting however I’m more focused on implementing this algorithm right now.
from qiskit import QuantumCircuit
value_1 = 42
value_2 = 69
binary_value_1 = bin(value_1)[2:][::-1]
binary_value_2 = bin(value_2)[2:][::-1]
qc = QuantumCircuit(16, 16)
for index, value in enumerate(binary_value_1):
if value == "1":
qc.x(index)
for index, value in enumerate(binary_value_2, start=8):
if value == "1":
qc.x(index)
for left_qubit, right_qubit in enumerate(range(8, 16)):
qc.swap(qubit1=left_qubit, qubit2=right_qubit)
qc.measure(range(16), range(16))
# the image is getting really big
# display(qc.draw(output="mpl"))
# the memory is arranged from high to low so the first 8 characters correspond to qubits 8-15
binary_result = run(qc, shots=1).get_memory()[0]
binary_result_1 = binary_result[8:]
print(f"input 1: {value_1}")
print(f"output 1: {binary_result_1}")
print(f"as int 1: {int(binary_result_1, base=2)}")
print()
binary_result_2 = binary_result[:8]
print(f"input 2: {value_2}")
print(f"output 2: {binary_result_2}")
print(f"as int 2: {int(binary_result_2, base=2)}")input 1: 42
output 1: 01000101
as int 1: 69
input 2: 69
output 2: 00101010
as int 2: 42There is a Control Swap Gate as well which will be vital for correctly using this. So with everything set up can we create a probabilistic int flipper?
It’s easy to put a qubit into the 50/50 superposition, as the Hadamard gate does that for us directly. If we use a qubit that has been put into that superposition as a control qubit then we can conditionally swap. Then roughly half of the time the two numbers will be swapped.
from qiskit import QuantumCircuit
value_1 = 42
value_2 = 69
binary_value_1 = bin(value_1)[2:][::-1]
binary_value_2 = bin(value_2)[2:][::-1]
qc = QuantumCircuit(17, 16) # one extra bit to determine if we swap
for index, value in enumerate(binary_value_1):
if value == "1":
qc.x(index)
for index, value in enumerate(binary_value_2, start=8):
if value == "1":
qc.x(index)
# the conditional qubit
qc.h(qubit=16)
for left_qubit, right_qubit in enumerate(range(8, 16)):
qc.cswap(control_qubit=16, target_qubit1=left_qubit, target_qubit2=right_qubit)
qc.measure(range(16), range(16))
# the memory is arranged from high to low so the first 8 characters correspond to qubits 8-15
# show_results(qc) ; None
output = run(qc).get_counts()
for key, count in output.items():
print(f"{int(key[8:], base=2)} is the first value {count:,} times")69 is the first value 509 times
42 is the first value 515 timesWith this we have our first number shuffler!
Extending this is tricky because there is no looping construct that I have found. The instructions that you are defining proceed forward in time, so every loop needs to be unrolled. I need to write this in a way where I don’t have to operate over previous “loop” results as that could get really messy.
With this in mind I want to be able to set one index at a time and then ignore it going forward. This leads to a simple algorithm, which was introduced at the start of this post:
To make it easy to start with I am going to implement this manually as a 3 way shuffle. This will lead to large blocks of code. After that we can look at how to make it N-way.
I will start by defining a few functions that will be useful going forward:
def to_rotation(probability: float) -> float:
return math.acos(-(probability * 2 - 1))
def initialize(qc: QuantumCircuit, value: int, offset: int) -> None:
binary_value = bin(value)[2:][::-1]
for index, value in enumerate(binary_value, start=offset):
if value == "1":
qc.x(index)
def conditional_swap(qc: QuantumCircuit, control_qubit: int, left_offset: int, right_offset: int) -> None:
for left_qubit, right_qubit in zip(
range(left_offset, left_offset+8),
range(right_offset, right_offset+8)
):
qc.cswap(
control_qubit=control_qubit,
target_qubit1=left_qubit,
target_qubit2=right_qubit
)And we can start by performing the first step, which will:
from qiskit import QuantumCircuit
from random import randint
import math
control_qubit_1 = 8*3
control_qubit_2 = 8*3 + 1
value_1 = randint(0, 255)
value_2 = randint(0, 255)
value_3 = randint(0, 255)
qc = QuantumCircuit((8 * 3) + 2, 8 * 3)
initialize(qc, value=value_1, offset=0)
initialize(qc, value=value_2, offset=8)
initialize(qc, value=value_3, offset=8*2)
# the conditional qubits
qc.rx(theta=to_rotation(1/3), qubit=control_qubit_1)
qc.h(qubit=control_qubit_2)
# turn off conditional qubit 2 if conditional qubit 1 is set
qc.ch(control_qubit=control_qubit_1, target_qubit=control_qubit_2)
conditional_swap(qc, control_qubit=control_qubit_1, left_offset=0, right_offset=8)
conditional_swap(qc, control_qubit=control_qubit_2, left_offset=0, right_offset=8*2)
qc.measure(range(8*3), range(8*3))
output = run(qc).get_counts()
print([value_1, value_2, value_3])
for key, count in output.items():
# remember that the high qubit values are in the _low_ string indexes,
# so we step through this backwards
values = [
int(key[index-8:index], base=2)
for index in range(8*3, 0, -8)
]
print(f"{values}: {count:,}")[0, 48, 158]
[48, 0, 158]: 331
[0, 48, 158]: 350
[158, 48, 0]: 343To complete this 3 way shuffler the second and third numbers need to be conditionally swapped.
from qiskit import QuantumCircuit
from random import randint
from itertools import permutations
control_qubit_1_1 = 8*3
control_qubit_1_2 = 8*3 + 1
control_qubit_2 = 8*3 + 2
# easier to see that it is working with these values
value_1 = 1
value_2 = 2
value_3 = 3
qc = QuantumCircuit((8 * 3) + 3, 8 * 3)
initialize(qc, value=value_1, offset=0)
initialize(qc, value=value_2, offset=8)
initialize(qc, value=value_3, offset=8*2)
# first swap: 0 with either 1 or 2 (or no swap)
# the conditional qubits
qc.rx(theta=to_rotation(1/3), qubit=control_qubit_1_1)
qc.h(qubit=control_qubit_1_2)
# turn off conditional qubit 2 if conditional qubit 1 is set
qc.ch(control_qubit=control_qubit_1_1, target_qubit=control_qubit_1_2)
conditional_swap(
qc,
control_qubit=control_qubit_1_1,
left_offset=0,
right_offset=8
)
conditional_swap(
qc,
control_qubit=control_qubit_1_2,
left_offset=0,
right_offset=8*2
)
# second swap: 1 with 2 (or no swap)
# the conditional qubit
qc.h(qubit=control_qubit_2)
conditional_swap(
qc,
control_qubit=control_qubit_2,
left_offset=8,
right_offset=8*2
)
qc.measure(range(8*3), range(8*3))
# at this point, with barriers, the circuit image wraps 3 times
# output = run(qc, shots=1_000_000).get_counts()
output = run(qc).get_counts()
print([value_1, value_2, value_3])
for key in permutations([value_1, value_2, value_3]):
index = "".join(
f"{bin(value)[2:]:>08}"
for value in key[::-1]
)
print(f"{key}: {output[index]:,}")[1, 2, 3]
(1, 2, 3): 151
(1, 3, 2): 178
(2, 1, 3): 182
(2, 3, 1): 170
(3, 1, 2): 182
(3, 2, 1): 161This can now produce all 6 permutations of a three element list.
One thing that is annoying is the fact that you can’t use a qubit as a control and as the target, which makes it hard to reliably clear a control qubit to reuse it later. For now I am going to use a new set of control qubits each time.
Given we have a way to pick a value out of the list and then swap it with the first value, we can implement the shuffler by repeatedly invoking this, moving one value forward each time. This is nice because it has a fixed shape - the list elements to ignore are known irrespective of the outcome of the random process.
from typing import List
from qiskit import QuantumCircuit
import math
def nway_shuffle(values: List[int]) -> List[int]:
qc = nway_shuffler(values)
output = run(qc, shots=1).get_memory()[0]
return [
int(output[index*8:(index+1)*8], base=2)
for index in range(len(output)//8)
][::-1]
def nway_shuffler(values: List[int], barriers: bool = False) -> QuantumCircuit:
value_count = len(values)
# one choice is skipped at each level
control_qubits = sum(range(1, value_count))
qc = QuantumCircuit(
(8 * value_count) + control_qubits,
8 * value_count,
)
for index, value in enumerate(values):
initialize(qc, value=value, offset=index*8)
if barriers:
qc.barrier()
control_qubit_offset = len(values) * 8
for iteration in range(value_count-1):
count = value_count - iteration
nway_shuffle_step(
qc=qc,
n=count,
value_qubit_offset=8*iteration,
control_qubit_offset=control_qubit_offset,
)
control_qubit_offset += count-1
if barriers:
qc.barrier()
qc.measure(range(8*value_count), range(8*value_count))
return qc
def nway_shuffle_step(
qc: QuantumCircuit,
n: int,
value_qubit_offset: int,
control_qubit_offset: int,
) -> None:
nway_chooser(
qc=qc,
n=n,
control_qubit_offset=control_qubit_offset
)
for index in range(n-1):
conditional_swap(
qc,
control_qubit=control_qubit_offset + index,
left_offset=value_qubit_offset,
right_offset=value_qubit_offset + (8*(index+1))
)
def nway_chooser(
qc: QuantumCircuit,
n: int,
control_qubit_offset: int,
) -> None:
# n less one because not doing something is a choice
n = n - 1
rotations = {
v: to_rotation(1/(1+n-v))
for v in range(n)
}
for qubit, rotation in rotations.items():
qc.rx(theta=rotation, qubit=control_qubit_offset + qubit)
for control_index in range(n - 1):
for target_index in range(control_index+1, n):
control_qubit = control_qubit_offset + control_index
target_qubit = control_qubit_offset + target_index
target_rotation = rotations[target_index]
qc.crx(
theta=-target_rotation,
control_qubit=control_qubit,
target_qubit=target_qubit
)nway_shuffle([1, 2, 3])[3, 2, 1]Running this once does not establish if it actually works though. We need to run it repeatedly to ensure that all of the possible values can be generated and that they are evenly distributed.
qc = nway_shuffler([1, 2, 3])
results = run(qc, shots=1_000_000)
results.get_counts(){'000000110000001000000001': 166426,
'000000010000001100000010': 166050,
'000000100000001100000001': 166576,
'000000010000001000000011': 166890,
'000000100000000100000011': 167434,
'000000110000000100000010': 166624}counts = list(results.get_counts().values())
variance = (max(counts) / min(counts)) - 1
print(f"The most frequent result occurs {variance*100:0.3f}% more than the least frequent result")The most frequent result occurs 0.833% more than the least frequent resultThere are 6 different values and they are within ~1% of each other. I’m reasonably happy with this.
How fast does it run?
import timeit
def time_quantum_shuffle(n: int, runs: int = 3) -> float:
total_time = timeit.timeit(
lambda: nway_shuffle(list(range(n))),
globals=globals(),
number=runs
)
return total_time / runs
quantum_shuffle_times = pd.DataFrame([
{
"number": n,
"time (seconds)": time_quantum_shuffle(n),
}
for n in range(2, 7)
]).set_index("number")
quantum_shuffle_times.plot()
quantum_shuffle_times| time (seconds) | |
|---|---|
| number | |
| 2 | 0.011618 |
| 3 | 1.721771 |
| 4 | 0.088402 |
| 5 | 2.040435 |
| 6 | 635.043775 |
While this algorithm works it scales hilariously badly. My CPU / RAM weren’t even taxed that heavily. I think that this has become too large to effectively simulate?
Ultimately the qubit simulator maintains a matrix of every qubit and their interactions. It may well be that this matrix became too big to operate efficiently.
While I think that the implementation that I have is fine the underlying problem is that the system ends up in an extremely entangled state and there are a large number of qubits to simulate. I can address this by returning to the N-bit chooser from the previous post and gluing the other parts together in plain python. While it’s not ideal I do hope that this would be significantly faster.
from qiskit import QuantumCircuit
import math
def quantum_choice(n: int) -> int:
qc = nway_chooser(n)
output = run(qc, shots=1).get_memory()[0]
return output.index("1")
def nway_chooser(n: int) -> QuantumCircuit:
qc = QuantumCircuit(n, n)
rotations = {
v: to_rotation(1/(n-v))
for v in range(n)
}
for qubit, rotation in rotations.items():
qc.rx(theta=rotation, qubit=qubit)
for control_qubit in range(n):
for target_qubit in range(control_qubit+1, n):
target_rotation = rotations[target_qubit]
qc.crx(
theta=-target_rotation,
control_qubit=control_qubit,
target_qubit=target_qubit
)
qc.measure(range(n), range(n))
return qcThe scaling here is approaching linear for the individual invocations. If we look at the cumulative time it has vastly improved compared to the full quantum implementation.
import timeit
def time_one_shuffle(n: int, runs: int = 10) -> float:
total_time = timeit.timeit(
lambda: quantum_choice(n),
globals=globals(),
number=runs
)
return total_time / runs
one_shuffle_times = pd.DataFrame([
{
"number": n,
"time (seconds)": time_one_shuffle(n),
}
for n in range(2, 20)
]).set_index("number")
one_shuffle_times["cumulative time (seconds)"] = one_shuffle_times["time (seconds)"].cumsum()
one_shuffle_times.plot()
one_shuffle_times| time (seconds) | cumulative time (seconds) | |
|---|---|---|
| number | ||
| 2 | 0.012632 | 0.012632 |
| 3 | 0.024178 | 0.036810 |
| 4 | 0.041090 | 0.077900 |
| 5 | 0.063083 | 0.140983 |
| 6 | 0.089879 | 0.230862 |
| 7 | 0.122087 | 0.352949 |
| 8 | 0.161558 | 0.514507 |
| 9 | 0.201130 | 0.715636 |
| 10 | 0.250159 | 0.965795 |
| 11 | 0.304193 | 1.269988 |
| 12 | 0.360929 | 1.630918 |
| 13 | 0.432081 | 2.062998 |
| 14 | 0.501075 | 2.564073 |
| 15 | 0.584941 | 3.149014 |
| 16 | 0.670605 | 3.819619 |
| 17 | 0.762924 | 4.582543 |
| 18 | 0.843015 | 5.425558 |
| 19 | 0.940590 | 6.366148 |
If we use this to implement the full shuffler:
def quantum_shuffle(values: List[int]) -> List[int]:
if not values:
return []
values = list(values)
shuffled = []
while len(values) > 1:
index = quantum_choice(len(values))
shuffled.append(values.pop(index))
shuffled.append(values[0])
return shuffledWe can see that this works as it completely reorders this list:
quantum_shuffle(list(range(10)))[0, 5, 2, 6, 9, 1, 3, 8, 7, 4]Does this address the performance problem?
import timeit
def time_full_shuffle(n: int, runs: int = 10) -> float:
total_time = timeit.timeit(
lambda: quantum_shuffle(list(range(n))),
globals=globals(),
number=runs
)
return total_time / runs
full_shuffle_times = pd.DataFrame([
{
"number": n,
"time (seconds)": time_one_shuffle(n),
}
for n in range(2, 15)
]).set_index("number")
full_shuffle_times.plot()
full_shuffle_times| time (seconds) | |
|---|---|
| number | |
| 2 | 0.012342 |
| 3 | 0.023739 |
| 4 | 0.040371 |
| 5 | 0.062664 |
| 6 | 0.088807 |
| 7 | 0.119778 |
| 8 | 0.156464 |
| 9 | 0.196033 |
| 10 | 0.242926 |
| 11 | 0.294877 |
| 12 | 0.354454 |
| 13 | 0.420764 |
| 14 | 0.487811 |
This is fast enough to be usable for non tiny lists.
These times are extremely similar to the times for the individual invocations which surprises me. Something to think about.
It’s a shame that I could not implement the full shuffler using quantum computing. I do think that this optimized algorithm is fast enough for my unstated goal…
"".join(
quantum_shuffle(list("puzzler"))
)'rpzzeul'