#Libraries needed to implement and simulate quantum circuits
from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, transpile, qasm2, qasm3
from qiskit_aer import Aer
from qiskit.primitives import BackendSampler
#Custem functions to simplify answers
import Our_Qiskit_Functions as oq #a part of the library presented in arXiv:1903.04359v1.

import numpy as np
import math as m
import scipy as sci

#Initialize backends simulators to visualize circuits
S_simulator = Aer.backends(name='statevector_simulator')[0]
M_simulator = Aer.backends(name='qasm_simulator')[0]

def XOR(a,b):
    if( (a+b)==1 ):
        return 1
    else:
        return 0

def OR(a,b):
    if( (a+b)>=1 ):
        return 1
    else:
        return 0

def AND(a,b):
    if( (a+b)==2 ):
        return 1
    else:
        return 0
    
#====================================
a = [1,1,0,0]
b = [0,1,1,1]
C = 0
Sum = []
for i in range( len(a) ):
    xor1 = XOR(a[0-i-1],b[0-i-1])
    S = XOR(xor1,C)
    and1 = AND(a[0-i-1],b[0-i-1])
    and2 = AND(xor1,C)
    C = OR(and1,and2)
    Sum.insert(0,S)
if( C == 1 ):
    Sum.insert(0,C)

print('Adding a = ',a,' + b = ',b,' ===> ',Sum)

q = QuantumRegister(3,name='q')
qc= QuantumCircuit(q,name='qc')
#------------------------------
qc.h( q[0] )
qc.h( q[1] )
qc.h( q[2] )
qc.p( 3*m.pi/2, q[0] )
qc.p( m.pi, q[1] )
qc.p( 0, q[2] )
oq.Wavefunction(qc);

q = QuantumRegister(3,name='q')
qc= QuantumCircuit(q,name='qc')
#------------------------------
qc.x( q[0] )
qc.x( q[1] )
oq.QFT(qc,q,3)
oq.Wavefunction(qc);

qa = QuantumRegister(2,name='a')
qb = QuantumRegister(2,name='b')
qc = QuantumCircuit(qa,qb,name='qc')
qc.x( qa[1] )
qc.x( qb[0] )
oq.Wavefunction(qc,systems=[2,2]);

qa = QuantumRegister(2,name='a')
qb = QuantumRegister(2,name='b')
qc = QuantumCircuit(qa,qb,name='qc')
#------------------------------------
qc.x( qa[1] )
qc.x( qb[0] )
print('_____ States to Add Together _____')
oq.Wavefunction(qc,systems=[2,2])
oq.QFT(qc,qa,2)
#----------------------------- phase contributions from |b>
qc.cp( m.pi, qb[0], qa[0] )
qc.cp( m.pi/2, qb[1], qa[0] )
qc.cp( m.pi, qb[1], qa[1] )
#-----------------------------
oq.QFT_dgr(qc,qa,2)
print('\n___ Sum Stored in |a> ___')
oq.Wavefunction(qc,systems=[2,2],show_systems=[True,False]);

A_States = [[0,1],[1,0],[1,1]]
B_States = [[1,0],[1,1]]
#------------------------------------
for a in range(len(A_States)):
    A = A_States[a]
    for b in range(len(B_States)):
        B = B_States[b]
        qa = QuantumRegister(2,name='a')
        qb = QuantumRegister(2,name='b')
        qc = QuantumCircuit(qa,qb,name='qc')
        #-----------------------------------
        if(A[0]==1):
            qc.x( qa[0] )
        if(A[1]==1):
            qc.x( qa[1] )
        if(B[0]==1):
            qc.x( qb[0] )
        if(B[1]==1):
            qc.x( qb[1] )
        oq.QFT(qc,qa,2)
        qc.cp( m.pi, qb[0], qa[0] )
        qc.cp( m.pi/2, qb[1], qa[0] )
        qc.cp( m.pi, qb[1], qa[1] )
        oq.QFT_dgr(qc,qa,2)
        print('\nA:',A,' B:',B)
        oq.Wavefunction(qc,systems=[2,2],show_systems=[True,False])


qa = QuantumRegister(3,name='a')
qb = QuantumRegister(2,name='b')
qc = QuantumCircuit(qa,qb,name='qc')
#-----------------------------------
qc.x( qa[1] )
qc.x( qb[0] )
qc.x( qb[1] )
print('_____ States to Add Together _____')
oq.Wavefunction(qc,systems=[3,2])
oq.QFT(qc,qa,3)
#------------------------------ phase contributions from |b>
qc.cp( m.pi/2, qb[0], qa[0] )
qc.cp( m.pi/4, qb[1], qa[0] )
qc.cp( m.pi, qb[0], qa[1] )
qc.cp( m.pi/2, qb[1], qa[1] )
qc.cp( m.pi, qb[1], qa[2] )
#------------------------------
oq.QFT_dgr(qc,qa,3)
print('\n___ Sum Stored in |a> ___')
oq.Wavefunction(qc,systems=[3,2],show_systems=[True,False]);

A = [0,1,1,0]
B = [1,1,0,1] #A and B need to be arrays of equal length (don't include the extra 0 qubit for A)
print('States to Add Together: ',A,' + ',B)
#=========================================
qa = QuantumRegister(len(A)+1,name='a')
qb = QuantumRegister(len(B),name='b')
qc = QuantumCircuit(qa,qb,name='qc')
#--------------------------------------
oq.Quantum_Adder(qc,qa,qb,A,B)
print('\n___ Sum Stored in |a> ___')
oq.Wavefunction(qc,systems=[len(A)+1,len(B)],show_systems=[True,False]);

qa = QuantumRegister(3,name='a')
qb = QuantumRegister(3,name='b')
qc = QuantumCircuit(qa,qb,name='qc')
#-----------------------------------
qc.x( qa[0] )
qc.x( qa[1] )
qc.x( qb[1] )
qc.x( qb[2] )
print('____ States to Subtract ____')
oq.Wavefunction(qc,systems=[3,3])
oq.QFT(qc,qa,3)
#------------------------------ phase contributions from |b>
qc.cp( -m.pi, qb[0], qa[0] )
qc.cp( -m.pi/2, qb[1], qa[0] )
qc.cp( -m.pi/4, qb[2], qa[0] )
qc.cp( -m.pi, qb[1], qa[1] )
qc.cp( -m.pi/2, qb[2], qa[1] )
qc.cp( -m.pi, qb[2], qa[2] )
#------------------------------
oq.QFT_dgr(qc,qa,3)
print('\n___ Difference Stored in |a> ___')
oq.Wavefunction(qc,systems=[3,3],show_systems=[True,False]);

