Simulation of anyonic statistics and its topological path independence using a 7qubit quantum simulator
Abstract
Anyons, quasiparticles living in twodimensional spaces with exotic exchange statistics, can serve as the fundamental units for faulttolerant quantum computation. However, experimentally demonstrating anyonic statistics is a challenge due to the technical limitations of current experimental platforms. Here, we take a state perpetration approach to mimic anyons in the Kitaev lattice model using a 7qubit nuclear magnetic resonance quantum simulator. Anyons are created by dynamically preparing the ground and excited states of the 7qubit Kitaev lattice model, and are subsequently braided along two distinct, but topologically equivalent, paths. We observe that the phase acquired by the anyons is independent of the path, and coincides with the ideal theoretical predictions when decoherence and implementation errors are taken into account. As the first demonstration of the topological path independence of anyons, our experiment helps to study and exploit the anyonic properties towards the goal of building a topological quantum computer.
pacs:
03.67.Ac, 03.67.Lx, 05.30.Pr1 Introduction
It is a fundamental question to investigate the physical properties of exchanging two identical particles. In onedimensional space, the exchange is trivial as the two particles inevitably collide. In threedimensional spaces, a wave function of a system acquires either +1 or 1 phase factor for bosons and fermions upon the exchange, respectively. Alternatively, in twodimensional spaces, exotic statistical properties emerge when exchanging two identical particles, leading to the theoretical existence of anyons. For instance, the wave function can acquire an arbitrary phase factor , ranging continuously from +1 to 1 for Abelian anyons, or undergo nontrivial unitary evolutions for nonAbelian anyons. This exotic feature of anyons called fractional statistics has attracted great interest over the past few decades Khare (2005); Tsui et al. (1982); Moore and Read (1991); Wen (1991); Nayak et al. (2008); Stern (2010).
As truly twodimensional systems do not exist in nature, anyons appear in an effective twodimensional system as quasiparticles, a collective behaviour of a group of fundamental particles behaving as a single entity. The experimental evidence of quasiparticle anyons was first discovered in the fractional quantum hall effect in the 1980s Tsui et al. (1982). Since then, an experimental quest for anyons has interested researchers not only for their fundamentally intriguing feature, but also for their application in performing protected quantum information processing (QIP). As the goal of QIP is to exploit quantum mechanical properties for computation, manipulating quantum properties in a precise manner is critical Ladd et al. (2010). Topological properties of anyons have potential to achieve such a goal, as these properties are resilient to small fluctuations. Therefore, the prospect of utilizing anyons’ fractional statistics to achieve robust control has gained much attention. In the past two decades, many powerful quantum computing schemes using anyons have been proposed. Such schemes constitute topological quantum computing (TQC) Sarma et al. (2005); Nayak et al. (2008); Freedman et al. (2003) and topological quantum error correction Raussendorf et al. (2007); Kitaev (2003). Of the many proposals, the Kitaev model (KM) of spin lattice Kitaev (2003) is one of the most renowned. By artificially designing a spin lattice model with highly nontrivial ground states, individual localized Abliean anyons can be created and manipulated, leading to the realization of TQC.
Despite the prospective applications of anyonic statistics, realizing these ideas in experiments remains a challenge, as such tasks typically require generating and manipulating complex manybody quantum systems. Nevertheless, significant progress has been made in smallscale systems theoretically Han et al. (2007); Aguado et al. (2008) and experimentally Pachos et al. (2009); Lu et al. (2009); Du et al. (2007); Feng et al. (2013). Through the quantum simulation approach Feynman (1982); Lloyd (1996); Georgescu et al. (2014); Friedenauer et al. (2008); Kim et al. (2010); Lanyon et al. (2010); Du et al. (2010); Lanyon et al. (2011); Lu et al. (2011); Islam et al. (2013); Li et al. (2014), in which the experimental setup acts as a processor to mimic the dynamics of anyonic systems, several experiments have been implemented to demonstrate the exotic properties of anyons in small systems Pachos et al. (2009); Lu et al. (2009); Du et al. (2007); Feng et al. (2013). These experiments provide better understanding of braiding operations in realistic noise, opening up the possibility of fully utilizing the advantages of anyonic fractional statistics.
However, the path independent nature of the anyons’ braiding statistics has not been demonstrated yet as it requires larger quantum simulators with highfidelity coherent control. In this paper, we study a 7qubit system with three different paths to braid anyons: two nontrivial paths where the wave function picks up the phase, and one trivial path where the wave function remains unchanged. The additional nontrivial loop allows the experimental proofofprinciple demonstration of path independence (i.e. the statistics do not depend on the shape of the path taken by the anyons as long as the exchange takes place) Han et al. (2007). This topological feature is one of the key advantages for utilizing anyons. In our experiments, we use a 7qubit NMR quantum simlulator to realize the three braiding paths through the state preparation approach, and observe that the two phases acquired during the two nontrivial loops agree within experimental uncertainty, although they are below the theoretical value of . Our primary source of error is decoherence that leads to deviations between the theoretically predicted phases and the experimental ones; however, we analyze the errors quantitatively and show that the two nontrivial phases are close to after accounting for such errors.
The remainder of the paper proceeds as follows. In Sec. 2, we review the original KM and the simplified 7qubit KM and describe in detail how anyons gain a phase after braiding regardless of the shapes of the nontrivial paths. In Sec. 3, we briefly introduce our experimental setup and the mapping between the KM model and the experimental system, followed by stepbystep description of the implementation in an NMR system. Finally, in Sec. 4, we show the experimental results and analyze the errors. We account for experimental deviations from the theoretical predictions using numerical simulations that take realistic error sources into account.
2 Theoretical Model
1 Kitaev Model
In this model, qubits are located on the edges of a twodimensional lattice as shown in Fig. 1(a). The Hamiltonian of the system consists of two different types of four neighbouringqubit interactions, and ( and are Pauli matrices) at a vertex and plaquette , respectively. Hence,
(1) 
where and , which are referred to as stabilizer operators. Here, star () is a set of four spins that share a link with the vertex , and bond () is a set of four spins placed at the edges of the plaquette . Since all the stabilizer operators commute with each other, the ground state of this Hamiltonian is a +1 eigenstate of the and operators (note the minus sign in the Hamiltonian in Eq. 1): and . The case with a periodic boundary condition on the lattice exemplifies a toric code Kitaev (2003) where the ground states are fourfold degenerate. The degenerate ground states form a protected subspace from possible noiseinduced excited states.
Excited states of this Hamiltonian are created by applying singlequbit operators and/or to the ground state. These operations create two types of quasiparticles, particle at a vertex or particles on a plaquette, as described in Fig. 1(b). Subsequently, as shown in Fig. 2(a) and (b), one can move a particle to a different plaquette by applying to a qubit nearby. Particles created on the same site annihilate each other so the X operation effectively moves the particle. Analogously, applying to a relevant nearby qubit moves an particle.
One can demonstrate anyonic statistics between and particles by moving one around the other, making a closed loop as shown in Fig. 2(c). This braiding operation is equivalent to the two successive particle exchanges. Note that it is not possible to exchange their positions once, since one is located at a vertex and the other at a plaquette. It can be shown that the wave function acquires a 1 phase factor (corresponding to a phase) after such braiding, indicating that a single exchange of Abelian anyon and particles would result in a phase. Therefore, the anyonic statistics of and particles is .
2 The 7qubit Kitaev Model
The 7qubit model used to demonstrate the path independent property of anyonic braiding is shown in Fig. 3(a). For the case of a periodic lattice, the Hamiltonian consists of the four body interactions (Eq. 1). However, for the 7qubit model, we consider a lattice with a rough boudary, which results in twobody interactions at the boundary. Therefore, the Hamiltonian of this system is
(2) 
where
Moreover, due to the absence of periodic boundary conditions, the ground state of this model is nondegenerate. Since the state is already a +1 eigenstate of , the ground state is given by projecting on to the +1 eigenstate of and :
(3)  
Due to the lattice structure, exciting qubit 1 with the operator creates a single particle at rather than creating a pair, whereas exciting qubit 5 with the operator still creates a pair of particles at the plaquettes associated with and . Refer to Fig. 3(b) for the particle locations. Starting with this excited state, there are three possible loops to braid particles as shown in Fig. 3(b,c): the trivial braiding operation where a particle braids around without an particle, and the two nontrivial braiding operations and where a particle braids around with an particle. The wave function remains the same when the braiding operation is trivial; however, if the operation is nontrivial, the wave function picks up a phase from the fractional statistics.
To experimentally demonstrate the path independence of anyonic braiding, we simulated anyonic physics manifested in the 7qubit KM using a liquidstate NMR quantum simulator. Since it is experimentally challenging to engineer the KM Hamiltonian which involves fourbody interactions, we took the state preparation approach: dynamically preparing the ground and excited states of the KM Hamiltonian in a NMR system, instead of generating the KM Hamiltonian and cooling the system.
The quantum circuit which simulates the anyonic physics is shown on the right in Fig. 3. The main idea is to prepare and then create a superposition , where is a state with the pair of particles, and is a state with both the particle and pair of particles. If braided along the nontrivial paths such as or in which circulates around , gains a phase due to the fractional statistics; otherwise, remains unchanged. By measuring the variation of the relative phase on before and after the braiding, one can deduce whether the braiding path is trivial or not and ,furthermore, demonstrate the path independence. The details are described as follows.
First, two Hadamard and six controlledNOT (CNOT) gates are applied to prepare the ground state of the 7qubit KM Hamiltonian from , as depicted in Fig. 3. Then, applying and on generates
(4) 
respectively. To create a superposition of the two, we apply since . When the anyons are braided along a nontrivial loop, the superposition picks up a relative phase on the component. Finally, anyons are annihilated by reversing the creation operator in order to measure this relative phase. The system ultimately evolves to either the ground state or the excited state depending on different braiding paths. Therefore, we can experimentally demonstrate the path independence nature of anyonic braiding if the two phases obtained under the two different nontrivial loops and are the same.
The states corresponding to each step of the circuit shown on the right in Fig. 3 are
(5)  
(6)  
(7)  
(8) 
where is the phase gained from the anyonic statistics for different loops . When the moves around the trivial loop , the final state ends up at the ground state (), whereas when the is moved around the nontrivial loops or , ends up at the excited state (). In order to demonstrate path independence of anyonic braiding experimentally, we need to implement the entire circuit and observe for different loops.
3 Experimental Implementation in NMR
Our 7qubit NMR processor is the perClabeled dichlorocyclobutanone derivative Johnson et al. (2008); Lu et al. (2015) dissolved in d6acetone. The molecule consists of seven C spins and the five H spins. We denoted the seven nuclear spins of C as qubits, while H nuclei were decoupled throughout all experiments except for the initialization step to boost polarization on C. The molecular structure is depicted in Fig. 4(a), where two nearestneighbouring Cs have stronger coupling strengths, implying the ability to implement a faster twoqubit gate. Therefore, by comparing the geometry of KM qubits and the structure of nuclear spins, we mapped each KM qubit to the nuclear spin in as shown Fig. 4(b). The natural Hamiltonian of this system is described as
(9) 
where is the chemical shift frequency of the spin, and is the coupling strength between the and spins (refer to Appendix A for values of the parameters). All experiments were conducted on a Bruker DRX 700 MHz spectrometer at room temperature. The experiment was divided into five main steps as shown in Fig. 4 (b), as follows:
PPS initialization. We first utilized the catstate method proposed in Knill et al. (2000) to initialize the system to a labeled pseudopure state (PPS) state. It can be represented by a deviation matrix of the form , where C is the labeling qubit. Two techniques were adopted before this initialization step to improve the signaltonoise (SNR) ratio. One is turning on C and H couplings temporarily at the very beginning and applying a SWAP gate between C and H, to achieve a 4 times higher polarization on C. The other one is performing RFselection Ryan et al. (2009) sequence to pick out a slice of the NMR sample which experiences much better radiofrequency (RF) homogeneity by randomizing the other part with worse RF homogeneity. Subsequently, the labeled PPS state was prepared using nonunitary transformations via gradient fields and phase cycles Knill et al. (2000). The total length of the initialization sequence is about 100 ms. Refer to Appendix A for more details about the PPS initialization step.
Ground state preparation. Unlike the theoretical circuit on the right in Fig. 3, the implemented circuit in NMR prepared the ground state from , rather than the required pure state . Nevertheless, since contains half of and half of , we can simply write the deviation matrix after the ground state preparation as:
(10) 
where results from . Under perfect unitary transformation, stays orthogonal to throughout the circuit, thus not interfering with the final result if it can be separated in the NMR spectra. In fact, the additional two CNOT gates in the beginning of the ground state preparation shown in Fig. 4(b) were specifically added to achieve this separation. However, in the presence of errors, did slightly modify with the final result, as analyzed in Sec. 4.
The entire ground state preparation step was optimized by a 60 ms GRadient Ascent Pulse Engineering (GRAPE) pulse Khaneja et al. (2005) based on a subspace approach Ryan et al. (2008). The simulated fidelity of this pulse is over 0.99. Additionally, a special rectification method was used in the experiment to ensure that all of the GRAPE pulses acting on the spins were very close to theoretical expectations Weinstein et al. (2004); Moussa et al. (2012). We performed modified stabilizer measurements after the ground state preparation step to verify the state. This step is explained in detail in Appendix Appendix C: Modified Stabilizer Measurements of the labeled PPS and Ground States.
Anyon creation, braiding and annihilation. These three parts shown in the emulation circuit (on the right of Fig. 3) are compressed together to simplify the circuit as they only involve singlequbit rotations. The trivial loop and nontrivial loops , are all depicted in Fig. 4(b), and in each experiment only one loop was implemented. The three braiding operators were realized by 1 ms GRAPE pulses, respectively. In principle, after this stage we can determine by measuring coefficients of and in Eq. 8, but it does require many measurements in a 7qubit system.
Measurement. This additional ‘measurement’ step is added to estimate with a few measurements, which allows us to measure diagonal elements of the final state and then extract the value of . It separates diagonal elements of and via basis transformation by evolving the state to
(11)  
Therefore, considering , the final density matrix is
(12)  
with
(13)  
(14) 
The coefficients , and , and originates from the neglected part . In this case,
(15) 
To evaluate , we estimated by measuring the diagonal elements of and similarly by measuring the diagonal elements of .
Diagonal elements readout. Since the diagonal elements cannot be directly observed in NMR, we indirectly measured them by applying the readout pulse which rotates C by /2 around the axis. This readout pulse generated single coherences from the diagonal elements, and thus a detectable signal with distinct frequencies depending on the state of the other qubits (see Appendix for detailed descriptions). In particular, the transitions relevant to and estimations are at four distinct frequencies centered around (resonant frequency of C): 61.25Hz, 24.09Hz, 32.24Hz, and 4.93Hz. Therefore, the real coefficients of the peaks at these specified frequencies can indirectly estimate the diagonal elements of interest.
There is one assumption in the measurement of diagonal elements in the above method. The peaks are actually generated by the subtraction of two relevant diagonal elements after rotating C by /2 around the axis (see Appendix). For example, the intensity of the peak at 61.25Hz corresponds to , but we only need the value of the first term. So we assume that the latter term is 0 in order to get the value of the first term. We simulated the contributions from such small elements and found that this assumption should be good enough for the accurate estimation of .
theory  experiment  theory  experiment  theory  experiment  

No BD  1  0.830.01  0  0.010.01  0  (12.19.5) 
BD0  1  0.830.01  0  0.020.01  0  (17.46.0) 
BD1  0  0.050.01  1  0.850.01  (180)  (153.93.8) 
BD2  0  0.050.01  1  0.810.02  (180)  (151.43.8) 
4 Result and Discussion
We measured the anyonic phases of the four different cases:

BD0: PPS GD BD0 MM Readout

BD1: PPS GD BD1 MM Readout

BD2: PPS GD BD2 MM Readout

noBD: PPS GD MM Readout
where noBD and BD0 ideally have , and BD1 and BD2 have . GD and MM refer to the ground state preparation and measurement steps, respectively. Fig. 5 shows the C spectra of the labeled PPS and the above four cases. The experimental spectra agree qualitatively with our theoretical predictions. First, in theory, we expect to observe the same spectra for the noBD and BD0 cases and the same spectra for the two nontrivial braiding cases (BD1 and BD2) due to the path independent nature. From Fig. 5, it is clear that the spectra of noBD and BD0 match well, and also that BD1 and BD2 match well. Second, our experimental spectra matched well with the simulated spectra. In theory, the spectra resulting from the four cases are expected to show four peaks with equal height (two generated from and the other two from ), which is a quarter of the labeledPPS peak. The spectra shown in Fig. 5 qualitatively illustrate the expected behaviours. Third, recalling Eq. 15, we expect to observe no peaks at Hz and Hz for noBD and BD, resulting in , and no peaks at Hz and Hz for BD1 and BD2, resulting in . It should be noted that the other large peaks located not at the four frequencies in the spectra result from and are neglected in the analysis.
We estimated and by evaluating the intensities of the peaks at the frequencies of and , and frequencies and , respectively. The intensities of peaks at and are averaged to estimate , and the peaks at and are averaged to estimate . To evaluate the numbers, we fitted the spectra with a Lorentzian function of 64 peaks (the maximum number of observable peaks on C) using the leastsquare method. The experimental results of , and are displayed in Table 1 for all noBD, BD0, BD1 and BD2 cases.
The experimental results show that the anyonic phases under the two different nontrivial braiding paths and agree within the errors: (153.93.8) and (151.43.8). These experimental values clearly demonstrate path independence and the phase gained under the nontrivial paths compared to the cases of the trivial and no braiding paths [(17.46.0) and (12.19.5), respectively]. However, the experimental have discrepancies with the theoretical values, which are 0 for the trivial and no braiding paths, and 180 for the two nontrivial paths. For the nontrivial cases, this deviation is mostly attributed to the tiny peaks at and (Fig. 5), which result in , because is highly sensitive to as it is small and in the denominator (Eq. 15). For instance, consider a theoretical case when is 0. In this case, regardless of a value of , is always . Similarly, for the trivial and no braiding cases, the deviation of is mostly caused by the tiny peaks at and (Fig. 5), resulting in rather than the theoretical value of 0.
To investigate how the unwanted small peaks arise, we numerically simulated the NMR circuit starting from the ideal labeled PPS state using 99% fidelity unitaries calculated from the GRAPE pulses in the presence of the decoherence effect. The assumptions that we used to simulate decoherence are shown in Appendix. The results of the simulation indicate that the errors increase the trivial loop phases to 20 whereas the nontrivial loop phases decrease to 160, blurring the difference between the two. It should be noted that most of the phase deviation comes from the decoherence effect; simulating only the gate imperfection from 99% fidelity unitaries results in the nontrivial phases of 177. Now we discuss the different sources of error in detail.
First, the error primarily comes from the decoherence effect, and the ground state and measurement pulses contribute the most in causing the biases in the determination. In particular, the ground state we prepared was the ground state of the 7qubit KM, not a ground state of our physical system. Therefore, the ground state preparation step is susceptible to decoherence, as there is no protection of the ground state by the energy gap in our NMR system.
Second, to a much lesser extent, Eq. 15 no longer accurately determines the anyonic phase in the presence of gate imperfections. Therefore, to estimate the anyonic phase independent of imperfections of ground state and the measurement pulses, a different equation is required. However, it is difficult to find such an equation that is accurate and whose variables can be easily measured. Since the braiding operation is 1 ms, whereas the ground state and measurement pulses are 60 ms, the ground state and measurement pulse imperfections contribute more significantly to the determination. Moreover, we expect that gate imperfections are worse in experiments than in simulation, which could explain the 10 discrepancy between the simulation and experimental values after accounting for the other sources of error.
Third, we also examined the effect of on the determination through numerical simulations. We simulated two scenarios with one started from and the other from the labeled PPS. As mentioned above, the one starting with the labeled PPS results in the nontrivial of 160, whereas the one started from the pure state results in 150. This signifies that the contribution from cannot be neglected completely when both gate imperfections and decoherence effects are present.
5 Conclusions
We have successfully demonstrated path independence of anyonic braiding statistics by braiding two anyons under two different nontrivial paths in a 7qubit NMR quantum simulator. The anyonic phases of the two nontrivial paths and agree within the errors: and for and , respectively. As references, the cases of no braiding and braiding along a trivial path are also implemented. We measured significantly smaller phases for these trivial cases compared to the nontrivial cases, confirming the extra phase acquired by the anyons in the nontrivial cases. The deviation of the anyonic phases from the theoretical value are well accounted for by the inherent errors of decoherence and imperfect gates. These contributions can be mostly attributed to the ground state preparation and measurement steps, as these steps are significantly longer than the braiding step. Other experimental schemes or setups where such a long preparation step can be prevented may be less prone to such errors. Moreover, the measurement step which is used to remarkably reduce the number of experiments in our NMR system may be eliminated in other settings.
As a step towards the realization of topological quantum computing, we do not simulate the manybody interactions in the KM Hamiltonian but alternatively use a state preparation approach to simulating the KM. This method is sufficient to simulate some particular anyonic properties such as the path independent nature shown in this paper; however, realizing faulttolerant topological quantum computation would ultimately require engineering such Hamiltonians with manybody interactions. Fortunately, quantum simulation provides exponential speedup, outperforming classical computers as well as highly controllable systems instead of the natural intractable solid state systems. Hence, quantum simulation is a promising solution for creating and engineering the full KM Hamiltonians Negrevergne et al. (2005); Zhang et al. (2011); Bloch et al. (2012); Franchini et al. (2014) in the near future, and it may shed light on the goal of building a topological quantum computer in a faulttolerant manner.
Acknowledgements
We thank Aharon Brodutch, Jonathan Baugh, Guanru Feng and Hemant Katiyar for helpful discussions and comments, and Anthony P. Krismanich, Ahmad Ghavami, and Gary I. Dmitrienko for synthesizing the NMR sample. This work is supported by Industry Canada, NSERC and CIFAR.
Appendix A: Sample and Initialization
Our NMR quantum processor is a racemic mixture of perC labeled (1S,4S,5S)7,7dichloro6oxo2thiabicyclo[3.2.0]heptane4carboxylic acid and its enantiomer. The unlabeled compound was synthesized previously by us and its structure was established unambiguously by a single crystal Xray diffraction study Johnson et al. (2008). By decoupling the H channel throughout the experiment, this sample can be regarded as a 7qubit quantum processor which involves seven C spins. The and values in Eq. 9 of the natural Hamiltonian, as well as the relaxation time scales T and T, are shown in Fig. 6.
We initialized the thermal equilibrium to the labeled PPS using the NMR circuit shown in Fig. 7(a), where the entire circuit can be divided into five sections ae. The input state of this 12qubit system is the thermal equilibrium state
(16) 
where is the gyromagnetic ratio of the nuclear spins, is the identity matrix, and represents the polarization of the system. Typically, and with a constant factor ignored. As the large identity matrix part does not evolve under unital operators (which is roughly the case in our experiment as the experimental time is far less than T) and it cannot be measured in NMR experiments, we can simply neglect the identity part and rewrite the input state as
(17) 
In the following calculations we only focus on this deviation density matrix assuming that the identity has no influence on the entire experiment.
a. Rotate C to by a 1 ms GRAPE pulse around axis on C channel, and then crush it by a 2 ms gradient pulse. The total length is 3 ms and the state at step is .
b. SWAP the signal of C and H by applying a 8 ms GRAPE pulse. This GRAPE pulse was designed via statetostate approach and hence not a universal SWAP gate. The reason of implementing this SWAP operation is to improve the C signal by four times in principle, which enables a much better signaltonoise ratio (SNR) in experiment. The state at step is .
c. Turn on the Waltz16 decoupling sequence on H channel. It averages out the signals of all H spins and their interactions with the C spins. In quantum information, this step is equivalent to reducing the 12qubit system to 7 qubits which only involve C spins. Hence, the state at step is . Compared to the input thermal equilibrium state of , the signal of C has been boosted by four times.
d. RFselection technique is used to pick out a subsample which has much better RF homogeneity. As the sample in NMR has some volume in centimeters, the RF pulse applied to the sample may have inhomogeneity. Some molecules located in the centre of the RF coil experience the ideal RF amplitude, while majority of molecules experience overrotation or lessrotation for the sake of RF amplitude inhomogeneity along the sample size. Since NMR readout is an ensemble average, the large portion with bad homogeneity contributes a lot to the final signal and causes accumulated errors when multiple pulses are implemented. RFselection sequence Ryan et al. (2009) is such a technique to randomize this inhomogenous portion to plane while keeping the homogenous portion in the thermal equilibrium state, followed by a gradient pulse in direction to destroy all plane signals. It is usually applied before the primary circuit, and the inhomogenous portion will stay at nosignal case during the following pulse sequence. A typical RFselection sequence with 64 loops is
(18) 
where . When the molecules feel perfect RF amplitude, their states remain as thermal equilibrium after this sequence. By contrast, when the molecules feel for example 4.5% error in RF amplitude, their states mostly evolve to plane and thus be killed by the following gradient field. Note that although RFselection enables a better SNR as the RF pulses are much more precise, the cost of this technique is the absolute loss of signal as many molecules have no contributions to the signal any longer.
In our experiment, we used a GRAPE pulse instead of the long sequence to realize this RFselection technique. This GRAPE pulse was designed on a singlequbit system via the statetostate approach, by setting two constraints: evolve to plane when the RF inhomogeneity is more than 1%, or else do nothing to . After applying this GRAPE pulse on our 7qubit system, we found the signal reduced to about 30% but the RF pulses were indeed much more homogeneous by running the Rabi oscillation experiment. The two gradients and rotations in step are used to kill the minor signal of multicoherence generated by the Jcoupling evolution during the RFselection sequence. The state at step is the same as step , but with some loss that . For convenience, we simply mark this state as . Compared to the original thermal equilibrium state, this new state gains signal boost from H and owns much better RF homogeneity.
e. The main body of catstate method Knill et al. (2000) is implemented which creates the labeled PPS from . It consists of three steps: encoding, phase cycling, and decoding. The detailed NMR sequence is shown in Fig. 7(b). Starting from , the system evolves to after the encoding step. The phase cycling step contains seven loops, and in each loop the axis of the rotation is chosen as (the rotating angle is always ). The state after the phase cycling is . The decoding step is just the inverse of the encoding part and simplified according to our molecular information. The final state after the decoding step is .
Till now the labeled PPS has been successfully prepared. Regarding the performance of this state see Fig. 5 for its NMR spectrum.
Appendix B: Assumptions used in the Simulation of Decoherence
The list below shows the assumptions we used when numerically simulating the decoherence effects.

The environment is Markovian.

The system and the environment are uncorrelated at t=0.

We only considered the effect of dephasing due to T effect and neglect the effect of amplitude damping, since T is much larger than the circuit time.

The dephasing noise is independent (or uncorrelated) between the qubits. The probability of an error happening on a given qubit does not affect the probability of an error happening on other qubits.

When solving the master equation, we assumed that the dissipator and the total Hamiltonian commute for short times. Therefore, the evolution of of the state was simulated in a sequence of two steps: evolution by and subsequently, dephasing for , where was chosen to match the pulse discretization. The dephasing channel implements exponential decay of offdiagonal elements according to relevant linear combinations of T values of C.
Appendix C: Modified Stabilizer Measurements of the labeled PPS and Ground States
If the pure state is prepared as an initial state, the stabilizer operators () of such a state are . When these stabilizer operators evolve under the ground state preparation circuit shown in Fig. 3, one can reconstruct the stabilizer operators of the ground state of the 7qubit KM. However, since our circuit starts from , the are modified to , , , , and . The expectation values of these operators are +1 as Tr() =1, where is one of the modified stabilizer operators. These operators transform to the following operators under the implemented ground state gate which is shown in Fig. 4(b):
(19)  
Therefore, the experimentally prepared ground state have +1 expectation values of the above transformed operators . We measured the expectation values of of the labeled PPS state and the expectation values of of the ground state. For the measurements, a single readout pulse which rotates C by around axis is sufficient to measure all six operators; whereas five different readout pulses are required (thus, five different measurements) to measure the operators. The readout pulses are composed of the single qubit rotations that transform the product operator components of a density matrix corresponding to the operators to the measurable product operators in C spectra, which are a combination of or and different , where indicates the th C. For instance, the readout pulse required to measure the expectation value of the second (#2 in the above list) operator is which rotates to , and thus produces observable peaks at C spectrum. The experimentally measured expectation values are shown in Fig. 8.
Figure 9 shows the experimental and simulated spectra of C after the ground state preparation which were measured by the five different readout pulses. These spectra were used to estimate the operators.
References
 Khare (2005) A. Khare, Fractional statistics and quantum theory (World Scientific, Singapore, 2005).
 Tsui et al. (1982) D. C. Tsui, H. L. Stormer, and A. C. Gossard, Phys. Rev. Lett. 48, 1559 (1982).
 Moore and Read (1991) G. Moore and N. Read, Nucl. Phys. B 360, 362 (1991).
 Wen (1991) X.G. Wen, Phys. Rev. Lett. 66, 802 (1991).
 Nayak et al. (2008) C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. D. Sarma, Rev. Mod. Phys. 80, 1083 (2008).
 Stern (2010) A. Stern, Nature 464, 187 (2010).
 Ladd et al. (2010) T. D. Ladd, F. Jelezko, R. Laflamme, Y. Nakamura, C. Monroe, and J. L. OBrien, Nature 464, 45 (2010).
 Sarma et al. (2005) S. D. Sarma, M. Freedman, and C. Nayak, Phys. Rev. Lett. 94, 166802 (2005).
 Freedman et al. (2003) M. Freedman, A. Kitaev, M. Larsen, and Z. Wang, Bull. Am. Math. Soc. 40, 31 (2003).
 Raussendorf et al. (2007) R. Raussendorf, J. Harrington, and K. Goyal, New J. Phys. 9, 199 (2007).
 Kitaev (2003) A. Y. Kitaev, Ann. Phys. 303, 2 (2003).
 Han et al. (2007) Y.J. Han, R. Raussendorf, and L.M. Duan, Phys. Rev. Lett. 98, 150404 (2007).
 Aguado et al. (2008) M. Aguado, G. Brennen, F. Verstraete, and J. I. Cirac, Phys. Rev. Lett. 101, 260501 (2008).
 Pachos et al. (2009) J. Pachos, W. Wieczorek, C. Schmid, N. Kiesel, R. Pohlner, and H. Weinfurter, New J. Phys. 11, 083010 (2009).
 Lu et al. (2009) C.Y. Lu, W.B. Gao, O. Gühne, X.Q. Zhou, Z.B. Chen, and J.W. Pan, Phys. Rev. Lett. 102, 030502 (2009).
 Du et al. (2007) J.F. Du, J. Zhu, M.G. Hu, and J.L. Chen, arXiv:0712.2694 (2007).
 Feng et al. (2013) G. Feng, G. Long, and R. Laflamme, Phys. Rev. A 88, 022305 (2013).
 Feynman (1982) R. P. Feynman, Int. J. Theor. Phys. 21, 467 (1982).
 Lloyd (1996) S. Lloyd, Science , 1073 (1996).
 Georgescu et al. (2014) I. Georgescu, S. Ashhab, and F. Nori, Rev. Mod. Phys. 86, 153 (2014).
 Friedenauer et al. (2008) A. Friedenauer, H. Schmitz, J. T. Glueckert, D. Porras, and T. Schätz, Nat. Phys. 4, 757 (2008).
 Kim et al. (2010) K. Kim, M.S. Chang, S. Korenblit, R. Islam, E. Edwards, J. Freericks, G.D. Lin, L.M. Duan, and C. Monroe, Nature 465, 590 (2010).
 Lanyon et al. (2010) B. P. Lanyon, J. D. Whitfield, G. Gillett, M. E. Goggin, M. P. Almeida, I. Kassal, J. D. Biamonte, M. Mohseni, B. J. Powell, M. Barbieri, et al., Nat. Chem. 2, 106 (2010).
 Du et al. (2010) J. Du, N. Xu, X. Peng, P. Wang, S. Wu, and D. Lu, Phys. Rev. Lett. 104, 030502 (2010).
 Lanyon et al. (2011) B. Lanyon, C. Hempel, D. Nigg, M. Müller, R. Gerritsma, F. Zähringer, P. Schindler, J. Barreiro, M. Rambach, G. Kirchmair, et al., Science 334, 57 (2011).
 Lu et al. (2011) D. Lu, N. Xu, R. Xu, H. Chen, J. Gong, X. Peng, and J. Du, Phys. Rev. Lett. 107, 020501 (2011).
 Islam et al. (2013) R. Islam, C. Senko, W. Campbell, S. Korenblit, J. Smith, A. Lee, E. Edwards, C.C. Wang, J. Freericks, and C. Monroe, Science 340, 583 (2013).
 Li et al. (2014) Z. Li, H. Zhou, C. Ju, H. Chen, W. Zheng, D. Lu, X. Rong, C. Duan, X. Peng, and J. Du, Phys. Rev. Lett. 112, 220501 (2014).
 Knill et al. (2000) E. Knill, R. Laflamme, R. Martinez, and C.H. Tseng, Nature 404, 368 (2000).
 Johnson et al. (2008) J. W. Johnson, D. P. Evanoff, M. E. Savard, G. Lange, T. R. Ramadhar, A. Assoud, N. J. Taylor, and G. I. Dmitrienko, J. Org. Chem. 73, 6970 (2008).
 Lu et al. (2015) D. Lu, H. Li, D.A. Trottier, J. Li, A. Brodutch, A. P. Krismanich, A. Ghavami, G. I. Dmitrienko, G. Long, J. Baugh, et al., Phys. Rev. Lett. 114, 140505 (2015).
 Ryan et al. (2009) C. Ryan, M. Laforest, and R. Laflamme, New J. Phys. 11, 013034 (2009).
 Khaneja et al. (2005) N. Khaneja, T. Reiss, C. Kehlet, T. SchulteHerbrüggen, and S. J. Glaser, J. Magn. Reson. 172, 296 (2005).
 Ryan et al. (2008) C. Ryan, C. Negrevergne, M. Laforest, E. Knill, and R. Laflamme, Phys. Rev. A 78, 012328 (2008).
 Weinstein et al. (2004) Y. S. Weinstein, T. F. Havel, J. Emerson, N. Boulant, M. Saraceno, S. Lloyd, and D. G. Cory, J. Chem. Phys. 121, 6117 (2004).
 Moussa et al. (2012) O. Moussa, M. P. da Silva, C. A. Ryan, and R. Laflamme, Phys. Rev. Lett. 109, 070504 (2012).
 Negrevergne et al. (2005) C. Negrevergne, R. Somma, G. Ortiz, E. Knill, and R. Laflamme, Phys. Rev. A 71, 032344 (2005).
 Zhang et al. (2011) J. Zhang, T.C. Wei, and R. Laflamme, Phys. Rev. Lett. 107, 010501 (2011).
 Bloch et al. (2012) I. Bloch, J. Dalibard, and S. Nascimbène, Nat. Phys. 8, 267 (2012).
 Franchini et al. (2014) F. Franchini, J. Cui, L. Amico, H. Fan, M. Gu, V. Korepin, L. C. Kwek, and V. Vedral, Phys. Rev. X 4, 041028 (2014).
 Laflamme et al. (2002) R. Laflamme, E. Knill, D. G. Cory, E. M. Fortunato, T. Havel, C. Miquel, R. Martinez, C. Negrevergne, G. Ortiz, M. A. Pravia, et al., arXiv:quantph/0207172 (2002).