`calculate_unitary` raises error for one-qubit circuits on non-`q0` register

[ryanhill1]

Describe the bug Calling calculate_unitary on a one-qubit circuit where the register is not q0 raises an error.

To reproduce

import numpy as np
from braket.circuits import Circuit, Instruction, gates
from braket.circuits.unitary_calculation import calculate_unitary

qubit_index = 1
circuit = Circuit()
pgates = [
    gates.Rx,
    gates.Ry,
    gates.Rz,
    gates.PhaseShift,
]
angles = np.random.RandomState(11).random(len(pgates))
instructions = [
    Instruction(rot(a), target=qubit_index)
    for rot, a in zip(pgates, angles)
]
for instr in instructions:
    circuit.add_instruction(instr)

unitary = calculate_unitary(circuit.qubit_count, circuit.instructions)

The above works fine for qubit_index=0 but raises an IndexError for any qubit_index!=0

Expected behavior The calculate_unitary function should account for the empty/unused registers and calculate the unitary of the circuit.

Screenshots or logs

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
/var/folders/nm/_31bsg9s76v7_cxbkkw9f20w0000gn/T/ipykernel_35500/857523477.py in <module>
----> 1 unitary = calculate_unitary(circuit.qubit_count, circuit.instructions)

~/opt/miniconda3/envs/braket-env/lib/python3.9/site-packages/braket/circuits/unitary_calculation.py in calculate_unitary(qubit_count, instructions)
     71         targets = instr.target
     72 
---> 73         gate_indexes, un_indexes, result_indexes = _einsum_subscripts(targets, qubit_count)
     74         gate_matrix = np.reshape(matrix, len(targets) * [2, 2])
     75 

~/opt/miniconda3/envs/braket-env/lib/python3.9/site-packages/braket/circuits/unitary_calculation.py in _einsum_subscripts(targets, qubit_count)
     28     un_right_indexes = list(range(target_count + qubit_count, target_count + 2 * qubit_count))
     29 
---> 30     gate_right_indexes = [un_left_indexes[-1 - target] for target in targets]
     31 
     32     result_left_indexes = un_left_indexes.copy()

~/opt/miniconda3/envs/braket-env/lib/python3.9/site-packages/braket/circuits/unitary_calculation.py in <listcomp>(.0)
     28     un_right_indexes = list(range(target_count + qubit_count, target_count + 2 * qubit_count))
     29 
---> 30     gate_right_indexes = [un_left_indexes[-1 - target] for target in targets]
     31 
     32     result_left_indexes = un_left_indexes.copy()

IndexError: list index out of range

System information

• Amazon Braket Python SDK version: 1.11.1
• Python version: 3.9.7

[ryanhill1]

Update: The issue seems to persist for any circuit with an empty/unused register indexed within the range of the maximally indexed register. For example,

• Circuits for which calculate_unitary works:

  
  T  : |0|                      T  : |0|                      T  : |0|                               

q0 : -X- q0 : -X- q0 : -X-

q1 : -Y- q1 : -Y- T : |0|

q2 : -Z- T : |0|

T : |0|

- Circuits for which `calculate_unitary` *does not* work
```python3
T  : |0|                      T  : |0|                      T  : |0|                      T  : |0|                                

q1 : -Y-                      q0 : -X-                      q1 : -Y-                      q2 : -Z-

q2 : -Z-                      q2 : -Z-                      T  : |0|                      T  : |0| 

T  : |0|                      T  : |0|

[krneta]

Hi @ryanhill1,

Thank you very much for reaching out!

We're aware of this limitation and our documentation states "When constructing a circuit with the simulators, Amazon Braket currently requires that you use contiguous qubits/indices."

Could you please tell us more about why this functionality (having empty/unused register index) is important to your use case?

[ryanhill1]

@krneta thank you for your quick response! I'm sorry I wasn't already aware of this limitation. I stumbled upon it while working on the unitaryfund/mitiq error mitigation toolkit. Their project interfaces between a few different front-end circuit-building packages. They verify each circuit conversion by comparing the circuit unitary before and after conversion, but for some packages, don't require qubit equality. I was investigating why that was, and now it makes sense! This functionality isn't overly important to me, but I'm curious, does it have to do with a limitation of the backends themselves? Or is it to do with how the circuits are currently processed through the SDK?

[krneta]

My pleasure.

It's just a simplification we made early on that we haven't had time to go back to. That's why I was curious about your use case - to see if we need to look into doing that work again.

[ryanhill1]

Gotcha. If you did decide to add support for circuits with non-contiguous qubits/indices, how do you think you would handle the unused registers? Suppose you started with the following circuit:

T  : |0|1|2|

q0 : -H-C---
        |   
q2 : ---X-C-
          | 
q4 : -----X-

T  : |0|1|2|

There seem to be two main approaches to handling the unused q1 and q3 registers in this example.

1. You can "vertically compress" the circuit, eliminating the unused registers. So at runtime, the circuit would act like:

T  : |0|1|2|

q0 : -H-C---
        |   
q1 : ---X-C-
          | 
q2 : -----X-

T  : |0|1|2|

2. You can maintain the dimensionality of the circuit, treating the unused registers as if they were acted on by the Identity gate. So at runtime, the circuit would act like:

   
   T  : |0|1|2|

q0 : -H-C---	q1 : -I-	---

q2 : ---X-C-	q3 : -I---	-

q4 : -----X-

T : |0|1|2|

[krneta]

Personally I would prefer to go with the first option (mostly because of performance implications), but I'd want to discuss it with a few more people on the team before we'd actually pick it up.

[ryanhill1]

Sounds good!

[speller26]

Hi @ryanhill1,

You can use the big-endian calculate_unitary_big_endian method:

import numpy as np
from braket.circuits import Circuit, Instruction, gates
from braket.circuits.unitary_calculation import calculate_unitary

qubit_index = 1
circuit = Circuit()
pgates = [
    gates.Rx,
    gates.Ry,
    gates.Rz,
    gates.PhaseShift,
]
angles = np.random.RandomState(11).random(len(pgates))
instructions = [
    Instruction(rot(a), target=qubit_index)
    for rot, a in zip(pgates, angles)
]
for instr in instructions:
    circuit.add_instruction(instr)

# Alternatively, and more idiomatically,
# unitary = circuit.to_unitary()
unitary = calculate_unitary_big_endian(circuit.instructions, circuit.qubits)

Is there anything else that you need for this issue?

[ryanhill1]

@speller26 Nothing else! Seems like the issue was fixed with Release v1.19.0