<jupyter_start><jupyter_text>TODO: look at pg 14 of: https://arxiv.org/pdf/quant-ph/0406176v5.pdfThis is what Qiskit useshttps://quantumcomputing.stackexchange.com/questions/6755/controlled-initialize-instruction<jupyter_code>from qiskit.extensions import Initialize
import numpy as np
def Single_qubit_rotation(zero_state_amp, one_state_amp):
    """
    Get rotation to create passed in qubit state
    
    |psi> = cos(theta/2) |0> + e^{-1j * phi/2} sin(theta/2) |1>
    
    See equation 8 and 9 of https://arxiv.org/pdf/quant-ph/0406176v5.pdf
    
    """
    zero_state_amp = complex(zero_state_amp)
    one_state_amp = complex(one_state_amp)
    
    norm = np.sqrt(zero_state_amp**2 + one_state_amp**2 )
    
    if np.isclose(norm, 0):
        theta = 0
        a_arg = 0
        b_arg = 0
        final_t = 0
        phi = 0
    else:
        theta = 2 * np.arccos(np.abs(zero_state_amp) / norm)
        a_arg = np.angle(zero_state_amp)
        b_arg = np.angle(one_state_amp)
        final_t = a_arg + b_arg
        phi = b_arg - a_arg
    
    return norm * np.exp(1.J * final_t / 2), theta, phi

state,theta, phi= Single_qubit_rotation(1/np.sqrt(2), 1/np.sqrt(2))
np.cos(theta/2) 
np.exp(-1j*phi/2)*np.sin(theta/2)
np.exp(-1j*phi/2)*np.sin(theta/2)
 Single_qubit_rotation(1/np.sqrt(2), 1/np.sqrt(2))
Initialize._bloch_angles((1/np.sqrt(2), 1/np.sqrt(2)))
def Rotations_to_disentangle(qubit_state_vector):
    """
    pg 11 of https://arxiv.org/pdf/quant-ph/0406176v5.pdf
    
    Method to work out Ry and Rz rotation angles to disentangle the least significant bit (LSB).
    These rotations make up a block diagonal matrix U == multiplexor
    
    
    
    ### futher background:
    
    Given |ψ> an (n+1) qubit state, seperate the state into a separable (unentanged) state by the following circuit:
    
               n  :──\\───(C)───────────────(C)────────────
   |ψ>                     │                 │             
               n+1:─────── Rz (-phi) ─────── Ry(-theta)──── |ψ''>


    the 2^{n+1} state vector is split into TWO 2^{n} states... This can be done by the circuit above.
    
    Overall |ψ> (2^{n+1} state vector) is split into 2^{n} contiguous 2-element blocks. This can be interpreted as a
    2D complex vector. We can label this |ψ_{c}>
    
    Then:
    
    Rz(-φ_{c}) Ry(-θ_{c}) |ψ> = r_{c} exp(1j*t_{c}) |0 >
    
    |ψ''> is the n-qubit state given by the 2^{n}-element row vector with c-th entry r_{c} exp(1j*t_{c}).
    
    If we let U be the block diagonal sum ⊕_{c} Ry(-θ_{c}) Rz(-φ_{c}). THEN:
    
                U |ψ> = |ψ''> |0>
    
    We can implement U as a multiplexed Rz gate followed by a multiplexed Ry gate!
    
    """
    
    remaining_vector = []
    theta_list = []
    phi_list = []
    
    param_len = len(qubit_state_vector)
    for state_ind in range(param_len//2):
            # Ry and Rz rotations to move bloch vector from 0 to "imaginary"
            # qubit
            # (imagine a qubit state signified by the amplitudes at index 2*i
            # and 2*(i+1), corresponding to the select qubits of the
            # multiplexor being in state |i>)
            
            amp_2i = qubit_state_vector[2*state_ind] # amp at 2i
            amp_2i_2= qubit_state_vector[(2*state_ind)+1] #  amp at 2(i+1)
            remaining_qubit_state_vector, theta, phi = Single_qubit_rotation(amp_2i, amp_2i_2)
            
            remaining_vector.append(remaining_qubit_state_vector)
            theta_list.append(-1*theta)
            phi_list.append(-1*phi)
    
    return remaining_vector, theta_list, phi_list
state = [np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25)]
Rotations_to_disentangle(state)

Initialize._rotations_to_disentangle(state)
import cirq
C = cirq.Circuit()
C.append(cirq.rx(np.pi/2).controlled(num_controls=1, control_values=[1]).on(cirq.LineQubit(1), cirq.LineQubit(9)))
C
def recursive_multiplex(target_gate, list_of_angles, start_qubit_num, end_qubit_num, last_cnot=True):
    """
    Args:
        target_gate (Gate): Ry or Rz gate to apply to target qubit,
                            multiplexed over all other "select" qubits
                            
        list_of_angles (list[float]): list of rotation angles to apply Ry and Rz
        
        last_cnot (bool): add the last cnot if last_cnot = True
    """
    number_angles = len(list_of_angles)
    local_num_qubits = int(np.log2(number_angles)) + 1 # +1 for n+1 qubits!
    
    qubits_list = cirq.LineQubit.range(start_qubit_num, end_qubit_num)
    
    LSB = qubits_list[0] # least significant bit
    MSB = qubits_list[local_num_qubits-1] # most significant bit
    
    circuit = cirq.Circuit()
    
    # case of no multiplexing: base case for recursion
    if local_num_qubits == 1:
        if target_gate == 'Ry':
            Ry_gate = cirq.ry(list_of_angles[0])
            circuit.append(Ry_gate.on(LSB))
        elif target_gate == 'Rz':
            Rz_gate = cirq.rz(list_of_angles[0])
            circuit.append(Rz_gate.on(LSB))
        else:
            raise ValueError(f'Incorrect gate specificed: {target_gate}')
        
        return circuit
    
    angle_weight = np.kron([[0.5, 0.5], [0.5, -0.5]],
                               np.identity(2 ** (local_num_qubits - 2)))
    
    # calc the combo angles
    list_of_angles = angle_weight.dot(np.array(list_of_angles)).tolist()
    
    # recursive step on half the angles fulfilling the above assumption
    multiplex_1 = recursive_multiplex(target_gate, list_of_angles[0:(number_angles // 2)],
                                      start_qubit_num,
                                      end_qubit_num-1,
                                      False)
    circuit = cirq.Circuit(
       [
           circuit.all_operations(),
           *multiplex_1.all_operations(),
       ]
    )
    
    circuit.append(cirq.CNOT(MSB, LSB))

    # implement extra efficiency from the paper of cancelling adjacent
    # CNOTs (by leaving out last CNOT and reversing (NOT inverting) the
    # second lower-level multiplex)
    multiplex_2 = recursive_multiplex(target_gate, list_of_angles[(number_angles // 2):],
                                      start_qubit_num,
                                      end_qubit_num-1,
                                      False)
    
    if number_angles > 1:
        circuit = cirq.Circuit(
                               [
                                   circuit.all_operations(),
                                   *list(multiplex_2.all_operations())[::-1],
                               ]
                            )
    else:
        circuit = cirq.Circuit(
                               [
                                   circuit.all_operations(),
                                   *multiplex_2.all_operations(),
                               ]
                            )
    # attach a final CNOT
    if last_cnot:
        circuit.append(cirq.CNOT(MSB, LSB))
    
    return circuit
recursive_multiplex( 'Ry', [np.pi/2, np.pi/2, np.pi/2, np.pi/4],0,3, last_cnot=True)
from qiskit.circuit.library.standard_gates.ry import RYGate
qiskit_test= Initialize(state)
qiskit_test._multiplex(RYGate, [np.pi/2, np.pi/2, np.pi/2, np.pi/4], last_cnot=True).decompose().decompose().draw()
# from copy import deepcopy
# def disentangle_circuit(qubit_state_vector, n_qubits):
#     """
#     """
#     circuit = cirq.Circuit()

#     remaining_vector = deepcopy(qubit_state_vector)
#     for qubit_ind in range(n_qubits):
#         # work out which rotations must be done to disentangle the LSB
#         # qubit (we peel away one qubit at a time)
#         remaining_vector, theta_list, phi_list = Rotations_to_disentangle(remaining_vector)
        
        
#         add_last_cnot = True
#         if np.linalg.norm(phi_list) != 0 and np.linalg.norm(theta_list) != 0:
#             add_last_cnot = False

#         if np.linalg.norm(phi_list) != 0:
#             rz_mult_circuit = recursive_multiplex('Rz',
#                                                   phi_list,
#                                                   qubit_ind,
#                                                   n_qubits,
#                                                   last_cnot=add_last_cnot)
#             circuit.append(rz_mult_circuit)

#         if np.linalg.norm(theta_list) != 0:
#             ry_mult_circuit = recursive_multiplex('Ry',
#                                                   theta_list,
#                                                   qubit_ind,
#                                                   n_qubits,
#                                                   last_cnot=add_last_cnot)
#             circuit = cirq.Circuit(
#                        [
#                            circuit.all_operations(),
#                            *list(ry_mult_circuit.all_operations())[::-1],
#                        ]
#                     )
#     return circuit
from copy import deepcopy
def disentangle_circuit(qubit_state_vector, start_qubit_ind, end_qubit_ind):
    """
    """

    circuit = cirq.Circuit()
    n_qubits = np.log2(len(qubit_state_vector))
    end_ind_corr = end_qubit_ind+1

    if n_qubits!= len(list(range(start_qubit_ind, end_ind_corr))):
        raise ValueError('incorrect qubit defined qubit indices!')

    remaining_vector = deepcopy(qubit_state_vector)
    for qubit_ind in range(start_qubit_ind, end_qubit_ind+1):
        # work out which rotations must be done to disentangle the LSB
        # qubit (we peel away one qubit at a time)
        remaining_vector, theta_list, phi_list = Rotations_to_disentangle(remaining_vector)
        
        
        add_last_cnot = True
        if np.linalg.norm(phi_list) != 0 and np.linalg.norm(theta_list) != 0:
            add_last_cnot = False

        if np.linalg.norm(phi_list) != 0:
            rz_mult_circuit = recursive_multiplex('Rz',
                                                  phi_list,
                                                  qubit_ind,
                                                  start_qubit_ind+end_ind_corr,
                                                  last_cnot=add_last_cnot)
            circuit.append(rz_mult_circuit)

        if np.linalg.norm(theta_list) != 0:
            ry_mult_circuit = recursive_multiplex('Ry',
                                                  theta_list,
                                                  qubit_ind,
                                                  start_qubit_ind+end_ind_corr,
                                                  last_cnot=add_last_cnot)
            circuit = cirq.Circuit(
                       [
                           circuit.all_operations(),
                           *list(ry_mult_circuit.all_operations())[::-1],
                       ]
                    )
    return circuit
state = [np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25)]
n_qubits=2

test = disentangle_circuit(state, 0, 1)
test
qiskit_test.gates_to_uncompute().decompose().decompose().draw()
test.unitary()[:,0]
# def intialization_circuit(qubit_state_vector):
    
#     n_qubits = np.log2(len(qubit_state_vector))
    
#     if np.ceil(n_qubits) != np.floor(n_qubits):
#         raise ValueError('state vector is not a qubit state')
    
#     if not np.isclose(sum(np.abs(qubit_state_vector)**2), 1):
#         raise ValueError('state vector is not normalized')
    
#     n_qubits = int(np.log2(len(qubit_state_vector)))
    
#     disentangling_circuit = disentangle_circuit(qubit_state_vector, n_qubits)
#     inverse_circuit = cirq.inverse(disentangling_circuit)
    
#     return inverse_circuit
    
def intialization_circuit(qubit_state_vector, start_qubit_ind, end_qubit_ind):
    """
    """

    n_qubits = np.log2(len(qubit_state_vector))
    end_ind_corr = end_qubit_ind+1
    
    if np.ceil(n_qubits) != np.floor(n_qubits):
        raise ValueError('state vector is not a qubit state')
    
    if not np.isclose(sum(np.abs(qubit_state_vector)**2), 1):
        raise ValueError('state vector is not normalized')
    
    
    disentangling_circuit = disentangle_circuit(qubit_state_vector, start_qubit_ind, end_qubit_ind)
    inverse_circuit = cirq.inverse(disentangling_circuit)
    
    return inverse_circuit
    
# state = [np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25),np.sqrt(0.25)]
# state = [np.sqrt(0.9),np.sqrt(0.1), 0 , 0, 0 , 0, 0 , 0]
state = [np.sqrt(0.9),np.sqrt(0.05), 0 , 0, 0 , 0, 0 , np.sqrt(0.05)]
circuit = intialization_circuit(state, 0, 2)
circuit
# Initialize Simulator
s=cirq.Simulator()

print('Simulate the circuit:')
results=s.simulate(circuit)
print(results)
print()

results.state_vector()
qiskit_test= Initialize(state)
qiskit_test._define_synthesis().decompose().decompose().decompose().decompose().draw()
circuit.unitary()[:,0]
np.around(circuit.unitary()[:,0], 3)
circuit.unitary()[0,:]
# remove reset part of circuit!
import qiskit

qiskit_circuit = qiskit_test._define_synthesis()


qcirc = qiskit_circuit.decompose().reverse_bits() #### <--- note reversing bits here!!! (flips order from top to bottom!)

# need to remove reset part of circuit
new_data = []
for index, tup in enumerate(qcirc.data):
    op_type, _, _ = tup
    if isinstance(op_type, qiskit.circuit.reset.Reset):
        continue
    else:
        new_data.append(tup)
qcirc.data = new_data
qcirc.decompose().draw()
from qiskit.quantum_info import Operator
matrix = Operator(qcirc).data
np.allclose(circuit.unitary(), matrix)
qiskit_circuit.decompose().decompose().decompose().decompose().draw()
circuit
cirq.inverse(circuit)
out = np.around(Operator(CIRC).data[:,0], 6)

from quchem.Qcircuit.misc_quantum_circuit_functions import Get_state_as_str

n_qubits = int(np.log2(len(out)))
print('N_qubits =', n_qubits)
for i, amp in enumerate(out):
    print(Get_state_as_str(n_qubits, i), amp)<jupyter_output>N_qubits = 3
000 (0.948683-0j)
001 (0.223607-0j)
010 (-0+0j)
011 -0j
100 -0j
101 -0j
110 0j
111 (0.223607-0j)
{'fraction_non_alphanumeric': 0.11119856887298747, 'fraction_numerical': 0.021323792486583185, 'mean_word_length': 2.020531661984007, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 20, 'https://': 4, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 44, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '576398', 'license': 'MIT', 'n_tokens_mistral': 4820, 'n_tokens_neox': 4514, 'n_words': 923, 'path': 'quchem_examples/(4.2  hacky notes on Prepare arbitrary State via Circuit SECOND APPROACH).ipynb', 'repo_name': 'AlexisRalli/VQE-code', 'size': '13975'}