from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi

# the angles of rotations
theta1 = pi/4
theta2 = pi/6

# the circuit with two qubits
qreg = QuantumRegister(2) 
creg = ClassicalRegister(2) 
mycircuit = QuantumCircuit(qreg,creg)


# when the second qubit is in |0>, the first qubit is rotated by theta1
mycircuit.x(qreg[1])
mycircuit.cu3(2*theta1,0,0,qreg[1],qreg[0])
mycircuit.x(qreg[1])

# when the second qubit is in |1>, the first qubit is rotated by theta2
mycircuit.cu3(2*theta2,0,0,qreg[1],qreg[0])
              
# we read the unitary matrix
job = execute(mycircuit,Aer.get_backend('unitary_simulator'),optimization_level=0)
u=job.result().get_unitary(mycircuit,decimals=3)
# we print the unitary matrix in nice format
for i in range(len(u)):
    s=""
    for j in range(len(u)):
        val = str(u[i][j].real)
        while(len(val)<8): val  = " "+val
        s = s + val
    print(s)
    



from math import pi, sin, cos

theta1 = pi/4
theta2 = pi/6

print(round(cos(theta1),3),-round(sin(theta1),3),0,0)
print(round(sin(theta1),3),-round(cos(theta1),3),0,0)
print(0,0,round(cos(theta2),3),-round(sin(theta2),3))
print(0,0,round(sin(theta2),3),-round(cos(theta2),3))

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi, sin, cos

# the angle of rotation
theta1 = pi/23
theta2 = 2*pi/23
theta3 = 4*pi/23


precision = 3

print("a1 = theta3 => sin(a1) = ",round(sin(theta3),precision))
print("a2 = theta2+theta3 => sin(a2) = ",round(sin(theta2+theta3),precision))
print("a3 = theta1 => sin(a3) = ",round(sin(theta1),precision))
print("a4 = theta1+theta2 => sin(a4) = ",round(sin(theta1+theta2),precision))
print()

qreg = QuantumRegister(3) 
creg = ClassicalRegister(3) 
circuit = QuantumCircuit(qreg,creg)

# controlled rotation when the third qubit is |1>
circuit.cu3(2*theta1,0,0,qreg[2],qreg[0])

# controlled rotation when the second qubit is |1>
circuit.cu3(2*theta2,0,0,qreg[1],qreg[0])

# controlled rotation when the third qubit is |0>
circuit.x(qreg[2])
circuit.cu3(2*theta3,0,0,qreg[2],qreg[0])
circuit.x(qreg[2])

# read the corresponding unitary matrix
job = execute(circuit,Aer.get_backend('unitary_simulator'),optimization_level=0)
unitary_matrix=job.result().get_unitary(circuit,decimals=precision)
for i in range(len(unitary_matrix)):
    s=""
    for j in range(len(unitary_matrix)):
        val = str(unitary_matrix[i][j].real)
        while(len(val)<precision+4): val  = " "+val
        s = s + val
    print(s)

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi, sin, cos
from random import randrange

# the angle of rotation
k1 = randrange(1,31)
theta1 = k1*2*pi/31
k2 = randrange(1,31)
theta2 = k2*2*pi/31
k3 = randrange(1,31)
theta3 = k3*2*pi/31

max_percentange = 0
    
# for each stream of length from 1 to 31
for i in range(1,32):    
    # initialize the circuit
    qreg = QuantumRegister(3) 
    creg = ClassicalRegister(3) 
    circuit = QuantumCircuit(qreg,creg)
    
    # Hadamard operators before reading the stream
    for m in range(3):
        circuit.h(qreg[m])   
        
    # read the stream of length i
    print("stream of length",i,"is being read")
    for j in range(i):    
        # controlled rotation when the third qubit is |1>
        circuit.cu3(2*theta1,0,0,qreg[2],qreg[0])

        # controlled rotation when the second qubit is |1>
        circuit.cu3(2*theta2,0,0,qreg[1],qreg[0])

        # controlled rotation when the third qubit is |0>
        circuit.x(qreg[2])
        circuit.cu3(2*theta3,0,0,qreg[2],qreg[0])
        circuit.x(qreg[2])
        
    # Hadamard operators after reading the stream
    for m in range(3):
        circuit.h(qreg[m])  
    # we measure after reading the whole stream
    circuit.measure(qreg,creg)
    # execute the circuit N times
    N = 1000
    job = execute(circuit,Aer.get_backend('qasm_simulator'),shots=N)
    counts = job.result().get_counts(circuit)
    print(counts)
    if '000' in counts.keys():
        c = counts['000']
    else:
        c = 0
    print('000 is observed',c,'times out of',N)
    percentange = round(c/N*100,1)
    if max_percentange < percentange and i != 31: max_percentange = percentange
    print("the ration of 000 is ",percentange,"%")
    print()  
print("maximum percentage of observing unwanted '000' is",max_percentange)

#
# your solution is here
#


from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi

# initialize the circuit
qreg = QuantumRegister(4)  
circuit = QuantumCircuit(qreg)

# we use the fourth qubit as the auxiliary

# apply a rotation to the first qubit when the third and second qubits are in states |0> and |1>
# change the state of the third qubit to |1>
circuit.x(qreg[2])
# if both the third and second qubits are in states |1>, the state of auxiliary qubit is changed to |1> 
circuit.ccx(qreg[2],qreg[1],qreg[3])
# the rotation is applied to the first qubit if the state of auxiliary qubit is |1>
circuit.cu3(2*pi/6,0,0,qreg[3],qreg[0])
# reverse the effects
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.x(qreg[2])

circuit.draw()

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi,sin

# the angles of rotations
theta1 = pi/10
theta2 = 2*pi/10
theta3 = 3*pi/10
theta4 = 4*pi/10

# for verification, print sin(theta)'s
print("sin(theta1) = ",round(sin(theta1),3))
print("sin(theta2) = ",round(sin(theta2),3))
print("sin(theta3) = ",round(sin(theta3),3))
print("sin(theta4) = ",round(sin(theta4),3))
print()

qreg =  QuantumRegister(4) 
circuit = QuantumCircuit(qreg)

# the third qubit is in |0>
# the second qubit is in |0>
circuit.x(qreg[2])
circuit.x(qreg[1])
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.cu3(2*theta1,0,0,qreg[3],qreg[0])
# reverse the effects
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.x(qreg[1])
circuit.x(qreg[2])


# the third qubit is in |0>
# the second qubit is in |1>
circuit.x(qreg[2])
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.cu3(2*theta2,0,0,qreg[3],qreg[0])
# reverse the effects
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.x(qreg[2])

# the third qubit is in |1>
# the second qubit is in |0>
circuit.x(qreg[1])
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.cu3(2*theta3,0,0,qreg[3],qreg[0])
# reverse the effects
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.x(qreg[1])

# the third qubit is in |1>
# the second qubit is in |1>
circuit.ccx(qreg[2],qreg[1],qreg[3])
circuit.cu3(2*theta4,0,0,qreg[3],qreg[0])
# reverse the effects
circuit.ccx(qreg[2],qreg[1],qreg[3])

# read the corresponding unitary matrix
job = execute(circuit,Aer.get_backend('unitary_simulator'),optimization_level=0)
unitary_matrix=job.result().get_unitary(circuit,decimals=3)
for i in range(len(unitary_matrix)):
    s=""
    for j in range(len(unitary_matrix)):
        val = str(unitary_matrix[i][j].real)
        while(len(val)<7): val  = " "+val
        s = s + val
    print(s)

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi,sin
from random import randrange

# the angle of rotation
k1 = randrange(1,61)
theta1 = k1*2*pi/61
k2 = randrange(1,61)
theta2 = k2*2*pi/61
k3 = randrange(1,61)
theta3 = k3*2*pi/61
k4 = randrange(1,61)
theta4 = k4*2*pi/61

max_percentange = 0

# for each stream of length of 1, 11, 21, 31, 41, 51, and 61
for i in [1,11,21,31,41,51,61]: 
#for i in range(1,62): 
    # initialize the circuit
    qreg =  QuantumRegister(4) 
    creg = ClassicalRegister(4)
    circuit = QuantumCircuit(qreg,creg)

    # Hadamard operators before reading the stream
    for m in range(3):
        circuit.h(qreg[m])   
        
    # read the stream of length i
    print("stream of length",i,"is being read")
    for j in range(i):         
        # the third qubit is in |0>
        # the second qubit is in |0>
        circuit.x(qreg[2])
        circuit.x(qreg[1])
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.cu3(2*theta1,0,0,qreg[3],qreg[0])
        # reverse the effects
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.x(qreg[1])
        circuit.x(qreg[2])


        # the third qubit is in |0>
        # the second qubit is in |1>
        circuit.x(qreg[2])
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.cu3(2*theta2,0,0,qreg[3],qreg[0])
        # reverse the effects
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.x(qreg[2])

        # the third qubit is in |1>
        # the second qubit is in |0>
        circuit.x(qreg[1])
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.cu3(2*theta3,0,0,qreg[3],qreg[0])
        # reverse the effects
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.x(qreg[1])

        # the third qubit is in |1>
        # the second qubit is in |1>
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        circuit.cu3(2*theta4,0,0,qreg[3],qreg[0])
        # reverse the effects
        circuit.ccx(qreg[2],qreg[1],qreg[3])
        
    # Hadamard operators after reading the stream
    for m in range(3):
        circuit.h(qreg[m])  
    # we measure after reading the whole stream
    circuit.measure(qreg,creg)
    # execute the circuit N times
    N = 1000
    job = execute(circuit,Aer.get_backend('qasm_simulator'),shots=N)
    counts = job.result().get_counts(circuit)
    print(counts)
    if '0000' in counts.keys():
        c = counts['0000']
    else:
        c = 0
    print('0000 is observed',c,'times out of',N)
    percentange = round(c/N*100,1)
    if max_percentange < percentange and i != 61: max_percentange = percentange
    print("the ration of 0000 is ",percentange,"%")
    print()  
print("maximum percentage of observing unwanted '0000' is",max_percentange)

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi

qreg11 =  QuantumRegister(4)
creg11 = ClassicalRegister(4)

theta = pi/4

# define our quantum circuit
mycircuit11 = QuantumCircuit(qreg11,creg11)

def ccc_ry(angle,q1,q2,q3,q4):
    mycircuit11.cu3(angle/2,0,0,q3,q4)
    mycircuit11.ccx(q1,q2,q4)
    mycircuit11.cu3(-angle/2,0,0,q3,q4)
    mycircuit11.ccx(q1,q2,q4)
    
ccc_ry(2*theta,qreg11[3],qreg11[2],qreg11[1],qreg11[0])

mycircuit11.draw(output='mpl')

from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, execute, Aer
from math import pi

qreg12 =  QuantumRegister(4)
creg12 = ClassicalRegister(4)

theta1 = pi/16
theta2 = 2*pi/16
theta3 = 3*pi/16
theta4 = 4*pi/16
theta5 = 5*pi/16
theta6 = 6*pi/16
theta7 = 7*pi/16
theta8 = 8*pi/16

# define our quantum circuit
mycircuit12 = QuantumCircuit(qreg12,creg12)

def ccc_ry(angle,q1,q2,q3,q4):
    mycircuit12.cu3(angle/2,0,0,q3,q4)
    mycircuit12.ccx(q1,q2,q4)
    mycircuit12.cu3(-angle/2,0,0,q3,q4)
    mycircuit12.ccx(q1,q2,q4)
    
mycircuit12.x(qreg12[3])
mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[1])
ccc_ry(2*theta1,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
mycircuit12.x(qreg12[1])
mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[3])

mycircuit12.x(qreg12[3])
mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[1])
ccc_ry(2*theta2,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
#mycircuit12.x(qreg12[1])
mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[3])

mycircuit12.x(qreg12[3])
#mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[1])
ccc_ry(2*theta3,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
mycircuit12.x(qreg12[1])
#mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[3])

mycircuit12.x(qreg12[3])
#mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[1])
ccc_ry(2*theta4,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
#mycircuit12.x(qreg12[1])
#mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[3])

#mycircuit12.x(qreg12[3])
mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[1])
ccc_ry(2*theta5,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
mycircuit12.x(qreg12[1])
mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[3])

#mycircuit12.x(qreg12[3])
mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[1])
ccc_ry(2*theta6,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
#mycircuit12.x(qreg12[1])
mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[3])

#mycircuit12.x(qreg12[3])
#mycircuit12.x(qreg12[2])
mycircuit12.x(qreg12[1])
ccc_ry(2*theta7,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
mycircuit12.x(qreg12[1])
#mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[3])

#mycircuit12.x(qreg12[3])
#mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[1])
ccc_ry(2*theta8,qreg12[3],qreg12[2],qreg12[1],qreg12[0])
#mycircuit12.x(qreg12[1])
#mycircuit12.x(qreg12[2])
#mycircuit12.x(qreg12[3])

job = execute(mycircuit12,Aer.get_backend('unitary_simulator'),optimization_level=0)
u=job.result().get_unitary(mycircuit12,decimals=3)
for i in range(len(u)):
    s=""
    for j in range(len(u)):
        val = str(u[i][j].real)
        while(len(val)<7): val  = " "+val
        s = s + val
    print(s)