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'))
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)

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

# 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)

