From b2658ea2e3ba535427684c06c3a2ad7e8a5273b0 Mon Sep 17 00:00:00 2001 From: Renato Mello Date: Fri, 4 Oct 2024 15:34:04 +0400 Subject: [PATCH] singular matrices --- src/qibojit/backends/gpu.py | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/src/qibojit/backends/gpu.py b/src/qibojit/backends/gpu.py index 54b3430..cf45251 100644 --- a/src/qibojit/backends/gpu.py +++ b/src/qibojit/backends/gpu.py @@ -543,9 +543,21 @@ def calculate_matrix_power( self, matrix, power: Union[float, int], precision_singularity: float = 1e-14 ): - if isinstance(power, int): + if isinstance(power, int) and power >= 0.0: return self.cp.linalg.matrix_power(matrix, power) + if power < 0.0: + # negative powers of singular matrices via SVD + determinant = self.cp.linalg.det(matrix) + if abs(determinant) < precision_singularity: + U, S, Vh = self.cp.linalg.svd(matrix) + S_inv = self.cp.where( + self.cp.abs(S) < precision_singularity, 0.0, S**power + ) + return ( + self.cp.linalg.inv(Vh) @ self.cp.diag(S_inv) @ self.cp.linalg.inv(U) + ) + copied = self.to_numpy(matrix) copied = super().calculate_matrix_power(copied, power, precision_singularity)