diff --git a/hamlib/qiskit/hamlib_simulation_benchmark.py b/hamlib/qiskit/hamlib_simulation_benchmark.py index bdf5d730..b7054e3d 100644 --- a/hamlib/qiskit/hamlib_simulation_benchmark.py +++ b/hamlib/qiskit/hamlib_simulation_benchmark.py @@ -38,6 +38,8 @@ verbose = False +# contains the correct bitstring for a random pauli circuit +bitstring_dict = {} # Creates a key for distribution of initial state for method = 3. def key_from_initial_state(num_qubits, num_shots, init_state, random_pauli_flag): @@ -150,19 +152,21 @@ def analyze_and_print_result( if verbose: print(f"... exact computation time = {round((time.time() - ts), 3)} sec") - - # for method 3, compute expected distribution from the initial state - elif method == 3: - - # check simple distribution if not inserting random Paulis - if not random_pauli_flag: - correct_dist = key_from_initial_state(num_qubits, num_shots, init_state, random_pauli_flag) - - # if using random paulis, distribution is based on all ones in initial state - # random_pauli_flag is (for now) set to always return bitstring of 1's - else: - correct_dist = {"1" * num_qubits: num_shots} + # for method 3, compute expected distribution from the initial state + elif method == 3: + + # check simple distribution if not inserting random Paulis + if not random_pauli_flag: + correct_dist = key_from_initial_state( + num_qubits, num_shots, init_state, random_pauli_flag + ) + + # if using random paulis, a potentially random bitstring is collected from circuit generation + else: + global bitstring_dict + correct_bitstring = bitstring_dict[qc.name] + correct_dist = {correct_bitstring: num_shots} else: raise ValueError("Method is not 1 or 2 or 3, or hamiltonian is not valid.") @@ -220,7 +224,9 @@ def print_top_measurements(label, counts, top_n): def run(min_qubits: int = 2, max_qubits: int = 8, max_circuits: int = 1, skip_qubits: int = 1, num_shots: int = 100, hamiltonian: str = "TFIM", method: int = 1, - random_pauli_flag: bool = False, init_state: str = None, + random_pauli_flag: bool = False, + random_init_flag: bool = False, + init_state: str = None, K: int = None, t: float = None, backend_id: str = None, provider_backend = None, hub: str = "ibm-q", group: str = "open", project: str = "main", exec_options = None, @@ -325,19 +331,26 @@ def execution_handler(qc, result, num_qubits, type, num_shots): ####################################################################### - # Loop over only the same circuit, executing it num_circuits times + # in the case of random paulis, method = 3: loop over multiple random pauli circuits + # otherwise, loop over the same circuit, executing it num_circuits times for circuit_id in range(num_circuits): + ts = time.time() - - # create the HamLibSimulation kernel and the associated Hamiltonian operator - qc, ham_op = HamiltonianSimulation( - num_qubits, - hamiltonian=hamiltonian, - K=K, t=t, - init_state=init_state, - method = method, - random_pauli_flag=random_pauli_flag - ) + + #used to store random pauli correct bitstrings + global bitstring_dict + + # create the HamLibSimulation kernel, random pauli bitstring, and the associated Hamiltonian operator + qc, bs, ham_op = HamiltonianSimulation( + num_qubits, + K=K, t=t, + hamiltonian=hamiltonian, + init_state=init_state, + method = method, + random_pauli_flag=random_pauli_flag, + random_init_flag=random_init_flag) + + bitstring_dict[qc.name] = bs metrics.store_metric(num_qubits, circuit_id, 'create_time', time.time() - ts) @@ -390,6 +403,7 @@ def get_args(): parser.add_argument("--global_rinst", "-param_rinst", default=None, help="paramater rinst") parser.add_argument("--num_steps", "-steps", default=None, help="Number of Trotter steps", type=int) parser.add_argument("--time", "-time", default=None, help="Time of evolution", type=float) + parser.add_argument("--random_init_flag", "-rani", action="store_true", help="random pauli flag") return parser.parse_args() # if main, execute method @@ -420,6 +434,7 @@ def get_args(): hamiltonian=args.hamiltonian, method=args.method, random_pauli_flag=args.random_pauli_flag, + random_init_flag=args.random_init_flag, init_state = args.init_state, K = args.num_steps, t = args.time, diff --git a/hamlib/qiskit/hamlib_simulation_kernel.py b/hamlib/qiskit/hamlib_simulation_kernel.py index 858357cb..ee6870a0 100644 --- a/hamlib/qiskit/hamlib_simulation_kernel.py +++ b/hamlib/qiskit/hamlib_simulation_kernel.py @@ -228,6 +228,7 @@ def create_circuit( method: int = 1, init_state: str = None, random_pauli_flag: bool = False, + random_init_flag: bool = False ): """ Create a quantum circuit based on the Hamiltonian data from an HDF5 file. @@ -292,6 +293,8 @@ def create_circuit( # Append K Trotter steps of inverse, if method 3 inv = None + bitstring = None + if method == 3: # if not adding random Paulis, just create simple inverse Trotter steps @@ -305,17 +308,23 @@ def create_circuit( # if adding Paulis, do that here, with code from pyGSTi else: from pygsti_mirror import convert_to_mirror_circuit - circuit, _ = convert_to_mirror_circuit(circuit_without_initial_state) - + + # if random init flag state is set, then discard the inital state input and use a completely (harr) random one + if random_init_flag: + circuit, bitstring = convert_to_mirror_circuit(circuit_without_initial_state, random_pauli = True, init_state=None) + else: + init_state = initial_state(n_spins, init_state) + circuit, bitstring = convert_to_mirror_circuit(circuit_without_initial_state, random_pauli = True, init_state=init_state) + # convert_to_mirror_circuit adds its own measurement gates if not random_pauli_flag: circuit.measure_all() - return circuit, ham_op, evo - + return circuit, bitstring, ham_op, evo + else: # print(f"Dataset not available for n_spins = {n_spins}.") - return None, None, None + return None, None, None, None ############### Initial Circuit Definition @@ -357,7 +366,8 @@ def HamiltonianSimulation( K: int = 5, t: float = 1.0, init_state = None, method: int = 1, - random_pauli_flag = False + random_pauli_flag = False, + random_init_flag = False ) -> QuantumCircuit: """ Construct a Qiskit circuit for Hamiltonian simulation. @@ -382,15 +392,15 @@ def HamiltonianSimulation( qc = QuantumCircuit(qr, cr, name=f"hamsim-{num_qubits}-{secret_int}") hamiltonian = hamiltonian.strip().lower() - - # create the quantum circuit for this Hamiltonian, along with the operator and trotter evolution circuit - qc, ham_op, evo = create_circuit( + # create the quantum circuit for this Hamiltonian, along with the correct pauli bstring, the operator and trotter evolution circuit + qc, bitstring, ham_op, evo = create_circuit( n_spins=n_spins, time=t, method=method, init_state=init_state, num_trotter_steps=K, random_pauli_flag=random_pauli_flag, + random_init_flag=random_init_flag ) # Save smaller circuit example for display @@ -403,10 +413,12 @@ def HamiltonianSimulation( # Collapse the sub-circuits used in this benchmark (for Qiskit) qc2 = qc.decompose().decompose() - - # return both the circuit created and the Hamiltonian operator - return qc2, ham_op - + + # return both the circuit created, the bitstring, and the Hamiltonian operator + # if random_pauli_flag is false or method isn't 3, bitstring will be None + return qc2, bitstring, ham_op + + ############### Circuit Drawer @@ -442,4 +454,3 @@ def kernel_draw(hamiltonian: str = "hamlib", method: int = 1): print(" ... circuit too large!") - \ No newline at end of file diff --git a/hamlib/qiskit/pygsti_mirror.py b/hamlib/qiskit/pygsti_mirror.py index 43ee9ded..88344af6 100644 --- a/hamlib/qiskit/pygsti_mirror.py +++ b/hamlib/qiskit/pygsti_mirror.py @@ -3,18 +3,20 @@ import numpy as np import scipy as sp from qiskit import transpile, QuantumCircuit +import re # Below code created by Timothy Proctor, July 16 2024, using code developed by various members of Sandia's QPL. +# Code has been edited to accept random pauli bitstrings by Sonny Rappaport, July 23 2024. +# Additional functions below also by Sonny based on Tim's code. -def sample_mirror_circuits(c, num_mcs=100, randomized_state_preparation=True, rand_state=None): +def sample_mirror_circuits(c, num_mcs=100, randomized_state_preparation=True, rand_state = None, random_pauli = True): if rand_state is None: rand_state = np.random.RandomState() - + mirror_circuits = [] mirror_circuits_bitstrings = [] for j in range(num_mcs): - mc, bs = central_pauli_mirror_circuit(c, randomized_state_preparation=randomized_state_preparation, - rand_state=rand_state) + mc, bs = central_pauli_mirror_circuit(c, randomized_state_preparation=randomized_state_preparation,random_pauli = True, rand_state = rand_state) mirror_circuits.append(mc) mirror_circuits_bitstrings.append(bs) @@ -36,7 +38,7 @@ def sample_reference_circuits(qubits, num_ref=100, randomized_state_preparation= return ref_circuits, ref_circuits_bitstrings -def central_pauli_mirror_circuit(circ, randomized_state_preparation=True, rand_state=None): +def central_pauli_mirror_circuit(circ, randomized_state_preparation=True, random_pauli = True, rand_state = None): if rand_state is None: rand_state = np.random.RandomState() @@ -50,10 +52,13 @@ def central_pauli_mirror_circuit(circ, randomized_state_preparation=True, rand_s n = circuit_to_mirror.width d = circuit_to_mirror.depth - # SONNY'S EDIT TO TIM CODE: MAKE THE CENTRAL PAULI (for now) NOT RANDOM - # Note that we could rework the code to pass in a desired bitstring here later, rather than randomly generating one here. - #central_pauli = 2 * rand_state.randint(0, 2, 2*n) - central_pauli = 2 * np.ones(2*n, dtype = np.int8 ) + # Sonny's edit here: Rather than generating central_pauli randomly, now we can use a non-random one if desired. + if random_pauli is True: + rand_state = np.random.RandomState() + central_pauli = 2 * np.random.randint(0, 2, 2*n) #old code + else: + central_pauli = 2 * np.ones(2*n, dtype = np.int8) + central_pauli_layer = pauli_vector_to_u3_layer(central_pauli, qubits) q = central_pauli.copy() @@ -125,8 +130,8 @@ def haar_random_u3(q, rand_state=None): if rand_state is None: rand_state = np.random.RandomState() - a, b = 2 * np.pi * rand_state.rand(2) - theta = mod_2pi(2 * np.arcsin(np.sqrt(rand_state.rand(1)))[0]) + a, b = 2 * np.pi * np.random.rand(2) + theta = mod_2pi(2 * np.arcsin(np.sqrt(np.random.rand(1)))[0]) phi = mod_2pi(a - b + np.pi) lamb = mod_2pi(-1 * (a + b + np.pi)) return pygsti.baseobjs.Label('Gu3', q, args=(theta, phi, lamb)) @@ -325,23 +330,36 @@ def pygsti_string(qc): return pygstr.strip() -def convert_to_mirror_circuit(qc): +def convert_to_mirror_circuit(qc, random_pauli, init_state): """ Takes an arbitrary quantum circuit, transpiles it to CX and u3, and eventually returns a mirror circuit version of the quantum circuit. The mirror circuit will have a random initial state in addition to random paulis. """ + # if no initial state is provided, use a ranmdom one + random_init_flag = True + + if init_state: + random_init_flag = False + init_state.append(qc, range(qc.num_qubits)) + qc = init_state + qc = transpile(qc, basis_gates=["cx","u3"]) pygstr = pygsti_string(qc) pygsti_circuit = pygsti.circuits.Circuit(pygstr) - mcs, bss = sample_mirror_circuits(pygsti_circuit, num_mcs=1) - print(bss) + mcs, bss = sample_mirror_circuits(pygsti_circuit, num_mcs=1, random_pauli = random_pauli, randomized_state_preparation=random_init_flag) + + qasm_string = mcs[0].convert_to_openqasm() + + # use regex magic to remove lines that begin with Delay gates + delays_removed_string = re.sub(r'^delay.*\n?', '', qasm_string, flags=re.MULTILINE) - qiskit_circuit = QuantumCircuit.from_qasm_str(mcs[0].convert_to_openqasm()) + qiskit_circuit = QuantumCircuit.from_qasm_str(delays_removed_string) # this circuit is made out of u3 and cx gates by pygsti default - return qiskit_circuit , bss[0] + # the bitstring is reversed to account for qiskit ordering + return qiskit_circuit , bss[0][::-1]