Quantum Computing#
Author:
Date:
Time Spent on this Assignment:
!pip install qiskit[visualization];
!pip install qutip
!pip install qiskit_aer;
import numpy as np
import pylab as plt
import qiskit
import qutip
from qiskit import QuantumCircuit, transpile
from qiskit_aer import AerSimulator
from qiskit.visualization import plot_histogram
def RunMe(circuit,a=None):
"""
Runs a quantum circuit and returns the state vector and measurement results.
Args:
circuit (QuantumCircuit): The quantum circuit to run.
a (list, optional): A list of initial state amplitudes. If provided, the circuit will be initialized with these amplitudes. Otherwise it will start from |0>.
Returns:
tuple: A tuple containing the state vector and measurement results.
"""
numQubits=circuit.num_qubits
if a!=None:
numQubits=max(numQubits,a.num_qubits)
initCircuit=QuantumCircuit(a.num_qubits,a.num_qubits)
initCircuit.initialize(a)
circuit=AddCircuits([initCircuit,circuit])
circuit.save_statevector(label='myStateVector')
compiled_circuit = transpile(circuit, simulator)
resultA = simulator.run(compiled_circuit).result()
forward=list(range(0,numQubits))
reverse=forward[::-1]
circuit.measure(forward,reverse)
compiled_circuit = transpile(circuit, simulator)
resultB = simulator.run(compiled_circuit).result()
return resultA.data()['myStateVector'],resultB.data()['counts']
def MakeState(v):
"""
Creates a qiskit state vector from a numpy array of amplitudes.
"""
numWires=int(round(np.log2(len(np.array(v)))))
qc1=QuantumCircuit(numWires,numWires)
qc1.initialize(v)
a,_=RunMe(qc1)
return a
def PrintDirac(out_state_a):
"""
Prints the state vector in Dirac notation.
"""
out_state=np.asarray(out_state_a)
l=len(out_state)
num_qubits=int(round(np.log(len(out_state))/np.log(2)))
for i in range(0,l):
if not np.isclose(out_state[i],0):
print(str(out_state[i])+'|'+bin(i)[2:].zfill(num_qubits)[::-1]+'>',end=' + ')
print()
def PlotBloch(v,bb):
"""
Plots the Bloch sphere representation of a quantum state.
"""
vv=np.asarray(v)
a=vv[0]*np.conj(vv[0])
b=vv[1]*np.conj(vv[0])
x = np.real(2.0 * b.real)
y = np.real(2.0 * b.imag)
z = np.real(2.0 * a - 1.0)
bb.add_vectors([x,y,z])
return
def PlotState(v,plotAll=True,myTitle=""):
"""
Plots the real and imaginary parts of a quantum state vector.
The x-axis is the index of the state in the list of 2^N, where N is the number of qubits.
The y-axis is the real and imaginary part of the state vector.
Args:
v (list): The quantum state vector to plot.
plotAll (bool): If True, plot all states. If False, plot only the first half.
myTitle (str): The title of the plot.
"""
if plotAll:
N=len(np.asarray(v))
else:
N=len(np.asarray(v))//2
vv=np.zeros(N,dtype=complex)
vv[:]=np.asarray(v)[0:N]
vv=vv*np.exp(-1.j*np.angle(vv[0]))
plt.axhline(0)
for idx,i in enumerate(vv):
plt.plot([int(idx),int(idx)],[0,np.real(i)],color='red')
plt.plot([int(idx)],[np.real(i)],'o',color='red')
for idx,i in enumerate(vv):
plt.plot([idx+0.1,idx+0.1],[0,np.imag(i)],color='blue')
plt.plot([idx+0.1],[np.imag(i)],'o',color='blue')
avg=np.average(vv)
plt.axhline(np.real(avg),linestyle='--')
plt.ylim(-1,1)
plt.title(myTitle)
plt.show()
def AddCircuits(theCircuits):
"""
Concatenates a list of circuits into a new circuit.
Args:
theCircuits (list): A list of QuantumCircuit objects to concatenate.
Returns:
QuantumCircuit: A new QuantumCircuit object that is the concatenation of the input circuits.
"""
numQubits=np.array([c.num_qubits for c in theCircuits])
numQubits=np.max(numQubits)
circuit=QuantumCircuit(numQubits,numQubits)
for i in range(0,len(theCircuits)):
circuit=circuit.compose(theCircuits[i],qubits=list(range(0,theCircuits[i].num_qubits)))
return circuit
def Mark(r,N):
"""
Creates a quantum circuit that marks a specific state |r> in the computational basis. (this is for Grover's algorithm)
Args:
r (int): The index of the state to mark in the computational basis.
N (int): The number of qubits in the circuit.
"""
circuit=QuantumCircuit(N,N)
circuit.barrier()
myString=np.binary_repr(r,width=N)[::-1]
for i in range(0,len(myString)):
if myString[i]=='0':
circuit.x(i)
circuit.barrier()
circuit.h(N-1)
circuit.mcx(list(range(0,N-1)), N-1,mode='noancilla')
circuit.h(N-1)
circuit.barrier()
for i in range(0,len(myString)):
if myString[i]=='0':
circuit.x(i)
circuit.barrier()
return circuit
def InitializeCircuitRandom(N):
"""
Creates a quantum circuit that initializes a random state vector of size 2^N.
Args:
N (int): The number of qubits in the circuit.
Returns:
QuantumCircuit: A QuantumCircuit object that initializes a random state vector.
"""
r=np.random.random(2**N)
r=r/np.linalg.norm(r)
circuit=QuantumCircuit(N,N)
circuit.initialize(list(r))
return circuit
def InitializeCircuit(my_string,N):
"""
Creates a quantum circuit that initializes a specific state vector.
Args:
my_string (str): A string of 0s and 1s representing the state vector. That is, all ups would be '00000' and all downs would be '11111'.
N (int): The number of qubits in the circuit.
"""
circuit = QuantumCircuit(N,N)
circuit.initialize(my_string[::-1])
return circuit
simulator = AerSimulator()
Exercise 1: A Single Qubit (33 points)#
a. A one qubit state#
A quantum state consists of a certain fraction of \(|0\rangle\) and a certain fraction of \(|1\rangle\) - i.e. \(\sqrt{0.3} |0\rangle - \sqrt{0.7} |1\rangle\).
We can plot this state in a variety of ways.
state=np.array([np.sqrt(0.3),-np.sqrt(0.7)])
state=MakeState(state)
PrintDirac(state)
PlotState(state)
b=qutip.Bloch()
PlotBloch(state,b)
b.render()
b.show()
Go ahead and try it out
### ANSWER HERE
Let’s go ahead and plot these three states
\(|0\rangle \equiv [1,0]\)
\(|1\rangle \equiv [0,1]\)
\(1/\sqrt{2} |0\rangle + 1/\sqrt{2} |1\rangle \equiv [1/\sqrt{2},1/\sqrt{2}]\)
Pay special attention to where those three states are on the Bloch Sphere
### ANSWER HERE
b. One qubit gates#
Quantum circuits are made out of quantum gates. Let us start by considering 1-qubit gates. The 1-qubit gate rotates states around the Bloch sphere. Three special gates are rx, ry, and rz which respectively rotate around the X, Y, and Z axis. These gates take an angle and a wire - i.e. rx(1/3*np.pi,0) will rotate \(1/3 \pi\) around the X axis.
To set up your circuit you can do
circuit=QuantumCircuit(1,1)
#Build your circuit here
state,measure = RunMe(circuit) #<--- This runs the circuit
# now you can plot the circuit
Go ahead and figure out how to use rx gate to move from a \(|0\rangle\) to a \(|1\rangle\) state.
### ANSWER HERE
Now we want to go ahead and use ry to go from \(|0\rangle \rightarrow 1/\sqrt{2}|0\rangle + 1/\sqrt{2}|1\rangle\). Write such a circuit. Then check to see what your circuit does to the state \(|1\rangle\) (you can get this using your other circuit from above).
### ANSWER HERE
So far we’ve seen we can get half \(|0\rangle\) and half \(|1\rangle\). We can also adjust states to get more of \(|1\rangle\) then \(|0\rangle\) (or visa versa). Generically we can get a state \(\cos(\theta) |0\rangle + \sin(\theta) e^{i\phi} |1\rangle\).
If we had such a state, after measurement we get “0” with probability \(\cos^2\theta\) and “1” with probability \(\sin^2\theta\).
Notice that \(\cos^2\theta + \sin^2\theta =1\) so we either get “0” or “1”.
To produce this state with \(\phi=0\), we can use the gate qc.ry(2*theta,wire) which takes \(|0\rangle \rightarrow \cos \theta |0\rangle + \sin \theta |1\rangle\)
Produce the state \(\cos(0.5)|0\rangle + \sin(0.5)|1\rangle\) and check to make sure it gives you the right amplitudes.
### ANSWER HERE
Now, make a circuit that rotates rx and rz both by \(0.3 \pi\).
### ANSWER HERE
Notice that by smartly choosing the angles around rx and rz and then rx again you can produce an arbitrary rotation on the sphere and therefore take \(|0\rangle\) and \(|1\rangle\) to any two orthogonal states.
c. Seeing the rotation#
We would like to plot the rotation around the bloch sphere. We can do this as follows:
b=qutip.Bloch()
for angle in np.arange(0,2*np.pi*0.8,0.1):
quantumWires=1
qc = QuantumCircuit(quantumWires,quantumWires)
# rotate around the x-axis by angle theta
vec,output=RunMe(qc)
PlotBloch(vec,b)
b.render()
b.show()
### ANSWER HERE
d. Measuring the state#
So far we’ve been cheating by looking at the quantum state. In the real world, you can’t do that. Instead you have to measure at the end of your circuit. The RunMe circuit returns two parameters. The second one is measurement outcomes - i.e.
vector,measure=RunMe(circuit)
can then plot the measurement outcomes by doing
plot_histgoram(measure)
Go ahead and measure the outcomes of a state rotated around the X axis by \(0.3\pi\). See that it gives you the expected raction of 0 and 1. Remember that the probability should be equal to the amplitude squared.
### ANSWER HERE
Exercise 2: Two qubits (34 points)#
In this exercise we will see how to build quantum states of two qubits.
a. Control-not gates#
The controlled-not gate: The key gate for two qubits is the control-not gate (circuit.cx(0,1)). The key gate “nots” the second wire (wire 1) if the first wire is “1” (wire 0). Let’s check it out. Apply the “control-not” to a state \(|00\rangle\) and to a state \(|10\rangle\) and print out the state (not the measurement outcome)
### ANSWER HERE
b. Build an EPR Pair#
The most interesting two qubit state is an EPR pair, \(\frac{1}{\sqrt{2}}\left(|00\rangle + |11\rangle\right)\). You can build it with a Hadamard and a CNOT. Go ahead and try different things and figure out how to build yourself an EPR pair. You can start out with
quantumWires=2
classicalWires=2
qc=QuantumCircuit(quantumWires,classicalWires)
to start out a quantum circuit with two wires. You can’t plot this on the bloch sphere but you can still plot it with PlotState. =
### ANSWER HERE
c. Measuring EPR Pairs#
Now that you’ve made an EPR pair, let’s go ahead and measure your EPR circuit above plotting the histogram of your measurements.
### ANSWER HERE
There is something very powerful here. If you get a “0” on the top wire, then you always get a “0” on the bottom wire. If you get a “1” on the top wire, then you get a “1” on bottom wire. This is even the case if these wires are taken miles apart before you measure.
Exercise 3: Grover’s Algorithm#
In this exercise, we are going to implement Grover’s algorithm. Grover’s algorithm solves the search problem
Grover’s algorithm takes a unitary which marks an element \(i\) and uses that as a subroutine with more circuit elements so that, after measurement, you sample the marked element \(i\).
To do this we are going to need to put together a lot of different circuit elements. Let’s start with the unitary which does the marking which we will give you. In this case it’s going to mark element 3 by taking it and changing its sign. Go ahead and run this to see explicitly what it does.
## Check Mark
circuit=InitializeCircuitRandom(5)
a,b=RunMe(circuit)
PlotState(a,True,myTitle="Initial Random State")
circuit=Mark(3,5)
a,b=RunMe(circuit,a)
PlotState(a,True,myTitle="Hopefully marked state")
### ANSWER HERE
a. Control-Control-Control Z#
Define a function controlZ(N) which generates the control-control-control-control Z on \(n\) qubits. The first \(n-1\) qubits have to be all 1 for the last Z qubit to happen.
To do this, implement the following circuit for general values of N. This is an example for controlZ(5).
You can get a control-control-control-control-x from qiskit by doing the following:
circuit.mcx(list(range(0,N-1)), N-1,mode='noancilla')
where \(N-1\) is the numer of target wires.
Test out your circuit by using the following code:
initCircuit=InitializeCircuit('11110',5)
#initCircuit=InitializeCircuit('11010',5)
initCircuit.h(4)
circuit=AddCircuits([initCircuit,ControlZ(5)])
print(circuit)
a,b=RunMe(circuit)
PrintDirac(a)
You should get something close to
(0.7071067811865475+8.659560562354932e-17j)|11110> + (-0.7071067811865476-8.659560562354934e-17j)|11111>
and then if you comment out the second initialization with 11010 you will get
(0.7071067811865476+0j)|11010> + (0.7071067811865475+0j)|11011>
### ANSWER HERE
b. All Hadamards#
Define a function AllHadamard(n) which builds \(n\) wires of Hadamards
Check it by the following test:
circuit=AllHadamard(5)
a,b=RunMe(circuit)
PlotState(a,True)
and seeing that it produces a state which is uniform. This is what we will use for the start of Grover’s algorithm.
### ANSWER HERE
c. FlipAllButZero#
Define FlipAllButZero(n). which should flip all but the zero configurations on \(n\) bits. It will us \(n+1\) wires (you need an ancilla, which is an extra bit that doesn’t get read out later. We will plot our state using PlotState(a,False) to ignore this extra bit). The example below is the output for FlipAllbutZero(4)
You might want to add some circuits together. To add these circuits you can do
AddCircuits([circuit1,circuit2,circuit3])
You can test it your result by doing
circuit=InitializeCircuitRandom(4)
a,b =RunMe(circuit)
PlotState(a,True,myTitle="Initial State")
circuit=FlipAllButZero(4)
a,b =RunMe(circuit,a)
PlotState(a,False,myTitle="Flipped Circuit")
You should find that everything but the zero flipped. It is ok if you are finding the zero flipped and nothing else. Quantum states aren’t defined up to a global phase and so the flipping of everything is equivalent to doing nothing.
### ANSWER HERE
d. Invert#
Define the Invert(N) function. This will also be on \(N+1\) wires. Here you need to add together the AllHadamard and FlipAllButZero and AllHadamard again.
You can test it in the following way:
circuit=InitializeCircuitRandom(4)
a,b=RunMe(circuit)
PlotState(a,True)
circuit=Invert(4)
a,b=RunMe(circuit,a)
PlotState(a,False)
This will generate a random state and then invert it around the mean (shown in a blue dotted line). Your mean before and after should be the same and everything else should flip around that dotted line. (Also everything totally flipping is still allowed)
### ANSWER HERE
e. Running Grover’s Algorithm#
Put everything together to run Grover’s algorithm.
### ANSWER HERE
Acknowledgement:
Bryan Clark (original)
Copyright: 2025