New algorithm for RB building Cliffords by layers

[merav-aharoni]

Summary

We implement a new algorithm for RB, 1 and 2 qubits.

Details and comments

Here are the main ideas behind the implementation:

1. For 1-qubits Cliffords, we transpile all the Cliffords in advance and store them in a list. When we add a Clifford to an RB circuit, we select it out of this list.
2. For 2-qubit Cliffords, we use the partition of the Cliffords into 'layers' as shown at the bottom of this summary. We store the transpiled Cliffords for each of the layers. When we add a Clifford to an RB circuit, we select one Clifford from each layer-list, and compose them together.
3. Every Clifford is mapped by a number. We store a list of the compositions of Cliffords represented as numbers. For example, if Clifford1.compose(Clifford2) == Clifford3, then we conceptually, we store {(1, 2) : 3}.
4. Similarly, we store for each number representing a Cliffords, the number representing the inverse Clifford. This data, along with the data in (2) is stored in the file qiskit_experiments/library/randomized_benchmarking/clifford_data.py. The data is generated by the file, in the same directory, create_clifford_map.py. clifford_data.py is constant. As a rule, create_clifford_map.py does not need to be executed unless some overall change is needed.
5. For (2) above, we don't actually store the map (as a dict), but only the results of the compose in a list. This is more efficient in performance. The index of the required Clifford is computed using the indices of the input Cliffords. Similarly, for (3) we also only store the list of inverse numbers.
6. For the compose map, we don't actually store the full compose table of all-cliffords X all-cliffords. Instead, we define an array of single-gate-cliffords. This comprises all Cliffords that consist of a single gate. These arrays are stored in clifford_data.py as well. There are 9 such Cliffords for 1-qubit, and 21 such Cliffords for 2 qubits. It is sufficient to store the compose table of all-cliffords X single-gate-cliffords, since for every Clifford on the right hand side, we can break it down into single gate Cliffords, and do the composition one at a time. This greatly reduces the storage space for the array of composition results (from O(n^2) to O(n)).
7. The only case where a user might want to change clifford_data.py is if they wish to add a new single-gate Clifford, which is not equivalent to any of the already supported gates. In this case, one needs to add the gate to one of the gate lists in create_clifford_map.py.
8. I changed the class CliffordUtils, so that it now has members for num_qubits, basis_gates and 'backend`. If several RB experiments are run in sequence, e.g., parallel, the class will be updated only if these parameters have changed. This improved performance for cases where the parameters have not changed.
9. I split the code for [1,2] qubits and >3 qubits on the highest level - i.e., in circuits(). The reason for this is that in addition to the changes listed above, I made a change in the high level flow. In the existing code - _sample_circuits (now for 3+ qubits), first we created sequences of Cliffords and then create circuits out of the sequences. In the code for 1,2 qubits - _build_circuits, I do not create the sequences. I immediately create the circuits, skipping the sequences step.

[merav-aharoni]

This PR is now ready for review - @nkanazawa1989 , @ShellyGarion , @gadial .

[merav-aharoni]

A question to the reviewers – @itoko (in his recent PR) and I each took a different approach in the split to <2 and >2 qubits. I made the split on the level of circuits (see point (9) above, because my entire algorithm is different than the existing one. @itoko made the split inside the low-level methods, assuming the main flow for the two cases would remain the same. As a result, in the current implementation, in all such methods (for example, _to_instruction), the number of qubits will always be >2.
I see two options to resolve this:

1. Clean up the unnecessary checks for the number of qubits where not needed.
2. Leave the code as is, and in a future PR, re-unite the two flows at high level. Then we can use Itoko’s approach. If we go for this approach, I would not like to include this in the current PR, because it would require significant changes.

I welcome your inputs on this.

[nkanazawa1989]

Thanks @merav-aharoni . My main concerns are the (1) CliffordUtil tied to basis gates and (2) enforced layer approach. To reduce variance of RB, it would be nice to decompose a single Clifford gate always with 3 CNOTs. However we are no longer able to obtain original Clifford circuit and thus I cannot transpile it by myself by talking .circuits(). Apart from this, I and @eliarbel are okey with moving forward. Let's see what we can get from the synergy of you and @itoko 's PRs.

Here are my comments:

For 1-qubits Cliffords, we transpile all the Cliffords in advance and store them in a list. When we add a Clifford to an RB circuit, we select it out of this list.

Since 1q Clifford contains small number of elements, perhaps it could be faster to generate them on the fly with cache. Please confirm in your PR @itoko

For 2-qubit Cliffords, we use the partition of the Cliffords into 'layers' as shown at the bottom of this summary. We store the transpiled Cliffords for each of the layers. When we add a Clifford to an RB circuit, we select one Clifford from each layer-list, and compose them together.

This is really nice technique to reduce cache data size. Great job!

Every Clifford is mapped by a number. We store a list of the compositions of Cliffords represented as numbers. For example, if Clifford1.compose(Clifford2) == Clifford3, then we conceptually, we store {(1, 2) : 3}.

This is also reasonable.

Similarly, we store for each number representing a Cliffords, the number representing the inverse Clifford. This data, along with the data in (2) is stored in the file qiskit_experiments/library/randomized_benchmarking/clifford_data.py. The data is generated by the file, in the same directory, create_clifford_map.py. clifford_data.py is constant. As a rule, create_clifford_map.py does not need to be executed unless some overall change is needed.

The file size is about 1.4 MB and this is quite heavy. You should switch to another compressed data format, otherwise you could include script to generate this data and run the script at the time of setup. Please do some followup @itoko

For (2) above, we don't actually store the map (as a dict), but only the results of the compose in a list. This is more efficient in performance. The index of the required Clifford is computed using the indices of the input Cliffords. Similarly, for (3) we also only store the list of inverse numbers.

Dict is O(1) lookup and sounds like fast enough and I'm curious how it's further boosted. I'm fine with your implementation as long as performance is guaranteed.

For the compose map, we don't actually store the full compose table of all-cliffords X all-cliffords. Instead, we define an array of single-gate-cliffords. This comprises all Cliffords that consist of a single gate. These arrays are stored in clifford_data.py as well. There are 9 such Cliffords for 1-qubit, and 21 such Cliffords for 2 qubits. It is sufficient to store the compose table of all-cliffords X single-gate-cliffords, since for every Clifford on the right hand side, we can break it down into single gate Cliffords, and do the composition one at a time. This greatly reduces the storage space for the array of composition results (from O(n^2) to O(n)).

This is also nice. We need to write user/developer guide at some point to describe this data structure for future reviewers.

The only case where a user might want to change clifford_data.py is if they wish to add a new single-gate Clifford, which is not equivalent to any of the already supported gates. In this case, one needs to add the gate to one of the gate lists in create_clifford_map.py.

User should be able to run RB experiment with new basis gates without regenerating file (i.e. try-except) as long as they don't care about the speed. They may want to choose to wait for another dozens of seconds with conventional approach, rather than spending time for learning how to use of this generator function and waiting another minutes to regenerate full Clifford table. This is why I prefer a single interface @itoko designed, and we can hide actual logic behind the interface as a black box. For example, forest benchmarking has a dedicated generator class that takes RB configuration with basis gates (likely they generates sequence through Cloud service). With this, users and casual developers don't need to know actual implementation. Proposed implementation exposes everything to the main class and researchers may need to read handles of lines of code to figure out how to make small update. Personally I prefer the implementation of forest benchmarking. I'm fine with promoting CliffordUtil to a generator class in this sense. It would be great if you could do some followup @itoko

I changed the class CliffordUtils, so that it now has members for num_qubits, basis_gates and 'backend`. If several RB experiments are run in sequence, e.g., parallel, the class will be updated only if these parameters have changed. This improved performance for cases where the parameters have not changed.

I feel this is no longer util. I was thinking to remove this class since its usage is just a name space. You can still explicitly add basis gates to the arguments of methods/functions that requires the basis information. On the other hand, I think we can promote this class to a generator as I mentioned above.

I split the code for [1,2] qubits and >3 qubits on the highest level - i.e., in circuits(). The reason for this is that in addition to the changes listed above, I made a change in the high level flow. In the existing code - _sample_circuits (now for 3+ qubits), first we created sequences of Cliffords and then create circuits out of the sequences. In the code for 1,2 qubits - _build_circuits, I do not create the sequences. I immediately create the circuits, skipping the sequences step.

This design may increases maintenance overhead. I feel we need extra abstraction layer for circuit generation. Please do followup @itoko

Overall the performance improvement is really nice and I'm happy with your awesome work! I'll just merge this PR now, and I hope @itoko will address my comments above with his rework without performance regression.

[itoko]

Great job, @merav-aharoni ! This PR has dramatically improved the performance of the RB experiments.
Actually, it achieves about 7.7x speed-up for _transpiled_circuits() of 1Q StandardRB
(5.7x for 2Q StandardRB, 3.8x for 1Q InterleavedRB and 5.3x for 2Q InterleavedRB),
measuring the performance with a benchmark: https://gist.github.com/itoko/f078218ea6458f32d1c0f9be827f614f
See tables below for the details of benchmark result.

Note: The benchmark measures run time without "setup time" (import time and caching time in circuits()), in that sense, the resulting speedup may be a bit optimistic. As far as I checked the setup time with another code, it seems negligible. (I observed a few seconds loading time sometimes, but they were due to slow loading of FakeBackend objects, so not relevant with the RB performance).

[itoko]

main: Before this PR (0359f4c)
new: After this PR (0b15ccb)

method = time_transpiled_circuits for both settings

Standard RB 1Q:

backend	max_length	num_samples	main[s]	new[s]	main->new
fake_manila	1000	3	6.858569	0.851461	0.124146
fake_manila	1000	6	12.919110	1.660500	0.128531
fake_manila	2000	3	12.545703	1.642936	0.130956
fake_manila	2000	6	25.421245	3.265286	0.128447
fake_mumbai	1000	3	6.696955	0.868488	0.129684
fake_mumbai	1000	6	13.016662	1.817998	0.139667
fake_mumbai	2000	3	12.516374	1.767419	0.141209
fake_mumbai	2000	6	24.813701	3.380879	0.136250
fake_washington	1000	3	10.513716	0.965256	0.091809
fake_washington	1000	6	13.176188	1.944821	0.147601
fake_washington	2000	3	12.700285	1.831261	0.144191
fake_washington	2000	6	25.937748	3.700574	0.142671

Standard RB 2Q:

backend	max_length	num_samples	main[s]	new[s]	main->new
fake_manila	1000	3	18.058308	4.673694	0.258811
fake_manila	1000	6	35.537921	9.215669	0.259319
fake_manila	2000	3	35.208998	9.280712	0.263589
fake_manila	2000	6	71.094951	18.541894	0.260805
fake_mumbai	1000	3	17.977917	4.591728	0.255409
fake_mumbai	1000	6	35.188215	9.284733	0.263859
fake_mumbai	2000	3	34.964976	8.912482	0.254897
fake_mumbai	2000	6	71.040783	18.087106	0.254602
fake_washington	1000	3	22.380718	5.125869	0.229031
fake_washington	1000	6	36.875958	10.234201	0.277530
fake_washington	2000	3	36.186682	10.308951	0.284882
fake_washington	2000	6	71.671835	19.716134	0.275089

Interleaved RB 1Q:

backend	max_length	num_samples	main[s]	new[s]	main->new
fake_manila	1000	3	18.058308	4.673694	0.258811
fake_manila	1000	6	35.537921	9.215669	0.259319
fake_manila	2000	3	35.208998	9.280712	0.263589
fake_manila	2000	6	71.094951	18.541894	0.260805
fake_mumbai	1000	3	17.977917	4.591728	0.255409
fake_mumbai	1000	6	35.188215	9.284733	0.263859
fake_mumbai	2000	3	34.964976	8.912482	0.254897
fake_mumbai	2000	6	71.040783	18.087106	0.254602
fake_washington	1000	3	22.380718	5.125869	0.229031
fake_washington	1000	6	36.875958	10.234201	0.277530
fake_washington	2000	3	36.186682	10.308951	0.284882
fake_washington	2000	6	71.671835	19.716134	0.275089

Interleaved RB 2Q:

backend	max_length	num_samples	main[s]	new[s]	main->new
fake_manila	100	3	3.804730	0.663383	0.174357
fake_manila	100	6	7.130793	1.335523	0.187290
fake_manila	200	3	6.238839	1.277523	0.204769
fake_manila	200	6	12.598263	2.683764	0.213026
fake_mumbai	100	3	3.748988	0.694320	0.185202
fake_mumbai	100	6	7.109178	1.358674	0.191115
fake_mumbai	200	3	6.365602	1.334885	0.209703
fake_mumbai	200	6	12.608916	2.633852	0.208888
fake_washington	100	3	7.953887	0.819290	0.103005
fake_washington	100	6	7.819094	1.608669	0.205736
fake_washington	200	3	6.671117	1.493850	0.223928
fake_washington	200	6	13.213095	2.938669	0.222406