From d8ec46674d9a0f7fdd4eea9eaf05dcaaf7b76e96 Mon Sep 17 00:00:00 2001
From: Jan Petykiewicz <jan@mpxd.net>
Date: Tue, 23 May 2023 12:56:38 -0700
Subject: [PATCH] sqrtm increases precision, so cast back to double

---
 meanas/fdfd/bloch.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/meanas/fdfd/bloch.py b/meanas/fdfd/bloch.py
index 83c45d8..9dacdee 100644
--- a/meanas/fdfd/bloch.py
+++ b/meanas/fdfd/bloch.py
@@ -578,7 +578,7 @@ def eigsolve(
 
         trace_U = real(trace(U))
         if trace_U > 1e8 * num_modes:
-            Z = Z @ scipy.linalg.sqrtm(U).conj().T
+            Z = Z @ scipy.linalg.sqrtm(U).astype(numpy.complex128).conj().T
             prev_traceGtKG = 0
             continue
         break
@@ -741,7 +741,7 @@ def eigsolve(
     # Recover eigenvectors from Z
     #
     U = numpy.linalg.inv(ZtZ)
-    Y = Z @ scipy.linalg.sqrtm(U)
+    Y = Z @ scipy.linalg.sqrtm(U).astype(numpy.complex128)
     W = Y.conj().T @ (scipy_op @ Y.copy())
 
     eigvals, W_eigvecs = numpy.linalg.eig(W)