diff --git a/qiskit/providers/aer/pulse/controllers/digest_pulse_qobj.py b/qiskit/providers/aer/pulse/controllers/digest_pulse_qobj.py index 71bdd5e5a9..35d6484db9 100644 --- a/qiskit/providers/aer/pulse/controllers/digest_pulse_qobj.py +++ b/qiskit/providers/aer/pulse/controllers/digest_pulse_qobj.py @@ -269,8 +269,7 @@ def format_pulse_samples(pulse_samples): return [[samp.real, samp.imag] for samp in new_samples] -def experiment_to_structs(experiment, ham_chans, pulse_inds, - pulse_to_int, dt, qubit_list=None): +def experiment_to_structs(experiment, ham_chans, pulse_inds, pulse_to_int, dt, qubit_list=None): """Converts an experiment to a better formatted structure Args: @@ -334,9 +333,22 @@ def experiment_to_structs(experiment, ham_chans, pulse_inds, structs['channels'][chan_name][0].extend([inst['t0'] * dt, None, index, cond]) pv_needs_tf[ham_chans[chan_name]] = 1 - # Frame changes + # ShiftPhase instructions elif inst['name'] == 'fc': - structs['channels'][chan_name][1].extend([inst['t0'] * dt, inst['phase'], cond]) + # get current phase value + current_phase = 0 + if len(structs['channels'][chan_name][1]) > 0: + current_phase = structs['channels'][chan_name][1][-2] + + structs['channels'][chan_name][1].extend([inst['t0'] * dt, + current_phase + inst['phase'], + cond]) + + # SetPhase instruction + elif inst['name'] == 'setp': + structs['channels'][chan_name][1].extend([inst['t0'] * dt, + inst['phase'], + cond]) # A standard pulse else: diff --git a/src/open_pulse/numeric_integrator.cpp b/src/open_pulse/numeric_integrator.cpp index 53fd204a8a..c0d398450d 100644 --- a/src/open_pulse/numeric_integrator.cpp +++ b/src/open_pulse/numeric_integrator.cpp @@ -69,24 +69,21 @@ complex_t chan_value( // TODO floating point comparsion with complex ?! // Seems like this is equivalent to: out != complex_t(0., 0.) if(out != 0.){ - double phase = 0.; num_times = floor_div(fc_array.shape[0], 3); - for(auto i = 0; i < num_times; ++i){ - // TODO floating point comparison - if(t >= fc_array[3 * i]){ - bool do_fc = true; - if(fc_array[3 * i + 2] >= 0){ - if(!reg[static_cast(fc_array[3 * i + 2])]){ - do_fc = false; - } - } - if(do_fc){ - phase += fc_array[3 * i + 1]; - } - }else{ - break; - } + + // get the index of the phase change + // this loop will result in finding the index of the phase to use +1 + auto phase_idx = 0; + while(phase_idx < num_times){ + if(t < fc_array[3 * phase_idx]) break; + phase_idx++; } + + double phase = 0.; + if(phase_idx > 0){ + phase = fc_array[3 * (phase_idx - 1) + 1]; + } + if(phase != 0.){ out *= std::exp(complex_t(0., 1.) * phase); } diff --git a/test/terra/backends/test_pulse_simulator.py b/test/terra/backends/test_pulse_simulator.py index b7e1e9e1d0..e21e597645 100644 --- a/test/terra/backends/test_pulse_simulator.py +++ b/test/terra/backends/test_pulse_simulator.py @@ -21,14 +21,13 @@ import numpy as np from scipy.linalg import expm from scipy.special import erf -from scipy.linalg.blas import get_blas_funcs from qiskit.providers.aer.backends import PulseSimulator from qiskit.compiler import assemble from qiskit.quantum_info import state_fidelity -from qiskit.pulse import (Schedule, Play, ShiftPhase, Acquire, SamplePulse, DriveChannel, - ControlChannel, AcquireChannel, MemorySlot) +from qiskit.pulse import (Schedule, Play, ShiftPhase, SetPhase, Delay, Acquire, SamplePulse, + DriveChannel, ControlChannel, AcquireChannel, MemorySlot) from qiskit.providers.aer.pulse.de.DE_Methods import ScipyODE from qiskit.providers.aer.pulse.de.DE_Options import DE_Options from qiskit.providers.aer.pulse.system_models.pulse_system_model import PulseSystemModel @@ -38,8 +37,6 @@ from .pulse_sim_independent import (simulate_1q_model, simulate_2q_exchange_model, simulate_3d_oscillator_model) -dznrm2 = get_blas_funcs("znrm2", dtype=np.float64) - class TestPulseSimulator(common.QiskitAerTestCase): r"""PulseSimulator tests.""" @@ -144,7 +141,7 @@ def test_x_gate_rwa(self): meas_map=[[0]], qubit_lo_freq=[omega_d], memory_slots=1, - shots=256) + shots=1) # set backend backend_options including initial state y0 = np.array([1.0, 0.0]) @@ -556,7 +553,7 @@ def test_2Q_interaction(self): meas_map=[[0]], qubit_lo_freq=[omega_d0, omega_d1], memory_slots=2, - shots=1.) + shots=1) y0 = np.kron(np.array([1., 0.]), np.array([0., 1.])) backend_options = {'seed' : 9000, 'initial_state': y0} @@ -603,7 +600,7 @@ def test_subsystem_restriction(self): meas_map=[[0]], qubit_lo_freq=[omega_d, omega_d, omega_d], memory_slots=2, - shots=1.) + shots=1) y0 = np.kron(np.array([1., 0.]), np.array([0., 1.])) backend_options = {'seed' : 9000, 'initial_state': y0} @@ -637,7 +634,7 @@ def test_subsystem_restriction(self): meas_map=[[0]], qubit_lo_freq=[omega_d, omega_d, omega_d], memory_slots=2, - shots=1.) + shots=1) y0 = np.kron(np.array([1., 0.]), np.array([0., 1.])) backend_options = {'seed' : 9000, 'initial_state': y0} @@ -792,7 +789,7 @@ def test_gaussian_drive(self): meas_map=[[0]], qubit_lo_freq=[omega_d], memory_slots=2, - shots=1000) + shots=1) y0 = np.array([1., 0.]) backend_options = {'seed' : 9000, 'initial_state' : y0} @@ -876,6 +873,250 @@ def test_2Q_exchange(self): # Check fidelity of statevectors self.assertGreaterEqual(state_fidelity(pulse_sim_yf, yf), 1-(10**-5)) + def test_delay_instruction(self): + """Test for delay instruction.""" + + # construct system model specifically for this + hamiltonian = {} + hamiltonian['h_str'] = ['0.5*r*X0||D0', '0.5*r*Y0||D1'] + hamiltonian['vars'] = {'r': np.pi} + hamiltonian['qub'] = {'0': 2} + ham_model = HamiltonianModel.from_dict(hamiltonian) + + u_channel_lo = [] + subsystem_list = [0] + dt = 1. + + system_model = PulseSystemModel(hamiltonian=ham_model, + u_channel_lo=u_channel_lo, + subsystem_list=subsystem_list, + dt=dt) + + # construct a schedule that should result in a unitary -Z if delays are correctly handled + # i.e. do a pi rotation about x, sandwiched by pi/2 rotations about y in opposite directions + # so that the x rotation is transformed into a z rotation. + # if delays are not handled correctly this process should fail + sched = Schedule() + sched += Play(SamplePulse([0.5]), DriveChannel(1)) + sched += Delay(1, DriveChannel(1)) + sched += Play(SamplePulse([-0.5]), DriveChannel(1)) + + sched += Delay(1, DriveChannel(0)) + sched += Play(SamplePulse([1.]), DriveChannel(0)) + + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[0., 0.], + memory_slots=2, + shots=1) + + + # Result of schedule should be the unitary -1j*Z, so check rotation of an X eigenstate + backend_options = {'initial_state': np.array([1., 1.]) / np.sqrt(2)} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + + statevector = results.get_statevector() + expected_vector = np.array([-1j, 1j]) / np.sqrt(2) + + self.assertGreaterEqual(state_fidelity(statevector, expected_vector), 1 - (10**-5)) + + # verify validity of simulation when no delays included + sched = Schedule() + sched += Play(SamplePulse([0.5]), DriveChannel(1)) + sched += Play(SamplePulse([-0.5]), DriveChannel(1)) + + sched += Play(SamplePulse([1.]), DriveChannel(0)) + + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[0., 0.], + memory_slots=2, + shots=1) + + backend_options = {'initial_state': np.array([1., 1.]) / np.sqrt(2)} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + + statevector = results.get_statevector() + U = expm(1j * np.pi * self.Y /4) @ expm(-1j * np.pi * (self.Y / 4 + self.X / 2)) + expected_vector = U @ np.array([1., 1.]) / np.sqrt(2) + + self.assertGreaterEqual(state_fidelity(statevector, expected_vector), 1 - (10**-5)) + + def test_shift_phase(self): + """Test ShiftPhase command.""" + + omega_0 = 1.123 + r = 1. + + system_model = self._system_model_1Q(omega_0, r) + + # run a schedule in which a shifted phase causes a pulse to cancel itself. + # Also do it in multiple phase shifts to test accumulation + sched = Schedule() + amp1 = 0.12 + sched += Play(SamplePulse([amp1]), DriveChannel(0)) + phi1 = 0.12374 * np.pi + sched += ShiftPhase(phi1, DriveChannel(0)) + amp2 = 0.492 + sched += Play(SamplePulse([amp2]), DriveChannel(0)) + phi2 = 0.5839 * np.pi + sched += ShiftPhase(phi2, DriveChannel(0)) + amp3 = 0.12 + 0.21 * 1j + sched += Play(SamplePulse([amp3]), DriveChannel(0)) + + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[omega_0], + memory_slots=2, + shots=1) + + y0 = np.array([1., 0]) + backend_options = {'initial_state': y0} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + pulse_sim_yf = results.get_statevector() + + #run independent simulation + samples = np.array([[amp1], + [amp2 * np.exp(1j * phi1)], + [amp3 * np.exp(1j * (phi1 + phi2))]]) + indep_yf = simulate_1q_model(y0, omega_0, r, np.array([omega_0]), samples, 1.) + + self.assertGreaterEqual(state_fidelity(pulse_sim_yf, indep_yf), 1 - (10**-5)) + + # run another schedule with only a single shift phase to verify + sched = Schedule() + amp1 = 0.12 + sched += Play(SamplePulse([amp1]), DriveChannel(0)) + phi1 = 0.12374 * np.pi + sched += ShiftPhase(phi1, DriveChannel(0)) + amp2 = 0.492 + sched += Play(SamplePulse([amp2]), DriveChannel(0)) + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[omega_0], + memory_slots=2, + shots=1) + + y0 = np.array([1., 0]) + backend_options = {'initial_state': y0} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + pulse_sim_yf = results.get_statevector() + + #run independent simulation + samples = np.array([[amp1], [amp2 * np.exp(1j * phi1)]]) + indep_yf = simulate_1q_model(y0, omega_0, r, np.array([omega_0]), samples, 1.) + + self.assertGreaterEqual(state_fidelity(pulse_sim_yf, indep_yf), 1 - (10**-5)) + + def test_set_phase(self): + """Test SetPhase command. Similar to the ShiftPhase test but includes a mixing of + ShiftPhase and SetPhase instructions to test relative vs absolute changes""" + + omega_0 = 1.3981 + r = 1. + + system_model = self._system_model_1Q(omega_0, r) + + # intermix shift and set phase instructions to verify absolute v.s. relative changes + sched = Schedule() + amp1 = 0.12 + sched += Play(SamplePulse([amp1]), DriveChannel(0)) + phi1 = 0.12374 * np.pi + sched += ShiftPhase(phi1, DriveChannel(0)) + amp2 = 0.492 + sched += Play(SamplePulse([amp2]), DriveChannel(0)) + phi2 = 0.5839 * np.pi + sched += SetPhase(phi2, DriveChannel(0)) + amp3 = 0.12 + 0.21 * 1j + sched += Play(SamplePulse([amp3]), DriveChannel(0)) + phi3 = 0.1 * np.pi + sched += ShiftPhase(phi3, DriveChannel(0)) + amp4 = 0.2 + 0.3 * 1j + sched += Play(SamplePulse([amp4]), DriveChannel(0)) + + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[omega_0], + memory_slots=2, + shots=1) + + y0 = np.array([1., 0.]) + backend_options = {'initial_state': y0} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + pulse_sim_yf = results.get_statevector() + + #run independent simulation + samples = np.array([[amp1], + [amp2 * np.exp(1j * phi1)], + [amp3 * np.exp(1j * phi2)], + [amp4 * np.exp(1j * (phi2 + phi3))]]) + indep_yf = simulate_1q_model(y0, omega_0, r, np.array([omega_0]), samples, 1.) + + self.assertGreaterEqual(state_fidelity(pulse_sim_yf, indep_yf), 1 - (10**-5)) + + def test_set_phase_rwa(self): + """Test SetPhase command using an RWA approximate solution.""" + omega_0 = 5.123 + r = 0.01 + + system_model = self._system_model_1Q(omega_0, r) + + sched = Schedule() + sched += SetPhase(np.pi / 2, DriveChannel(0)) + sched += Play(SamplePulse(np.ones(100)), DriveChannel(0)) + + sched |= Acquire(1, AcquireChannel(0), MemorySlot(0)) << sched.duration + + qobj = assemble([sched], + backend=self.backend_sim, + meas_level=2, + meas_return='single', + meas_map=[[0]], + qubit_lo_freq=[omega_0], + memory_slots=2, + shots=1) + + y0 = np.array([1., 1.]) / np.sqrt(2) + backend_options = {'initial_state': y0} + + results = self.backend_sim.run(qobj, system_model, backend_options).result() + pulse_sim_yf = results.get_statevector() + + #run independent simulation + phases = np.exp((-1j * 2 * np.pi * omega_0 * np.array([1, -1]) / 2) * 100) + approx_yf = phases * (expm(-1j * (np.pi / 2) * self.Y) @ y0) + + self.assertGreaterEqual(state_fidelity(pulse_sim_yf, approx_yf), 0.99) def _system_model_1Q(self, omega_0, r): """Constructs a standard model for a 1 qubit system.