diff --git a/filter_functions/basis.py b/filter_functions/basis.py index 64d8d78..a69fefd 100644 --- a/filter_functions/basis.py +++ b/filter_functions/basis.py @@ -318,11 +318,8 @@ def H(self) -> 'Basis': @property def T(self) -> 'Basis': """Return the basis transposed element-wise.""" - if self.ndim == 3: - return self.transpose(0, 2, 1) - - if self.ndim == 2: - return self.transpose(1, 0) + if self.ndim >= 2: + return self.swapaxes(-1, -2) return self diff --git a/filter_functions/numeric.py b/filter_functions/numeric.py index ff71548..1c15498 100644 --- a/filter_functions/numeric.py +++ b/filter_functions/numeric.py @@ -1108,9 +1108,20 @@ def calculate_cumulant_function( cumulant_function[..., 1:, 1:] -= frequency_shifts[..., 1:, 1:] cumulant_function[..., 1:, 1:] += frequency_shifts[..., 1:, 1:].swapaxes(-1, -2) else: - # Multi qubit case. Use general expression. Drop imaginary part since - # result is guaranteed to be real (if we didn't do anything wrong) + # Multi qubit case. Use general expression. traces = pulse.basis.four_element_traces + # g_iklj = ( + # + traces.transpose(0, 1, 2, 3) + # - traces.transpose(0, 1, 3, 2) + # - traces.transpose(0, 3, 1, 2) + # + traces.transpose(0, 3, 2, 1) + # ) + # g_jikl = ( + # + traces.transpose(3, 0, 1, 2) + # - traces.transpose(2, 0, 1, 3) + # - traces.transpose(2, 0, 3, 1) + # + traces.transpose(1, 0, 3, 2) + # ) cumulant_function = - ( + oe.contract('...kl,klji->...ij', decay_amplitudes, traces, backend='sparse') - oe.contract('...kl,kjli->...ij', decay_amplitudes, traces, backend='sparse') @@ -1118,6 +1129,22 @@ def calculate_cumulant_function( + oe.contract('...kl,kijl->...ij', decay_amplitudes, traces, backend='sparse') ) / 2 if second_order: + # f_iklj = ( + # + traces.transpose(0, 1, 2, 3) + # - traces.transpose(0, 2, 1, 3) + # - traces.transpose(0, 2, 3, 1) + # + traces.transpose(0, 3, 2, 1) + # ) + # + # f_iklj = -g_iljk: + # f = -g.transpose(0, 2, 3, 1) + # + # f_jikl = ( + # + traces.transpose(3, 0, 1, 2) + # - traces.transpose(3, 0, 2, 1) + # - traces.transpose(1, 0, 2, 3) + # + traces.transpose(1, 0, 3, 2) + # ) cumulant_function -= ( + oe.contract('...kl,klji->...ij', frequency_shifts, traces, backend='sparse') - oe.contract('...kl,lkji->...ij', frequency_shifts, traces, backend='sparse') @@ -1125,6 +1152,7 @@ def calculate_cumulant_function( + oe.contract('...kl,lkij->...ij', frequency_shifts, traces, backend='sparse') ) / 2 + # Drop imaginary part since result is guaranteed to be real (if we didn't do anything wrong) return cumulant_function.real @@ -1851,17 +1879,20 @@ def infidelity( The ``PulseSequence`` instance for which to calculate the infidelity for. spectrum: array_like, shape ([[n_nops,] n_nops,] omega) or callable - The two-sided noise power spectral density in units of inverse - frequencies as an array of shape (n_omega,), (n_nops, n_omega), - or (n_nops, n_nops, n_omega). In the first case, the same - spectrum is taken for all noise operators, in the second, it is - assumed that there are no correlations between different noise - sources and thus there is one spectrum for each noise operator. + The noise power spectral density in units of inverse frequencies + as an array of shape (n_omega,), (n_nops, n_omega), or + (n_nops, n_nops, n_omega). In the first case, the same spectrum + is taken for all noise operators, in the second, it is assumed + that there are no correlations between different noise sources + and thus there is one spectrum for each noise operator. In the third and most general case, there may be a spectrum for each pair of noise operators corresponding to the correlations between them. n_nops is the number of noise operators considered and should be equal to ``len(n_oper_identifiers)``. + See :ref:`Notes ` for a discussion on one- and two-sided + power spectral densities. + If *test_convergence* is ``True``, a function handle to compute the power spectral density from a sequence of frequencies is expected. @@ -1946,15 +1977,34 @@ def infidelity( infidelities that can be computed by setting ``which='correlations'``. + **One- and two-sided spectral densities** - To convert to the average gate infidelity, use the - following relation given by Horodecki et al. [Hor99]_ and - Nielsen [Nie02]_: + Since the real (imaginary) part of filter function :math:`F(\omega)` + is even (odd), it does not matter for integral whether + :math:`S(\omega)` is taken to be the one- or two-sided spectral + density. However, care should be taken that, if it is one or the + other, the frequencies :math:`\omega` are positive or symmetric + about zero, respectively. + + To convert between one- and two-sided PSDs, use the following + relationship: + + .. math:: + + S_\mathrm{onesided}(\omega) = 2 S_\mathrm{twosided}(\omega). + + **Conversion to the Average Gate Infidelity (AGI)** + + To convert the entanglement infidelity to the average gate + infidelity, use the following relation given by Horodecki et al. + [Hor99]_ and Nielsen [Nie02]_: .. math:: \mathcal{I}_\mathrm{avg} = \frac{d}{d+1}\mathcal{I}. + **Goodness of approximation** + The smallness parameter is given by .. math:: diff --git a/tests/test_basis.py b/tests/test_basis.py index 5fc0e90..9f5d3e4 100644 --- a/tests/test_basis.py +++ b/tests/test_basis.py @@ -156,6 +156,25 @@ def test_basis_properties(self): self.assertTrue(basis.isorthonorm) self.assertArrayEqual(basis.T, basis.view(np.ndarray).T) + def test_transpose(self): + arr = rng.normal(size=(2, 3, 3)) + b = arr.view(ff.Basis) + self.assertArrayEqual(b.T, arr.transpose(0, 2, 1)) + + arr = rng.normal(size=(2, 2)) + b = arr.view(ff.Basis) + self.assertArrayEqual(b.T, arr.transpose(1, 0)) + + # Edge cases for not officially intended usage + arr = rng.normal(size=2) + b = arr.view(ff.Basis) + self.assertArrayEqual(b.T, arr) + + arr = rng.normal(size=(2, 3, 4, 4)) + b = arr.view(ff.Basis) + self.assertArrayEqual(b.T, arr.transpose(0, 1, 3, 2)) + + def test_basis_expansion_and_normalization(self): """Correct expansion of operators and normalization of bases""" # dtype