Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
raise ValueError("The covariance matrix does not correspond to a pure state")
rpt = i
beta = Beta(mu, hbar=hbar)
Q = Qmat(cov, hbar=hbar)
A = Amat(cov, hbar=hbar)
(n, _) = cov.shape
N = n // 2
B = A[0:N, 0:N].conj()
alpha = beta[0:N]
if np.linalg.norm(alpha) < tol:
# no displacement
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
B_rpt = reduction(B, rpt)
haf = hafnian(B_rpt)
else:
haf = hafnian_repeated(B, rpt)
else:
gamma = alpha - B @ np.conj(alpha)
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
B_rpt = reduction(B, rpt)
np.fill_diagonal(B_rpt, reduction(gamma, rpt))
haf = hafnian(B_rpt, loop=True)
else:
haf = hafnian_repeated(B, rpt, mu=gamma, loop=True)
if include_prefactor:
pref = np.exp(-0.5 * (np.linalg.norm(alpha) ** 2 - alpha @ B @ alpha))
haf *= pref
return haf / np.sqrt(np.prod(fac(rpt)) * np.sqrt(np.linalg.det(Q)))
probs1 = np.zeros([cutoff + 1], dtype=np.float64)
kk = np.arange(k + 1)
mu_red, V_red = reduced_gaussian(local_mu, cov, kk)
if approx:
Q = Qmat(V_red, hbar=hbar)
A = Amat(Q, hbar=hbar, cov_is_qmat=True)
for i in range(cutoff):
indices = result + [i]
ind2 = indices + indices
if approx:
factpref = np.prod(fac(indices))
mat = reduction(A, ind2)
probs1[i] = (
hafnian(np.abs(mat.real), approx=True, num_samples=approx_samples) / factpref
)
else:
probs1[i] = density_matrix_element(
mu_red, V_red, indices, indices, include_prefactor=True, hbar=hbar
).real
if approx:
probs1 = probs1 / np.sqrt(np.linalg.det(Q).real)
probs2 = probs1 / prev_prob
probs3 = np.maximum(
probs2, np.zeros_like(probs2)
) # pylint: disable=assignment-from-no-return
ssum = np.sum(probs3)
if ssum < 1.0:
probs3[-1] = 1.0 - ssum
B = A[0:N, 0:N].conj()
alpha = beta[0:N]
if np.linalg.norm(alpha) < tol:
# no displacement
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
B_rpt = reduction(B, rpt)
haf = hafnian(B_rpt)
else:
haf = hafnian_repeated(B, rpt)
else:
gamma = alpha - B @ np.conj(alpha)
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
B_rpt = reduction(B, rpt)
np.fill_diagonal(B_rpt, reduction(gamma, rpt))
haf = hafnian(B_rpt, loop=True)
else:
haf = hafnian_repeated(B, rpt, mu=gamma, loop=True)
if include_prefactor:
pref = np.exp(-0.5 * (np.linalg.norm(alpha) ** 2 - alpha @ B @ alpha))
haf *= pref
return haf / np.sqrt(np.prod(fac(rpt)) * np.sqrt(np.linalg.det(Q)))
beta = Beta(mu, hbar=hbar)
A = Amat(cov, hbar=hbar)
if np.linalg.norm(beta) < tol:
# no displacement
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
A_rpt = reduction(A, rpt)
haf = hafnian(A_rpt)
else:
haf = hafnian_repeated(A, rpt)
else:
# replace the diagonal of A with gamma
gamma = beta.conj() - A @ beta
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
A_rpt = reduction(A, rpt)
np.fill_diagonal(A_rpt, reduction(gamma, rpt))
haf = hafnian(A_rpt, loop=True)
else:
haf = hafnian_repeated(A, rpt, mu=gamma, loop=True)
if include_prefactor:
haf *= prefactor(mu, cov, hbar=hbar)
return haf / np.sqrt(np.prod(fac(rpt)))
probs1 = np.zeros([cutoff + 1], dtype=np.float64)
kk = np.arange(k + 1)
mu_red, V_red = reduced_gaussian(local_mu, cov, kk)
if approx:
Q = Qmat(V_red, hbar=hbar)
A = Amat(Q, hbar=hbar, cov_is_qmat=True)
for i in range(cutoff):
indices = result + [i]
ind2 = indices + indices
if approx:
factpref = np.prod(fac(indices))
mat = reduction(A, ind2)
probs1[i] = (
hafnian(np.abs(mat.real), approx=True, num_samples=approx_samples) / factpref
)
else:
probs1[i] = density_matrix_element(
mu_red, V_red, indices, indices, include_prefactor=True, hbar=hbar
).real
if approx:
probs1 = probs1 / np.sqrt(np.linalg.det(Q).real)
probs2 = probs1 / prev_prob
probs3 = np.maximum(
probs2, np.zeros_like(probs2)
) # pylint: disable=assignment-from-no-return
ssum = np.sum(probs3)
if ssum < 1.0:
probs3[-1] = 1.0 - ssum
used to scale the result.
tol (float): tolerance for determining if displacement is negligible
hbar (float): the value of :math:`\hbar` in the commutation
relation :math:`[\x,\p]=i\hbar`.
Returns:
complex: the density matrix element
"""
rpt = i + j
beta = Beta(mu, hbar=hbar)
A = Amat(cov, hbar=hbar)
if np.linalg.norm(beta) < tol:
# no displacement
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
A_rpt = reduction(A, rpt)
haf = hafnian(A_rpt)
else:
haf = hafnian_repeated(A, rpt)
else:
# replace the diagonal of A with gamma
gamma = beta.conj() - A @ beta
if np.prod([k + 1 for k in rpt]) ** (1 / len(rpt)) < 3:
A_rpt = reduction(A, rpt)
np.fill_diagonal(A_rpt, reduction(gamma, rpt))
haf = hafnian(A_rpt, loop=True)
else:
haf = hafnian_repeated(A, rpt, mu=gamma, loop=True)
if include_prefactor:
haf *= prefactor(mu, cov, hbar=hbar)
return haf / np.sqrt(np.prod(fac(rpt)))