Photon counting probabilities of the output field for a single-photon input

The scattering of light on a quantum system is one of core issues in quantum optics. Many efforts were put recently to describe the scattering of light prepared in $N$-photon state in a one-dimensional (1D) waveguide. The reason for this interest is the rapid development of techniques of generating and manipulating single-photon and multi-photon states of light Banaszek05 ; Cooper13 ; Yukawa13 ; Rybarczyk13 ; Scarani2013 ; Rempe15 ; Lodahl15 ; Ogawa16 ; Lodahl17 ; Sun18 ; Leong2016 . There are several theoretical techniques allowing to tackle this problem. The scattering in waveguide can be studied by using the Heisenberg picture approach Domokos02 ; WMSS11 ; WMS12 ; Stolyarov2013 ; Scarani2016a , methods based on the Lippmann-Schwinger equation ShenFan2005 ; Fan2005b ; ShenFan2007a ; ShenFan2007b ; Fan2009 ; Fan2012 ; Gritsev2012 ; Shen2015 , input-output formalism Fan2010 ; Fan2013 , generalized master equation GEPZ98 ; Baragiola2012 ; Cirac2015 , pure-state analysis Konyk16 ; Nysteen2015 or a stochastic approach Gough12a ; Gough12b ; Gough13 ; Dong15 ; Song16 ; Pan16 ; Baragiola17 ; Dabrowska17 ; Dabrowska18 ; Dabrowska19 ; Zhang19 ; Dong19 ; Dabrowska2020a . In Scarani2016b a method based on the operational translation of the system nonlinear response is described. In Fischer2018a ; Fischer2018b the scattering into 1D is tackled by a discrete approximation of bath Hilbert space. A comprehensive review on this subject was given in Zubairy2016 ; Roy2017 .

In this paper, we analyze the problem of scattering of a single-photon wave packet on a quantum system, within the Wesskopf-Wigner approximation Scully1997 , by making use of quantum filtering theory BarBel91 ; Car93 ; BGM04 ; B05 ; GZ10 ; WM10 . We describe the properties of the output field (the field after an interaction with quantum system) in time domain. The input-output formalism and the quantum filtering theory were formulated in the framework of the quantum stochastic Itô calculus (QSC) HP84 ; Par92 . To find the exclusive counting probability densities BarBel91 ; Car93 ; Srinivas81 ; Srinivas1982 ; Mollow68 , which fully characterize the statistics of photons in the output field, we determine quantum trajectories Car93 — the conditional states of the system depending on the results of continuous in time measurement performed on the output field. These states are also called a posteriori states BarBel91 . Due to the temporal correlations of the input field, the evolution of open quantum system is non-Markovian in this case Ciccarello2018 ; Dabrowska2020 . It is therefore not governed by Gorini-Kossakowski-Linblad-Sudarshan master equation GKS76 ; Lin76 but by a set of coupled master equations. Instead of studying an averaged evolution of the quantum system, given by the set of master equations, we analyze the set of filtering equations — the set of stochastic equations for the conditional operators of the system.

We formulate the problem of scattering by making use of a discrete model of repeated interactions and measurements (collision model) B02 ; BH08 ; BHJ09 ; Kretschmer2016 ; Breuer2016 ; Ciccarello2017 ; Gross18 ; M95 ; G04 ; GS04 ; Attal06 ; P08 ; PP09 ; P10 ; Ciccarello2021 with a bidirectional field defined as two collections of “ancillas” which are two-level systems (qubits). One of these collections is taken in an entangled state — a discrete analogue of continuous-mode single photon state L00 ; RMS07 ; Milburn08 , and the second collection is prepared in the vacuum, which means that all its qubits are in the ground state. We assume that bath qubits do not interact with each other and there is no initial correlation between the quantum system and its environment. Successive interactions (“collisions”) of the environment with the quantum system are described by unitary operators involving the system and two bath qubits. A schematic sketch of the discrete dynamics of the composed system, consisting of the travelling environment and the quantum system, is shown in Fig. 1. To characterize the properties of the output field, we assume that the measurements are performed on all bath elements after their interactions with the system. The results of these measurements give rise to the conditional dynamics of the quantum system, which is described by a set of difference filtering equations in the discrete case, and by a set of differential filtering equations in the continuous-time limit. Instead of solving the set differential filtering equations numerically, as it is usually done, we derive analytical formulae for the corresponding quantum trajectories. This allows us to determine general expressions for exclusive counting probability densities. We describe in detail the results for the case when a single-photon field interacts with a two-level atom. We fully characterize the statistics of photons of the output field in this case, and we show how to use the photon counting probability densities to determine the probability densities of successive counts and the mean number of photons registered up to a given time. In this paper, we present a direct generalization of the results for a single-photon unidirectional field described in Dabrowska17 ; Dabrowska19 ; Dabrowska2020a .

We would like to emphasize that the collision model provides the results which are consistent with the studies based on quantum stochastic calculus Gough12a ; Gough12b ; Baragiola2012 ; Baragiola17 , but the discrete approach allows one to simplify calculations and to provide an intuitive and rigorous interpretation to quantum trajectories. It allows the reader unfamiliar with QSC to understand the stochastic approach. The decomposition of the reduced system dynamics into quantum trajectories helps also one to better understand the nontrivial evolution of a quantum system interacting with a temporally correlated field.

The derivation of the Markovian case of conditional and unconditional evolution of an open quantum system in the discrete model of repeated interactions and measurements, with the bath defined as a sequence of qubits, one can find, for instance, in B02 ; GS04 ; Gross18 . For a discussions on physical assumptions leading to collision models in quantum optics see, for instance, Gross18 ; Ciccarello2017 ; Fischer2018a .

The paper is organized as follows. We define the collision model for a bidirectional field in Sec. II. In Sec. III, we provide a description of repeated measurements and we derive the formula for the conditional state of a system for discrete stochastic evolution. In Sec. IV, we determine the set of difference filtering equations and we obtain its continuous-time limit. Next, in Sec. V, we give analytical formulae for conditional vectors defining the quantum trajectories and we provide the expressions for photon counting probability densities in the continuous-time limit. In Sec. VI, we derive the statistics of counts for the output field for a two-level atom.

We consider a quantum system $\mathcal{S}$ interacting with an environment $\mathcal{E}$ being a bidirectional field modeled by two chains of qubits. One chain describes the field going to the right and the second one refers to the field going to the left. We assume that the qubits do not interact with each other but qubits of each chain interact subsequently with the system $\mathcal{S}$. At a given moment $\mathcal{S}$ interacts with only two qubits: one coming from the left and the other one coming from the right. Any interaction (“collision”) has the same duration $\tau$ and each of the bath qubits interacts with the system only once. The schematic of the collision dynamics is shown in Fig. 1. We describe the dynamics of the composed system $\mathcal{E}+\mathcal{S}$ up to time $T=N\tau$, where $N$ denotes the number of qubits in each chain. The Hilbert space of the environment is then defined as

$$\mathcal{H}_{\mathcal{E}}=\mathcal{H}_{\mathcal{E}_{1}}\otimes\mathcal{H}_{\mathcal{E}_{2}},$$

(1)

where

$$\mathcal{H}_{\mathcal{E}_{l}}=\bigotimes_{k=0}^{N-1}\mathcal{H}_{\mathcal{E}_{l,k}},\;l=1,2$$

(2)

and $\mathcal{H}_{\mathcal{E}_{l,k}}=\mathbb{C}^{2}$ is the Hilbert space of the qubit of the $l$-th part of the environment which interacts with $\mathcal{S}$ in the time interval $[k\tau,(k+1)\tau)$. Note that the Hilbert spaces $\mathcal{H}_{\mathcal{E}_{l}},\;l=1,2$ can be split as

$$\mathcal{H}_{\mathcal{E}_{l}}=\mathcal{H}_{\mathcal{E}_{l}}^{j-1]}\otimes\mathcal{H}_{\mathcal{E}_{l}}^{[j},$$

(3)

where

$$\mathcal{H}_{\mathcal{E}_{l}}^{j-1]}=\bigotimes_{k=0}^{j-1}\mathcal{H}_{\mathcal{E}_{l,k}},\;\;\;\mathcal{H}_{\mathcal{E}_{l}}^{[j}=\bigotimes_{k=j}^{N-1}\mathcal{H}_{\mathcal{E}_{l,k}}.$$

(4)

Thus, if $j\tau$ is a current moment, then

$$\mathcal{H}_{\mathcal{E}}^{j-1]}=\mathcal{H}_{\mathcal{E}_{1}}^{j-1]}\otimes\mathcal{H}_{\mathcal{E}_{2}}^{j-1]}$$

(5)

refers to the part of the environment which has already interacted with $\mathcal{S}$, constituting the output field, and

$$\mathcal{H}_{\mathcal{E}}^{[j}=\mathcal{H}_{\mathcal{E}_{1}}^{[j}\otimes\mathcal{H}_{\mathcal{E}_{2}}^{[j}$$

(6)

refers to the part which has not interacted with $\mathcal{S}$ yet—the input field. We shall call $\mathcal{H}_{\mathcal{E}}^{j-1]}$ and $\mathcal{H}_{\mathcal{E}}^{[j}$ respectively the past and future environment spaces.

The evolution of the composed system is described by the sequence of unitary operators, $U_{j}$ for $0\leq j\leq N-1$, defined by

$$U_{j\tau}={V}_{j-1}{V}_{j-2}\ldots{V}_{0},\;\;\;\;\;U_{0}=\mathbbm{1},$$

(7)

where $V_{k}$ for $0\leq k\leq N-1$ describes the interaction between $\mathcal{S}$ and $\mathcal{E}$ in the time-interval $[k\tau,(k+1)\tau)$. The operator ${V}_{k}$ acts non-trivially only in the space $\mathcal{H}_{\mathcal{E}_{1},k}\otimes\mathcal{H}_{\mathcal{E}_{2},k}\otimes\mathcal{H}_{\mathcal{S}}$, where $\mathcal{H}_{\mathcal{S}}$ is the Hilbert space of $\mathcal{S}$, and $V_{k}$ has the form

$${V}_{k}=\exp\left(-i\tau H_{k}\right)$$

(8)

with the Hamiltonian G04 ; Fischer2018a ; Gross18 ; Ciccarello2017

$\displaystyle\ H_{k}=H_{\mathcal{S}}+\sum_{l=1}^{2}\frac{i}{\sqrt{\tau}}\left(\sigma_{l,k}^{+}\otimes L_{l}-\sigma_{l,k}^{-}\otimes L^{\dagger}_{l}\right).$

(9)

We set $\hbar=1$ throughout the paper and to simplify the notation we omit a multiplication by identity operators. The model is formulated in the framework of some standard assumptions made in quantum optics: rotating wave-approximation, a flat coupling constant, and the extension of the lower limit of integration over frequency to minus infinity Scully1997 . The bandwidth of the spectrum is assumed to be much smaller that the central frequency of the pulse. The Hamiltonian $H_{k}$ is written in the interaction picture eliminating the free evolution of the field. Here $H_{\mathcal{S}}$ stands for the Hamiltonian of $\mathcal{S}$, $\sigma^{+}_{l,k}$ and $\sigma^{-}_{l,k}$ denote respectively the raising and lowering operators acting in $\mathcal{H}_{\mathcal{E}_{l,k}}$. From the mathematical point of view, $L_{l}$ for $l=1,2$ are arbitrary bounded operators on $\mathcal{H}_{S}$. They are called jump operators or the Itô coefficients. In Sec. V we consider $\mathcal{S}$ which is a two-level atom and we define $L_{l}$ as $\sqrt{\Gamma_{l}}\sigma_{-}$, where $\Gamma_{l}$ is a positive coupling constant and $\sigma_{-}$ is the atom lowering operator. If $\mathcal{S}$ is a two-sided cavity, then $L_{l}=\sqrt{\Gamma_{l}}a$ and $a$ is the annihilation operator of a cavity mode. A discussion about the structure of the collision model and physical assumptions leading to it one can find, for instance, in Gross18 ; Fischer2018a ; Ciccarello2017 . We use the representation of (8) in the basis $\{|00\rangle_{k},|01\rangle_{k},|10\rangle_{k},|11\rangle_{k}\}$, such that

$$\exp\left(-i\tau H_{k}\right)=\sum_{i_{1},i_{2},i_{3},i_{4}=0,1}|i_{1}i_{2}\rangle_{k}{}_{k}\langle i_{3}i_{4}|\otimes V_{i_{1}i_{2},i_{3}i_{4}}$$

(10)

where $|i_{l_{1}}i_{l_{2}}\rangle_{k}=|i_{l_{1}}\rangle_{k}\otimes|i_{l_{2}}\rangle_{k}$, $V_{i_{1}i_{2},i_{3}i_{4}}$ are operators on $\mathcal{S}$, and their explicit forms are given in A. Clearly, in order to approximate the continuous evolution, we take $\tau$ that satisfies the conditions

$$\tau\ll\Gamma_{l}^{-1},\;\;l=1,2.$$

(11)

Note that we are interested in terms up to order $\tau$ for determining first-order differential equations, it means that we have to expand $V_{k}$ to the second order. To obtain the continuous dynamics, we will take finally the limit of $\tau\to 0$ and $N\to\infty$ such that $T=N\tau$ is fixed.

The initial state of the composed system is assumed to be the product state vector of the form

$$|\Psi_{0}\rangle=|1_{\xi}\rangle\otimes|vac\rangle\otimes|\psi_{0}\rangle,$$

(12)

where $|\psi_{0}\rangle$ is the initial state of $\mathcal{S}$ and

$\displaystyle|1_{\xi}\rangle=\sum_{k=0}^{N-1}\sqrt{\tau}\xi_{k}\sigma^{+}_{1,k}|vac\rangle$

(13)

with the vacuum vector $|vac\rangle=|0\rangle_{0}\otimes|0\rangle_{1}\otimes\ldots\otimes|0\rangle_{N-1}$, where $|0\rangle_{k}$ is the ground state in $\mathbb{C}^{2}$, and $\displaystyle{\sum_{k=0}^{N-1}}\tau|\xi_{k}|^{2}=1$. Note that vector $|1_{\xi}\rangle$ poses the additive decomposition property

$$|1_{\xi}\rangle=\sum_{k=0}^{j}\sqrt{\tau}\xi_{k}\sigma^{+}_{1,k}|vac\rangle+\sum_{k=j+1}^{N-1}\sqrt{\tau}\xi_{k}\sigma^{+}_{1,k}|vac\rangle$$

(14)

and it can be written in the form

$$|1_{\xi}\rangle=\sum_{k=0}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle,$$

(15)

where $|1_{k}\rangle=|0\rangle_{0}\otimes|0\rangle_{1}\otimes\ldots|0\rangle_{k-1}\otimes|1\rangle_{k}\otimes|0\rangle_{k+1}\otimes\ldots|0\rangle_{N-1}$. Thus $|1_{\xi}\rangle$ is a superposition of vectors with one qubit prepared in the excited state and all the others in the ground state. Clearly, $|1_{\xi}\rangle$ an entangled state and $|\xi_{k}|^{2}\tau$ is the probability that the qubit of the $k$ number in the first chain is in the excited state and all the others qubits in this chain are in the ground state. One can easily check the identities

$$\sigma^{-}_{1,k}|1_{\xi}\rangle=\sqrt{\tau}\xi_{k}|vac\rangle,$$

(16)

$$\sigma^{+}_{1,k}\sigma^{-}_{1,k}|1_{\xi}\rangle=\sqrt{\tau}\xi_{k}|1_{k}\rangle.$$

(17)

The state $|1_{\xi}\rangle$ is a discrete version of a continuous-mode single-photon state L00 ; RMS07 ; Milburn08 .

After $j$ interactions the state of the composed system is given as

$$U_{j-1}|\Psi_{0}\rangle\langle\Psi_{0}|U_{j-1}^{\dagger}.$$

(18)

And taking the partial trace over the environment, we obtain the reduced state of $\mathcal{S}$ at time $j\tau$:

$$\varrho_{j}=\mathrm{Tr}_{\mathcal{E}}\left[U_{j-1}|\Psi_{0}\rangle\langle\Psi_{0}|U_{j-1}^{\dagger}\right].$$

(19)

We describe now the setup of repeated measurements performed on the environment chains. We assume that after each step of interaction, the measurements are performed on two qubits which have just interacted with $\mathcal{S}$. Clearly, the first chain is going to the right and its output is measured on the right side of $\mathcal{S}$, the second chain is going to the left and its output is measured on the left side of $\mathcal{S}$. We consider measurements of the input observables:

$$\sigma_{l,k}^{+}\sigma_{l,k}^{-},\;\;l=1,2,\;k=0,1,\ldots N-1.$$

(20)

It means that we study the measurements in the orthogonal basis: $\{|00\rangle_{k},|01\rangle_{k},|10\rangle_{k},$ $|11\rangle_{k}\}$.

One can check that after the first interaction, the state of the composed system has the following form

	$\displaystyle{V}_{0}|\Psi_{0}\rangle=|vac\rangle\otimes|vac\rangle\otimes(-\tau\xi_{0}L_{1}^{\dagger})|\psi_{0}\rangle$
	$\displaystyle+|0\rangle_{0}\otimes\sum_{k=1}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle_{[1}\otimes|vac\rangle\otimes\left(1-i\tau H_{\mathcal{S}}-\frac{\tau}{2}\left(L_{1}^{\dagger}L_{1}+L_{2}^{\dagger}L_{2}\right)\right)|\psi_{0}\rangle$
	$\displaystyle+|1\rangle_{0}\otimes\sum_{k=1}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle_{[1}\otimes|vac\rangle\otimes\sqrt{\tau}L_{1}|\psi_{0}\rangle$
	$\displaystyle+|1\rangle_{0}\otimes|vac\rangle_{[1}\otimes|vac\rangle\otimes\sqrt{\tau}\xi_{0}|\psi_{0}\rangle$
	$\displaystyle+|0\rangle_{0}\otimes\sum_{k=1}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle_{[1}\otimes|1\rangle_{0}\otimes|vac\rangle_{[1,}\otimes\sqrt{\tau}L_{2}|\psi_{0}\rangle+$
	$\displaystyle+|1\rangle_{0}\otimes|vac\rangle_{[1}\otimes|1\rangle_{0}\otimes|vac\rangle_{[1}\otimes\tau\xi_{0}L_{2}|\psi_{0}\rangle+$
	$\displaystyle+|1\rangle_{0}\otimes\sum_{k=1}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle_{[1}\otimes|1\rangle_{0}\otimes|vac\rangle_{[1}\otimes\frac{\tau}{2}\left(L_{1}L_{2}+L_{2}L_{1}\right)|\psi_{0}\rangle+O(\tau^{3/2}),$		(21)

where

$$|vac\rangle_{[1}=|0\rangle_{1}\otimes|0\rangle_{2}\otimes\ldots\otimes|0\rangle_{N-1},$$

(22)

$$|1_{k}\rangle_{[1}=|0\rangle_{1}\otimes|0\rangle_{2}\otimes\ldots|0\rangle_{k-1}\otimes|1\rangle_{k}\otimes|0\rangle_{k+1}\otimes\ldots|0\rangle_{N-1},$$

(23)

and $O(.)$ is the Landau symbol. We assume that now the measurement is performed in the basis $\{|00\rangle_{0},|01\rangle_{0},|10\rangle_{0},|11\rangle_{0}\}$. Thus, if we have observed, for instance, the result $(0,0)$ (zero counts for the both detectors), then the conditional state of the composed system is given by

$$\frac{|00\rangle_{0}\langle 00|{V}_{0}|\Psi_{0}\rangle\langle\Psi_{0}|{V}_{0}^{\dagger}|00\rangle_{0}\langle 00|}{\mathrm{Tr}_{\mathcal{E+S}}\left[{V}_{0}|\Psi_{0}\rangle\langle\Psi_{0}|{V}_{0}^{\dagger}|00\rangle_{0}\langle 00|\right]}.$$

(24)

Note that it is a pure state. By eliminating the degrees of freedom of the past environment, which will not interact with $\mathcal{S}$ in the future, we obtain the recipe for the a posteriori state of $\mathcal{S}$ and the input part of the environment. Thus the conditional state of $\mathcal{S}$ and the input part of the environment at time $\tau$ can be written in the form

$$|\tilde{\Psi}_{1|(\eta_{1,1},\eta_{2,1})}\rangle=\frac{|\Psi_{1|(\eta_{1,1},\eta_{2,1})}\rangle}{\sqrt{\langle\Psi_{1|(\eta_{1,1},\eta_{2,1})}|\Psi_{1|(\eta_{1,1},\eta_{2,1})}\rangle}},$$

(25)

where $(\eta_{1,1},\eta_{2,1})$ stands for the results of the first measurement performed on the qubits of two chains at time $\tau$. Clearly, $|\Psi_{1|(\eta_{1,1},\eta_{2,1})}\rangle$ is a random unnormalized vector from the Hilbert space $\mathcal{H}_{\mathcal{E}}^{[1}\otimes\mathcal{H_{S}}$ and its form depends on the observed outcomes. For instance, if the result is $(0,0)$, then we get

	$\displaystyle|\Psi_{1|(0,0)}\rangle$	$\displaystyle=$	$\displaystyle\sum_{k=1}^{N-1}\sqrt{\tau}\xi_{k}|1_{k}\rangle_{[1}\otimes|vac\rangle_{[1}\otimes\left(1-i\tau H_{\mathcal{S}}-\frac{\tau}{2}\left(L_{1}^{\dagger}L_{1}+L_{2}^{\dagger}L_{2}\right)\right)|\psi_{0}\rangle$		(26)
			$\displaystyle-|vac\rangle_{[1}\otimes|vac\rangle_{[1}\otimes\tau\xi_{0}L_{1}^{\dagger}|\psi_{0}\rangle.$

One can easily check that the results of the first measurement appear with the following probabilities

$$p_{1}\left((0,0)\right)=1-\tau\left(\langle\psi_{0}|(L_{1}^{\dagger}L_{1}+L_{2}^{\dagger}L_{2})|\psi_{0}\rangle+|\xi_{0}|^{2}\right)+O(\tau^{2}),$$

(27)

$$p_{1}\left((1,0)\right)=\tau\left(\langle\psi_{0}|L_{1}^{\dagger}L_{1}|\psi_{0}\rangle+|\xi_{0}|^{2}\right)+O(\tau^{2}),$$

(28)

$$p_{1}\left((0,1)\right)=\tau\langle\psi_{0}|L_{2}^{\dagger}L_{2}|\psi_{0}\rangle+O(\tau^{2}),$$

(29)

$$p_{1}\left((1,1)\right)=O(\tau^{2}).$$

(30)

Since the probability of the result $(1,1)$ is $O(\tau^{2})$, we ignore such detection. In the next step the system $\mathcal{S}$ interacts with the second pair of the qubits and after this interaction we perform the next measurement. Note that the conditional state of $\mathcal{S}$ and the input field in time $2\tau$ will depend on the results of the two past measurements. We formulate our result for time $j\tau$ $(1\leq j\leq N-1)$ in the form of a theorem.

The proof of Theorem 1 one can find in B.

Note that the vectors $|\alpha_{j|\boldsymbol{\eta}_{j}}\rangle$, $|\beta_{j|\boldsymbol{\eta}_{j}}\rangle$ depend on all results of measurements up to the time $j\tau$ and this is indicated by their subscript $\boldsymbol{\eta}_{j}$. It is seen from the formula (32) that the system $\mathcal{S}$ becomes entangled with the left input part of the environment. This property distinguishes the studied case from the Markovian one within after each measurement the state of the system and the part of environment which has not interacted with the system is a product state (see, for instance, G04 ). Clearly, the non-Markovian character of the evolution comes here from the initial entanglement of the bath qubits. Let us notice that the conditional state (31) has the following physical interpretation. The first term (with the conditional vector $|\alpha_{j|\boldsymbol{\eta}_{j}}\rangle$) describes the scenario that $\mathcal{S}$ has not met the qubit prepared in the upper state yet and it appears in the future. The second term (with the vector $|\beta_{j|\boldsymbol{\eta}_{j}}\rangle$) is associated with the scenario that $\mathcal{S}$ has already interacted with the qubit prepared in the upper state and it meets in the future only qubits in the ground states. We treat $\tau$ as a small time step in the sense that all elements of the recurrence equations (35)-(42) of order more than $\tau$ can be ignored. To understand why these terms disappear one needs to consider the continuous-time limits of the solutions to the derived equations.

The probability of a given sequence, $\boldsymbol{\eta}_{j}$, is

$$p(\boldsymbol{\eta}_{j})=\langle{\Psi}_{j|\boldsymbol{\eta}_{j}}|{\Psi}_{j|\boldsymbol{\eta}_{j}}\rangle$$

(43)

and one can easily check that it can be expressed by the conditional vectors as

$$p(\boldsymbol{\eta}_{j})=\langle\alpha_{j|\boldsymbol{\eta}_{j}}|\alpha_{j|\boldsymbol{\eta}_{j}}\rangle\;\displaystyle{\sum_{k=j}^{N-1}}\;\tau|\xi_{k}|^{2}+\langle\beta_{j|\boldsymbol{\eta}_{j}}|\beta_{j|\boldsymbol{\eta}_{j}}\rangle.$$

(44)

Let us note that the conditional probability of the outcome $\eta_{j+1}=(0,0)$ at $(j+1)\tau$ given the sequence $\boldsymbol{\eta}_{j}$ is defined as

$$p_{j+1}\left((0,0)|\boldsymbol{\eta}_{j}\right)=\frac{\langle{\Psi}_{j|\left((0,0),\boldsymbol{\eta}_{j}\right)}|{\Psi}_{j|\left((0,0),\boldsymbol{\eta}_{j}\right)}\rangle}{\langle{\Psi}_{j|\boldsymbol{\eta}_{j}}|{\Psi}_{j|\boldsymbol{\eta}_{j}}\rangle}.$$

(45)

In a similar way we define the conditional probabilities of the other outcomes. By using of the difference equations (35)-(42), one can check the following important properties

$$p_{j+1}\left((0,0)|\boldsymbol{\eta}_{j}\right)=1+O(\tau),$$

(46)

$$p_{j+1}\left((1,0)|\boldsymbol{\eta}_{j}\right)=O(\tau),$$

(47)

$$p_{j+1}\left((0,1)|\boldsymbol{\eta}_{j}\right)=O(\tau),$$

(48)

$$p_{j+1}\left((1,1)|\boldsymbol{\eta}_{j}\right)=O(\tau^{2}).$$

(49)

They mean that most of the time we observe two vacuums and from time to time we observe a count on the left or a count on the right. However, because the simultaneous counts at the left and at the right detectors appear with the probability of order $O(\tau^{2})$, we do not observe such result.

Note that the $\eta_{j}$, describing the $j$-th results, is a two-dimensional random variable which is statistically dependent on the sequence $\eta_{1},\ldots,\eta_{j-1}$. Let us introduce the two-dimensional discrete stochastic process

$$n_{j}=\left(n_{1,j},n_{2,j}\right)$$

(50)

where

$$n_{1,j}=\sum_{k=1}^{j}\eta_{1,j},\;\;n_{2,j}=\sum_{k=1}^{j}\eta_{2,j}$$

(51)

are the stochastic processes referring to counts registered respectively by the right and the left detector.

Taking the partial trace of $|\tilde{\Psi}_{j|\boldsymbol{\eta}_{j}}\rangle\langle\tilde{\Psi}_{j|\boldsymbol{\eta}_{j}}|$ over the environment, one obtains the a posteriori state of $\mathcal{S}$ at the time $j\tau$,

$$\tilde{\rho}_{j|\boldsymbol{\eta}_{j}}=\frac{\rho_{j|\boldsymbol{\eta}_{j}}}{\mathrm{Tr}\rho_{j|\boldsymbol{\eta}_{j}}},$$

(52)

where

$$\rho_{j|\boldsymbol{\eta}_{j}}=|\alpha_{j|\boldsymbol{\eta}_{j}}\rangle\langle\alpha_{j|\boldsymbol{\eta}_{j}}|\;\displaystyle{\sum_{k=j}^{N-1}}\;\tau|\xi_{k}|^{2}+|\beta_{j|\boldsymbol{\eta}_{j}}\rangle\langle\beta_{j|\boldsymbol{\eta}_{j}}|.$$

(53)

Clearly, $\tilde{\rho}_{j|\boldsymbol{\eta}_{j}}$ is the state of $\mathcal{S}$ depending on the results of all measurements performed on the output field up to $j\tau$. Thus the proceeding of the repeated interactions and measurements allows us to derive discrete quantum trajectories for $\mathcal{S}$ associated with the two-dimensional counting stochastic process (50). Note that a different realization of the stochastic process (50) means a different quantum trajectory given by (52).

Considering the counting stochastic process (50) it is convenient to replace the notation $\boldsymbol{\eta}_{j}$ of all results from $0$ to $j\tau$ by writing only the location of counts in $(\eta_{j},\ldots,\eta_{1})$. So the sequence

$$\left(D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}\right)$$

(54)

means the following scenario: the first photon was counted at time $l_{1}\tau$ at the detector $D_{1}$, the second photon at time $l_{2}\tau$ at the detector $D_{2}$ and so on, where $l_{1}<l_{2}<\ldots<l_{m}$, and no other photons were detected from $0$ to $j\tau$. We shall use the notation $R$ and $L$ to indicate respectively the right and left detector, thus $D_{1},D_{2},\ldots,D_{m}=R,L$.

The discrete conditional vectors, being the solutions to Eqs. (35)-(40) for some chosen sequences of the outcomes are given in C.

By taking the average over all possible outcomes one gets from (52) the formula for the a priori state of the system $\mathcal{S}$:

$$\varrho_{j}=\rho_{j|\boldsymbol{0}_{j}}+\sum_{m=1}^{j}\sum_{l_{m}=m}^{j}\sum_{l_{m-1}=m-1}^{j-1}\ldots\sum_{l_{1}=1}^{l_{2}-1}\sum_{D_{m},\ldots,D_{2},D_{1}=R,L}\rho_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}},$$

(55)

where

$$\rho_{j|\boldsymbol{0}_{j}}=|\alpha_{j|\boldsymbol{0}_{j}}\rangle\langle\alpha_{j|\boldsymbol{0}_{j}}|\sum_{k=j}^{N-1}\tau|\xi_{k}|^{2}+|\beta_{j|\boldsymbol{0}_{j}}\rangle\langle\beta_{j|\boldsymbol{0}_{j}}|$$

(56)

and

	$\displaystyle\rho_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}}$	$\displaystyle=$	$\displaystyle|\alpha_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}}\rangle\langle\alpha_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}}|\sum_{k=j}^{N-1}\tau|\xi_{k}|^{2}$		(57)
			$\displaystyle+|\beta_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}}\rangle\langle\beta_{j|D_{m},l_{m};D_{m-1},l_{m-1};\ldots;D_{1},l_{1}}|.$

Clearly, for the moment $j\tau$ there are $4^{j}$ of different scenarios (trajectories), but some of them occur with zero probability. The a priori state specifies the unconditional evolution of $\mathcal{S}$. We constructed a decomposition (unraveling) of the reduced density operator of $\mathcal{S}$ defined by (19). Let us stress that for different choices of the measured observables one gets different quantum trajectories (conditional evolution of $\mathcal{S}$) but all of them lead to the same average dynamics.

In this section, we shall determine the set stochastic equations defining the a posteriori evolution of $\mathcal{S}$. To simplify the notation, we skip in this section the conditional subscript.

First, we derive the recurrence formulae for a a posteriori state (52) to the order of $\tau$. To determine the conditional formulae, we make use of the recurrence equations (35)-(40). If the results of measurements at time $(j+1)\tau$ is $\eta_{j+1}=(0,0)$, then we get

	$\displaystyle\rho_{j+1}$	$\displaystyle=$	$\displaystyle\rho_{j}-i[H_{\mathcal{S}},\rho_{j}]\tau-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},\rho_{j}\}\tau$		(58)
			$\displaystyle-|\beta_{j}\rangle\langle\alpha_{j}|L_{1}\xi_{j}^{\ast}\tau-L_{1}^{\dagger}|\alpha_{j}\rangle\langle\beta_{j}|\xi_{j}\tau$
			$\displaystyle-|\alpha_{j}\rangle\langle\alpha_{j}||\xi_{j}|^{2}\tau+O(\tau^{2}),$

	$\displaystyle|\alpha_{j+1}\rangle\langle\beta_{j+1}|$	$\displaystyle=$	$\displaystyle|\alpha_{j}\rangle\langle\beta_{j}|-i[H_{\mathcal{S}},|\alpha_{j}\rangle\langle\beta_{j}|]\tau$		(59)
			$\displaystyle-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},|\alpha_{j}\rangle\langle\beta_{j}|\}\tau$
			$\displaystyle-|\alpha_{j}\rangle\langle\alpha_{j}|L_{1}\xi_{j}^{\ast}\tau+O(\tau^{2}),$

	$\displaystyle|\alpha_{j+1}\rangle\langle\alpha_{j+1}|$	$\displaystyle=$	$\displaystyle|\alpha_{j}\rangle\langle\alpha_{j}|-i[H_{\mathcal{S}},|\alpha_{j}\rangle\langle\alpha_{j}|]\tau$		(60)
			$\displaystyle-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},|\alpha_{j}\rangle\langle\alpha_{j}|\}\tau+O(\tau^{2}),$

where $\{a,b\}=ab+ba$. Let us notice that the unnormalized a posteriori operator $\rho_{j+1}$ depends on the results of the $(j+1)$-st measurements and the conditional operators: $\rho_{j}$, $|\alpha_{j}\rangle\langle\beta_{j}|$, $|\beta_{j}\rangle\langle\alpha_{j}|$, and $|\alpha_{j}\rangle\langle\alpha_{j}|$. The conditional probability of the outcome $(0,0)$ at time $(j+1)\tau$ given the a posteriori state $\tilde{\rho}_{j}$ is defined as

$$p_{j+1}\left((0,0)|\tilde{\rho}_{j}\right)=\frac{\mathrm{Tr}{\rho}_{j+1}}{\mathrm{Tr}{\rho}_{j}}.$$

(61)

One can check that

$$p_{j+1}\left((0,0)|\tilde{\rho}_{j}\right)=1-k_{j}\tau+O(\tau^{2}),$$

(62)

where

	$\displaystyle k_{j}$	$\displaystyle=$	$\displaystyle\mathrm{Tr}\left[\left(L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2}\right)\tilde{\rho}_{j}+\xi_{j}^{\ast}L_{1}|\tilde{\beta}_{j}\rangle\langle\tilde{\alpha}_{j}|\right.$		(63)
			$\displaystyle\left.+\xi_{j}\,|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|L^{\dagger}_{1}+|\xi_{j}|^{2}|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|\right]$

and $|\tilde{\alpha}_{j}\rangle=\frac{|{\alpha}_{j}\rangle}{\sqrt{\mathrm{Tr}\rho_{j}}}$, $|\tilde{\beta}_{j}\rangle=\frac{|{\beta}_{j}\rangle}{\sqrt{\mathrm{Tr}\rho_{j}}}$. Now, by making use of the property

$$\frac{1}{\mathrm{Tr}\rho_{j+1}}=\frac{1}{\mathrm{Tr}\rho_{j}}\left(1+k_{j}\tau\right)+O(\tau^{2}),$$

(64)

we obtain the difference equation for the normalized a posteriori density operator,

	$\displaystyle\tilde{\rho}_{j+1}$	$\displaystyle=$	$\displaystyle\tilde{\rho}_{j}+\tilde{\rho}_{j}k_{j}\tau-i[H_{\mathcal{S}},\tilde{\rho}_{j}]\tau-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},\tilde{\rho}_{j}\}\tau$		(65)
			$\displaystyle-|\tilde{\beta}_{j}\rangle\langle\tilde{\alpha}_{j}|L_{1}\xi_{j}^{\ast}\tau-L^{\dagger}_{1}|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|\xi_{j}\tau-|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}||\xi_{j}|^{2}\tau.$

In a similar way, one can find the following coupled difference equations

	$\displaystyle|\tilde{\alpha}_{j+1}\rangle\langle\tilde{\beta}_{j+1}|$	$\displaystyle=$	$\displaystyle|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|+|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|k_{j}\tau-i[H_{\mathcal{S}},|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|]\tau$		(66)
			$\displaystyle-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},|\tilde{\alpha}_{j}\rangle\langle\tilde{\beta}_{j}|\}\tau-|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|L_{1}\xi_{j}^{\ast}\tau,$

	$\displaystyle|\tilde{\alpha}_{j+1}\rangle\langle\tilde{\alpha}_{j+1}|$	$\displaystyle=$	$\displaystyle|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|+|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|k_{j}\tau-i[H_{\mathcal{S}},|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|]\tau$		(67)
			$\displaystyle-\frac{1}{2}\{L^{\dagger}_{1}L_{1}+L^{\dagger}_{2}L_{2},|\tilde{\alpha}_{j}\rangle\langle\tilde{\alpha}_{j}|\}\tau.$