QubitUnitary accepts sparse matrices as inputs

[JerryChen97]

Context:
As a follow-up of #6863, in this PR we would like to introduce the overloaded QubitUnitary capable of taking sparse matrices and processing along.

Currently, if a user create a new QubitUnitary with sparse matrix as input, it will succeed with unitary_check=False but fail with True. This is because lots of the internal processing make use of incompatible functionalities e.g. einsum for csr_matrix, for which we will create more sparse-specific routines to help build a reasonable workflow for the sparse cases.

Description of the Change:
QubitUnitary has quite some API that we need to consider, some of which might need further attention to decide whether or not to implement for sparse matrices right now.

• init: inserted extra conditional unitary check, and separated the previous logic as static _unitary_check
• matrix(): the method to generate QubitUnitary's matrix should be able to return a sparse as well
• adjoint: conditional treatment for sparse matrices, since our default dispatch method does not work for csr_matrix
• pow: similar as above, although this could be easily merged into our multi-dispatch in the future
• compute_decomposition: to avoid unnecessary extra work, one- and two-qubit sparse ops are still converted back to desnse matrices, although in practice it actually does not even make sense to use sparse matrices at this level (recalling our need for sparse matrices, it's for systems of size ~20)
• support for wires of size 20: should be able to init a sparse QubitUnitary obj with 20 wires.

Benefits:
Necessary extension for QubitUnitary implemented.

Possible Drawbacks:
Implementations are quite scattered in the source code. In the future we can think about adding overloading to qml.math instead to help unify the matrix manipulations of sparse matrices.

Also, to minimize the amount of conditional branches specifically for csr_matrix, we neglected the following functionalities:

• _controlled: this method actually involves the sparse implementation of ControlledQubitUnitary, which is beyond current scope; should be included in a future story of implementing sparse for generic controlled ops.
• Any other hidden interfaces that were inherited from parent classes: same as above; require further discussions and extra stories/epics if needed. Current interfaces already provide enough features for QubitUnitary to be called in the pipeline. For example, prod, sum etc.

P.S. please note that the multiplication between operators and states is not yet involved here since it is a matter of apply_operation which should be resolved in some follow-ups (e.g. #6883 (comment), or some other refined smaller PR).

Related GitHub Issues:
[sc-82068]

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 99.59%. Comparing base (be3ae76) to head (3e3b31e).
Report is 1 commits behind head on master.

Additional details and impacted files

@@           Coverage Diff           @@
##           master    #6889   +/-   ##
=======================================
  Coverage   99.59%   99.59%           
=======================================
  Files         480      480           
  Lines       45495    45514   +19     
=======================================
+ Hits        45310    45329   +19     
  Misses        185      185           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[PietropaoloFrisoni]

This seems flawless to me! 🚀

Just one question. Except for the _controlled method and any other interfaces inherited from parent classes (as well as multiplication between operators and states), is there something else that is supported without sparse matrices but not supported with sparse matrices?

[PietropaoloFrisoni]

Suggested change

	# Note: it is necessary to explicitly cast back to csr, or it will be come csc
	# Note: it is necessary to cast back to csr explicitly, or it will become csc

Curious, why does this happen?

[JerryChen97]

@PietropaoloFrisoni SciPy has its own optimization over matrix manipulation and storage. Whenever a CSR got transposed it automatically convert the obj to CSC no matter the actual matrix structure.

m1 = sparse.csr_matrix([[1, 0], [0, 1]])
print(m1) # CSR
m2 = m1.T
print(m2) # CSC

Ref: https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.csr_matrix.transpose.html

[AmintorDusko]

Nice work so far.
@JerryChen97, can you think of a way to accept any scipy sparse representation, and internally convert it to CSR? If you do that in a helper function it will be easier to use it everywhere. This might help you do the trick.

[AmintorDusko]

Sometimes the matrices users choose to run in a sparse mode are very difficult to handle in dense form.
Converting it to a dense matrix here might be a problem.
Do we need that? If yes, can we have a default where this will not happen but raise an exception, but the user can ask for this conversion to happen?

[AmintorDusko]

Same as before. Do we want to just convert such matrices to dense and use the already-in-place functionality?

[AmintorDusko]

If we are doing everything in CSR format, let's convert it as soon as possible.

Suggested change

	def adjoint(self) -> "QubitUnitary":
	U = self.matrix()
	if isinstance(U, csr_matrix):
	adjoint_sp_mat = U.conjugate().transpose()
	# Note: it is necessary to explicitly cast back to csr, or it will be come csc
	return QubitUnitary(csr_matrix(adjoint_sp_mat), wires=self.wires)
	return QubitUnitary(qml.math.moveaxis(qml.math.conj(U), -2, -1), wires=self.wires)
	def adjoint(self) -> "QubitUnitary":
	U = self.matrix()
	if isinstance(U, csr_matrix):
	adjoint_sp_mat = csr_matrix(U.conjugate().transpose())
	# Note: it is necessary to explicitly cast back to csr, or it will become csc
	return QubitUnitary(adjoint_sp_mat, wires=self.wires)
	return QubitUnitary(qml.math.moveaxis(qml.math.conj(U), -2, -1), wires=self.wires)

[JerryChen97]

Nice work so far. @JerryChen97, can you think of a way to accept any scipy sparse representation, and internally convert it to CSR? If you do that in a helper function it will be easier to use it everywhere. This might help you do the trick.

@AmintorDusko My plan is to put issparse in the init of both StatePrep and QubitUnitary, and then convert to csr_matrix specifically.
It might be also helpful if we know whether or not this is necessary. My own experience based on past practice is SciPy.sparse works fastest with csr or csc (e.g. one can try freely scipy.sparse.kron over 20 sites without enforcing everyone to be csr and it will be super slow) but not sure if this is still the generic case right now.