Create a gate to represent a PauliRotation

[NoureldinYosri]

Is your feature request related to a use case or problem? Please describe.
The operation $e^{i \theta P}$ where $P$ is a pauli string, appears a lot in discussions around rotations and magic states. Cirq has a a general way for representing exponentials of pauli strings $e^{i \theta \sum_k P_k}$ called PauliSumExponential however the pauli strings are stored as sparse strings (i.e. dropping identitiy operations) which results in the wrong unitary.

Example P=XI $\theta = \frac{\pi}{4}$ correct unitary should be

$$ \frac{1}{\sqrt{2}} (I_4 + i X \otimes I_2) $$

but cirq gives

>>> cirq.unitary(cirq.PauliSumExponential(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)), exponent=np.pi/4))
array([[0.70710678+0.j        , 0.        +0.70710678j],
       [0.        +0.70710678j, 0.70710678+0.j        ]])

Describe the solution you'd like
ideally a fix to the PauliSumExponential operation or a new gate for representing the exponential of one pauli string

What is the urgency from your perspective for this issue? Is it blocking important work?
P2 - we should do it in the next couple of quarters

[NoureldinYosri]

cc: @pavoljuhas

[roseopicta]

I dived into this as part of Unitaryhack.

A fix of PauliSumExponential seems unlikely, because the faulty behavior is not specific to this op. For example:

>>>cirq.unitary(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)))
[[0.+0.j 1.+0.j]
 [1.+0.j 0.+0.j]]

While one could expect a 4x4 unitary describing $X \otimes I$.

The general behavior is that DensePauliString "reduces" to a PauliString whenever applied to qubits, and then non-identity elements are dropped from the tensor product. This happens when passed as an argument to PauliSumExponential in particular but of course anywhere else!

The good news is that PauliStringPhasorGate, which serves as a backend for PauliSumExponential accepts DensePauliString as an argument and appears to do the job:

>>>cirq.unitary(cirq.PauliStringPhasorGate(
      cirq.DensePauliString('X'),
      exponent_neg=-1/4,
      exponent_pos=1/4).on(*cirq.LineQubit.range(2)))

[[0.70710678+0.j         0.        +0.j         0.        +0.70710678j
  0.        +0.j        ]
 [0.        +0.j         0.70710678+0.j         0.        +0.j
  0.        -0.70710678j]
 [0.        +0.70710678j 0.        +0.j         0.70710678+0.j
  0.        +0.j        ]
 [0.        +0.j         0.        -0.70710678j 0.        +0.j
  0.70710678+0.j        ]]

At this point three options to close this:

• PauliStringPhasorGate (with the little annoyance of a division by pi and repeated arguments) is a good enough solution.
• Syntaxic sugar returning such a gate is good enough.
• A proper Gate object dedicated to this task is preferred.

[NoureldinYosri]

@burlemarxiste thanks for the invisitigation. I assigned the task to you. can you try to modify the methods in DensePauliString to behave correctly (i.e. not drop identity terms). hepefully nothing breaks. if lots of things break then we may accept a patch solution that converts DensePauliString to PauliStringPhasorGate when needed. but we should first try to make DensePauliString work

[roseopicta]

Thanks! I actually found a flaw in my reasoning and PauliStringPhasorGate does not always produce the expected result and needs to be modified to handle identity correctly.

The documentation states that "When applied on qubits, a DensePauliString results in cirq.PauliString as an operation". Any insight on the original reason for this behavior, since it now needs to be fixed?

[NoureldinYosri]

@burlemarxiste no but I presume there was a good reason. anyway that can't be changed now since that will most likely break lots of things, but we can add an optional parameter to change the behaviour

Cirq/cirq-core/cirq/ops/dense_pauli_string.py

Lines 326 to 327 in 7df83d7

	def on(self, *qubits: 'cirq.Qid') -> 'cirq.PauliString':
	return self.sparse(qubits)

you can change the definition of on to be

    def on(self, *qubits: 'cirq.Qid', sparse: bool=True) :

so that the current behaviour is preserved and the same time we can use on(..., sparse=False) to get the new behaviour

[roseopicta]

After several red herrings in the internals of PauliString, PauliSum and PauliStringPhasor, I wonder if there is actually much to do at all!

The current implementation of PauliSumExponential can already be used to compute the exponential of a Pauli string, up to:

• a global phase shift. This is because I terms are not handled.
• unitary() ignoring the subspace on which the exponential acts trivially. The unitary displayed for the resulting object might thus be misleading because it only shows the subspace on which the unitary is non-trivial, the tensor products by identities on the rest of the space are implicit.

I would thus suggest the following contributions:

• Minor bug fixes: make PauliSumExponential well-behaved on I(fix the diagram display that is broken at the moment, and return a GlobalPhaseGate if the argument evaluates to just I).
• More serious bug fix: extend PauliStringPhasor to behave correctly on the extra qubits given to it (according to the documentation, they should be treated as identities, but for now they are handled as if they were extra 'Z's). For exponentiating single Pauli strings, that object is a neater alternative to calling PauliSumExponential, with the benefit of being a Gate object, and thus controllable. (cirq.PauliStringPhasor returns wrong unitary #6612)

I can create issues targeting these specific changes.

Some examples of "solutions" I tried to the proposed challenge of getting the right unitary to display:

1. Add an extra qubits argument to "PauliStringExponential" that could eventually be used to add extra qubits treated as I, with a syntax like: cirq.PauliSumExponential(cirq.DensePauliString('XI')(*qubits), exponent=np.pi/4), qubits=qubits). However, the semantics of substituting qubits to the ones originally present in the argument were very unclear for more complex expressions like PauliSums of several terms.
2. Through extra class members, let a PauliString keep track of the qubits it acted on during the dense->sparse reduction, accumulate those in the PauliSum, and later append them as extra Is in the output of PauliSumExponential.

Both feel like bloat for a very minor benefit.

Overall the entire infrastructure around PauliString is really engineered for observables and not very I friendly, and I have learnt to respect it for what it is.

[daxfohl]

Is it worth considering adding I as a proper Pauli gate? The way they silently get removed when working with Pauli strings has come up before, and the asymmetry can be hard to work with. #5715. We have optimizers now that can explicitly remove them when desired, so the implicit behavior isn't needed. It would be a fairly large change and there may be good reasons not to do it regardless (idk). But figured I'd throw it out there.

[roseopicta]

Current status:

• Fixed the extra qubits bug in PauliStringPhasor
• Made sure PauliStringPhasor returned the right unitaries by adding extra I gates whenever necessary.
• Fixed PauliStringPhasor to ensure correct phase (useful when PauliStringPhasorGate is controlled!).
• Fixed miscellaneous circuit display crash when the argument to PauliStringPhasor or PauliSumExponential reduces to I.
• Tests.

main...burlemarxiste:Cirq:main

At this stage, the only remaining quirk is that cirq.unitary on a PauliSumExponential produces a result correct up to extras $\otimes I_2$. The unit test checks the unitary of PauliSumExponential followed by Is on the remaining qubits.

One should decide if displaying a unitary for the whole space (as opposed to only the non-trivial subspace) is worth bloating PauliString and PauliSum.

A refactoring of PauliString to handle I might make sense, but I am not familiar with the core design principles (should we aim for phase-correctness? unitary-level correctness? which rewrites can be done silently? what does it mean for two Pauli strings to be equal?). It might make sense to refactor to get rid of other odd behaviors I have seen, but it is probably out of the scope of a hackathon challenge task.

[NoureldinYosri]

@burlemarxiste I think you are getting side tracked ... try to follow the logic that happens when creating a PauliSum from a DensePauliString ... find the place where the identity gets dropped and propose a solution for that

[roseopicta]

Yes, I tried to follow that logic, see my first message. To summarize:

• The DensePauliString is converted into a PauliString as soon as on() is applied to it, and there is no remaining trace of the original qubit list.
• PauliSum accumulates PauliStrings so effectively, there is no trace of that original list.
• The terms of the PauliSum are used to create PauliStringPhasors that are composed together.

The superficial solution for all this is to allow PauliSum and PauliString to keep track of the original qubit list they were created with, to ultimately pass them as extra qubits to PauliStringPhasor, which generates the final decomposition of the gate. I think this "genealogy" of qubits is what you are expecting me to work on since the beginning.

My point is that even with this change, we'll still get incorrect results or odd corner cases because PauliStringPhasor, the "backend" for PauliSumExponential has flaws: it treats qubits non-referenced in the PauliString as a Z() (the issue you closed), it doesn't behave well when the string only contains identities (including an exception when attempting to print the circuit), and it is not phase-accurate (thus not controllable).

To me, having a correctly working PauliStringPhasor is a necessary condition for a good PauliSumExponential. This is why I am focusing on it.

Could you please answer the following questions:

1. Do you agree with my point that PauliStringPhasor treats extra qubits as if they had a Z applied to them in the Pauli string (cirq.PauliStringPhasor returns wrong unitary #6612)? There is both practical evidence (wrong output) and a clear explanation for it (they are involved in the parity-checking CNOT but shouldn't).

Why are the second and third unitary identical to the first?

>>> cirq.unitary(cirq.PauliStringPhasor(cirq.DensePauliString('XZ').on(*cirq.LineQubit.range(2)), qubits=cirq.LineQubit.range(2), exponent_pos=1/4, exponent_neg=-1/4))
array([[0.70710678+0.j        , 0.        +0.j        ,
        0.        +0.70710678j, 0.        +0.j        ],
       [0.        +0.j        , 0.70710678+0.j        ,
        0.        +0.j        , 0.        -0.70710678j],
       [0.        +0.70710678j, 0.        +0.j        ,
        0.70710678+0.j        , 0.        +0.j        ],
       [0.        +0.j        , 0.        -0.70710678j,
        0.        +0.j        , 0.70710678+0.j        ]])
>>> cirq.unitary(cirq.PauliStringPhasor(cirq.DensePauliString('XI').on(*cirq.LineQubit.range(2)), qubits=cirq.LineQubit.range(2), exponent_pos=1/4, exponent_neg=-1/4))
array([[0.70710678+0.j        , 0.        +0.j        ,
        0.        +0.70710678j, 0.        +0.j        ],
       [0.        +0.j        , 0.70710678+0.j        ,
        0.        +0.j        , 0.        -0.70710678j],
       [0.        +0.70710678j, 0.        +0.j        ,
        0.70710678+0.j        , 0.        +0.j        ],
       [0.        +0.j        , 0.        -0.70710678j,
        0.        +0.j        , 0.70710678+0.j        ]])
>>> cirq.unitary(cirq.PauliStringPhasor(cirq.DensePauliString('X').on(*cirq.LineQubit.range(1)), qubits=cirq.LineQubit.range(2), exponent_pos=1/4, exponent_neg=-1/4))
array([[0.70710678+0.j        , 0.        +0.j        ,
        0.        +0.70710678j, 0.        +0.j        ],
       [0.        +0.j        , 0.70710678+0.j        ,
        0.        +0.j        , 0.        -0.70710678j],
       [0.        +0.70710678j, 0.        +0.j        ,
        0.70710678+0.j        , 0.        +0.j        ],
       [0.        +0.j        , 0.        -0.70710678j,
        0.        +0.j        , 0.70710678+0.j        ]])

The docstring for PauliStringPhasor mentions "The pauli_string contains only the non-identity component of the phasor, while the qubits supplied here and not in pauli_string are acted upon by identity. " This is clearly not the case here, the extra qubit is acted upon by a Z.

Whether it helps or not with issue 6598, shouldn't we fix this?

2. What is the end goal of this issue? Get correct results when simulating a circuit with exponentials of Pauli sums, or get the missing $\otimes I_2$ in the unitary display? Does phase / controllability matter at all? To get controllable Pauli rotations behaving well on the string I, phase accuracy is important and it's happening in PauliStringPhasor.

3. Is getting the whole space (including the trivial subspace) a requirement "by design" for cirq.unitary()? If so, in order to get the correct subspace, one should yield extra I()s on the qubits on which the operation acts trivially. Then, should it be done only in PauliSumExponential or also at the root of the issue, in PauliStringPhasor?

[roseopicta]

This is the CL that introduced the extra qubits to PauliStringPhasor:
#5565

I'm really sorry to resort to this, but could maybe @dabacon or @Strilanc take a look at the unitaries in my previous message and confirm that the extra qubits are acted upon an unwanted Z, not I.

[NoureldinYosri]

I did a hacky solution that only affect cirq-core/cirq/ops/dense_pauli_string.py and cirq-core/cirq/ops/pauli_string.py to preserve the identity (by casting it as X**0)

>>> cirq.PauliSum.from_pauli_strings(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)))
cirq.PauliSum(cirq.LinearDict({frozenset({(cirq.LineQubit(0), cirq.X), (cirq.LineQubit(1), (cirq.X**0.0))}): (1+0j)}))
>>> cirq.DensePauliString('XI')(*cirq.LineQubit.range(2))
((1+0j)*cirq.X(cirq.LineQubit(0))*(cirq.X**0.0).on(cirq.LineQubit(1)))
>>> cirq.PauliSum.from_pauli_strings(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)))
cirq.PauliSum(cirq.LinearDict({frozenset({(cirq.LineQubit(0), cirq.X), (cirq.LineQubit(1), (cirq.X**0.0))}): (1+0j)}))
>>> cirq.unitary(cirq.PauliSum.from_pauli_strings(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)))).round(3)
array([[0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j],
       [0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j],
       [1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],
       [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j]])
>>> cirq.PauliSumExponential(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)), exponent=np.pi/4)
cirq.PauliSumExponential(cirq.PauliSum(cirq.LinearDict({frozenset({(cirq.LineQubit(0), cirq.X), (cirq.LineQubit(1), (cirq.X**0.0))}): (1+0j)})), 0.7853981633974483)
>>> cirq.unitary(cirq.PauliSumExponential(cirq.DensePauliString('XI')(*cirq.LineQubit.range(2)), exponent=np.pi/4)).round(3)
array([[0.707+0.j   , 0.   +0.j   , 0.   +0.j   , 0.   +0.707j],
       [0.   +0.j   , 0.707+0.j   , 0.   +0.707j, 0.   +0.j   ],
       [0.   +0.j   , 0.   +0.707j, 0.707+0.j   , 0.   +0.j   ],
       [0.   +0.707j, 0.   +0.j   , 0.   +0.j   , 0.707+0.j   ]])

---

whether there is an issue with PauliStringPhasor is irrelevant to this task

---

note casting identity as X**0 is not the desired solution it's just a hacky way to show that densepaulistring is all that is needed for this issue.

---

Is it worth considering adding I as a proper Pauli gate?

I think: yes