Skip to content

Performance improvement of SparsePauliOp.to_matrix #9620

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 9 commits into from
Mar 23, 2023
Merged

Performance improvement of SparsePauliOp.to_matrix #9620

Changes from 2 commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
9b4adc0
optimize SparsePauliOp.to_matrix
t-imamichi Feb 17, 2023
e0f541f
optimize dense
t-imamichi Feb 18, 2023
fb1428c
proper vectorize
t-imamichi Feb 22, 2023
7cfb91c
Merge branch 'main' into opt-tomatrix2
t-imamichi Feb 22, 2023
d962498
update copyright year
t-imamichi Feb 22, 2023
cea780e
small opt
t-imamichi Feb 22, 2023
972b406
update a comment
t-imamichi Feb 22, 2023
8094903
Merge branch 'main' into opt-tomatrix2
t-imamichi Feb 23, 2023
6915b1e
Merge branch 'main' into opt-tomatrix2
t-imamichi Mar 21, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 8 additions & 5 deletions qiskit/quantum_info/operators/symplectic/base_pauli.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
# This code is part of Qiskit.
#
# (C) Copyright IBM 2017, 2020
# (C) Copyright IBM 2017, 2023
#
# This code is licensed under the Apache License, Version 2.0. You may
# obtain a copy of this license in the LICENSE.txt file in the root directory
Expand All @@ -26,6 +26,10 @@
from qiskit.quantum_info.operators.mixins import AdjointMixin, MultiplyMixin


# utility for _to_matrix
_PARITY = np.array([-1 if bin(i).count("1") % 2 else 1 for i in range(256)], dtype=complex)


class BasePauli(BaseOperator, AdjointMixin, MultiplyMixin):
r"""Symplectic representation of a list of N-qubit Paulis.

Expand Down Expand Up @@ -394,7 +398,7 @@ def _from_array(z, x, phase=0):

@staticmethod
def _to_matrix(z, x, phase=0, group_phase=False, sparse=False):
"""Return the matrix matrix from symplectic representation.
"""Return the matrix from symplectic representation.

The Pauli is defined as :math:`P = (-i)^{phase + z.x} * Z^z.x^x`
where ``array = [x, z]``.
Expand Down Expand Up @@ -432,7 +436,7 @@ def _to_matrix(z, x, phase=0, group_phase=False, sparse=False):
coeff = (-1j) ** phase
else:
coeff = 1
data = np.array([coeff * (-1) ** (bin(i).count("1") % 2) for i in z_indices & indptr])
data = coeff * _PARITY[(z_indices & indptr).astype(np.uint8)]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why is uint8 enough? Is it possible to have a larger number?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The original code computes the parity of each value of z_indices & indptr.
My code first applies astype(np.uint8) to take mod 256 (parity does not change) and then applies a lookup table _PARITY[256] to compute the parity. Does it make sense?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for your explanation. I misunderstood the definition of parity.

Copy link
Member Author

@t-imamichi t-imamichi Feb 22, 2023 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry. I noticed that I misunderstood the party of the code. Your concern is correct. This code does not work with many qubits. I will rethink the optimization.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

fixed at fb1428c

if sparse:
# Return sparse matrix
from scipy.sparse import csr_matrix
Expand All @@ -441,8 +445,7 @@ def _to_matrix(z, x, phase=0, group_phase=False, sparse=False):

# Build dense matrix using csr format
mat = np.zeros((dim, dim), dtype=complex)
for i in range(dim):
mat[i][indices[indptr[i] : indptr[i + 1]]] = data[indptr[i] : indptr[i + 1]]
mat[range(dim), indices[:dim]] = data[:dim]
return mat

@staticmethod
Expand Down