# Dimension-independent weak value estimation via controlled SWAP operations

Giulio Chiribella<EMAIL_ADDRESS>QICI Quantum Information and Computation
Initiative, Department of Computer Science, The University of Hong Kong,
Pokfulam Road, Hong Kong  Quantum Group, Department of Computer Science,
University of Oxford, Wolfson Building, Parks Road, Oxford, OX1 3QD, United
Kingdom  Perimeter Institute for Theoretical Physics, 31 Caroline Street
North, Waterloo, Ontario, Canada  Kyrylo Simonov<EMAIL_ADDRESS>Fakultät für Mathematik, Universität Wien, Oskar-Morgenstern-Platz 1, 1090
Vienna, Austria  Xuanqiang Zhao<EMAIL_ADDRESS>QICI Quantum
Information and Computation Initiative, Department of Computer Science, The
University of Hong Kong, Pokfulam Road, Hong Kong

###### Abstract

Weak values of quantum observables are a powerful tool for investigating a
broad spectrum of quantum phenomena. For this reason, several methods to
measure them in the laboratory have been proposed. Some of these methods
require weak interactions and postselection, while others are deterministic,
but require statistics over a number of experiments growing exponentially with
the number of measured particles. Here we propose a deterministic dimension-
independent scheme for estimating weak values of arbitrary observables. The
scheme, based on coherently controlled SWAP operations, does not require prior
knowledge of the initial and final states, nor of the measured observables,
and therefore can work with uncharacterized preparation and measurement
devices. As a byproduct, our scheme provides an alternative expression for
two-time states, that is, states describing quantum systems subject to pre-
and post-selections. Using this expression, we show that the controlled-SWAP
scheme can be used to estimate weak values for a class of two-time states
associated to bipartite quantum states with positive partial transpose.

##  I Introduction

In a seminal 1988 paper Aharonov  _et al._ (1988), Aharonov, Albert, and
Vaidman introduced the notion of weak values and showed that they can be
experimentally accessed by letting the measured system interact weakly with a
pointer in the time interval between a pre-selection and a post-selection.
Since then, weak values have proven a powerful tool for analyzing a broad
spectrum of quantum phenomena Aharonov  _et al._ (2010); Shikano (2012);
Kofman  _et al._ (2012); Dressel  _et al._ (2014). On the foundational side,
they provide a lens for experimentally investigating quantum paradoxes
Aharonov  _et al._ (2002); Lundeen and Steinberg (2009); Kocsis  _et al._
(2011) and Leggett-Garg inequalities Williams and Jordan (2008); Dressel  _et
al._ (2011); Goggin  _et al._ (2011a); Groen  _et al._ (2013), as well as a
quantitative indicator of non-classicality of quantum states Pusey (2014);
Dressel  _et al._ (2014); Dressel (2015); Kunjwal  _et al._ (2019). On the
applied side, they provide probabilistic amplification techniques for quantum
metrology Hosten and Kwiat (2008); Dixon  _et al._ (2009); Brunner and Simon
(2010); Starling  _et al._ (2010a, b); Hofmann (2011); Hofmann  _et al._
(2012); Strübi and Bruder (2013) and for the direct measurement of quantum
states Lundeen  _et al._ (2011); Lundeen and Bamber (2012); Salvail  _et al._
(2013); Malik  _et al._ (2014).

Several experimental schemes for measuring weak values have been proposed and
demonstrated in the laboratory. Some schemes, based on the original definition
of weak values, involve weak measurement interactions and post-selection (see
Dressel  _et al._ (2014) for a review). More recently, there has been a
growing interest in methods for estimating weak values without post-selection
and weak couplings. An ingenious method was developed by Hoffman Hofmann
(2012), who showed that weak values can be estimated by performing standard
measurements on the outputs of the universal quantum cloning machine. An
alternative method was proposed by Wagner et al Wagner  _et al._ (2023), based
on the so-called cyclic shift test Oszmaniec  _et al._ (2021), which allows
one to estimate the trace of the product of multiple density matrices.

A limitation of most of the existing deterministic estimation protocols,
however, is that their number of samples generally grows linearly in the
system dimension $d$. For example, the cloning method Hofmann (2012) provides
an estimate of the weak value multiplied by a terms of order $1/d$. Hence,
obtaining a reliable estimate requires a number of repetitions of the
experiment growing at least as $d$. For a quantum system composed of $n$
particles, this scaling results into an exponential increase of the sample
complexity with $n$. A similar issue arises in the cyclic test approach Wagner
_et al._ (2023), whose sample complexity is also growing with $d$ in the worst
case over all possible observables.

In this paper, we provide a deterministic method for estimating weak values
with dimension-independent sample complexity. In our method, two identical
quantum systems are prepared in two quantum states associated to the pre- and
post-selections appearing in the definition of weak value. Then, the two
systems undergo a controlled SWAP operation, which exchanges them or leaves
them unchanged depending on the quantum state of a control system. After the
controlled-SWAP operation, the two systems undergo a possibly noisy
measurement, whose measurement outcomes are used to estimate the weak value by
classical post-processing.

The controlled-SWAP scheme studied in this paper also provides insights into
the theory of two-time states, a generalized notion of states that describe
quantum systems subject to both pre- and post-selections Aharonov  _et al._
(1964); Aharonov and Vaidman (1991); Silva  _et al._ (2014). We show that the
expectation values associated to two-time states can be obtained from the
expectation values associated to ordinary density matrices, by applying a
linear fractional transformation that involves a SWAP operation and a partial
transpose. Using this result, we show that the controlled-SWAP method can be
used to estimate the expectation values of all two-time states associated to
density matrices with positive partial transpose (PPT) Peres (1996); Horodecki
_et al._ (1996). Finally, we discuss several adaptations of the controlled-
SWAP method that allow to access the expectation values for two-time states
associated to non-PPT states.

The rest of the paper is structured as follows. In Section II, we briefly
review the notion of weak value and put forward a dimension-independent
estimation scheme based on controlled-SWAP operations. In Section III, we
develop several generalizations of the proposed estimation scheme and
introduce a new type of quantities dubbed “double weak values”. In Section IV,
we provide the reader with an overview of the theory of two-time states. In
Section V, we discuss the matrix representation of two-time states. In Section
VI, we unravel a fundamental connection between two-time states and the
developed estimation scheme. In Section VII, we draw the conclusions and
outlook.

##  II Weak value estimation

Weak values were introduced by Aharonov, Albert and Vaidman, as part of a
framework describing pre- and post-selected ensembles Aharonov  _et al._
(1988). Their original definition referred to the scenario where a quantum
system is pre-selected and post-selected in two pure states, described by rays
in the system’s Hilbert space $\mathcal{H}$: given two (generally
unnormalized) vectors $|\psi\rangle\in\mathcal{H}$ and
$|\psi^{\prime}\rangle\in\mathcal{H}$ satisfying the condition
$\langle\psi|\psi^{\prime}\rangle\not=0$, the weak value of an observable $A$
is

$\displaystyle
W(A\,|\,\psi,\psi^{\prime}):=\frac{\langle\psi^{\prime}|A|\psi\rangle}{\langle\psi^{\prime}|\psi\rangle}\,.$
(1)

Here the vector $|\psi\rangle$ ($|\psi^{\prime}\rangle$) represents the
initial (final) state of the system. The observable $A$ is typically taken to
be a self-adjoint operator, but more generally could be any linear operator
acting on $\mathcal{H}$. Hereafter, the algebra of all linear operators on
$\mathcal{H}$ will be denoted by $\mathsf{Lin}(\mathcal{H})$.

The notion of weak value was later generalized to mixed states (see e.g.
Vaidman  _et al._ (2017)): for a pair of density matrices $\rho$ and
$\rho^{\prime}$ satisfying the condition
$\operatorname{Tr}[\rho\rho^{\prime}]\not=0$, the weak value is defined as

$\displaystyle
W(A\,|\,\rho,\rho^{\prime})=\frac{\operatorname{Tr}[\rho^{\prime}A\rho]}{\operatorname{Tr}[\rho^{\prime}\rho]}\,.$
(2)

Eq. (1) can be obtained as a special case of Eq. (2) by setting
$\rho=|\psi\rangle\langle\psi|$ and
$\rho^{\prime}=|\psi^{\prime}\rangle\langle\psi^{\prime}|$.

An important question is how to measure weak values. In their seminal paper
Aharonov  _et al._ (1988), Aharonov, Albert, and Vaidman provided an
experimental scheme using weak interactions and postselection. A number of
other schemes was subsequently devised by other authors Ritchie  _et al._
(1991); Pryde  _et al._ (2005); Hosten and Kwiat (2008); Dixon  _et al._
(2009); Brunner and Simon (2010); Goggin  _et al._ (2011b); Lundeen  _et al._
(2011); Li  _et al._ (2011); Dressel and Jordan (2012); Xu  _et al._ (2013);
Magaña-Loaiza  _et al._ (2014); Piacentini  _et al._ (2016); Vaidman (2017);
Cohen and Pollak (2018); Rebufello  _et al._ (2021); Adam and Diósi (2023);
Wiseman  _et al._ (2023). The most recent scheme was presented in Ref. Wagner
_et al._ (2023). In its simplest version, the scheme is defined for rank-one
observables, of the form $A=|\alpha\rangle\langle\alpha|$ for some vector
$|\alpha\rangle\in\mathcal{H}$. It consists in two subprotocols: a cyclic test
Oszmaniec  _et al._ (2021) for estimating the trace
$\operatorname{Tr}[\rho^{\prime}|\alpha\rangle\langle\alpha|\rho]]$, and a
swap test Barenco  _et al._ (1997); Buhrman  _et al._ (2001) for estimating
the trace $\operatorname{Tr}[\rho^{\prime}\rho]$. The extension to general
self-adjoint observables $A$ is done by using the spectral decomposition
$A=\sum_{i}\,a_{i}\,|\alpha_{i}\rangle\langle\alpha_{i}|$ and computing the
weak value of $A$ as a linear combination of the weak values of the rank-one
observables $A_{i}:=|\alpha_{i}\rangle\langle\alpha_{i}|$. A limitation of
this approach, however, is that in general an observable $A$ can have up to
$d$ distinct eigenvalues, and therefore the number of experimental settings
needed to estimate the weak value of $A$ grows with the system’s dimension in
the worst case. For a system of $n$ particles, the number of settings for the
estimation of a generic weak value grows exponentially with $n$.

We now provide a way to estimate weak values of arbitrary observables in a
dimension-independent way. Our scheme applies also to infinite-dimensional
systems, and does not require the use of ideal projective measurements.
Consider a general quantum measurement, described by a positive operator-
valued measure (POVM) $(P_{j})_{j=1}^{N}$, satisfying the conditions
$P_{j}\geq 0\,\forall j$ and $\sum_{j}P_{j}=I$. We say that this POVM allows
estimation of the observable $A$ if $A$ is a linear combination of the POVM
operators, namely

$\displaystyle A=\sum_{j}\,x_{j}\,P_{j}\,,$ (3)

for suitable coefficients $\\{x_{j}\\}_{j=1}^{N}\in\mathbb{R}^{N}$ (here and
in the following we consider the estimation of self-adjoint observables). An
example of POVM that satisfies condition (3) is a noisy measurement of $A$,
corresponding to $P_{j}=(1-p)\,|\alpha_{j}\rangle\langle\alpha_{j}|+p\,I/d$
for $j\in\\{1,\dots,d\\}$ and $p\in[0,1]$. In this case, Eq. (3) is satisfied
with $x_{j}=(a_{j}-p\,\operatorname{Tr}[A]/d)/(1-p)$.

The operators $P_{j}$ give rise to a complex measure

$\displaystyle
q(j|\rho,\rho^{\prime}):=\operatorname{Tr}[P_{j}\,\rho\,\rho^{\prime}]\,,$ (4)

normalized as
$\sum_{j}q(j|\rho,\rho^{\prime})=\operatorname{Tr}[\rho\,\rho^{\prime}]$. We
call this measure the weak value (WV) measure.

When condition Eq. (3) is satisfied, the weak value (2) can be written as the
ratio of two expectation values with respect to the WV measure: the
expectation value
$\mathbb{E}_{q}(X):=\sum_{j}\,x_{j}\,q(j|\rho,\rho^{\prime})$ of the random
variable $X$ and the expectation value
$\mathbb{E}_{q}(Y):=\sum_{j}\,q(j|\rho,\rho^{\prime})$ of the constant random
variable $Y$ with values $y_{j}=1\,\forall j$; in formula,

$\displaystyle
W(A\,|\,\rho,\rho^{\prime})=\frac{\mathbb{E}_{q}(X)}{\mathbb{E}_{q}(Y)}\,.$
(5)

In general, the WV measure $q(j|\rho,\rho^{\prime})$ is not a probability
measure, and therefore $\mathbb{E}_{q}(\cdot)$ are not proper expectation
values. Nevertheless, we now show that every expectation value with respect to
the WV measure can be evaluated by computing expectations of suitable random
variables with respect to a proper probability distribution. To this purpose,
we use the following protocol.

Protocol 1 Sampling from WV measure
---
Inputs. Two copies $S_{1}$, $S_{2}$ of the system, on which the weak value of
(3) is intended to be measured, prepared in the states $\rho$ and
$\rho^{\prime}$, respectively, and an auxiliary qubit $C$ prepared in the
state $|+\rangle\langle+|$ with $|+\rangle=(|0\rangle+|1\rangle)/\sqrt{2}$. 1.
Application of the controlled SWAP gate $U=I_{S_{1}}\otimes
I_{S_{2}}\otimes|0\rangle\langle 0|_{C}+{\tt
SWAP}_{S_{1}S_{2}}\otimes|1\rangle\langle 1|_{C}$ (6) to all three systems. 2.
Measurement of $S_{1}$ with the POVM $(P_{j})_{j}$. 3. Measurement of $C$ with
the four-outcome POVM $(R_{c})_{c=0}^{3}$ with $\displaystyle R_{0}$
$\displaystyle=$ $\displaystyle 1/2|+\rangle\langle+|,$ (7) $\displaystyle
R_{1}$ $\displaystyle=$ $\displaystyle 1/2|-\rangle\langle-|,$ (8)
$\displaystyle R_{2}$ $\displaystyle=$ $\displaystyle
1/2|+i\rangle\langle+i|,$ (9) $\displaystyle R_{3}$ $\displaystyle=$
$\displaystyle 1/2|-i\rangle\langle-i|,$ (10) where
$|\pm\rangle=(|0\rangle\pm|1\rangle)/\sqrt{2}$, and $|\pm
i\rangle=(|0\rangle\pm i|1\rangle)/\sqrt{2}$. Physically, the POVM
$(R_{c})_{c}$ can be realized as a random measurement, by measuring either on
the basis $\\{|+\rangle,|-\rangle\\}$ or on the basis
$\\{|+i\rangle,|-i\rangle\\}$. Output. The protocol produces a pair of
outcomes $(j,c)$, distributed with probability $\displaystyle
p(j,c|\rho,\rho^{\prime})$ $\displaystyle:=$
$\displaystyle\operatorname{Tr}\left[(P_{j}\otimes I\otimes
R_{c})\,U(\rho\otimes\rho^{\prime}\otimes|+\rangle\langle+|)U^{\dagger}\right]$
(11) $\displaystyle=$
$\displaystyle\frac{1}{8}\Big{\\{}\operatorname{Tr}[P_{j}\,\rho]+\operatorname{Tr}[P_{j}\,\rho^{\prime}]$
$\displaystyle+\;2\theta(1-c)\,(-1)^{c}\,{\rm
Re}\left(\operatorname{Tr}[P_{j}\rho\rho^{\prime}]\right)$
$\displaystyle-\;2\theta(c-2)\,(-1)^{c}\,{\rm
Im}\left(\operatorname{Tr}[P_{j}\rho\rho^{\prime}]\right)\Big{\\}}\,,$ where
$\theta(t)$ is the Heaviside step function, defined as $\theta(t)=1$ for
$t\geq 0$ and $\theta(t)=0$ for $t<0$ (see Appendix A for a derivation of Eq.
(11)).

The probability distribution (11) generated by the proposed protocol leads to
the following:

######  Lemma 1.

For every random variable $Z:j\mapsto z_{j}$, the expectation value of $Z$
with respect to the WV measure $q(j|\rho,\rho^{\prime})$ in Eq. (4) coincides
with the expectation value of the random variable
$\widetilde{Z}:(j,c)\mapsto\widetilde{z}_{j,c}$ defined by

$\displaystyle\widetilde{z}_{j,c}:=2z_{j}\,(-1)^{c}\,\Big{[}\,\theta(1-c)-i\theta(c-2)\Big{]}$
(12)

with respect to the probability distribution $p(j,c|\rho,\rho^{\prime})$.

###### Proof.

See Appendix B. $\blacksquare$

The above lemma, combined with Eq. (5), implies that the weak value
$W(A|\rho,\rho^{\prime})$ can be estimated by computing the empirical average
of the random variables $X$ and $Y$ with respect to the experimental
frequencies generated by Protocol II. The sample complexity of the protocol is
provided by the following:

######  Theorem 1.

With $K$ copies of the input state pair $\rho,\rho^{\prime}$, the weak value
$W(A\,|\,\rho,\rho^{\prime})$ of an observable $A=\sum_{j}x_{j}P_{j}$ with a
POVM $(P_{j})_{j}$ and real coefficients $\\{x_{j}\\}_{j}$ can be estimated up
to a small additive error $\epsilon$ with a probability no less than
$1-\delta$ given

$\displaystyle K$ $\displaystyle=$
$\displaystyle\frac{8\ln\Bigl{(}\frac{6}{\delta}\Bigr{)}}{\epsilon^{2}}\Biggl{(}\frac{x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|}{\operatorname{Tr}[\rho\rho^{\prime}]}\Biggr{)}^{2}$
(13) $\displaystyle+\;O\left(\frac{\ln\frac{1}{\delta}}{\epsilon}\frac{x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|}{(\operatorname{Tr}[\rho\rho^{\prime}])^{2}}\right)$

where $x_{\rm max}\coloneqq\max_{j}|x_{j}|$.

###### Proof.

See Appendix C. $\blacksquare$

When the POVM is the spectral decomposition of the observable $A$, the maximum
value of $x_{j}$ is simply the norm
$\|A\|=\sup_{\||\psi\rangle\|=1}\,\langle\psi|A|\psi\rangle$. Note that the
sample complexity depends on the overlap between the states $\rho$ and
$\rho^{\prime}$, and can become large when the states are nearly orthogonal.
However, there is no dependence of the sample complexity on the dimension of
the system. As a result, our method provides a reduction of the sample
complexity with respect to the method proposed in Ref. Wagner  _et al._ (2023)
in all the cases where the observable $A$ has rank larger than one. Another
appealing feature from the experimental point of view is that our method only
requires a controlled-SWAP operation, which is generally easier to implement
than the controlled shift operations in Ref. Wagner  _et al._ (2023).

Beside the dimension-independent sample complexity, an interesting feature of
the above protocol is that it works with uncharacterized preparation devices:
to estimate the weak value $W(A|\,\rho,\rho^{\prime})$, the experimenter does
not need to know in advance the initial state $\rho$ and final state
$\rho^{\prime}$. Similarly, the experimenter does not need to know in advance
the POVM $(P_{j})_{j}$ describing the measurement apparatus: the weak value of
any observable of the form $A=\sum_{j}\,x_{j}\,P_{j}$ can be computed by
averaging the random variable $x_{j}$ with respect to the frequency
distribution of the measurement outcome. These features make the above
protocol suitable for the tomographic characterization of preparation devices
(or measurement devices) using weak values.

##  III Double weak values

The controlled-SWAP protocol presented in the previous section lends itself to
several generalizations, discussed in this section. These generalizations
allow to estimate another type of quantities, which we call “double weak
values,” and facilitate the connection with the theory of two-time states,
discussed in Section VI.

###  III.1 Two local measurements

A first generalization of Protocol II is to measure both systems $S_{1}$ and
$S_{2}$, instead of measuring only $S_{1}$. Measuring $S_{1}$ and $S_{2}$ with
two POVMs $(P_{j})_{j}$ and $(Q_{k})_{k}$, respectively, gives rise to the
probability distribution

$\displaystyle
p(j,k,c|\rho,\rho^{\prime}):=\operatorname{Tr}\left[(P_{j}\otimes Q_{k}\otimes
R_{c})\,U(\rho\otimes\rho^{\prime}\otimes|+\rangle\langle+|)U^{\dagger}\right]\,.$
(14)

By sampling over this probability distribution, one can estimate the averages
of arbitrary random variables with respect to the complex measure

$\displaystyle q(j,k\,|\,\rho,\rho^{\prime}):=\operatorname{Tr}[P_{j}\rho
Q_{k}\rho^{\prime}]\,,$ (15)

thereby estimating all quantities of the form

$\displaystyle\operatorname{Tr}[A\rho B\rho^{\prime}]\,,$ (16)

for arbitrary observables $A$ and $B$ in the linear span of $\\{P_{i}\\}_{i}$
and $\\{Q_{j}\\}_{j}$, respectively.

These quantities give rise to a generalization of the notion of weak value:

######  Definition 1.

The double weak value of a pair of observables $A\in\mathsf{Lin}(\mathcal{H})$
and $B\in\mathsf{Lin}(\mathcal{H})$ with respect to the initial state $\rho$
and the final state $\rho^{\prime}$ is the quantity

$\displaystyle
W_{2}(A,B\,|\,\rho,\rho^{\prime}):=\frac{\operatorname{Tr}[A\,\rho\,B\,\rho^{\prime}]}{\operatorname{Tr}[\rho\rho^{\prime}]}\,.$
(17)

The double weak value generalizes the standard weak value (2), whose
expression can be retrieved from Eq. (17) by setting $B=I$; in short, one has

$\displaystyle
W(A\,|\,\rho,\rho^{\prime})=W_{2}(A,I\,|\,\rho,\rho^{\prime})\qquad\forall
A,\rho,\rho^{\prime}\in\mathsf{Lin}(\mathcal{H})\,.$ (18)

The double weak value has an interesting physical interpretation. For pure
states $\rho=|\psi\rangle\langle\psi|:=\psi$ and
$\rho^{\prime}=|\psi^{\prime}\rangle\langle\psi^{\prime}|:=\psi^{\prime}$, one
has

$\displaystyle
W_{2}(A,B\,|\,\psi,\psi^{\prime})=W(A|\psi,\psi^{\prime})~{}W(B|\psi^{\prime},\psi)\,,$
(19)

meaning that the double weak value is the product between the weak value of
$A$ in the forward time direction, with initial state $\psi$ and final state
$\psi^{\prime}$, and the weak value of $B$ in the backward time direction,
with initial state $\psi^{\prime}$ and final state $\psi$.

For general mixed states
$\rho=\sum_{m}\,r_{m}\,|\psi_{m}\rangle\langle\psi_{m}|$ and
$\rho^{\prime}=\sum_{n}\,r^{\prime}_{n}\,|\psi^{\prime}_{n}\rangle\langle\psi^{\prime}_{n}|$,
the double weak value quantifies the correlations between the weak values of
$A$ in the forward time direction and the weak values of $B$ in the backward
time direction; in formula:

$\displaystyle
W_{2}(A,B\,|\,\rho,\rho^{\prime})=\sum_{m,n}\,p(m,n)\,W(A|\psi_{m},\psi_{n}^{\prime})~{}W(B|\psi_{n}^{\prime},\psi_{m})\,,$
(20)

where $p(m,n)$ is the probability distribution

$\displaystyle
p(m,n):=\frac{r_{m}r^{\prime}_{n}\,|\langle\psi_{m}|\psi_{n}^{\prime}\rangle|^{2}}{\sum_{i,j}r_{i}r^{\prime}_{j}\,|\langle\psi_{i}|\psi_{j}^{\prime}\rangle|^{2}}\,.$
(21)

Recalling that the weak values were originally defined in terms of weak
measurement processes, the double weak value in Eq. (20) can be interpreted as
the correlation between the quantities observed in two weak measurement
processes where the roles of the pre- and post-selections are exchanged. These
correlations would appear in an exotic scenario where two different agents
operate in two opposite time directions Chiribella and Liu (2022), with one
agent preparing inputs in the past and selecting outputs in the future, and
the other agent preparing inputs in the future and selecting outputs in the
past. In this setting, the double weak value is the correlation between the
values observed by the two agents. It is quite remarkable that such exotic
correlations can be experimentally observed through the controlled-SWAP
protocol.

###  III.2 Local measurements and general bipartite states

A further generalization of Protocol II is to replace the two uncorrelated
states $\rho$ and $\rho^{\prime}$ with a single bipartite state $\rho_{\rm
in,fin}$. If systems $S_{1}$ and $S_{2}$ are measured with POVMs $(P_{j})_{j}$
and $(Q_{k})_{k}$, respectively, the controlled-SWAP protocol produces a
triple of outcomes $(j,k,c)$ distributed with probability

$\displaystyle
p(j,k,c|\rho,\rho^{\prime}):=\operatorname{Tr}\left[(P_{j}\otimes Q_{k}\otimes
R_{c})\,U(\rho_{\rm in,fin}\otimes|+\rangle\langle+|)U^{\dagger}\right]\,.$
(22)

By sampling over this probability distribution, one can simulate the averages
of arbitrary random variables with respect to the complex measure

$\displaystyle q(j,k\,|\,\rho_{\rm in,fin}):=\operatorname{Tr}[(P_{j}\otimes
Q_{k})\,\rho_{\rm in,fin}\,\tt SWAP]\,,$ (23)

thereby estimating quantities of the form

$\displaystyle\operatorname{Tr}[(A\otimes B)\,\rho_{\rm in,fin}\,\tt SWAP]$
(24)

for arbitrary observables $A$ and $B$ in the linear span of $\\{P_{i}\\}_{i}$
and $\\{Q_{j}\\}_{j}$, respectively.

As a biproduct, one can also estimate the ratios

$\displaystyle\frac{\operatorname{Tr}[(A\otimes B)\,\rho_{\rm in,fin}\,\tt
SWAP]}{\operatorname{Tr}[\rho_{\rm in,fin}\,\tt SWAP]}\,.$ (25)

Note that we do not use the term “(double) weak values” for the ratios in Eq.
(25). This omission is intentional: as we will see later in the paper, the
proper notion of (double) weak value does not, in general, coincide with the
ratios in Eq. (25), but rather with a variant of Eq. (25) where the bipartite
state $\rho_{\rm in,fin}$ is subject to a partial transposition.

###  III.3 General bipartite states and joint measurements

A third generalization of the controlled-SWAP protocol consists in performing
a joint measurement on the two copies of the system. If the two copies are
initially in the bipartite state $\rho_{\rm in,fin}$ and are measured with the
joint POVM $(\Pi_{j})_{j}$ after the controlled-SWAP operation, the protocol
produces a pair of outcomes $(j,c)$ distributed with probability

$\displaystyle
p(j,c|\rho,\rho^{\prime}):=\operatorname{Tr}\left[(\Pi_{j}\otimes
R_{c})\,U(\rho_{\rm in,fin}\otimes|+\rangle\langle+|)U^{\dagger}\right]\,.$
(26)

By sampling from this probability distribution, an experimenter can estimate
the expectation value of arbitrary random variables with respect to the
complex measure

$\displaystyle q(j|\rho)=\operatorname{Tr}[\Pi_{j}\,\rho_{\rm in,fin}\,{\tt
SWAP}]\,,$ (27)

and therefore the value of every quantity of the form

$\displaystyle\operatorname{Tr}[C\,\rho_{\rm in,fin}\,{\tt SWAP}]\,,$ (28)

where $C\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$ is an arbitrary
observable in the linear span of the POVM operators $\\{\Pi_{j}\\}$. Again,
this allows one to estimate the ratios

$\displaystyle\frac{\operatorname{Tr}[C\,\rho_{\rm in,fin}\,\tt
SWAP]}{\operatorname{Tr}[\rho_{\rm in,fin}\,\tt SWAP]}.$ (29)

In the following sections, we will compare the experimentally accessible
quantities (28) and (29) with the weak values associated to general two-time
states.

##  IV Two-time states

Weak values can be interpreted as expectation values with respect to a
generalized type of quantum states, known as two-time states Aharonov  _et
al._ (1988); Aharonov and Vaidman (1991); Silva  _et al._ (2014). In the
following, we provide a brief review of the notion of two-time state.

The prototype of a two-time state, introduced in the seminal work of Aharonov,
Bergmann, and Lebowitz Aharonov  _et al._ (1964), is the linear functional
$\lambda:\mathsf{Lin}(\mathcal{H})\to\mathbb{C}$ defined by the relation

$\displaystyle\lambda(A):=\langle\psi^{\prime}|A|\psi\rangle\,,\qquad\forall
A\in\mathsf{Lin}(A)\,,$ (30)

where $|\psi\rangle\in\mathcal{H}$ and $|\psi^{\prime}\rangle\in\mathcal{H}$
are two arbitrary vectors. The complex number $\lambda(A)$ is then interpreted
as the expectation value of the observable $A$ on the two-time state
$\lambda$. The functional $\lambda$ plays the role of an unnormalized state
vector in textbook quantum mechanics. In the following we will call $\lambda$
a two-time vector.

The notion of two-time vector was later generalized by Aharonov and Vaidman
Aharonov and Vaidman (1991), who considered general linear combinations of the
form

$\displaystyle\lambda(A):=\sum_{i=1}^{N}\langle\psi^{\prime}_{i}|\,A\,|\psi_{i}\rangle\,,$
(31)

where $N$ is a positive integer and
$\\{|\psi_{i}\rangle\\}_{i=1}^{N}\subset\mathcal{H}$ and
$\\{|\psi^{\prime}_{i}\rangle\\}_{i=1}^{N}\subset\mathcal{H}$ are arbitrary
sets of vectors.

Mathematically, the set of two-time vectors (31) is the set of all linear
functionals on the observables of the system.

######  Proposition 2.

The set of all two-time vectors (31) is the vector space consisting of all
linear functionals from $\mathsf{Lin}(\mathcal{H})$ to $\mathbb{C}$.

Proof. For every linear functional
$\lambda:\mathsf{Lin}(\mathcal{H})\to\mathbb{C}$, there exists one and only
one matrix $L\in\mathsf{Lin}(\mathcal{H})$ such that
$\lambda(A)=\operatorname{Tr}[L\,A]$. By the singular value decomposition
Datta (1995), the matrix $L$ can be rewritten as
$L=\sum_{i=}^{d}\,|\psi_{i}\rangle\langle\psi_{i}^{\prime}|$, where
$\\{|\psi_{i}\rangle\\}$ and $\\{|\psi_{i}^{\prime}\rangle\\}$ are two sets of
orthogonal states. Hence, the functional $\lambda$ has the Aharonov-Vaidman
form (31). This proves that every linear functional is a valid two-time
vector. The converse is trivial, since every two-time vector is, by
definition, a linear functional. $\blacksquare$

The general notion of two-time vector (31) has been extended from pure to
mixed states in the work by Silva et al Silva  _et al._ (2014), who developed
a framework for describing pre- and post-selected ensembles of quantum states.
We now review this framework, using a slightly different notation that
facilitates the connection with the notion of double weak value introduced in
this paper.

The transition from pure to mixed two-time states is similar to the transition
from state vectors to density matrices.

######  Definition 2.

The two-time matrix corresponding to a two-time vector
$\lambda:\mathsf{Lin}(\mathcal{H})\to\mathbb{C}$ is the linear functional
$E_{\lambda}:\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})\to\mathbb{C}$
uniquely defined by the relation

$\displaystyle E_{\lambda}(A\otimes
B):=\lambda(A)\,\lambda^{\dagger}(B)\,,\qquad\forall
A,B\in\mathsf{Lin}(\mathcal{H})\,,$ (32)

where $\lambda^{\dagger}$ is the functional defined by the relation
$\lambda^{\dagger}(A):=\overline{\lambda(A^{\dagger})}\,,\forall
A\in\mathsf{Lin}(\mathcal{H})$.

######  Definition 3.

A (generally unnormalized) two-time density matrix is linear functional
$\omega:\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})\to\mathbb{C}$ of the form

$\displaystyle\omega=\sum_{n}\,E_{\lambda_{n}}\,,$ (33)

where $(\lambda_{n})_{n}$ are arbitrary linear functionals and
$(E_{\lambda_{n}})_{n}$ are two-time matrices defined as in Eq. (32).

In the following, two-time density matrices will also be called two-time
states. Mathematically, the set of all unnormalized two-time states is a
convex cone, i.e. it contains all convex combinations and all positive
multiples of its elements. Silva et al Silva  _et al._ (2014) showed that the
(generally unnormalized) two-time states are in one-to-one correspondence with
(generally unnormalized) bipartite density matrices. In the next section we
will make this correspondence explicit, providing an expression that connects
two-time states to the controlled-SWAP protocol introduced earlier in the
paper.

##  V Matrix representation of two-time states

We now provide an explicit matrix representation of two-time states. To this
purpose, we recall that every linear functional
$\omega:\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})\to\mathbb{C}$ is of the
form

$\displaystyle\omega(C):=\operatorname{Tr}[O\,C]$ (34)

for a suitable matrix $O\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$. Our
goal is to characterize the constraints that the matrix $O$ has to satisfy in
order for the functional $\omega$ to be a valid two-time state. This goal is
achieved by the following theorem:

######  Theorem 3.

The functional $\omega$ in Eq. (34) is an unnormalized two-time state if and
only if the matrix $O$ is of the form

$\displaystyle O=P^{T_{2}}\,{\tt SWAP}\,,$ (35)

where $P\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$ is a positive
semidefinite matrix and $P^{T_{2}}$ is the partial transpose of $P$ over the
second Hilbert space.

###### Proof.

See Appendix D. $\blacksquare$

Theorem 3 gives an explicit matrix representation of the two-time states: in
short, a two-time state must be a functional of the form

$\displaystyle\omega(C)=\operatorname{Tr}\left[C\,P^{T_{2}}\,{\tt
SWAP}\right]\,,$ (36)

for some positive matrix $P$.

Eq. (36) shows in an explicit way the correspondence between the convex cone
of two-time states and the convex cone of positive bipartite matrices. We now
analyze the correspondence further, by characterizing the structure of the
normalized two-time states, defined as follows

######  Definition 4.

A normalized two-time state is a two-time state $\omega$ such that
$\omega(I)=1$, where $I$ is the identity operator. The set of normalized two-
time states will be denoted by ${\sf T}_{*}(\mathcal{H}\otimes\mathcal{H})$.

We now provide a characterization of the normalized states, showing that they
are in one-to-one correspondence with a convex subset of bipartite density
matrices.

######  Theorem 4.

A linear functional
$\omega:\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})\to\mathbb{C}$ is a
normalized two-time state if and only if it is of the form

$\displaystyle\omega(C)=\frac{\operatorname{Tr}[C\,\rho^{T_{2}}\,{\tt
SWAP}]}{\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]}:=\omega_{\rho}(C)\,,$
(37)

where $\rho\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$ is a normalized
density matrix such that $\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]\not=0$.

Proof. Suppose that $\omega$ is a functional of the form (37). Clearly,
$\omega$ satisfies the normalization condition $\omega(I)=1$. Moreover,
$\omega$ is of the form $\omega(C)=\operatorname{Tr}[C\,P^{T_{2}}\,{\tt
SWAP}]$ with $P:=\rho/\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]$. Note that
$P$ is positive, because $\rho$ is positive and

$\displaystyle\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]=\sum_{i,j}\langle
i|\langle i|\rho|j\rangle|j\rangle\,,$ (38)

which implies $\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]>0$. Hence, $\omega$
is of the form (36), and therefore it is a (normalized) two-time state.

Conversely, suppose that $\omega$ is a normalized two-time state. Since
$\omega$ is a two-time state, it must be of the form (36) for some positive
matrix $P$. Since $\omega$ is normalized, one has
$\operatorname{Tr}[P^{T_{2}}\,{\tt SWAP}]=\omega(I)=1$, which implies in
particular $P\not=0$, and, since $P$ is positive,
$\operatorname{Tr}[P]\not=0$. We can then define a normalized density matrix
$\rho:=P/\operatorname{Tr}[P]$ satisfying the condition
$\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]=\operatorname{Tr}[P^{T_{2}}\,{\tt
SWAP}]/\operatorname{Tr}[P]=1/\operatorname{Tr}[P]\not=0$. Hence, we have

$\displaystyle\omega(C)$ $\displaystyle=\operatorname{Tr}[C\,P^{T_{2}}\,{\tt
SWAP}]$
$\displaystyle=\operatorname{Tr}[P]\,\operatorname{Tr}[C\,\rho^{T_{2}}\,{\tt
SWAP}]$ $\displaystyle=\operatorname{Tr}[C\,\rho^{T_{2}}\,{\tt
SWAP}]/\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]\,.$ (39)

$\blacksquare$

Note that the set of bipartite density matrices satisfying the condition
$\operatorname{Tr}[\rho^{T_{2}}\,{\tt SWAP}]\not=0$ can be equivalently
characterized as the set of bipartite density matrices that have non-zero
overlap with the canonical unnormalized maximally entangled state
$|\Phi\rangle:=\sum_{i=1}^{d}|i\rangle\otimes|i\rangle$, as one can see from
Eq. (38). We denote this set by

$\displaystyle{\sf D}_{*}(\mathcal{H}\otimes\mathcal{H}):=\Big{\\{}\rho\geq
0\,,\operatorname{Tr}[\rho]=1\,,\langle\Phi|\rho|\Phi\rangle>0\Big{\\}}\,.$
(40)

Note also that the correspondence between two-time states and bipartite
matrices is a homeomorphism (that is, it is invertible and continuous), as it
is given by the linear fractional transformation $\rho\mapsto\omega_{\rho}$ in
Eq. (37). Summarizing, we have the following

######  Corollary 1.

The sets of two-time states and bipartite density matrices with nonzero
overlap with the maximally entangled state are homeomorphic; in formula, ${\sf
T}_{*}(\mathcal{H}\otimes\mathcal{H})\simeq{\sf
D}_{*}(\mathcal{H}\otimes\mathcal{H})$.

Note that the mapping $\rho\mapsto\omega_{\rho}$ in Eq. (37) is non-linear and
it does not preserve convex combinations. Nevertheless, it is linear
fractional, and therefore it maps convex combinations into convex
combinations, although with generally different weights. Explicitly, a convex
combination of density matrices, say $\rho=\sum_{i}\,p_{i}\,\rho_{i}$, is
mapped into a convex combination of two-time states
$\omega_{\rho}=\sum_{i}q_{i}\,\omega_{\rho_{i}}$, with
$q_{i}=p_{i}\,\langle\Phi|\rho_{i}\,|\Phi\rangle/(\sum_{j}p_{j}\,\langle\Phi|\rho_{j}\,|\Phi\rangle)$.
This condition implies that the mapping $\rho\mapsto\omega_{\rho}$ maps pure
bipartite states into extreme points of the set of two-time states, and vice-
versa.

##  VI Two-time states and the control-SWAP protocol

Theorem 4 reveals a fundamental connection between two-time states and the
control-SWAP protocol presented earlier in the paper. As we saw in Eq. (29),
the control-SWAP protocol allows an experimenter to estimate any quantity of
the form

$\displaystyle\widetilde{\omega}_{\rho}(C):=\frac{\operatorname{Tr}[C\,\rho\,{\tt
SWAP}]}{\operatorname{Tr}[\rho\,{\tt SWAP}]}\,,$ (41)

for every bipartite observable $C$ and every bipartite quantum state $\rho$.

The difference between the quantities (41) and the expectation values (37) is
the presence of the partial transpose on the second Hilbert space. Physically,
one can interpret the partial transpose as the signature of the difference
between the spatial correlations accessible with the control-SWAP protocol and
the time correlations associated to two-time states.

By comparing Eqs. (41) and (37) we can also have a clear view of the strengths
and limitations of the control-SWAP protocol. First, if the density matrix
$\rho$ is invariant under partial transpose, an experimentalist who has access
to two systems in the joint state $\rho$ can directly estimate the expectation
values on the two-time state $\omega_{\rho}$, by using the control-SWAP
protocol.

More generally if the density matrix is positive under partial transpose (PPT)
Peres (1996); Horodecki  _et al._ (1996), the control-SWAP protocol can
provide an estimate of the expectation values on the two time-state
$\omega_{\rho}$, if the experimenter is given access to the quantum state
$\rho^{T_{2}}$. This is the case of the protocols shown earlier in the paper,
where we saw how to estimate the weak values (1) and (2), which correspond to
the expectation values of the observable $A\otimes I$ on the two-time states
$\omega_{\psi\otimes(\psi^{\prime})^{T}}$ and
$\omega_{\rho\otimes(\rho^{\prime})^{T}}$, respectively. Similarly, the double
weak value $\operatorname{Tr}[\rho^{\prime}A\rho
B]/\operatorname{Tr}[\rho^{\prime}\rho]$ is the expectation of the product
observable $A\otimes B$ with respect to the two-time state
$\omega_{\rho\otimes(\rho^{\prime})^{T}}$.

In constrast, the control-SWAP protocol does not provide, in general, an
estimate of the expectation values on a two-time state when the density matrix
$\rho$ is not PPT. One way to circumvent this problem would be to approximate
the partial transpose operation with a physical process, such as the optimal
universal transpose map Horodecki (2001); Buscemi  _et al._ (2003). However,
the estimation protocol resulting from this approach would be dimension-
dependent: indeed, optimal universal transpose is the completely positive
trace-preserving map given by

$\displaystyle\mathscr{T}(\rho)=\frac{\rho^{\rm T}+I}{d+1}\,,$ (42)

and the $1/(d+1)$ factor implies that the sample complexity for the estimation
of the expectation values $\omega_{\rho}(C)$ grows linearly with $d$.

Another approach is to introduce an auxiliary system, and to reduce two-time
states associated to non-PPT density matrices to extended two-time states
associated to PPT density matrices. This approach can be easily illustrated in
the pure state case. In this case, Aharonov and Vaidman Aharonov and Vaidman
(1991) showed that every two-time vector
$\lambda(A):=\sum_{i=1}^{N}\langle\psi^{\prime}_{i}|\,A\,|\psi_{i}\rangle$ can
be obtained by introducing an auxiliary system of dimension $N$, and a joint
two-time vector $\Lambda(\cdot):=\langle\Psi^{\prime}|\,\cdot\,|\Psi\rangle$
corresponding to the initial vector
$|\Psi\rangle=\sum_{i}\,|\psi_{i}\rangle\otimes|i\rangle\in\mathcal{H}\otimes\mathcal{H}_{\rm
aux}$ and final vector
$|\Psi^{\prime}\rangle=\sum_{i}\,|\psi_{i}^{\prime}\rangle\otimes|i\rangle\in\mathcal{H}\otimes\mathcal{H}_{\rm
aux}$, where $\mathcal{H}_{\rm aux}$ is the Hilbert space of the auxiliary
system and $\\{|i\rangle\\}_{i=1}^{N}$ is an orthonormal basis for
$\mathcal{H}_{\rm aux}$. With this definition, one has the relation

$\displaystyle\lambda(A)=\Lambda(A\otimes I_{\rm aux})\,,\qquad\forall
A\in\mathsf{Lin}(\mathcal{H})\,,$ (43)

meaning that the expectation value of the observable $A$ on the two-time
vector $\lambda$ coincides with the expectation value of the observable
$A\otimes I_{\rm aux}$ on the the joint two-time vector $\Lambda$. In turn,
the two-time vector $\Lambda$ consists just in a pre-selection and a post-
selection to pure states, and therefore it can be reproduced by a product (and
therefore PPT) state. Explicitly, we have

$\displaystyle\frac{\Lambda(A\otimes I_{\rm aux})}{\Lambda(I\otimes I_{\rm
aux})}=\omega_{\Psi\otimes(\Psi^{\prime})^{T}}(A\otimes I_{\rm aux}\otimes
I\otimes I_{\rm aux})\,,$ (44)

where $\omega_{\Psi\otimes(\Psi^{\prime})^{T}}$, defined as in Eq. (37), is a
normalized two-time state in ${\sf T}_{*}(\mathcal{H}\otimes\mathcal{H}_{\rm
aux}\otimes\mathcal{H}\otimes\mathcal{H}_{\rm aux})$. Combining Eqs. (43) and
(44), we can see that the weak values $\lambda(A)/\lambda(I)$ can be estimated
by the control-SWAP protocol, if the experimenter has access to the pure
quantum states states proportional to the vectors $|\Psi\rangle$ and
$|\overline{\Psi}^{\prime}\rangle$.

The above argument can be extended from two-time vectors to two-time matrices.
For example, the expectation values with respect to the two-time matrix
$E_{\lambda}$ associated to the two-time vector $\lambda$ in the previous
paragraph can be computed as

$\displaystyle\frac{E_{\lambda}(C)}{E_{\lambda}(I)}=\omega_{\Psi\otimes(\Psi^{\prime})^{T}}\Big{(}W\,(C\otimes
I_{\rm aux}\otimes I_{\rm aux})\,W^{\dagger}\Big{)}\,,$ (45)

where $W$ is the unitary operator that permutes the Hilbert spaces
$\mathcal{H}\otimes\mathcal{H}\otimes\mathcal{H}_{\rm
aux}\otimes\mathcal{H}_{\rm aux}$ into $\mathcal{H}\otimes\mathcal{H}_{\rm
aux}\otimes\mathcal{H}\otimes\mathcal{H}_{\rm aux}$. Eq. (45) means that the
expectation value on the l.h.s. can be estimated by the controlled-SWAP
protocol, if the experimenter has access to the pure quantum states
proportional to the vectors $|\Psi\rangle$ and
$|\overline{\Psi}^{\prime}\rangle$.

Overall, the extension approach has the benefit of being dimension-
independent, but it requires either the two-time state $\omega$ to be known,
or the appropriate PPT state of the system and the auxiliary to be provided to
the experimenter.

##  VII Conclusions and outlook

In this paper we have designed a dimension-independent scheme for estimating
weak values of arbitrary observables. The scheme is based on the controlled-
SWAP gate, and generates a probability distribution that can be used to sample
from a complex measure underlying the weak values of interest. Crucially, the
initial and final states of the system as well as the observable itself, whose
weak value is estimated, do not have to be a priori known to experimenter.

As a byproduct, the structure of the controlled-SWAP scheme provides several
insights to theory of two-time states. In particular, we have derived an
alternative expression for two-time states, which provides an explicit
characterization of the correspondence between two-time states and (a subset
of) bipartite density matrices. Using this expression, we characterized the
domain of applicability of the control-SWAP protocol, showing that it can be
used to estimate the expectation values of all two-time states corresponding
to PPT density matrices (including, of course, the product density matrices
corresponding to the usual definition of weak values). For two-time states
corresponding to non-PPT density matrices, the controlled-SWAP protocol can
still be used if the experimenter is given access to an extendend quantum
states involving a pair of auxiliary systems.

An interesting open question for future research is whether dimension-
independent sampling schemes like ours could be constructed for the estimation
of Kirkwood-Dirac (KD) quasiprobability distributions, an important type of
complex measures that often arise in quantum information and foundations
Dressel (2015); Yunger Halpern  _et al._ (2018); Lostaglio  _et al._ (2023).
KD distributions have a close connection with weak values, and their
negativity provides a signature of quantum contextuality in Spekkens’
formulation Spekkens (2005). In turn, KD distributions have found numerous
applications beyond quantum foundations, including quantum metrology
Arvidsson-Shukur  _et al._ (2020); Jenne and Arvidsson-Shukur (2021); Lupu-
Gladstein  _et al._ (2022), condensed matter physics Yunger Halpern  _et al._
(2018); González Alonso  _et al._ (2019, 2022), and thermodynamics Lostaglio
_et al._ (2023). The approach of Wagner et al Wagner  _et al._ (2023) provides
a way to estimate the value of the KD distribution at every fixed point. The
open question is whether there exist ways to simulate sampling from the KD
distribution.

Another interesting direction concerns the theory of two-time states and its
relation of the study of causality in quantum theory Chiribella  _et al._
(2009, 2013); Oreshkov  _et al._ (2012); Č. Brukner (2014); Barrett  _et al._
(2019, 2021). Recent findings have suggested that anomalous weak values of
observables shared between several parties can assist in witnessing the causal
relationships between the parties’ laboratories Abbott  _et al._ (2019).
Moreover, two-time states (and multiple-time states in general) themselves can
carry information about the underlying causal structure Silva  _et al._ (2017)
and can be used to witness indefinite causal order of operations Liu  _et al._
(2023). Finally, protocols exploiting controlled causal order of operations in
a protocol known as the quantum SWITCH Chiribella  _et al._ (2009, 2013) have
been shown to be useful for efficient estimation of quantities that can be
connected with weak values and KD distributions, such as out-of-time-
correlators Swingle  _et al._ (2016) and incompatibility of quantum
observables Gao  _et al._ (2023). Therefore, an interplay between weak values,
KD distributions, and indefinite causal order in quantum SWITCH appears as a
promising direction for future investigations.

###### Acknowledgements.

We thank Nicole Yunger Halpern and Billy F. Braasch, Jr. (Quantum Steampunk
Laboratory, University of Maryland) for many fruitful discussions and useful
comments that helped us improve the manuscript. This work was supported by
funding from the Hong Kong Research Grant Council through grant no. 17307520
and through the Senior Research Fellowship Scheme SRFS2021-7S02, and from the
John Templeton Foundation through grant 62312, The Quantum Information
Structure of Spacetime (qiss.fr). The opinions expressed in this publication
are those of the authors and do not necessarily reflect the views of the John
Templeton Foundation. Research at the Perimeter Institute is supported by the
Government of Canada through the Department of Innovation, Science and
Economic Development Canada and by the Province of Ontario through the
Ministry of Research, Innovation and Science.

## References

  * Aharonov  _et al._ (1988) Y. Aharonov, D. Z. Albert, and L. Vaidman, “How the result of a measurement of a component of the spin of a spin-$\frac{1}{2}$ particle can turn out to be $100$,” Phys. Rev. Lett. 60, 1351 (1988).
  * Aharonov  _et al._ (2010) Y. Aharonov, S. Popescu, and J. Tollaksen, “A time-symmetric formulation of quantum mechanics,” Phys. Today 63, 27 (2010).
  * Shikano (2012) Y. Shikano,  _Theory of “weak value” and quantum mechanical measurements_ (InTech Shanghai and Rijeka, 2012).
  * Kofman  _et al._ (2012) A. G. Kofman, S. Ashhab, and F. Nori, “Nonperturbative theory of weak pre-and post-selected measurements,” Phys. Rep. 520, 43 (2012).
  * Dressel  _et al._ (2014) J. Dressel, M. Malik, F. M. Miatto, A. N. Jordan, and R. W. Boyd, “Colloquium: Understanding quantum weak values: Basics and applications,” Rev. Mod. Phys. 86, 307–316 (2014).
  * Aharonov  _et al._ (2002) Y. Aharonov, A. Botero, S. Popescu, B. Reznik, and J. Tollaksen, “Revisiting Hardy’s paradox: counterfactual statements, real measurements, entanglement and weak values,” Phys. Lett. A 301, 130 (2002).
  * Lundeen and Steinberg (2009) J. S. Lundeen and A. M. Steinberg, “Experimental joint weak measurement on a photon pair as a probe of Hardy’s paradox,” Phys. Rev. Lett. 102, 020404 (2009).
  * Kocsis  _et al._ (2011) S. Kocsis, B. Braverman, S. Ravets, M. J. Stevens, R. P. Mirin, L. K. Shalm, and A. M. Steinberg, “Observing the Average Trajectories of Single Photons in a Two-Slit Interferometer,” Science 332, 1170 (2011).
  * Williams and Jordan (2008) N. S. Williams and A. N. Jordan, “Weak Values and the Leggett-Garg Inequality in Solid-State Qubits,” Phys. Rev. Lett. 100, 026804 (2008).
  * Dressel  _et al._ (2011) J. Dressel, C. J. Broadbent, J. C. Howell, and A. N. Jordan, “Experimental violation of two-party Leggett-Garg inequalities with semiweak measurements,” Phys. Rev. Lett. 106, 040402 (2011).
  * Goggin  _et al._ (2011a) M. E. Goggin, M. P. Almeida, M. Barbieri, B. P. Lanyon, J. L. O’Brien, A. G. White, and G. J. Pryde, “Violation of the Leggett–Garg inequality with weak measurements of photons,” Proc. Natl. Acad. Sci. 108, 1256 (2011a).
  * Groen  _et al._ (2013) J. P. Groen, D. Ristè, L. Tornberg, J. Cramer, P. C. de Groot, T. Picot, G. Johansson, and L. DiCarlo, “Partial-measurement backaction and nonclassical weak values in a superconducting circuit,” Phys. Rev. Lett. 111, 090506 (2013).
  * Pusey (2014) M. F. Pusey, “Anomalous Weak Values Are Proofs of Contextuality,” Phys. Rev. Lett. 113, 200401 (2014).
  * Dressel (2015) J. Dressel, “Weak values as interference phenomena,” Phys. Rev. A 91, 032116 (2015).
  * Kunjwal  _et al._ (2019) R. Kunjwal, M. Lostaglio, and M. F. Pusey, “Anomalous weak values and contextuality: Robustness, tightness, and imaginary parts,” Phys. Rev. A 100, 042116 (2019).
  * Hosten and Kwiat (2008) O. Hosten and P. Kwiat, “Observation of the Spin Hall Effect of Light via Weak Measurements,” Science 319, 787 (2008).
  * Dixon  _et al._ (2009) P. B. Dixon, D. J. Starling, A. N. Jordan, and J. C. Howell, “Ultrasensitive Beam Deflection Measurement via Interferometric Weak Value Amplification,” Phys. Rev. Lett. 102, 173601 (2009).
  * Brunner and Simon (2010) N. Brunner and C. Simon, “Measuring Small Longitudinal Phase Shifts: Weak Measurements or Standard Interferometry?” Phys. Rev. Lett. 105, 010405 (2010).
  * Starling  _et al._ (2010a) D. J. Starling, P. B. Dixon, A. N. Jordan, and J. C. Howell, “Precision frequency measurements with interferometric weak values,” Phys. Rev. A 82, 063822 (2010a).
  * Starling  _et al._ (2010b) D. J. Starling, P. B. Dixon, N. S. Williams, A. N. Jordan, and J. C. Howell, “Continuous phase amplification with a Sagnac interferometer,” Phys. Rev. A 82, 011802(R) (2010b).
  * Hofmann (2011) H. F. Hofmann, “Uncertainty limits for quantum metrology obtained from the statistics of weak measurements,” Phys. Rev. A 83, 022106 (2011).
  * Hofmann  _et al._ (2012) H. F. Hofmann, M. E. Goggin, M. P. Almeida, and M. Barbieri, “Estimation of a quantum interaction parameter using weak measurements: Theory and experiment,” Phys. Rev. A 86, 040102(R) (2012).
  * Strübi and Bruder (2013) G. Strübi and C. Bruder, “Measuring ultrasmall time delays of light by joint weak measurements,” Phys. Rev. Lett. 110, 083605 (2013).
  * Lundeen  _et al._ (2011) J. S. Lundeen, B. Sutherland, A. Patel, C. Stewart, and C. Bamber, “Direct measurement of the quantum wavefunction,” Nature 474, 188 (2011).
  * Lundeen and Bamber (2012) J. S. Lundeen and C. Bamber, “Procedure for direct measurement of general quantum states using weak measurement,” Phys. Rev. Lett. 108, 070402 (2012).
  * Salvail  _et al._ (2013) J. Z. Salvail, M. Agnew, A. S. Johnson, E. Bolduc, J. Leach, and R. W. Boyd, “Full characterization of polarization states of light via direct measurement,” Nature Photon. 7, 316 (2013).
  * Malik  _et al._ (2014) M. Malik, M. Mirhosseini, M. P. J. Lavery, J. Leach, M. J. Padgett, and R. W. Boyd, “Direct measurement of a 27-dimensional orbital-angular-momentum state vector,” Nat. Commun. 5, 3115 (2014).
  * Hofmann (2012) H. F. Hofmann, “How Weak Values Emerge in Joint Measurements on Cloned Quantum Systems,” Phys. Rev. Lett. 109, 020408 (2012).
  * Wagner  _et al._ (2023) R. Wagner, Z. Schwartzman-Nowik, I. L. Paiva, A. Te’eni, A. Ruiz-Molero, R. Soares Barbosa, E. Cohen, and E. F. Galvão, “Quantum circuits measuring weak values and Kirkwood-Dirac quasiprobability distributions, with applications,” (2023), arXiv:2302.00705 .
  * Oszmaniec  _et al._ (2021) M. Oszmaniec, D. J. Brod, and E. F. Galvão, “Measuring relational information between quantum states, and applications,” (2021), arXiv:2109.10006 .
  * Aharonov  _et al._ (1964) Y. Aharonov, P. G. Bergmann, and J. L. Lebowitz, “Time Symmetry in the Quantum Process of Measurement,” Phys. Rev. 134, B1410 (1964).
  * Aharonov and Vaidman (1991) Y. Aharonov and L. Vaidman, “Complete description of a quantum system at a given time,” J. Phys. A: Math. Gen. 24, 10 (1991).
  * Silva  _et al._ (2014) R. Silva, Y. Guryanova, N. Brunner, N. Linden, A. J. Short, and S. Popescu, “Pre- and postselected quantum states: Density matrices, tomography, and Kraus operators,” Phys. Rev. A. 89, 012121 (2014).
  * Peres (1996) A. Peres, “Separability Criterion for Density Matrices,” Phys. Rev. Lett. 77, 1413 (1996).
  * Horodecki  _et al._ (1996) M. Horodecki, P. Horodecki, and R. Horodecki, “Separability of mixed states: necessary and sufficient conditions,” Phys. Lett. A 223, 1 (1996).
  * Vaidman  _et al._ (2017) L. Vaidman, A. Ben-Israel, J. Dziewior, L. Knips, M. Weissl, J. Meinecke, C. Schwemmer, R. Ber, and H. Weinfurter, “Weak value beyond conditional expectation value of the pointer readings,” Phys. Rev. A 96, 032114 (2017).
  * Ritchie  _et al._ (1991) N. W. M. Ritchie, J. G. Story, and R. G. Hulet, “Realization of a measurement of a “weak value”,” Phys. Rev. Lett. 66, 1107 (1991).
  * Pryde  _et al._ (2005) G. J. Pryde, J. L. O’Brien, A. G. White, T. C. Ralph, and H. M. Wiseman, “Measurement of quantum weak values of photon polarization,” Phys. Rev. Lett. 94, 220405 (2005).
  * Goggin  _et al._ (2011b) M. E. Goggin, M. P. Almeida, M. Barbieri, B. P. Lanyon, J. L. O’Brien, A. G. White, and G. J. Pryde, “Violation of the Leggett-Garg inequality with weak measurements of photons,” Proc. Natl. Acad. Sci. USA 108, 1256 (2011b).
  * Li  _et al._ (2011) C.-F. Li, X.-Y. Xu, J.-S. Tang, J.-S. Xu, and G.-C. Guo, “Ultrasensitive phase estimation with white light,” Phys. Rev. A 83, 044102 (2011).
  * Dressel and Jordan (2012) J. Dressel and A. N. Jordan, “Significance of the imaginary part of the weak value,” Phys. Rev. A 85, 012107 (2012).
  * Xu  _et al._ (2013) X.-Y. Xu, Y. Kedem, K. Sun, L. Vaidman, C.-F. Li, and G.-C. Guo, “Phase Estimation with Weak Measurement Using a White Light Source,” Phys. Rev. Lett. 111, 033604 (2013).
  * Magaña-Loaiza  _et al._ (2014) O. S. Magaña-Loaiza, M. Mirhosseini, B. Rodenburg, and R. W. Boyd, “Amplification of angular rotations using weak measurements,” Phys. Rev. Lett. 112, 200401 (2014).
  * Piacentini  _et al._ (2016) F. Piacentini, A. Avella, M. P. Levi, R. Lussana, F. Villa, A. Tosi, F. Zappa, M. Gramegna, G. Brida, I. P. Degiovanni, and M. Genovese, “Experiment Investigating the Connection between Weak Values and Contextuality,” Phys. Rev. Lett. 116, 180401 (2016).
  * Vaidman (2017) L. Vaidman, “Weak value controversy,” Phil. Trans. R. Soc. A 375, 20160395 (2017).
  * Cohen and Pollak (2018) E. Cohen and E. Pollak, “Determination of weak values of quantum operators using only strong measurements,” Phys. Rev. A 98, 042112 (2018).
  * Rebufello  _et al._ (2021) E. Rebufello, F. Piacentini, A. Avella, M. A. de Souza, M. Gramegna, J. Dziewior, E. Cohen, L. Vaidman, I. P. Degiovanni, and M. Genovese, “Anomalous weak values via a single photon detection,” Light Sci. Appl. 10, 106 (2021).
  * Adam and Diósi (2023) P. Adam and L. Diósi, “Sequential unsharp measurement of photon polarization,” Phys. Rev. A 107, 063706 (2023).
  * Wiseman  _et al._ (2023) H. M. Wiseman, A. M. Steinberg, and M. Hallaji, “Obtaining a single-photon weak value from experiments using a strong (many-photon) coherent state,” AVS Quantum Sci. 5, 024401 (2023).
  * Barenco  _et al._ (1997) A. Barenco, A. Berthiaume, D. Deutsch, A. Ekert, R. Jozsa, and C. Macchiavello, “Stabilization of Quantum Computations by Symmetrization,” SIAM J. Comput. 26, 1541 (1997).
  * Buhrman  _et al._ (2001) H. Buhrman, R. Cleve, J. Watrous, and R. de Wolf, “Quantum Fingerprinting,” Phys. Rev. Lett. 87, 167902 (2001).
  * Chiribella and Liu (2022) G. Chiribella and Z. Liu, “Quantum operations with indefinite time direction,” Commun. Phys. 5, 190 (2022).
  * Datta (1995) B. Datta,  _Numerical Linear Algebra and Application_ (Brooks/Cole Publishing Company, 1995).
  * Horodecki (2001) P. Horodecki, “From limits of quantum nonlinear operations to multicopy entanglement witnesses and state spectrum estimation,” (2001), arXiv:quant-ph/0111036 .
  * Buscemi  _et al._ (2003) F. Buscemi, G. M. D’Ariano, P. Perinotti, and M. F. Sacchi, “Optimal realization of the transposition maps,” Phys. Lett. A 314, 374 (2003).
  * Yunger Halpern  _et al._ (2018) N. Yunger Halpern, B. Swingle, and J. Dressel, “Quasiprobability behind the out-of-time-ordered correlator,” Phys. Rev. A 97, 042105 (2018).
  * Lostaglio  _et al._ (2023) M. Lostaglio, A. Belenchia, A. Levy, S. Hernández-Gómez, N. Fabbri, and S. Gherardini, “Kirkwood-Dirac quasiprobability approach to quantum fluctuations: Theoretical and experimental perspectives,” Quantum 7, 1128 (2023).
  * Spekkens (2005) R. W. Spekkens, “Contextuality for preparations, transformations, and unsharp measurements,” Phys. Rev. A 71, 052108 (2005).
  * Arvidsson-Shukur  _et al._ (2020) D. R. M. Arvidsson-Shukur, N. Yunger Halpern, H. V. Lepage, A. A. Lasek, C. H. W. Barnes, and S. Lloyd, “Quantum advantage in postselected metrology,” Nat. Commun. 11, 3775 (2020).
  * Jenne and Arvidsson-Shukur (2021) J. H. Jenne and D. R. M. Arvidsson-Shukur, “Quantum Learnability is Arbitrarily Distillable,” (2021), arXiv:2104.09520 .
  * Lupu-Gladstein  _et al._ (2022) N. Lupu-Gladstein, Y. B. Yilmaz, D. R. M. Arvidsson-Shukur, A. Brodutch, A. O. T. Pang, A. M. Steinberg, and N. Yunger Halpern, “Negative Quasiprobabilities Enhance Phase Estimation in Quantum-Optics Experiment,” Phys. Rev. Lett. 128, 220504 (2022).
  * González Alonso  _et al._ (2019) J. R. González Alonso, N. Yunger Halpern, and J. Dressel, “Out-of-Time-Ordered-Correlator Quasiprobabilities Robustly Witness Scrambling,” Phys. Rev. Lett. 122, 040404 (2019).
  * González Alonso  _et al._ (2022) J. R. González Alonso, N. Shammah, S. Ahmed, F. Nori, and J. Dressel, “Diagnosing quantum chaos with out-of-time-ordered-correlator quasiprobability in the kicked-top model,” (2022), arXiv:2201.08175 .
  * Chiribella  _et al._ (2009) G. Chiribella, G. M. D’Ariano, P. Perinotti, and B. Valiron, “Beyond quantum computers,” (2009), arXiv:0912.0195 .
  * Chiribella  _et al._ (2013) G. Chiribella, G. M. D’Ariano, P. Perinotti, and B. Valiron, “Quantum computations without definite causal structure,” Phys. Rev. A 88, 022318 (2013).
  * Oreshkov  _et al._ (2012) O. Oreshkov, F. Costa, and Č. Brukner, “Quantum correlations with no causal order,” Nat. Commun. 3, 1092 (2012).
  * Č. Brukner (2014) Č. Brukner, “Quantum causality,” Nature Phys. 10, 259 (2014).
  * Barrett  _et al._ (2019) J. Barrett, R. Lorenz, and O. Oreshkov, “Quantum causal models,” (2019), arXiv:1906.10726 .
  * Barrett  _et al._ (2021) J. Barrett, R. Lorenz, and O. Oreshkov, “Cyclic quantum causal models,” Nat. Commun. 12, 885 (2021).
  * Abbott  _et al._ (2019) A. A. Abbott, R. Silva, J. Wechs, N. Brunner, and C. Branciard, “Anomalous Weak Values Without Post-Selection,” Quantum 3, 194 (2019).
  * Silva  _et al._ (2017) R. Silva, Y. Guryanova, A. J. Short, P. Skrzypczyk, N. Brunner, and S. Popescu, “Connecting processes with indefinite causal order and multi-time quantum states,” New J. Phys. 19, 103022 (2017).
  * Liu  _et al._ (2023) X. Liu, Z. Jia, Y. Qiu, F. Li, and O. Dahlsten, “Unification of spatiotemporal quantum formalisms: mapping between process and pseudo-density matrices via multiple-time states,” (2023), arXiv:2306.05958 .
  * Swingle  _et al._ (2016) B. Swingle, G. Bentsen, M. Schleier-Smith, and P. Hayden, “Measuring the scrambling of quantum information,” Phys. Rev. A 94, 040302(R) (2016).
  * Gao  _et al._ (2023) N. Gao, D. Li, A. Mishra, J. Yan, K. Simonov, and G. Chiribella, “Measuring incompatibility and clustering quantum observables with a quantum switch,” Phys. Rev. Lett. 130, 170201 (2023).
  * Hoeffding (1963) W. Hoeffding, “Probability Inequalities for Sums of Bounded Random Variables,” J. Am. Stat. Assoc. 58, 13–30 (1963).

##  Appendix A Derivation of Eq. (11)

At the beginning of the protocol, the two copies of the system ($S_{1}$ and
$S_{2}$) and the control qubit are in the product state

$\omega=\rho\otimes\rho^{\prime}\otimes|+\rangle\langle+|.$ (46)

Application of the controlled SWAP gate $U$ produces the new state

$\displaystyle U\omega U^{\dagger}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\Bigl{\\{}\rho\otimes\rho^{\prime}\otimes|0\rangle\langle
0|+\rho^{\prime}\otimes\rho\otimes|1\rangle\langle 1|$ (47)
$\displaystyle+\;(\rho\otimes\rho^{\prime}){\tt SWAP}\otimes|0\rangle\langle
1|$ $\displaystyle+\;{\tt
SWAP}(\rho\otimes\rho^{\prime})\otimes|1\rangle\langle 0|\Bigr{\\}}\,.$

A measurement of the first system and auxiliary qubit with the POVM
$(P_{j})_{j}$ and $(R_{c})_{c}$, respectively, produces a pair of outcomes
$(j,c)$ distributed with probability

$\displaystyle p(j,c|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{2}\operatorname{Tr}\Bigl{[}P_{j}\rho\otimes\rho^{\prime}\otimes
R_{c}|0\rangle\langle 0|$
$\displaystyle+\;P_{j}\rho^{\prime}\otimes\rho\otimes R_{c}|1\rangle\langle
1|$ $\displaystyle+\;(P_{j}\otimes I)(\rho\otimes\rho^{\prime}){\tt
SWAP}\otimes R_{c}|0\rangle\langle 1|\Bigr{]}$ $\displaystyle+\;(P_{j}\otimes
I){\tt SWAP}(\rho\otimes\rho^{\prime})\otimes R_{c}|1\rangle\langle 0|\,.$

Using the relations

$\displaystyle\langle 0|R_{c}|0\rangle=\langle
1|R_{c}|1\rangle=\frac{1}{4}\qquad\forall c\in\\{0,1,2,3\\}$ (48)

and

$\displaystyle\langle
0|R_{c}|1\rangle=\frac{(-1)^{c}}{4}\,\Bigl{(}\theta(1-c)+i\theta(c-2)\Bigr{)}\qquad\forall
c\in\\{0,1,2,3\\}$ (49)

we then obtain

$\displaystyle p(j,c|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{8}\Bigl{\\{}\operatorname{Tr}[P_{j}\rho]+\operatorname{Tr}[P_{j}\rho^{\prime}]$
(50) $\displaystyle+\;(-1)^{c}\Bigl{(}\theta(1-c)+i\theta(c-2)\Bigr{)}$
$\displaystyle\cdot\;\operatorname{Tr}\Bigl{[}(P_{j}\otimes
I)(\rho\otimes\rho^{\prime}){\tt SWAP}\Bigr{]}$
$\displaystyle+\;(-1)^{c}\Bigl{(}\theta(1-c)-i\theta(c-2)\Bigr{)}$
$\displaystyle\cdot\;\operatorname{Tr}\Bigl{[}(P_{j}\otimes I){\tt
SWAP}(\rho\otimes\rho^{\prime})\Bigr{]}\Bigr{\\}}.$

Finally, using the relation $\operatorname{Tr}[(A\otimes B)\,{\tt
SWAP}]=\operatorname{Tr}[AB]\,,\forall A,B\in\mathsf{Lin}(\mathcal{H})$, we
obtain

$\displaystyle p(j,c|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{8}\Bigl{\\{}\operatorname{Tr}[P_{j}\rho]+\operatorname{Tr}[P_{j}\rho^{\prime}]$
(51)
$\displaystyle+\;(-1)^{c}\theta(1-c)\Bigl{(}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]+\operatorname{Tr}[\rho
P_{j}\rho^{\prime}]\Bigr{)}$
$\displaystyle+\;i(-1)^{c}\theta(c-2)\Bigl{(}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]-\operatorname{Tr}[\rho
P_{j}\rho^{\prime}]\Bigr{)}\Bigr{\\}}$ $\displaystyle=$
$\displaystyle\frac{1}{8}\Bigl{\\{}\operatorname{Tr}[P_{j}\rho]+\operatorname{Tr}[P_{j}\rho^{\prime}]$
$\displaystyle+\;2(-1)^{c}\theta(1-c)\operatorname{Re}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]$
$\displaystyle-\;2(-1)^{c}\theta(c-2)\operatorname{Im}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]\Bigr{\\}},$

where the second equality follows from the identity $\operatorname{Tr}[\rho
P_{j}\rho^{\prime}]=\operatorname{Tr}[\rho^{\prime}\rho
P_{j}]=\operatorname{Tr}[(P_{j}\rho\rho^{\prime})^{\dagger}]=\overline{\operatorname{Tr}[P_{j}\rho\rho^{\prime}]}$.
This concludes the proof of Eq. (11).

##  Appendix B Proof of Lemma 1

The proof of Lemma 1 is a straightforward calculation: we only need to
calculate the expectation value of the random variable $\widetilde{Z}$ in Eq.
(12) with respect to the probability distribution $p(j,c|\rho,\rho^{\prime})$.
Explicitly, the expectation value is

$\displaystyle\mathbb{E}_{p}[\widetilde{Z}]$ $\displaystyle=$
$\displaystyle\sum_{j,c}\,\widetilde{z}_{j,c}\,p(j,c|\rho,\rho^{\prime})$ (52)
$\displaystyle=$ $\displaystyle
2\sum_{j,c}z_{j}(-1)^{c}\Bigl{[}\theta(1-c)-i\theta(c-2)\Bigr{]}p(j,c|\rho,\rho^{\prime})$
$\displaystyle=$
$\displaystyle\frac{1}{4}\sum_{j}z_{j}\Bigl{(}\operatorname{Tr}[P_{j}\rho]+\operatorname{Tr}[P_{j}\rho^{\prime}]\Bigr{)}$
$\displaystyle\cdot\;\sum_{c}(-1)^{c}\Bigl{[}\theta(1-c)-i\theta(c-2)\Bigr{]}$
$\displaystyle+\;\frac{1}{2}\sum_{c}(-1)^{2c}\theta(1-c)\Bigl{[}\theta(1-c)-i\theta(c-2)\Bigr{]}$
$\displaystyle\cdot\;\sum_{j}z_{j}\operatorname{Re}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]$
$\displaystyle-\;\frac{1}{2}\sum_{c}(-1)^{2c}\theta(c-2)\Bigl{[}\theta(1-c)-i\theta(c-2)\Bigr{]}$
$\displaystyle\cdot\;\sum_{j}z_{j}\operatorname{Im}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]$
$\displaystyle=$
$\displaystyle\sum_{j}z_{j}\Bigl{(}\operatorname{Re}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]+i\operatorname{Im}\operatorname{Tr}[P_{j}\rho\rho^{\prime}]\Bigr{)}$
$\displaystyle=$ $\displaystyle\sum_{j}z_{j}q(j|\rho,\rho^{\prime})$
$\displaystyle=$ $\displaystyle\mathbb{E}_{q}[Z]\,.$

In summary, the expectation value of the random variable $\widetilde{Z}$ with
respect to the probability distribution $p(j,c|\rho,\rho^{\prime})$ is equal
to the expectation value of the random variable $Z$ with respect to the
complex distribution $q(j|\rho,\rho^{\prime})$.

##  Appendix C Proof of Theorem 1

The proof is based on two lemmas, provided in the following.

######  Lemma 2.

Let $\rho$ and $\rho^{\prime}$ be a pair of states, let $(P_{j})_{j}$ be a
POVM, and let $\\{x_{j}\\}_{j}$ be a set of real numbers, with
$x_{\max}:=\max_{j}|x_{j}|$. Let $\nu_{\rm Re}:={\rm
Re}(\operatorname{Tr}[\rho^{\prime}A\rho])$ and $\nu_{\rm Im}:={\rm
Im}(\operatorname{Tr}[\rho^{\prime}A\rho])$ be the real and imaginary parts of
the weak value of the observable $A:=\sum_{j}\,x_{j}\,P_{j}$, respectively.
The estimate of $\nu_{\rm Re}$ and $\nu_{\rm Im}$ obtained from $K$ runs of
Protocol II has error at most

$\epsilon_{\nu}=2x_{\rm max}\sqrt{\frac{1}{K}\ln\frac{2}{\delta}}$ (53)

with probability at least $1-\delta$.

###### Proof.

In accordance with (5),
$\operatorname{Tr}[\rho^{\prime}A\rho]=\mathbb{E}_{q}[X]$, where $q$ is the WV
measure defined in (4), and $X$ is a random variable taking values in
$\\{x_{j}\\}_{j}$. Lemma 1 guarantees that the probability distribution (11)
generated by Protocol II can be equivalently used to estimate
$\mathbb{E}_{q}[X]$. In turn, it can be rewritten as

$\displaystyle p(j,c|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{2}\Bigl{(}\theta(1-c)p_{\rm
Re}(j,c|\rho,\rho^{\prime})$ (54) $\displaystyle+\;\theta(c-2)p_{\rm
Im}(j,c-2|\rho,\rho^{\prime})\Bigr{)},$

with probability distributions

$\displaystyle p_{\rm Re}(j,\bar{c}|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{4}\Bigl{\\{}\operatorname{Tr}[P_{j}(\rho+\rho^{\prime})]+2(-1)^{\bar{c}}\operatorname{Re}\operatorname{Tr}[\rho^{\prime}P_{j}\rho]\Bigr{\\}}$
$\displaystyle p_{\rm Im}(j,\bar{c}|\rho,\rho^{\prime})$ $\displaystyle=$
$\displaystyle\frac{1}{4}\Bigl{\\{}\operatorname{Tr}[P_{j}(\rho+\rho^{\prime})]-2(-1)^{\bar{c}}\operatorname{Im}\operatorname{Tr}[\rho^{\prime}P_{j}\rho]\Bigr{\\}},$

where $\bar{c}\in\\{0,1\\}$.

Now, recall that the POVM $(R_{c})_{c}$ in Protocol II can be obtained by
randomly choosing between the projective measurements $(2R_{0},2R_{1})$ and
$(2R_{2},2R_{3})$, which give rise to the probability distributions $p_{\rm
Re}(j,\bar{c}|\rho,\rho^{\prime})$ and $p_{\rm
Im}(j,\bar{c}|\rho,\rho^{\prime})$, respectively. Then, a straightforward
calculation demonstrates that $\mathbb{E}_{p_{\rm Re}}[\widetilde{X}]=\nu_{\rm
Re}$ and $\mathbb{E}_{p_{\rm Im}}[-\widetilde{X}]=\nu_{\rm Im}$, so that

$\mathbb{E}_{q}[X]=\mathbb{E}_{p_{\rm Re}}[\widetilde{X}]+i\mathbb{E}_{p_{\rm
Im}}[-\widetilde{X}].$ (55)

In turn, $\widetilde{X}$ is a random variable associated with the
corresponding two-outcome projective measurement of auxiliary qubit and taking
values in $\\{\widetilde{x}_{j,\bar{c}}\\}_{j,\bar{c}}$, where
$\widetilde{x}_{j,\bar{c}}=(-1)^{\bar{c}}x_{j}$, and $j$ and $\bar{c}$ are the
measurement outcomes. Therefore, for $K$ runs of Protocol II, we can consider
a projective measurement $\\{2R_{0},2R_{1}\\}$ of the auxiliary qubit
performed in $K/2$ runs and the measurement $\\{2R_{2},2R_{3}\\}$ performed in
other $K/2$ runs.

In order to estimate $\nu_{\rm Re}$, in the $k$-th run of the protocol, we
associate the measurement $\\{2R_{0},2R_{1}\\}$ with a random variable
$X_{k}\coloneqq\widetilde{X}/(K/2)$ taking values in
$\\{2(-1)^{\bar{c}}x_{j}/K\\}_{j,\bar{c}}$. Hence, the estimator $X_{\rm
Re}\coloneqq\sum_{k=1}^{K/2}X_{k}$ is an unbiased estimator for $\nu_{\rm
Re}$. According to the Hoeffding’s inequality Hoeffding (1963), the
probability of estimation with an additive error at least $\epsilon_{\nu}$ is
upper-bounded as

$\Pr(|X_{\rm Re}-\nu_{\rm Re}|\geq\epsilon_{\nu})\leq
2\exp\left(-\frac{2\epsilon_{\nu}^{2}}{\sum_{k}(\max(X_{k})-\min(X_{k}))^{2}}\right).$
(56)

As each of $K/2$ estimators is bounded as $-2x_{\rm max}/K\leq X_{k}\leq
2x_{\rm max}/K$, where $x_{\rm max}=\operatorname{max}_{j}|x_{j}|$, we obtain

$\sum_{k}(\max(X_{k})-\min(X_{k}))^{2}=\frac{8x_{\rm max}^{2}}{K}.$ (57)

Therefore, the probability (56) is upper-bounded by

$\Pr(|X_{\rm Re}-\nu_{\rm Re}|\geq\epsilon_{\nu})\leq
2\exp\left(-\frac{\epsilon_{\nu}^{2}K}{4x_{\max}^{2}}\right).$ (58)

Requiring that it does not exceed $\delta$ is equivalent to the upper bound

$2\exp\left(-\epsilon_{\nu}^{2}K/4x_{\max}^{2}\right)\leq\delta,$ (59)

bounding hence the estimation error as

$\epsilon_{\nu}\geq 2x_{\rm max}\sqrt{\frac{1}{K}\ln\frac{2}{\delta}}.$ (60)

This means that, for $K/2$ runs of Protocol II, we can make sure that the
estimation of $\nu_{\rm Re}$ is within an error $\epsilon_{\nu}$ with a
probability no less than $1-\delta$. This relation between the error and the
number of copies also holds for the estimation of the imaginary part $\nu_{\rm
Im}$, since, in accordance with (55), it is associated with the same estimator
$X_{\rm Re}$ up to an overall minus sign, which does not change the upper
bound (56). hence the proof. $\blacksquare$

######  Corollary 2.

Given a pair of states $\rho$, $\rho^{\prime}$, for $K$ runs of Protocol II,
the estimate of $\mu:=\operatorname{Tr}[\rho^{\prime}\rho]$ can be guaranteed
to have an error at most

$\epsilon_{\mu}=2\sqrt{\frac{1}{K}\ln\frac{2}{\delta}},$ (61)

with probability at least $1-\delta$.

###### Proof.

First, we note that $\mu=\nu_{\rm Re}$ if $A=I$. Therefore, from the
$(K/2)$-round measurement outcomes for estimating $\nu_{\rm Re}$, we can
simultaneously construct an unbiased estimator
$Y\coloneqq\sum_{k=1}^{K/2}Y_{k}$ for $\mu$, where $Y_{k}\coloneqq
2(-1)^{c}/K$. Since $\mu$ is a nonnegative number, it is not necessary to
estimate its imaginary part. Applying the Hoeffding’s inequality (56) and
taking into account that $x_{\rm max}=1$ for $\mu$, we find that the
estimation error is bounded by

$\epsilon_{\mu}\leq 2\sqrt{\frac{1}{K}\ln\frac{2}{\delta}}.$ (62)

Therefore, for $K/2$ runs of the Protocol II, we can make sure that the
estimation of $\mu$ is within a error $\epsilon_{\mu}$ with a probability
$1-\delta$. $\blacksquare$

######  Lemma 3.

Given an observable $A=\sum_{j}x_{j}P_{j}$, where $(P_{j})_{j}$ is a set of
POVM elements, and a pair of states $\rho$, $\rho^{\prime}$, for $K$ runs of
Protocol II, the estimate of its weak value $W(A\,|\,\rho,\rho^{\prime})$ can
be guaranteed to have an error at most:

$\epsilon=\frac{\sqrt{2}(x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|)}{\mu/\epsilon_{\mu}-1},$ (63)

where $x_{\rm max}:=\operatorname{max}_{j}|x_{j}|$,
$\mu=\operatorname{Tr}[\rho^{\prime}\rho]$, and
$\epsilon_{\mu}=2\sqrt{\frac{1}{K}\ln\frac{2}{\delta}}$, with probability at
least $1-3\delta$.

###### Proof.

First, we we recall that $W(A\,|\,\rho,\rho^{\prime})=W_{\rm Re}+iW_{\rm Im}$,
where $W_{\rm Re}=\nu_{\rm Re}/\mu$ and $W_{\rm Im}=\nu_{\rm Im}/\mu$.
Therefore, the corresponding estimation error is given by

$\epsilon=\sqrt{\epsilon_{W_{\rm Re}}^{2}+\epsilon_{W_{\rm Im}}^{2}},$ (64)

where $\epsilon_{W_{\rm Re}}$ and $\epsilon_{W_{\rm Re}}$ are errors in
estimation of the real and imaginary parts of the weak value, respectively.
The error in estimation of the real part $W_{\rm Re}$ can be upper-bounded as

$\displaystyle\epsilon_{W_{\rm Re}}$ $\displaystyle=\left|\frac{X_{\rm
Re}}{Y}-\frac{\nu_{\rm Re}}{\mu}\right|$ $\displaystyle=\left|\frac{X_{\rm
Re}\mu-\nu_{\rm Re}\mu+\mu\nu_{\rm Re}-Y\nu_{\rm Re}}{Y\mu}\right|$
$\displaystyle\leq\frac{|X_{\rm Re}-\nu_{\rm Re}|\mu+|\mu-Y||\nu_{\rm
Re}|}{|Y|\mu}$ $\displaystyle\leq\frac{\epsilon_{\nu}+\epsilon_{\mu}|W_{\rm
Re}|}{|Y|},$ (65)

where the first inequality follows from the triangle inequality
$|a+b|\leq|a|+|b|$, while the second inequality, in accordance with Lemma 2
and Corollary 2, follows from the upper bounds (53) and (61), respectively. As
the latter can be given as $\mu-\epsilon_{\mu}\leq Y\leq\mu+\epsilon_{\mu}$,
we can assume that the corresponding estimation error is small enough to
fulfill the condition $\epsilon_{\mu}<\mu$, so that
$|Y|\geq\mu-\epsilon_{\mu}\geq 0$ and, thus,

$\displaystyle\epsilon_{{\rm
Re}(W)}\leq\frac{\epsilon_{\nu}+\epsilon_{\mu}|W_{\rm
Re}|}{\mu-\epsilon_{\mu}}.$ (66)

Comparing (53) and (61), it is easy to note that $\epsilon_{\nu}=x_{\rm
max}\epsilon_{\mu}$. Therefore,

$\displaystyle\epsilon_{W_{\rm Re}}\leq\frac{x_{\rm max}+|W_{\rm
Re}|}{\mu/\epsilon_{\mu}-1}.$ (67)

Similarly, the error of estimating ${\rm Im}(W(A\,|\,\rho,\rho^{\prime}))$ is
bounded from above as

$\displaystyle\epsilon_{{\rm Im}(W)}\leq\frac{x_{\rm max}+|W_{\rm
Im}|}{\mu/\epsilon_{\mu}-1}.$ (68)

Therefore, an upper bound on the total error in estimation of the weak value
$W(A|\rho,\rho^{\prime})$

$\displaystyle\epsilon$ $\displaystyle\leq\frac{\sqrt{(x_{\rm max}+|W_{\rm
Re}|)^{2}+(x_{\rm max}+|W_{\rm Im}|)^{2}}}{\mu/\epsilon_{\mu}-1}$
$\displaystyle=\frac{\sqrt{2x_{\rm max}^{2}+2x_{\rm max}\left(|W_{\rm
Re}|+|W_{\rm
Im}|\right)+|W(A|\rho,\rho^{\prime})|^{2}}}{\mu/\epsilon_{\mu}-1}$
$\displaystyle\leq\frac{\sqrt{2x_{\rm max}^{2}+4x_{\rm
max}|W(A|\rho,\rho^{\prime})|+|W(A|\rho,\rho^{\prime})|^{2}}}{\mu/\epsilon_{\mu}-1}$
$\displaystyle\leq\frac{\sqrt{2}(x_{\rm
max}+|W(A|\rho,\rho^{\prime})|)}{\mu/\epsilon_{\mu}-1},$ (69)

is valid if the bounds $|X_{\rm Re}-\nu_{\rm Re}|\leq\epsilon_{\nu}$, $|X_{\rm
Im}-\nu_{\rm Im}|\leq\epsilon_{\nu}$, and $|Y-\mu|\leq\epsilon_{\mu}$ are
satisfied. Since each bound is violated with a probability not higher than
$\delta$, the violation probability of at least one of them is upper-bounded
by $3\delta$. Hence, the estimation error of $W(A\,|\,\rho,\rho^{\prime})$
does not exceed the upper bound (69) with probability at least $1-3\delta$.
$\blacksquare$

Proof of Theorem 1. Applying Lemma 3 and taking into account (61), we can
reverse (63) in order to obtain an expression for the number $K$ of runs of
the Protocol II in terms of the estimation error $\epsilon$,

$\displaystyle
K=\frac{4\ln\frac{2}{\delta}}{\mu^{2}}\left(\frac{\sqrt{2}\left(x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|\right)}{\epsilon}+1\right)^{2}.$ (70)

Denoting $\Lambda\coloneqq\sqrt{2}\left(x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|\right)/\epsilon$ for the sake of
simplicity, we obtain

$\displaystyle K$ $\displaystyle=$
$\displaystyle\frac{4\ln\frac{2}{\delta}\Lambda^{2}}{\mu^{2}}\left(1+\frac{1}{\Lambda}\right)^{2}$
(71) $\displaystyle=$
$\displaystyle\frac{4\ln\frac{2}{\delta}\Lambda^{2}}{\mu^{2}}\left(1+\frac{2}{\Lambda}+\frac{1}{\Lambda^{2}}\right)$
$\displaystyle=$
$\displaystyle\frac{4\ln\frac{2}{\delta}\Lambda^{2}}{\mu^{2}}+O\left(\frac{\ln\frac{1}{\delta}\Lambda}{\mu^{2}}\right),$

where the third equality follows from the assumption that $\Lambda\gg 1$.
Finally, relabelling $3\delta$ to $\delta$ for the sake of simplicity, we
conclude that

$\displaystyle K$ $\displaystyle=$
$\displaystyle\frac{8\ln\Bigl{(}\frac{6}{\delta}\Bigr{)}}{\epsilon^{2}}\Biggl{(}\frac{x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|}{\operatorname{Tr}[\rho\rho^{\prime}]}\Biggr{)}^{2}$
(72) $\displaystyle+\;O\left(\frac{\ln\frac{1}{\delta}}{\epsilon}\frac{x_{\rm
max}+|W(A\,|\,\rho,\rho^{\prime})|}{(\operatorname{Tr}[\rho\rho^{\prime}])^{2}}\right)$

runs of Protocol II are enough to ensure that the estimation of weak value
$W(A\,|\,\rho,\rho^{\prime})$ has an additive error at most $\epsilon$ with a
probability at least $1-\delta$. Hence the proof.

##  Appendix D Proof of Theorem 3

Here we present the proof of Theorem 3. The proof consists of a few lemmas,
which are of independent interest as they clarify the mathematical details of
the matrix representation of two-time states.

######  Lemma 4.

Let $\lambda:\mathsf{Lin}(\mathcal{H})\to\mathbb{C}$ be a two-time vector,
written as $\lambda(A)=\operatorname{Tr}[L\,A]$ for a suitable matrix
$L\in\mathsf{Lin}(\mathcal{H})$. Then, one has

$\displaystyle\lambda^{\dagger}(A)=\operatorname{Tr}[L^{\dagger}\,A]\qquad\forall
A\in\mathsf{Lin}(\mathcal{H})$ (73)

and

$\displaystyle E_{\lambda}(C)=\operatorname{Tr}\left[\left(L\otimes
L^{\dagger}\right)\,C\right]\qquad\forall
C\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})\,.$ (74)

Proof. By definition,
$\lambda^{\dagger}(A)=\overline{\lambda(A^{\dagger})}=\overline{\operatorname{Tr}[LA^{\dagger}]}=\operatorname{Tr}[L^{\dagger}A]$,
whence Eq. (73). Now, one has

$\displaystyle E_{\lambda}(A\otimes B)$
$\displaystyle=\operatorname{Tr}\big{[}L\,A\big{]}\,\operatorname{Tr}\big{[}L^{\dagger}\,B\big{]}$
$\displaystyle=\operatorname{Tr}\Big{[}\big{(}L\otimes
L^{\dagger}\big{)}\,\big{(}A\otimes B\big{)}\Big{]}\qquad\forall
A,B\in\mathsf{Lin}(\mathcal{H})\,.$ (75)

Since the product matrices $A\otimes B$ are a spanning set for
$\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$, Eq. (74) follows.
$\blacksquare$

######  Lemma 5.

The set of unnormalized two-time density matrices consists of functionals
$\omega$ of the form $\omega(C)=\operatorname{Tr}[O\,C]$, where

$\displaystyle O=\sum_{n}\,L_{n}\otimes L_{n}^{\dagger}\,,$ (76)

and $(L_{n})_{n}\subset\mathsf{Lin}(\mathcal{H})$ are arbitrary matrices.

Proof. Immediate from Eqs. (33) and (74). $\blacksquare$

######  Lemma 6.

For every bipartite vector $|\Psi\rangle\in\mathcal{H}\otimes\mathcal{H}$
there exists a matrix $M\in\mathsf{Lin}(\mathcal{H})$ such that

$\displaystyle{\tt
SWAP}\,\Big{(}|\Psi\rangle\langle\Psi|\Big{)}^{T_{2}}=M^{\dagger}\otimes M.$
(77)

Proof. Let us expand $|\Psi\rangle$ as
$|\Psi\rangle=\sum_{i}\,|\alpha_{i}\rangle\otimes|\beta_{i}\rangle$ for
suitable vectors $(|\alpha_{i}\rangle)_{i}\subset\mathcal{H}$ and
$(|\beta_{i}\rangle)_{i}\subset\mathcal{H}$. Then, we have

$\displaystyle{\tt SWAP}\,\left(|\Psi\rangle\langle\Psi|\right)^{T_{2}}$
$\displaystyle=\sum_{i,j}{\tt
SWAP}\,\left(|\alpha_{i}\rangle\langle\alpha_{j}|\otimes|\beta_{i}\rangle\langle\beta_{j}|\right)^{T_{2}}$
$\displaystyle=\sum_{i,j}{\tt
SWAP}\,|\alpha_{i}\rangle\langle\alpha_{j}|\otimes|\overline{\beta}_{j}\rangle\langle\overline{\beta}_{i}|$
$\displaystyle=\sum_{i,j}|\overline{\beta}_{j}\rangle\langle\alpha_{j}|\otimes|\alpha_{i}\rangle\langle\overline{\beta}_{i}|$
$\displaystyle=M^{\dagger}\otimes M\,,$ (78)

with $M:=\sum_{i}|\alpha_{i}\rangle\langle\overline{\beta}_{i}|$.
$\blacksquare$

######  Lemma 7.

A matrix $O\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$ is of the form (76)
if and only if it is of the form (35).

Proof. Suppose that $O$ is of the form (76). Let
$|\Psi\rangle=\sum_{i}\,|\alpha_{i}\rangle\otimes|\beta_{i}\rangle$ be an
arbitrary vector in $\mathcal{H}\otimes\mathcal{H}$. Then, one has

$\displaystyle\langle\Psi|\,\left(O\,{\tt SWAP}\right)^{T_{2}}\,|\Psi\rangle$
$\displaystyle=\operatorname{Tr}\left[\left(O\,{\tt
SWAP}\right)^{T_{2}}\,|\Psi\rangle\langle\Psi|\,\right]$
$\displaystyle=\operatorname{Tr}\left[O\,{\tt
SWAP}\,\left(|\Psi\rangle\langle\Psi|\right)^{T_{2}}\,\right]$
$\displaystyle=\operatorname{Tr}\left[O\,(M^{\dagger}\otimes M)\right]$
$\displaystyle=\sum_{n}\operatorname{Tr}\left[(L_{n}\otimes
L_{n}^{\dagger})\,(M^{\dagger}\otimes M)\right]$
$\displaystyle=\sum_{n}\left|\operatorname{Tr}[L_{n}^{\dagger}M]\right|^{2}\geq
0\,,$ (79)

the third equality following from Lemma 6 and the fourth from Eq. (76). Since
the vector $|\Psi\rangle$ is arbitrary, the matrix $P:=\left(O\,{\tt
SWAP}\right)^{T_{2}}$ is positive semidefinite. Moreover, one has
$O=P^{T_{2}}\,{\tt SWAP}$. Hence, $O$ is of the form (35).

Conversely, suppose that $O$ has the form (35), namely $O=P^{T_{2}}\,{\tt
SWAP}$ for some positive matrix $P$. Decomposing $P$ as
$P=\sum_{n}|\Psi_{n}\rangle\langle\Psi_{n}|$, we obtain

$\displaystyle P^{T_{2}}\,{\tt SWAP}$ $\displaystyle={\tt SWAP}\,{\tt
SWAP}\,P^{T_{2}}\,{\tt SWAP}$ $\displaystyle=\sum_{n}{\tt SWAP}\left({\tt
SWAP}\,|\Psi_{n}\rangle\langle\Psi_{n}|^{T_{2}}\right)\,{\tt SWAP}$
$\displaystyle=\sum_{n}{\tt SWAP}\left(M_{n}^{\dagger}\otimes
M_{n}\right)\,{\tt SWAP}$ $\displaystyle=\sum_{n}M_{n}\otimes M_{n}^{\dagger}$
(80)

where the third equality follows from Lemma 6. Hence, $P^{T_{2}}\,{\tt SWAP}$
is of the form (76) with $L_{n}:=M_{n}$. $\blacksquare$

Proof of Theorem 3. By Lemma 5, the functional $\omega$ defined by
$\omega(C)=\operatorname{Tr}[O\,C]$ is a two-time state if and only if
$O=\sum_{n}\,L_{n}\otimes L_{n}^{\dagger}$. Then, Lemma 7 guarantees that
there exists a positive semidefinite matrix
$P\in\mathsf{Lin}(\mathcal{H}\otimes\mathcal{H})$ such that $O=P^{T_{2}}{\tt
SWAP}$ if and only if $O=\sum_{n}\,L_{n}\otimes L_{n}^{\dagger}$. Hence the
proof.