Tomography of a multimode quantum black box
Abstract
We report a technique for experimental characterization of an mode quantum optical process, generalizing the singlemode coherentstate quantumprocess tomography method [M. Lobino et al., Science 322, 563 (2008); A. Anis and A.I. Lvovsky, New J. Phys. 14, 105021 (2012)]. By measuring effect of the process on multimode coherent states via balanced homodyne tomography, we obtain the process tensor in the Fock basis. This rank tensor, which predicts the effect of the process on an arbitrary density matrix, is iteratively reconstructed directly from the experimental data via the maximumlikelihood method. We demonstrate the capabilities of our method using the example of a beam splitter, reconstructing its process tensor within the subspace spanned by the first three Fock states. In spite of using purely classical probe states, we recover quantum properties of this optical element, in particular the HongOuMandel effect.
 PACS numbers

03.65.Wj, 03.65.Ta, 42.50.Xa
pacs:
03.65.Wj, 03.65.Ta, 42.50.XaI Introduction
Precise understanding of the performance of individual quantum systems is a key requirement for the development of compound devices, e.g. quantum computers or secure communication networks. This requirement gives rise to the problem of experimentally characterizing quantum systems as ‘black boxes: learning to predict their effect on arbitrary quantum states by measuring their effect on a limited number of “probe” states. The art of solving this problem is referred to as quantum process tomography (QPT).
A straightforward approach to QPT consists of measuring the action of the black box on a set of states whose density operators form a spanning set in the space of all operators over a particular Hilbert space. Because any quantum process is a linear map with respect to density operators, this information is sufficient to fully characterize the process Poyatos1997 (). However, such a direct method typically requires a large set of difficulttoprepare probe states, and is consequently restricted to systems of very low dimension. Another possibility is the ancillaassisted method D'Ariano2001 () utilizing an input state that is a part of a fully entangled state in a larger Hilbert space. Although in this case only a single input is necessary thanks to the Jamiolkowski isomorphism Jamiolkowski1972 (), both preparation of this state and tomography of the output state is, again, complicated, which dramatically limits the practicality of the method.
In application to optics, the coherentstate quantum process tomography (csQPT) Lobino2008 (); Lobino2009 (); Keshari2011 () offers a practical solution. While being a member of the direct methods family described above, this technique uses only coherent states for probing the process , relying on the fact that these states span the space of operators over the optical Hilbert space (the optical equivalence theorem) Glauber1963b (); Sudarshan1963 (). The prediction for the output of the black box in response to to an arbitrary input state then involves integration of the measured output states , weighted by corresponding GlauberSudarshan function , over the phase space. A similar coherentstate based approach can also be used for the tomography of quantum measurements Lundeen2009 (); Zhang2012 (). While being a case of the direct method described above, csQPT is relatively easy to implement in an experiment, since coherent states are readily obtained from lasers, and their amplitudes and phases are easy to control.
On the other hand, is a generalized function, typically highly singular. Therefore the process reconstruction involving that function may either suffer from inaccuracies or involve an unreasonably large number of required probe states. Moreover, the procedures proposed in Refs. Lobino2008 (); Keshari2011 () evaluate each element of the process tensor individually, and can hence lead to unphysical (nontrace preserving or nonpositive) process tensors.
The above shortcomings are absent in a method known as MaxLik csQPT, which exploits the Jamiolkowski isomorphism to reduce the QPT problem to the wellstudied problem of the quantum state estimation, and applies the likelihood maximization technique to estimate the process tensor Hradil2004 (). In this way, one can perform the reconstruction without leaving the physically plausible space. MaxLik csQPT has been proposed in Ref. Anis2012 () and successfully realized for nondeterministic singemode processes Kumar2013 (); Cooper2014 ().
In this work, we expand csQPT beyond the “single input — single output” case, which covers only a few of practically relevant quantum optical black boxes. The need for our study is dictated by the growing fields of quantum optical communication and logic, which are impossible without multimode processing. Examples include multimode quantum memories Gisin2009 (); Moiseev2010 () and logic gates for processing photonic qubits Monroe (); Pittman2003 (), to name a few. Although our experiment is in the optical domain, the theory and methodology of csQPT can be employed on a much broader scale. It applies to any physical system whose Hamiltonian is equivalent to that of the harmonic oscillator — such as superconducting cavities, atomic spin ensembles and nanomechanical systems. In all of these, coherent states are the simplest to prepare and are hence most suitable as probe states in QPT.
Ii Multimode MaxLik csQPT
Our method generalizes the singlemode MaxLik csQPT approach Anis2012 (), which we briefly outline below. We work in the Fock basis and represent a general mode quantum process by a tensor of rank which maps the input density matrix into the output one:
(1) 
where underlined symbols refer to multimode Fock states. In practice, the infinite dimensions of both input and output optical Hilbert spaces are truncated to the lowest Fock states, so that .
In the experiment, the black box is tested with a set of mode coherent probe states . For every probe state, the output channels are examined by homodyne measurements, which gives a set of quadrature data , where is the set of local oscillator (LO) phases associated with the th measurement.
To provide enough information about the process, the probe states should cover the volume of interest in the multimode phase space corresponding to the energies up to the chosen photon truncation number . Because the mean quadrature variance of the photon state equals , this volume corresponds to a dimensional hypersphere of radius . On the other hand, a single multimode set of coherent states corresponds to a hypersphere of radius . Therefore the number of the necessary probe states can be estimated as .
Our process reconstruction method relies on the Jamiolkowski isomorphism, relates the superoperator to an operator on the product of and spaces:
(2) 
In this way, the process reconstruction is reduced to a more familiar problem of state reconstruction. The physicality of the process requires it to be completely positive and trace preserving. These conditions are equivalent to the requirement that the corresponding Jamiolkowski operator be positive semidefinite and that , where is identity operator. The latter condition is readily extended to conditional (tracereducing) processes as discussed in Refs. Anis2012 (); Kumar2013 ().
The maximum likelihood reconstruction consists of finding an operator which maximizes the probability of obtaining the harvested data set . Mathematically, this is equivalent to maximization of the functional
(3) 
where is Hermitian matrix of Lagrange multipliers incorporating the tracepreservation condition, and
(4) 
is probability of registering th outcome for the probe state and is the projector corresponding to the th measurement outcome. For deterministic processes, operator maximizing the likelihood (3) satisfies the extremal condition Anis2012 (); Hradil2004 ()
(5) 
where
(6)  
(7) 
Equations (5)–(7) can be solved iteratively, starting from an unbiased . Due to the Hermitian nature of operators and , remains positive semidefinite at each iteration. Together with the trace preservation constraint, this assures physicality of the reconstructed process. The likelihood functional (3) is convex over the space of positive semidefinite operators, which eliminates the possibility of the iteration process stopping at a local maximum.
Iii Tomography of beam splitting
The process of choice for testing the capability of our method is beam splitting. Its paramount importance in quantum optics needs no proof: all linear optical devices (interferometers, waveguide couplers, loss channels, etc.) are equivalent to single beam splitters (BSs) or sets thereof. Any single BS was recently shown to be generator of universal linear optics Bouland2014 (). Accompanied by single photon sources and photon detectors, BSs enable quantum computation Knill2001 (). In some form, a BS is present in any imaginable optical setting. In addition, the BS Hamiltonian is paramount in interfacing quantum information between harmonic oscillator systems of different nature, e.g. between an electromagnetic mode and either an atomic ensemble Hammerer2010 (), or an electromagnetic mode and a nanomechanical oscillator Aspelmeyer2013 ().
Although the operation of the BS is consistent with classical physics (coherent state inputs lead to coherent state outputs, and vice versa), its response to nonclassical input gives rise to quantum phenomena. A striking example is the HongOuMandel effect: when two photons impinge upon a symmetric BS, they appear only in pairs at one of its outputs Hong1987 (). Our technique reveals this quantum effect in spite of using only classical states in measurements.
The BS has previously been characterized by QPT in the role of a Bellstate filter Mitchell2003 () and an amplitude damping channel Bongioanni2010 (). In both these studies, tomography of the BS as a process on a multimode Hilbert space has been incomplete: limited to a specific photon number subspace of that space. Our technique is free of this shortcoming. It allows one to predict output of the process for any arbitrary Fock states and their superpositions in the input, up to a certain cutoff photon number.
Our technique is different from a recently developed methods for characterizing linear optical networks RahimiKeshari2013 () and Gaussian processes Wang2013 () in that it makes no assumptions about the content of the black box, in particular about its Gaussianity or linearoptical character. Although for the demonstration we do use a device which is both linear and Gaussian, our approach can be successfully applied to a multimode process of any nature.
The light source in our experiment is a modelocked Ti:Sapphire laser (Coherent Mira 900), which emits pulses at nm with a repetition rate of 76 MHz and a pulse width of ps. In order to stabilize and control the relative phases of the inputs and outputs, we realize symmetric beam splitting with respect to the horizontal and vertical polarization modes in the same spatial channel, marked 1 in Fig. 1. The polarizations are mixed using an electrooptical modulator (EOM) with its optical axis oriented at to horizontal and a voltage applied to it. A polarizing beam splitter (PBS) subsequently separates the output modes spatially for detection. Our black box is thus implemented by combination EOM + PBS. In the Heisenberg picture, this process has the form
(8) 
where are photon annihilation operators of the input and output modes. The relative amplitudes and phases of the input coherent states are set using a waveplate pair.
The measurement of the output is performed using balanced homodyne detectors (BHDs) Kumar2012a () in both output channels. To this end, we introduce two LOs in orthogonal polarizations in spatial mode 2, so the central PBS directs them into the two output spatial channels (Fig. 1). In each output channel of the PBS, we then find the signal and LO in orthogonal polarizations. For homodyne detection, these polarizations are mixed in each channel using a combination of a plate oriented at to the horizontal and an additional PBS. The relative phases the LOs can be controlled by two wave plates, while their common phase is slowly scanned using a piezomounted mirror in the signal channel.
Iv Evaluating the process tensor
The process reconstruction is simplified by its invariance with respect to the global phase shift. That is, if both input phases are shifted by some phase , so will be the output state. This invariance is a consequence of the symmetric nature of time: a global phase shift by is equivalent to a shift in time by , where is the optical frequency. If the “black box” is not connected to any external clock (such as in our case), it will respond to a signal that is shifted in time by the same amount in the output. The effect of phase invariance on the process tensor can be determined from the fact that a phase shift of both modes will transform density matrix elements according to
Reconciling this with Eq. (1), we find that only elements such that can be nonzero in tensor .
The process reconstruction requires knowledge of the LO phase vector , at each moment in time both for the input and output of the black box. For a general phaseinvariant process, this is equivalent to unknown phase relations. This requirement makes a marked difference between the reconstruction of singlemode and multimode processes. In the singlemode case, many relevant processes exhibit intensityindependent phase behavior, which, in combination with the phase invariance, allows one to disregard phase relations between the input and output modes altogether. In multimode processes, however, this is almost always not the case: even in the relatively simple case of the present work, total phase control is essential for successful reconstruction.
We acquire the phase vector by periodically setting the EOM voltage to zero, so the black box becomes the identity process and the quadrature measurements correspond to the input states. This allows us to monitor all three required phase relations in real time. The inverse sine of the mean quadrature value for each set yields the differences of the LO and input state phases for both modes.
The switching between the BS and identity processes is performed with a period of 0.1 s, which is much faster than the characteristic time of phase fluctuations caused by air movements in the two interferometer channels. In this way, the evaluated LO phases can be translated to the process output measurements by taking into account the linear motion of the piezo.
We acquire a total of 48 sets of quadrature samples for three different relative phases of the LOs: , and rad and, in addition to the vacuum, pairs of input coherent states, obtained by setting each waveplate at and . Each pair of the input states has the same total energy corresponding to photons. This set of probe states is sufficient to reliably reconstruct the process up to a cutoff photon number of .
We implement a twostep reconstruction process as prescribed by Ref. Anis2012 (). In the first step, we artificially inflate the reconstruction Hilbert space by choosing the cutoff point at . This is necessary to ensure that both the input probe states and the output states are well accommodated in that space, which is required for the proper function of the reconstruction algorithm. However, the fraction of 3 and 4photon terms in the Fock decomposition of the probe coherent states is relatively low, and so is their contribution to the loglikelihood functional. As a result, the corresponding terms of the process tensor are not estimated accurately. To eliminate these inaccuracies, we truncate the reconstructed tensor to a lower maximum photon number after the iterations have been completed Anis2012 ().
The phase invariance property of the process kills about of tensor elements. The resulting dimensionality of the optimization space is close to that in the ion tomography done in work Haffner2005 () and is computationally intensive. The iterative algorithm runs on an Intel Core i7 processor. Paralleled onto 4 of 8 computing cores, each iteration takes about 2 hours. The maximumlikelihood reconstruction algorithm appears to converge at around iterations.
V Results
Fig. 2 shows the result of the process reconstruction with in comparison with the theoretical expectation according to Eq. (8) with an additional common phase delay of 0.8 rad. The elements of the process tensor associated with the diagonal elements of the input and output density matrices [Fig. 2(a)] have transparent physical meaning as probabilities of the corresponding transitions. In particular, the HongOuMandel effect is represented by the probability of transition, which is zero for ideally symmetrical BS and amounts to in the reconstructed tensor.
The data in Fig. 2(a) are only a small fraction of the full tensor shown in Fig. 2(b), which has nonzero, generally complex elements. These elements determine the phase behavior of the black box, and are equally important in the description of the process. The left and right columns of the grid present, respectively, the real and imaginary parts of the tensor, while the top and bottom rows correspond to the reconstruction result and the theoretical expectation. The insets in each panel shows the response of our black box to the HongOuMandel query, the input state. The diagonal of the left (real) panel in Fig. 2(b) corresponds to the full panel in Fig. 2(a).
To characterize the quality of the reconstructed tensor shown in Fig. 2, we calculate the fidelity between the ideal and reconstructed processes in the Jamiolkowski state representation:
(9) 
We perform a few tests to find the source of this nonideality. First, we quantify the physical imperfections of our black box by fitting the observed phasedependent mean quadrature data by the theoretical prediction corresponding to an arbitrary BS. We obtain that the power transmittance corresponding to the best fit is 0.502. The fidelity between the processes associated with that slightly asymmetric BS and a symmetric one is 0.998, which shows that the physical errors (at least those which manifest in change of splitting ratio) are insignificant. Second, we evaluate the statistical and systematic errors of the reconstruction using bootstrapping. Specifically, we simulate the quadrature data expected from a model BS and apply the MaxLik reconstruction algorithm to them to calculate a set of tensors . The numbers of simulated data points, the dimensionality of the reconstruction space and the number of algorithm iterations were taken the same as in the real reconstruction procedure. We find for all . Similar values are observed for the pairwise fidelities as well as for the fidelity between the mean of the bootstrapping tensors and the theoretical one. These statistics show that the experimental fidelity of 0.95 results from both the inaccuracy of the numerical reconstruction algorithm and the statistical error conditioned by the limited amount of experimental data.
Vi Summary
We presented experimental csQPT reconstruction of the most common multimode optical process, the beam splitter. Our technique can be readily generalized to other processes, other physical systems and scaled up to a higher number of channels and larger state spaces thanks to the simplicity of the required optical measurements and probe state preparation.
Acknowledgment
We thank Russian Quantum Center for support as well as A. Masalov, D. Mylnikov, A. Bozhenko, and S. Snigirev for helping with the experiment. AKF is a Fellow of the Dynasty foundation. AIL is a CIFAR Fellow.
References
 (1) J.F. Poyatos, J.I. Cirac, and P. Zoller, Phys. Rev. Lett. 78, 390 (1997).
 (2) G.M. D’Ariano and P. Lo Presti, Phys. Rev. Lett. 86, 4195 (2001).
 (3) A. Jamiolkowski, Rep. Math. Phys. 3, 275 (1972).
 (4) M. Lobino, D. Korystov, C. Kupchak, E. Figueroa, B.C. Sanders, and A.I. Lvovsky, Science, 322, 563 (2008).
 (5) M. Lobino, C. Kupchak, E. Figueroa, and A.I. Lvovsky, Phys. Rev. Lett. 102, 203601 (2009).
 (6) S. RahimiKeshari, A. Scherer, A. Mann, A.T. Rezakhani, A.I. Lvovsky, and B.C. Sanders, New J. Phys. 13, 013006 (2011).
 (7) R. Glauber, Phys. Rev. Lett. 10, 84 (1963).
 (8) E.C.G. Sudarshan, Phys. Rev. Lett. 10, 277 (1963).
 (9) J.S. Lundeen, A. Feito, H. ColdenstrodtRonge, K.L. Pregnell, Ch. Silberhorn, T.C. Ralph, J. Eisert, M.B. Plenio, and I.A. Walmsley, Nat. Phys. 5, 27 (2009).
 (10) L. Zhang, H.B. ColdenstrodtRonge, A. Datta, G. Puentes, J.S. Lundeen, X. M. Jin, B. J. Smith, M.B. Plenio, and I.A. Walmsley, Nat. Photonics 6, 364 (2012).
 (11) Z. Hradil, J. Řeháček, J. Fiurášek, M. Ježek, Lecture Notes in Physics 649, 59 (2004).
 (12) A. Anis and A.I. Lvovsky, New J. Phys. 14, 105021 (2012).
 (13) R. Kumar, E. Barrios, C. Kupchak, and A. I. Lvovsky, Phys. Rev. Lett. 110, 130403 (2013).
 (14) M. Cooper, E. Slade, M. Karpinski, and B.J. Smith, New J. Phys. 17, 033041 (2015).
 (15) M. Afzelius, C. Simon, H. de Riedmatten, and N. Gisin, Phys. Rev. A 79, 052329 (2009).
 (16) S.A. Moiseev, S.N. Andrianov, and F.F. Gubaidullin, Phys. Rev. A 82, 022311 (2010).
 (17) C. Monroe, Nature 416, 238 (2002).
 (18) T.B. Pittman, M.J. Fitch, B.C. Jacobs, and J.D. Franson, Phys. Rev. A 68, 032316 (2003)
 (19) A. Bouland and S. Aaronson, Phys. Rev. A 89, 062316 (2014).
 (20) E. Knill, R. Laflamme, and G.J. Milburn, Nature 409, 46 (2001).
 (21) K. Hammerer, A.S. Sørensen, and E.S. Polzik, Rev. Mod. Phys. 82, 1041 (2010).
 (22) M. Aspelmeyer, T.J. Kippenberg, and F. Marquardt, Cavity Optomechanics (Springer Berlin Heidelberg, Berlin, 2014).
 (23) C.K. Hong, Z.Y. Ou, and L. Mandel, Phys. Rev. Lett. 59, 2044 (1987).
 (24) M.W. Mitchell, C.W. Ellenor, S. Schneider, and A.M. Steinberg, Phys. Rev. Lett. 91, 120402 (2003).
 (25) I. Bongioanni, L. Sansoni, F. Sciarrino, G. Vallone, and P. Mataloni, Phys. Rev. A 82, 042307 (2010).
 (26) S. RahimiKeshari, M. Broome, R. Fickler, A. Fedrizzi, T. Ralph, and A. White, Opt. Express 21, 13450 (2013).
 (27) X.B. Wang, Z.W. Yu, J.Z. Hu, A. Miranowicz, and F. Nori, Phys. Rev. A 88, 022101 (2013).
 (28) R. Kumar, E. Barrios, A. MacRae, E. Cairns, E.H. Huntington, and A.I. Lvovsky, Opt. Commun. 285, 5259 (2012).
 (29) H. Häffner, W. Hänsel, C. F. Roos, J. Benhelm, D. Chekalkar, M. Chwalla, T. Körber, U.D. Rapol, M. Riebe, P.O. Schmidt, C. Becher, O. Gühne, W. Dür, and R. Blatt, Nature 438, 643 (2005).
 (30) U. Leonhardt, Measuring the quantum state of light, (Cambridge University Press, Cambridge, 1997).