Concurrence, first introduced in this paper is defined as $C[\rho] = \max(0,\lambda_1 - \lambda_2 - \lambda_3 - \lambda_4)$, where ${\lambdai}_{i=1}^{4}$ are the eigenvalues of matrix $\sqrt{\sqrt{\rho}\widetilde{\rho}\sqrt{\rho}}$ in descending order, $\rho = \frac{1}{4}(1\otimes1 + B_i\sigma_i\otimes1 + \bar{B}_i1\otimes\sigma_i + C\{ij}\sigma_i\otimes\sigma_j)$ is the spin density matrix and $\widetilde{\rho} = (\sigma_2\otimes\sigma_2)\rho^*(\sigma_2\otimes\sigma2)$ its conjugate matrix (for lack of a better term). Alternatively, one can compute concurrence from the square roots of the eigenvalues of matrix $\rho\widetilde{\rho}$, which avoids the inconvenience of taking square roots of 4x4 complex matrices. While $\rho$ and $\widetilde{\rho}$ are guaranteed to be Hermitian (and with unit trace) by construction, their product is not. This could mean that the eigenvalues of $\rho\widetilde{\rho}$ can be imaginary, even though it is widely claimed in the literature that the matrix should have real eigenvalues. Mathematically this is not guaranteed, which could mean that either the claim is wrong (or based on some physics assumption), or that the spin correlations $C{ij}$ and polarizations $B_i$ and $\bar{B}_i$ must adhere to some additional constraint(s), which would guarantee that $\rho\widetilde{\rho}$ indeed has real eigenvalues.
In practice, the eigenvalues can sometimes pick up a small imaginary component, which spoils the computation of concurrence. Given that we're dealing with complex-valued matrices, it is expected that some eigenvalues have imaginary component that is on the level of floating point precision. To exclude those spurious instances, we set the imaginary part to zero if it's below a certain threshold:
https://github.com/HEP-KBFI/tautau-Entanglement/blob/233882c374108480acfbffc1bd685a9356872710/src/comp_concurrence.cc#L96-L99
However, some measurements still produce eigenvalues that have much larger imaginary component. This can happen a lot, especially if the statistics of the MC sample is very small due to cuts on $|\cos\theta^*|$. In those cases we set the concurrence to some unphysical value like -1 and exclude this data point from the bootstrap estimation:
This workaround is not ideal, because we have to rely on fewer bootstrap samples to estimate the error on concurrence.
I get the impression from reading our log files that if we end up with an unphysical value of concurrence, it's because the smallest two eigenvalues have substantial imaginary components, on the level 1e-4. Fortunately, the real part of those same eigenvalues is similarly very small. Given that our estimated error on concurrence is on the level of 1e-2, we can ignore all eigenvalues with real parts smaller than 1e-4 (remember, we're taking square roots of eigenvalues to compute the concurrence) regardless of their imaginary parts.
In a very small number of cases, the real part of some eigenvalues can sometimes be negative, but in those cases their absolute value is also very small, O(1e-7). Since those eigenvalues would have negligible impact on the final results anyways, we can set them to zero.
I expect that the two workarounds greatly reduce the number of instances where we get unphysical values for concurrence, but it doesn't mean that the unphysical values would be completely eliminated. I plan to implement the workarounds very soon; the issue here is more of a formal record to motivate future changes in the code.
Concurrence, first introduced in this paper is defined as $C[\rho] = \max(0,\lambda_1 - \lambda_2 - \lambda_3 - \lambda_4)$, where ${\lambdai}_{i=1}^{4}$ are the eigenvalues of matrix $\sqrt{\sqrt{\rho}\widetilde{\rho}\sqrt{\rho}}$ in descending order, $\rho = \frac{1}{4}(1\otimes1 + B_i\sigma_i\otimes1 + \bar{B}_i1\otimes\sigma_i + C\{ij}\sigma_i\otimes\sigma_j)$ is the spin density matrix and $\widetilde{\rho} = (\sigma_2\otimes\sigma_2)\rho^*(\sigma_2\otimes\sigma2)$ its conjugate matrix (for lack of a better term). Alternatively, one can compute concurrence from the square roots of the eigenvalues of matrix $\rho\widetilde{\rho}$, which avoids the inconvenience of taking square roots of 4x4 complex matrices. While $\rho$ and $\widetilde{\rho}$ are guaranteed to be Hermitian (and with unit trace) by construction, their product is not. This could mean that the eigenvalues of $\rho\widetilde{\rho}$ can be imaginary, even though it is widely claimed in the literature that the matrix should have real eigenvalues. Mathematically this is not guaranteed, which could mean that either the claim is wrong (or based on some physics assumption), or that the spin correlations $C{ij}$ and polarizations $B_i$ and $\bar{B}_i$ must adhere to some additional constraint(s), which would guarantee that $\rho\widetilde{\rho}$ indeed has real eigenvalues.
In practice, the eigenvalues can sometimes pick up a small imaginary component, which spoils the computation of concurrence. Given that we're dealing with complex-valued matrices, it is expected that some eigenvalues have imaginary component that is on the level of floating point precision. To exclude those spurious instances, we set the imaginary part to zero if it's below a certain threshold: https://github.com/HEP-KBFI/tautau-Entanglement/blob/233882c374108480acfbffc1bd685a9356872710/src/comp_concurrence.cc#L96-L99
However, some measurements still produce eigenvalues that have much larger imaginary component. This can happen a lot, especially if the statistics of the MC sample is very small due to cuts on $|\cos\theta^*|$. In those cases we set the concurrence to some unphysical value like -1 and exclude this data point from the bootstrap estimation:
https://github.com/HEP-KBFI/tautau-Entanglement/blob/233882c374108480acfbffc1bd685a9356872710/src/SpinAnalyzer.cc#L168-L180
This workaround is not ideal, because we have to rely on fewer bootstrap samples to estimate the error on concurrence.
I get the impression from reading our log files that if we end up with an unphysical value of concurrence, it's because the smallest two eigenvalues have substantial imaginary components, on the level 1e-4. Fortunately, the real part of those same eigenvalues is similarly very small. Given that our estimated error on concurrence is on the level of 1e-2, we can ignore all eigenvalues with real parts smaller than 1e-4 (remember, we're taking square roots of eigenvalues to compute the concurrence) regardless of their imaginary parts.
In a very small number of cases, the real part of some eigenvalues can sometimes be negative, but in those cases their absolute value is also very small, O(1e-7). Since those eigenvalues would have negligible impact on the final results anyways, we can set them to zero.
I expect that the two workarounds greatly reduce the number of instances where we get unphysical values for concurrence, but it doesn't mean that the unphysical values would be completely eliminated. I plan to implement the workarounds very soon; the issue here is more of a formal record to motivate future changes in the code.