## Abstract

Unphysical solutions are ruled out in physical equations, as they lead to behavior that violates fundamental physical laws. One of the celebrated equations that allows unphysical solutions is the relativistic Majorana equation, thought to describe neutrinos and other exotic particles predicted in theories beyond the standard model. The neutrally charged Majorana fermion is the equation’s physical solution, whereas the charged version is, due to charge nonconservation, unphysical and cannot exist. Here, we present an experimental scheme simulating the dynamics of a charged Majorana particle by light propagation in a tailored waveguide chip. Specifically, we simulate the free-particle evolution as well as the unphysical operation of charge conjugation. We do this by exploiting the fact that the wave function is not a directly observable physical quantity and by decomposing the unphysical solution to observable entities. Our results illustrate the potential of investigating theories beyond the standard model in a compact laboratory setting.

© 2015 Optical Society of America

## 1. INTRODUCTION

When Ettore Majorana wrote down his famous equation in 1937 [1,2], he explicitly suggested describing the characteristics of neutrinos on its basis. He noted that Lorentz invariance allowed not only the Dirac equation, but also the expression ($\hslash \equiv c\equiv 1$)

for the wave function $\psi $ of a particle with (Majorana) mass $m$ and its charge conjugate ${\psi}_{\mathrm{c}}$. The appearance of the so-called Majorana mass term points to the violation of charge conservation, suggesting that a particle obeying the Majorana equation must be its own antiparticle. For this physical reason, $\psi $ is commonly taken to be charge-neutral; i.e., the Majorana equation is frequently supplemented by the condition $\psi ={\psi}_{\mathrm{c}}$ (the resulting particle is called the Majorana fermion) [3]. To date, no elementary particle has been identified as a Majorana fermion. However, as Majorana originally envisioned, there is the possibility that the neutrino is a Majorana fermion. In this case, the corresponding lepton number would not be conserved and the nature of the neutrino can therefore be tested by lepton number nonconserving processes such as neutrinoless double-beta decay [4]. The concept of Majorana fermions has also found use in condensed-matter physics, where quasi-particle excitations can be their own antiparticle. This happens, for example, in superconducting systems where the Bogoliubov quasi-particles, whose evolution is described by the Bogoliubov–De Gennes equations, are Majorana fermions [5]. Such quasi-particles can be non-Abelian anyons whose non-Abelian braiding statistics form the building blocks of topological quantum computation [6,7].On the other hand, the fact that the charged version of a Majorana fermion, the Majoranon [8,9], violates charge conservation may provide access to physics beyond the standard model. In many theories, a potential violation of charge conservation, for example, associated with higher spacetime dimensions [10] or a nonvanishing photon mass [11], is considered. In addition, simulating unphysical effects may yield unexpected benefits in other areas, as recently shown for the case of complex conjugation that provides an efficient method to measure entanglement [12].

In this work, we break new ground and devise an experimental scheme to simulate the dynamics of a Majoranon, thereby implementing a classical optical simulator of an unphysical particle. To this end, we consider the Majorana equation in $1+1$-dimensional spacetime, which reads for the two-component spinor $\psi =\left(\genfrac{}{}{0ex}{}{{\psi}_{1}}{{\psi}_{2}}\right)$ as

Here, ${p}_{x}$ is the momentum along the spatial coordinate and we have used the representation such that ${\psi}_{\mathrm{c}}=-i{\sigma}_{z}{\sigma}_{y}{\psi}^{*}=-\left(\genfrac{}{}{0ex}{}{\begin{array}{c}{\psi}_{2}^{*}\end{array}}{\begin{array}{c}{\psi}_{1}^{*}\end{array}}\right)$, where ${\sigma}_{x}$, ${\sigma}_{y}$, ${\sigma}_{z}$ are the Pauli matrices. One cannot directly simulate this equation due to the fact that it contains a complex conjugation, which renders its Hamiltonian formulation impossible [8]. However, somewhat surprisingly, this equation can be decomposed into two physical equations [9]: which are the Dirac equations with positive and negative mass terms. Here, ${\psi}_{\pm}$ are two charge-neutral Majorana fields, i.e., ${\psi}_{\mathrm{c};\pm}=-i{\sigma}_{z}{\sigma}_{y}{\psi}_{\pm}^{*}={\psi}_{\pm}$, and the original Majoranon field can be reconstructed from them via the simple relationshipImportantly, the Dirac equation itself is a physical equation and can be presented in the Hamiltonian form. As such, it can be simulated using various systems such as trapped ions [13,14] or light [15]. Physical operations in this decomposed Hilbert space of two independent Majorana fermions can be used to simulate unphysical operations acting on the Majoranon, i.e., complex conjugation and charge conjugation, to which the evolution is intrinsically linked [cf. Eq. (2)]. It is important to note, however, that while the two simulations of the Dirac equations can be performed in parallel, they need to be mutually coherent, such that the Majoranon field $\psi $ can be reconstructed. Using a photonic chip setup we implement such a coherent, parallel simulation of the free evolution of a Majoranon. On top of demonstrating the unphysical Majoranon dynamics directly by measuring the absolute values of the spinor components, we also compare the dynamics of a Majoranon with its Dirac cousin—the same initial spinor following the Dirac evolution. Note that discrepancies between the two arise from the difference in the term proportional to the mass that renders the Majorana equation unphysical. To illustrate these discrepancies, we evaluate the quantity $\u3008{\sigma}_{z}\u3009=\sum _{n}{|{\psi}_{1,n}|}^{2}-{|{\psi}_{2,n}|}^{2}$ [16]. For a Dirac particle at rest (${p}_{x}=0$, or equivalently $m\to \infty $), it measures the population difference between the positive and negative energy branches and is a conserved quantity. On the contrary, it is not conserved for the Majoranon at rest, but oscillates due to the unphysical mass term that continuously forces exchanges between the spinor components. Borrowing from the physics of the Dirac particle, we will hereafter call this quantity a pseudo-energy for convenience.

## 2. OPTICAL SIMULATION IN WAVEGUIDE LATTICES

Our system consists of two $1+1$-dimensional photonic lattices, each composed of a periodic array of waveguides that are evanescently coupled to one another. Such waveguide lattices have attracted considerable interest and have been used in the exploration of a number of fundamental wave-transport phenomena, including Anderson localization [17], discrete solitons [18], and photonic topological insulators [19]. In order to describe the light evolution along the longitudinal spatial axis $Z$ in a waveguide array, one commonly employs a coupled-mode approach [20], which yields

where ${\psi}_{k}$ is the field amplitude in the $k$th lattice site, $\kappa $ is the coupling between the waveguides, and ${\beta}_{k}$ is a position-dependent detuning.#### A. Dirac Dynamics

When a broad input beam with an initial phase shift of $\pi /2$ between adjacent guides (lateral waveguide spacing ${d}_{0}$) is launched into a binary waveguide array composed of two interleaved sublattices A and B with different refractive indices amounting to detunings $\pm \beta $, the light evolution can be approximated by [15,21]

*Zitterbewegung*of a relativistic electron [23]. Figure 1(a) shows an experimentally observed photonic

*Zitterbewegung*in a photonic lattice using a tailored input phase distribution (see Sections 2.B.1 and 2.B.2 for details on the experimental implementation). A numerical simulation of the

*Zitterbewegung*, based on Eq. (5), is shown in Fig. 1(b). The close correspondence proves the ability to simulate the Dirac equation in a waveguide lattice.

#### B. Majorana Dynamics

In our setting, we make use of the capacity for coherent Dirac simulation and let two light beams propagate along two parallel planar waveguide lattices with masses of opposite sign, such that the two Dirac equations for ${\psi}_{+}$ and ${\psi}_{-}$ [see Eq. (3)] are simulated in parallel, leading to *Zitterbewegung* in opposite directions (see Fig. 2). After the desired propagation distance (corresponding to a specific evolution time), the amplitude distributions are coherently combined in order to retrieve the Majoranon wave function according to Eq. (4). In the following, the structure of the simulator, its fabrication, and the encoding and readout of input and output states, respectively, will be discussed in detail.

### 1. Design and Fabrication of the Simulator

The experimental platform for the simulation of the Majorana equation consists of two binary waveguide lattices, which only differ in the ordering of the two sites A and B forming a unit cell. The first part of the sample (right-hand side in Fig. 2) is occupied by the encoding stage (see Section 2.B.2). In the central part, the Dirac equation, Eq. (6), with positive (negative) mass is simulated over the evolution length ${L}_{\mathrm{e}}$ in the upper (lower) lattice. In this discrete setting, each spinor amplitude ${\psi}_{+,n}$ in unit cell $n$ of the upper plane has its counterpart ${\psi}_{-,n}$ in the lower plane. By construction, the first spinor components ${\psi}_{\pm ,n,1}$ are distributed over the odd lattice sites, whereas the second components ${\psi}_{\pm ,n,2}$ are found on the even sites.

The evolution is terminated by a fan-out section of length ${L}_{\mathrm{f}}$, in which the waveguide separation is increased until no more significant evanescent coupling takes place. This fan-out trajectory follows a harmonic curve, and ${L}_{\mathrm{f}}$ is sufficiently long to ensure that bending losses are negligibly small. Due to the gradual reduction of the coupling strength in this section, some residual evolution takes place, which effectively extends the evolution length to some value ${L}_{\mathrm{e},\mathrm{eff}}>{L}_{\mathrm{e}}$ [24].

Finally, all waveguide pairs of the two planes are mutually connected by vertical directional couplers of length ${L}_{\mathrm{c}}$. For balanced couplers [20], the output amplitudes in the upper ports are proportional to the discrete Majorana spinor ${\psi}_{n}={\psi}_{+,n}+i{\psi}_{-,n}$. Thus the desired recombination of the two spinors, Eq. (4), is performed in an integrated and spatially resolved fashion.

Similarly, the amplitude in the lower port is proportional to the charge-conjugated spinor ${\psi}_{\mathrm{c};n}=-i{\sigma}_{z}{\sigma}_{y}{\psi}_{n}^{*}={\psi}_{+,n}-i{\psi}_{-,n}$. Two simulators with two different masses $\beta $ have been realized in the course of this work. In the configuration with the larger mass ($\beta =1.2\kappa $; presented in Figs. 4 and 5) these lower output ports reach the output facet of the sample and are experimentally accessible (not shown in Fig. 2).

The waveguides are inscribed in bulk fused silica by nonlinear absorption of focused (numerical aperture 0.35) pulsed laser radiation (wavelength 800 nm, pulse duration $\tau $, pulse energy ${E}_{\mathrm{p}}$, repetition rate 100 kHz). These nonlinear absorption processes lead to a permanent increase of the refractive index of the material. By translating the material with velocity ${v}_{0}$ on a certain path through the focus, a waveguide channel is written [15,19,22]. The fabrication parameters are $\tau =150\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{fs}$, ${E}_{\mathrm{p}}=300\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nJ}$, and ${v}_{0}=100\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}/\mathrm{min}$ for the Dirac lattice of Fig. 1 as well as the low-mass Majoranon simulator ($\beta =0.65\kappa $) of Fig. 3 and $\tau =120\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{fs}$, ${E}_{\mathrm{p}}=260\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nJ}$, and ${v}_{0}=90\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}/\mathrm{min}$ for the high-mass lattice shown in Figs. 4 and 5, respectively. The waveguide separation in the evolution section is ${d}_{0}=18.5(19.5)\pm 0.3\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{\mu m}$ for the low (high)-mass lattice, and the refractive index difference between sublattices A and B is realized by modulating the inscription velocity by $\pm 6(14)\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}/\mathrm{min}$.

### 2. Encoding of the Input State and Experimental Observation Technique

We investigate an initial Majoranon wave packet of width $\sigma $, centered around position ${n}_{0}$, with zero average momentum and occupation of only the first spinor component, i.e., ${\psi}_{n}(z=0)\propto \mathrm{exp}(-{(n-{n}_{0})}^{2}/2{\sigma}^{2})(\genfrac{}{}{0ex}{}{1}{0})$. The corresponding decomposed spinors are then given by ${\psi}_{+[-],n}\propto \mathrm{exp}(-{(n-{n}_{0})}^{2}/2{\sigma}^{2})(\genfrac{}{}{0ex}{}{1}{-1})[-(\genfrac{}{}{0ex}{}{i}{i})]$ [9]. In order to ensure equal amplitude distributions in the two planes simulating ${\psi}_{+}$ and ${\psi}_{-}$, balanced directional couplers, each with only a single input port, are introduced at the front end of the device, which is then illuminated by a spatially extended beam in the experiment (see Fig. 2). The beam has a flat-phased Gaussian profile with a waist radius ($1/\mathrm{e}$ intensity) of 40(50) μm for the low (high)-mass device, corresponding to $\sigma =1.1(1.3)$, and a wavelength of $\lambda =633\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}$.

Due to the mapping from the spinors to light amplitudes [21], the two Dirac lattices with opposing masses require a phase shift of $\pi /2$ between adjacent waveguides at the start of the evolution, but with opposite directions of the phase gradient. This is implemented by a tailored phase segmentation of the waveguides, i.e., an intentional periodic omission of waveguide sections [22,25]. The period of this segmentation is 40 μm, and the filling factor is $1/2$. For $\lambda =633\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}$, a phase shift of $-j\pi /2$ is introduced by a segmented section of length $js$, with $j=0,\dots ,3$ and $s=1.76(1.85)\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$ for the low (high)-mass lattice (see inset of Fig. 2).

The intensity evolution in a single Dirac lattice is observed directly by the fluorescence of color centers in the waveguides [26], whereas the evolution in the Majoranon simulator is inferred from the measured output intensity distributions after the recombination step. For each value of the mass, two samples with two different evolution lengths have been fabricated.

The measurement uncertainty of these classical light intensities is negligibly small due to a high signal-to-noise ratio on the camera. Yet, random errors are introduced to the simulator by fabrication tolerances of the waveguide system, triggering random deviations from the target parameters. The positioning precision of the waveguides relative to each other is $\pm 0.3\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{\mu m}$, yielding a relative uncertainty of $\frac{\mathrm{\Delta}\kappa}{\kappa}\approx 0.06$ for all in-plane coupling strengths and causing errors in the intensity splitting ratio of the vertical couplers by $\pm 4\%$ around the target value of 50%. It was determined from independent measurements that temporal fluctuations in the inscription parameters cause random variations of each waveguide’s propagation constant by $\frac{\mathrm{\Delta}{\beta}_{k}}{\kappa}\approx 0.04$. The error bars on the following experimental results (shown in Figs. 3, 4, and 5) are derived from these uncertainties via simulations of the light propagation through a Gaussian-distributed ensemble of 10,000 devices for each setting.

## 3. EXPERIMENTAL RESULTS

#### A. Observed Majorana Evolution

Figure 3 shows our experimental results in the low-mass lattice ($\beta =0.65\kappa $; $\kappa =0.064\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$) consisting of 26 waveguides, i.e., discretization points $n=1\dots ,13$ for the spinors. In Figs. 3(a) and 3(b), the computed parallel evolution of both components of the Majoranon spinor is presented. We observe that although initially all intensity is concentrated in ${\psi}_{1}$, it immediately starts to oscillate between the two spinor components and, at the same time, to spread along the transverse space coordinate. Using our photonic structure, we observe the population of both spinor components at two different propagation distances. For a small effective evolution distance of $Z={L}_{\mathrm{e},\mathrm{eff}}=0.55{\kappa}^{-1}=8.6\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$, the light mostly remains in odd waveguide sites, which heralds the prevalent occupation of ${\psi}_{1}$ (Fig. 3(c)). For a larger distance of $Z=4.4{\kappa}^{-1}$, one expects another minimum of spinor 2 accompanied by extensive spreading of the wave packet (cf. Figs. 3(a) and 3(b)). Indeed, most of the light is again trapped in the odd channels and the entire wave packet is spread over a much larger spatial region (Fig. 3(d)). The individual spinor intensities, which are equivalent to the light intensities on the odd/even sites, are shown in Figs. 3(e) and 3(f), together with the theoretical data. Differences in the spinor distribution between theory and experimental data can be attributed to the fabrication precision of the simulator (see error bars) as well as systematic errors arising at the phase preparation stage (mainly losses and the gradual accumulation of the phases in the segmented channels). At both lengths, the population of ${\psi}_{1}$ predominates ${\psi}_{2}$.

#### B. Pseudo-Energy Oscillation

In Fig. 3(g) we show the expected unphysical oscillations in the pseudo-energy of the Majoranon as discussed earlier. The measured values of $\u3008{\sigma}_{z}\u3009$ at the two evolution lengths agree with the expected values within one standard deviation, while displaying a significant difference from the calculated pseudo-energy of the same initial spinor subjected to the Dirac equation (6). Note that the oscillations in pseudo-energy for the Dirac particle and the Majoranon occur for entirely different reasons: the oscillation for the Dirac particle occurs due to nonzero momentum components in the initial wave packet, while the oscillation for the Majoranon is mainly due to the unphysical mass term.

To elaborate on this difference further, we investigate the evolution of a Majoranon with a larger mass, where the nonzero momentum contribution plays a smaller role. For this purpose, we measure the light evolution in the second system with $\beta =1.2\kappa $, where $\kappa =0.072\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{-1}$ and 30 lattice sites were used. The results are summarized in Fig. 4. Due to the reduced momentum contribution in the evolution, the amplitude of the oscillation in pseudo-energy has gotten smaller for the Dirac particle, resulting in larger discrepancies from the Majoranon, whose oscillation amplitude is not affected by the increase in mass (see Fig. 4(g)). The oscillation frequency, however, has increased, such that already at small distances $Z=0.9{\kappa}^{-1}$, mostly ${\psi}_{2}$ is populated (see Figs. 4(a)–4(c) and 4(e)). After a distance of $Z=3.5{\kappa}^{-1}$, a further oscillation period has occurred, leading again to a strong population of ${\psi}_{2}$. However, the transverse spreading of the wave packet is much less pronounced than for the smaller mass of $\beta =0.65\kappa $, as is clearly visible from Figs. 4(d) and 4(f). This is consistent with the fact that the amplitude of the *Zitterbewegung* of ${\psi}_{\pm}$ decreases for larger masses, whereas the frequency is increased [14,15].

#### C. Implementation of Charge Conjugation

Finally, we demonstrate how the unphysical operation of charge conjugation can be simulated in our system. To this end, one can make use of the fact that the final waveguide couplers possess two output ports. So far, only the upper ports have been used to produce the Majoranon wave function. In the lower ports, the two Dirac wave packets are recombined with the opposite phase, which yields the charge-conjugated spinor ${\psi}_{\mathrm{c}}(Z)={\psi}_{+}(Z)-i{\psi}_{-}(Z)$, up to a phase, after the propagation distance $Z$. The measured and numerically simulated output intensity distributions and reconstructed spinors are displayed in Figs. 5(a) and 5(b). It is evident that for both evolution lengths, the occupation of the spinor components is reversed: whereas the second component is prevalent in the Majorana spinor (cf. Fig. 4), its charge-conjugated brother has mostly its first component occupied. That is exactly what one would expect from the definition of charge conjugation: ${\psi}_{\mathrm{c}}=\left(\genfrac{}{}{0ex}{}{\begin{array}{c}-{\psi}_{2}^{*}\end{array}}{\begin{array}{c}{\psi}_{1}^{*}\end{array}}\right)$. Accordingly, the pseudo-energy $\u3008{\sigma}_{z}\u3009$ after charge conjugation oscillates in antiphase with that of the unconjugated particle. The experimental observation of the pseudo-energy at the two evolution lengths lies within one standard deviation of this theoretical prediction (symbols and solid lines in Fig. 5(c)). Dirac particles of either charge, on the other hand, would preserve their pseudo-energy, and merely momentum-induced oscillations would occur (dashed lines). This illustrates the striking difference between the physical Dirac equation, where charge is preserved, and the unphysical Majorana equation, which violates charge conservation.

## 4. DISCUSSION

In our work, we coherently simulated the dynamics of a Majoranon wave packet by classical optics in a compact integrated waveguide architecture. The evolution of the free particle is unphysical due to the fact that it involves charge conjugation and complex conjugation. Moreover, the simulator permits us to directly implement the operation of charge conjugation on the evolving wave packet. Here, evolution and charge conjugation are performed simultaneously for arbitrary, but predefined, evolution lengths. This differs from the proposed trapped ion scheme [8], where such operations are performed sequentially and independently from one another. A specific advantage of our classical scheme lies in the possibility that both the unconjugated and the charge-conjugated Majorana spinors can be accessed simultaneously.

Simulating such unphysical operations provides an entirely new approach for probing and understanding exotic phenomena and particles that cannot exist in nature, such as the Majoranon. Our scheme uses the fact that even for real particles the wave function itself is not a physical entity, but only its square modulus is. Hence, the superposition of such wave functions can result in an unphysical phenomenon, which means, conversely, that the latter can be reproduced by simulating the individual wave functions. Many interesting questions are prompted, concerning, e.g., possible decay mechanisms of the Majoranon, the impact of many-body effects and interactions, and their scattering characteristics. For example, the dynamics of a pair of Majoranons, which remain to be explored, will be influenced by their exotic exchange symmetry as well as their Coulomb interaction. An arbitrary Abelian exchange symmetry can be investigated directly in the same setup by resorting to photon pairs, which are entangled in another degree of freedom, such as polarization [27,28]. Interactions, on the other hand, are most effectively simulated by resorting to two-dimensional waveguide architectures [29,30]. We would also like to note that superconducting systems may provide an alternative platform to simulate the Majorana equation with the potential for studying the effects of quantization [compare Eq. (3) with Eq. (20) of Ref. [5], for example]. Finally, we anticipate that this first (to the best of our knowledge) explicit demonstration of unphysical operations in the laboratory will stimulate many exciting proposals that utilize the freedom of going beyond the ‘physical’ operations in areas such as exotic particle physics and quantum information processing. We note that after the completion of this work, an experimental simulation of the Majorana equation in a trapped ion system has appeared [31].

## FUNDING INFORMATION

German Ministry of Education and Research (03Z1HN31); German–Israeli Foundation for Scientific Research and Development (GIF) (1157-127.14/2011); Ministry of Education—Singapore (MOE) (MOE2012-T3-1-009); Singapore National Research Foundation; Thuringian Ministry for Education, Science and Culture (11027-514).

## REFERENCES

**1. **E. Majorana, “Teoria simmetrica dell’elettrone e del positrone,” Nuovo Cimento **14**, 171–184 (1937). [CrossRef]

**2. **A. Zee, *Quantum Field Theory in a Nutshell* (Princeton University, 2003).

**3. **C. Giunti and C. W. Kim, *Fundamentals of Neutrino Physics and Astrophysics* (Oxford University, 2007).

**4. **R. N. Mohapatra and P. B. Pal, *Massive Neutrinos in Physics and Astrophysics*, Vol. 72 of Lecture Notes in Physics (World Scientific, 2004).

**5. **C. W. J. Beenakker, “Random-matrix theory of Majorana fermions and topological superconductors,” arXiv:1407.2131 (2015).

**6. **A. Y. Kitaev, “Anyons in an exactly solved model and beyond,” Ann. Phys. (N.Y.) **321**, 2–111 (2006). [CrossRef]

**7. **C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. Das Sarma, “Non-Abelian anyons and topological quantum computation,” Rev. Mod. Phys. **80**, 1083–1159 (2008). [CrossRef]

**8. **J. Casanova, C. Sabin, J. León, I. L. Egusquiza, R. Gerritsma, C. F. Roos, J. J. Garcya-Ripoll, and E. Solano, “Quantum simulation of the Majorana equation and unphysical operations,” Phys. Rev. X **1**, 021018 (2011). [CrossRef]

**9. **C. Noh, B. M. Rodríguez-Lara, and D. G. Angelakis, “Proposal for realization of the Majorana equation in a tabletop experiment,” Phys. Rev. A **87**, 040102(R) (2013). [CrossRef]

**10. **S. N. Gninenko, “Limit on the electric charge-nonconserving *μ*^{+} → invisible decay,” Phys. Rev. D **76**, 055004 (2007). [CrossRef]

**11. **A. S. Goldhaber and M. M. Nieto, “Photon and graviton mass limits,” Rev. Mod. Phys. **82**, 939–979 (2010). [CrossRef]

**12. **R. Di Candia, B. Mejia, H. Castillo, J. S. Pedernales, J. Casanova, and E. Solano, “Embedding quantum simulators for quantum computation of entanglement,” Phys. Rev. Lett. **111**, 240502 (2013). [CrossRef]

**13. **L. Lamata, J. León, T. Schätz, and E. Solano, “Dirac equation and quantum relativistic effects in a single trapped ion,” Phys. Rev. Lett. **98**, 253005 (2007). [CrossRef]

**14. **R. Gerritsma, G. Kirchmair, F. Zähringer, E. Solano, R. Blatt, and C. F. Roos, “Quantum simulation of the Dirac equation,” Nature **463**, 68–71 (2010). [CrossRef]

**15. **F. Dreisow, M. Heinrich, R. Keil, A. Tünnermann, S. Nolte, S. Longhi, and A. Szameit, “Classical simulation of relativistic Zitterbewegung in photonic lattices,” Phys. Rev. Lett. **105**, 143902 (2010). [CrossRef]

**16. **L. Lamata, J. Casanova, I. L. Egusquiza, and E. Solano, “The nonrelativistic limit of the Majorana equation and its simulation in trapped ions,” Phys. Scr. **T147**, 014017 (2012). [CrossRef]

**17. **T. Schwartz, G. Bartal, S. Fishman, and M. Segev, “Transport and Anderson localization in disordered two-dimensional photonic lattices,” Nature **446**, 52–55 (2007). [CrossRef]

**18. **J. W. Fleischer, M. Segev, N. K. Efremidis, and D. N. Christodoulides, “Observation of two-dimensional discrete solitons in optically induced nonlinear photonic lattices,” Nature **422**, 147–150 (2003). [CrossRef]

**19. **M. C. Rechtsman, J. M. Zeuner, Y. Plotnik, Y. Lumer, D. Podolsky, F. Dreisow, S. Nolte, M. Segev, and A. Szameit, “Photonic Floquet topological insulators,” Nature **496**, 196–200 (2013). [CrossRef]

**20. **A. Yariv, “Coupled-mode theory for guided-wave optics,” IEEE J. Quantum Electron. **9**, 919–933 (1973). [CrossRef]

**21. **S. Longhi, “Photonic analog of Zitterbewegung in binary waveguide arrays,” Opt. Lett. **35**, 235–237 (2010). [CrossRef]

**22. **R. Keil, J. M. Zeuner, F. Dreisow, M. Heinrich, A. Tünnermann, S. Nolte, and A. Szameit, “The random mass Dirac model and long-range correlations on an integrated optical platform,” Nat. Commun. **4**, 1368 (2013).

**23. **E. Schrödinger, “Über die kräftefreie Bewegung in der relativistischen Quantenmechanik,” Sitz. Preuss. Akad. Wiss. Phys.-Math. Kl. **24**, 418–428 (1930).

**24. **A. Peruzzo, M. Lobino, J. C. F. Matthews, N. Matsuda, A. Politi, K. Poulios, X.-Q. Zhou, Y. Lahini, N. Ismail, K. Wörhoff, Y. Bromberg, Y. Silberberg, M. G. Thompson, and J. L. O’Brien, “Quantum walks of correlated photons,” Science **329**, 1500–1503 (2010). [CrossRef]

**25. **A. Szameit, F. Dreisow, M. Heinrich, T. Pertsch, S. Nolte, A. Tünnermann, E. Suran, F. Louradour, A. Barthélémy, and S. Longhi, “Image reconstruction in segmented femtosecond laser-written waveguide arrays,” Appl. Phys. Lett. **93**, 181109 (2008). [CrossRef]

**26. **F. Dreisow, M. Heinrich, A. Szameit, S. Döring, S. Nolte, A. Tünnermann, S. Fahr, and F. Lederer, “Spectral resolved dynamic localization in curved fs laser written waveguide arrays,” Opt. Express **16**, 3474–3483 (2008). [CrossRef]

**27. **L. Sansoni, F. Sciarrino, G. Vallone, P. Mataloni, A. Crespi, R. Ramponi, and R. Osellame, “Two-particle bosonic-fermionic quantum walk via integrated photonics,” Phys. Rev. Lett. **108**, 010502 (2012). [CrossRef]

**28. **J. C. F. Matthews, K. Poulios, J. D. A. Meinecke, A. Politi, A. Peruzzo, N. Ismail, K. Wörhoff, M. G. Thompson, and J. L. O’Brien, “Observing fermionic statistics with photons in arbitrary processes,” Sci. Rep. **3**, 1539 (2013).

**29. **G. Corrielli, A. Crespi, G. Della Valle, S. Longhi, and R. Osellame, “Fractional Bloch oscillations in photonic lattices,” Nat. Commun. **4**, 1555 (2013). [CrossRef]

**30. **X. Qin, Y. Ke, X. Guan, Z. Li, N. Andrei, and C. Lee, “Statistics-dependent quantum co-walking of two particles in one-dimensional lattices with nearest-neighbor interactions,” Phys. Rev. A **90**, 062301 (2014). [CrossRef]

**31. **X. Zhang, Y. Shen, J. Zhang, J. Casanova, L. Lamata, E. Solano, M.-H. Yung, J.-N. Zhang, and K. Kim, “Time reversal and charge conjugation in an embedding quantum simulator,” arXiv:1409.3681 (2014).