{'fraction_non_alphanumeric': 0.025438663703930733, 'fraction_numerical': 0.007279135698413914, 'mean_word_length': 4.868480215827338, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantum computation of silicon electronic band structure

Development of quantum architectures during the last decade has inspired hybrid classical-quantum algorithms in physics and quantum chemistry that promise simulations of fermionic systems beyond the capability of modern classical computers, even before the era of quantum computing fully arrives. Strong research efforts have been recently made to obtain minimal depth quantum circuits which could accurately represent chemical systems. Here, we show that unprecedented methods used in quantum chemistry, designed to simulate molecules on quantum processors, can be extended to calculate properties of periodic solids. In particular, we present minimal depth circuits implementing the variational quantum eigensolver algorithm and successfully use it to compute the band structure of silicon on a quantum machine for the first time. We are convinced that the presented quantum experiments performed on cloud-based platforms will stimulate more intense studies towards scalable electronic structure computation of advanced quantum materials.


I. INTRODUCTION
Quantum computing aims to leverage superposition, entanglement and interference of quantum bits in order to tackle computational tasks that scale exponentially on classical computers. 1,2 While renowned quantum algorithms, such as unsorted database search or integer factorization require resources that remain out of reach, 3,4 quantum chemistry calculations are gaining steam as a key application performed on available quantum architectures. 5,6 The idea of so-called quantum simulations originally proposed by Feynman, 7-9 relies on a mapping between the fermionic system and the set of qubits, so that the dynamics of the former is directly followed by the latter. Therefore, wave functions of complex many-body systems could be effectively reproduced in quantum measurements performed on qubits, providing a tool to compute desired quantities with an unprecedented accuracy. Even though available quantum computers contain merely few tens of qubits, 10 they have been employed to solve quantum chemistry problems, such as the estimation of nuclear binding energies or molecular ground states. [11][12][13][14] Remarkably, these successful quantum experiments relied on variational approaches that greatly reduced the required hardware resources, inspiring more active research in order to solve elusive condensed matter systems beyond quantum chemistry. [15][16][17][18][19] Here, we put forward an approach to calculate the electronic structure of the periodic crystal on a quantum computer. While developments of quantum computation for molecules were primarily focused on the ground state energies, to evaluate a band structure one needs to determine the excited states. We have shown that a standard hybrid quantum/classical algorithm, variational quantum eigensolver (VQE) can be easily adapted to provide an accurate estimation of the electronic bands in the solid. In particular, by casting a Si tight-binding (TB) Hamiltonian in terms of fermionic operators, we have designed a low-depth quantum circuit, enough robust to capture the electronic properties of a crystal in the reciprocal space. The quantum measurements have been performed on sets of qubits available remotely via cloudbased platforms provided by IBM and Rigetti Computing. Importantly, we have tested different classical optimization routines that minimize expectation values, corrected beforehand against the readout errors. Comparison between bands computed on the quantum processors, the quantum virtual machine and by classical diagonalization revealed a satisfactory agreement, confirming validity of the algorithm which could be generalized to explore materials more complex than crystalline silicon.

II. HAMILTONIAN REPRESENTATION
Let us consider a silicon lattice in the diamond cubic structure. The Hamiltonian describing the electronic system can be approximated, in atomic units, aŝ where r i (R i ) are the positions of electrons (nuclei) and Z j denotes the nuclear charge, respectively. We have assumed the Born-Oppenheimer approximation and considered the nuclei as stationary charges, thus neglecting their kinetic energy and treating the ion repulsion as a constant. The last term of Eq.1 represents the electronelectron interaction, whose correct estimation is one of the long-term goals of quantum simulation. However, we are now primarily focused on the proof-of-principle band structure calculations, and have disregarded the electronic correlations for the purpose of the present study. In order to convert the Hamiltonian into a computational problem, a suitable basis set needs to be selected. While different representations were proposed for quantum computation, 15 we introduce here a simple basis of atomic orbitals at each lattice site arising from the tightbinding (TB) approximation. The unit cell of silicon contains two tetrahedrally coordinated ions and is well described in terms of s, p x , p y and p z orbitals centered at each atom. Because magnetic order is absent, the spin degrees of freedom can be omitted in the analysis. Using the second quantization formalism, we can express the TB Hamiltonian via creation and annihilation operators (a † in and a in ) acting at the orbital n and the site R i : In this expression, E n correspond to the atomic energies and t in,jm denote the hopping integrals whose numerical values have been reported elsewhere. 20 Only the tunneling between pairs of nearest neighbors, denoted by the < i, j > summation, have been considered. The Hamiltonian can be then easily converted to the momentum space via standard Fourier transform applied to the raising and lowering operators. Last, such a representation (Ĥ k ) needs to be mapped onto the system of qubits.
In practice, qubits are manipulated on a quantum processor by operating on a set of Pauli matrices X, Y , Z and I, the latter denoting 2×2 identity matrix. Any Hermitian matrix can be decomposed using a complete Pauli basis for matrices of dimension N = 2 n with n = log 2 N terms, that can be generated by taking a tensor product: Thus, TB Hamiltonian can be decomposed as follows: where the set {σ} n is the set of 4 n possible basis matrices, and {c k } n is a set of complex coefficients. {c k } n is known as the spectral decomposition and can be determined easily. In particular, we can exploit the orthogonality of Pauli matrices and the trace inner product between two of them: By taking the inner product Tr(Ĥ k †σ i ), we can eliminate all terms but one from the sum, yielding: Therefore, the Hamiltonian is represented by a list of coefficients corresponding to each of the 4 n Pauli basis matrices suitable for simulation on a QPU. 21

III. VARIATIONAL QUANTUM EIGENSOLVER
We have computed the energy spectrum using the variational quantum eigensolver in conjunction with overlapbased techniques. VQE is a standard hybrid quantumclassical algorithm capable to determine the lowest or highest eigenvalue of an operator using minimal quantum resources, implemented by combining measurements on a quantum computer with classical routines. [22][23][24] The ground state wave function and energy can be found based on Rayleigh-Ritz variational principle, whereby the energy expectation value can be minimized by a specific set of parameters. In practice, the state preparation and the expectation value measurements are implemented on a quantum machine, while the optimization of the parameters is performed classically. The whole algorithm used for the ground state calculation can be summarized in three following steps: 1. We create a quantum circuitV (θ) depending on a set of parameters θ, known as a variational form.
Then, we prepare a trial wave function (or ansatz ) |ψ(θ) =V |0 , where |0 denotes an initial state ensuring the measurement of each qubit.
2. We measure the expectation value ofĤ k , which depend on the parameters θ, The Hamiltonian is represented by series of operators. The wave function |ψ is measured in the Pauli basis, yielding each σ i . We can then reconstruct Ĥ k with the spectrum {c k }: The measurement should be treated as a probabilistic element of the algorithm and needs to be performed several times. An arbitrary precision can be achieved with a sufficient number of repetitions.
3. We apply a classical optimization routine to explore the parameter space and minimize E(θ). We define 0 = ψ(θ min )|Ĥ k |ψ(θ min ) as a ground state energy, where θ min denotes the set of parameters minimizing the expectation value ofĤ k .

IV. ENERGIES BEYOND THE GROUND STATE
After having determined the ground state, we can calculate excited states using a procedure similar to the quantum deflation algorithm that exploits orthogonality of the Hamiltonian eigenvectors. [25][26][27] In particular, we define an effective Hamiltonian (Ĥ ) whose lowest eigenstate is the excited state of the original one (Ĥ). By subtracting from the latter a corresponding ground state projector weighted by the ground state energy, we obtain: We observe that the last equality provides the following spectral decomposition of the excited Hamiltonian: The procedure is used iteratively to determine as many eigenvalues as desired. Updating the spectral decomposition c i → c i − 0 σi 2 n effectively removes all ground state contributions from the Hamiltonian.
We note that the effect of subtracting the ground state density matrix weighted by its corresponding eigenvalue is to project that eigenstate onto the zero value. Because an arbitrary Hermitian matrix can have both positive and negative eigenvalues, special care must be taken to ensure that the zero is not erroneously computed as a ground state after all negative eigenvalues are determined. One reconciliation is to subtract a value greater than the maximum eigenvalue from the diagonal elements of the Hamiltonian, ensuring that all eigenvalues are lower than zero. Therefore, projecting an eigenstate to zero would not affect the remaining eigenvalues that need to be determined. Such a shift requires the modification of only one coefficient of the spectral decomposition, which stands before the identity matrix.

V. DATA ACQUISITION
Before discussing the results of quantum experiments, let us remark on the various techniques that we have employed to compute the band structure of silicon. A careful distinction must be made between the use of quantum processor, quantum virtual machine and quantum state simulation. In particular, simulated qubits helped us analyze the performance of variational forms and the effect of measurement uncertainty on a noiseless machine. Three independent techniques will be further referenced: 1. Quantum Processor Unit (QPU) is prepared for measurements under subsequent sets of parameters. The measurements are performed in real time.
The available APIs compile quantum programs and directly manipulate qubits, providing measured expectation values in the form of bitstrings.
2. Quantum Virtual Machine (QVM) chooses one of the possible outcomes to be "measured", weighted by its respective probability computed with the quantum state simulator (see below). The quantum processor is mimicked, providing a noiseless (unless noise is simulated) simulation of the measurement process. This method helps to analyze the effects in the band structure determined by discrete measurements of the energy expectations values.
3. Quantum State Simulator (QSS) carries out linear algebra to obtain an exact wave function which would represent the simulated state of a qubit on a quantum processor after the application of specified gates. It can serve as an analytical guideline for quantum measurements. Optimization can be easily performed with the quantum state simulator, providing a convenient framework to test the performance of variational forms.

VI. QUANTUM EXPERIMENTS
Quantum computations of the band structure have been performed following two different techniques, both yielding a correct spectrum while compared with the classical diagonalization of the TB Hamiltonian. The first approach relies on a true quantum measurement, employing one qubit that we access on remote quantum machines Rigetti Aspen and IBMQ Armonk. Although these cloud platforms permit the use of larger resources, the practical realization of the VQE algorithm for diagonalization of the 8 × 8 Hamiltonian of Si required a substantial amount of time. Therefore, we have started with a reduced Hamiltonian, considering only the interactions between s-states which give rise to the lowest bands of silicon. After neglecting s − p hopping parameters in the originalĤ k , a smaller 2×2 matrix block can be decoupled and diagonalized using VQE on the QPU. Figure 1 shows the two-gate circuit acting on a single qubit, often referred to as the mean field ansatz, 28 which has been used in the experiment. In principle, to ensure that finding the true minimum is possible, circuits must be designed to span every state allowed by the operating qubits, unless the space is restricted by physical arguments, such as fermionic commutation relations in the UCC strategies. 29 The ansatz below takes a pure state |0 and applies two rotations described by the angles θ = (θ, φ). A polar rotation brings the qubit into a superposition of |0 and |1 states, while an azimuthal rotation scans the sphere's latitude. The two rotations produce a state represented by the following wave function: FIG. 1. Mean-field circuit acting on a single qubit has been employed to determine the lowest bands of silicon. It consists of a polar rotation (Ry) followed by an azimuthal rotation (Rz). In last step, the expectation value ofĤ k is measured. The band structure has been computed along a highsymmetry line X − Γ − L by repeating the whole algorithm for each of the k-points. Figures 2(a-b) report the two-band electronic structure evaluated on the quantum machines of IBM (red squares) and Rigetti (green circles), complemented by the data from the classical diagonalization (black solid line). In addition, we present the results obtained via quantum-classical algorithm performed on QSS (blue squares) and QVM (yellow circles). While the latter directly follow the bands calculated classically, the quantum data reveal tiny deviations that can be noticed around the high-symmetry points Γ and L for Rigetti and IBM, respectively. The sources of errors in the experiment can be manifold. The probabilistic aspect can obviously play a role, despite a large number of measurements (8192) taken for each parameterization. Importantly, simulation of noise on QVM have revealed that any gate noise or readout error tends to increase the measured energy, shifting the expectation value toward different eigenstates. As described in the next sections, we have attempted to characterize and reduce the effects of errors arising from the qubit manipulation.
We note that the standard optimization routines have not been here applied. Instead, we have used the meanfield circuit to measure a dense grid of parameter angles in order to find the minimum expectation value. Sampling the entire parameter space provides a visual tool for analyzing the structure of parameter space. Figure 2(c) shows examples of the expectation value surfaces computed for one selected point k = π 4a 1, 1, 1 . The three subsequent panels report the surfaces obtained analyt-ically on QSS (blue) and experimentally on IBM (red) and Rigetti (green). The two latter have been smoothed by minimizing the root-mean-square error across all data points. Again, the data collected on IBM reveals largest irregularities in the energy contour lines, especially compared with the analytical surface evaluated on QSS.
The second approach, employed to diagonalize full 8×8 Hamiltonian, relies on QSS. Figure 3 presents a robust three-qubit circuit that we have designed to variationally minimize the expectation value ofĤ k at any k-point and each level of excitation. The set of twelve parameters θ = (θ 1 , θ 2 , ..., θ 12 ) in this ansatz, measured in the Pauli word basis from the Hamiltonian decomposition defined in Eq.(4), are varied to minimize the energy expectation values. Figure (4) displays the electronic structure computed using this circuit, demonstrating that it is indeed capable of representing the silicon Hamiltonian anywhere along the k-line. Although small discrepancies are again visible, the overall agreement with the bands calculated classically seems to be sufficient. We note that now the results do not depend on external factors that can perturb the behavior of qubits. The deviations are related to the optimization procedures whose proper choice is essential to correctly determine the energy spectrum.
Several classical optimization routines have been tested in conjunction with the three-qubit circuit used for the evaluation of full electronic structure. Minimizing a function in parameter space of twelve dimensions is rather challenging and requires a compromise between the number of measurements and the smoothness of the space being optimized. We have found that the Broyden-Fletcher- Goldfarb-Shanno (BFGS) and Constrained Optimization BY Linear Approximation (COBYLA) routines 30 yielded the most accurate results. The former requires fewer function evaluations to reach a minimum, but it suffers from instability due to the rough surface in parameter space. The latter, being a direct search method, entirely omits the idea of gradient decent which makes it more robust against becoming trapped in a local minimum. Even though it may provide more reliable global minima, 31 it occasionally fails to settle on the correct set of parameters. Figure 4 clearly shows that especially the excited energy levels are sensitive to fluctuations in the determined parameters. The comparison of both routines, BFGS and COBYLA, eventually indicates the superior performance of the former, at least in the present case.

VII. ADDITIONAL REMARKS ON MEASURING EXPECTATION VALUES
While the previous section was entirely focused on the realization and results of quantum experiments, the measurements of expectation values need a more detailed dis-cussion. The quantities we have measured on the quantum computer are the expectation values σ k , where the operator σ k is an n-length Pauli word consisting of an I, X, Y , or Z for each qubit. They depend on the state |ψ of the qubits, and could be written as the integral ψ|σ k |ψ . Because we do not know |ψ , we must measure the state of each qubit in the computational basis, resulting in a single bitstring (eg. |00101 ). Repeating the measurement a large number of times M , we construct the expectation value σ k from the ensemble of bitstrings. In the following paragraphs, we will first consider the simple single-qubit case σ k = Z, then the multiqubit case where σ k consists only of I and Z operators and last, the general case including X and Y operators.
The Pauli operator Z can be written in a matrix form: It is a diagonal matrix with eigenvalues +1, corresponding to the state |0 , and −1, corresponding to the state |1 . The expectation value Z is the average of these two eigenvalues, weighted by the number of measurements in each state. If p is the probability that we measure |0 rather than |1 , the expectation value Z is given by: Now, consider an operatorÂ defined as a Kronecker product of I and Z operators, each acting on their own qubit. It is a degenerate operator with half the eigenvalues +1 and half −1. Because its matrix form is diagonal, each bitstring we measure corresponds exactly to an eigenstate. The parity (±1) of a given bitstring z is precisely the parity of the substring z which omits any index corresponding to an I operator in A. For example, if A = I 5 Z 4 Z 3 I 2 Z 1 and z = |00101 , the substring z leaves off the second and fifth indices: z = |011 . This string has a weight of two, which is an even parity and therefore corresponds to the eigenvalue +1. The expectation value A is once again an average of +1 and −1, weighted by the frequency of bitstrings corresponding to each of the two states.
Last, let us consider a general Pauli word σ k . Half its eigenvalues are again +1 and half −1, but bitstrings in the computational basis do not correspond exactly to the eigenstates. We therefore need to diagonalize σ k . Let A k be the Pauli word which replaces all X and Y in σ k by Z, and the operator U k changes the basis so that σ k = U † k A k U k . Then, for each expectation value we have σ k = U † k A k U k . This is equivalent to measuring the expectation value A k in a new state |ψ = U k |ψ . Thus, we may apply at the end of the variational circuit the sequence of gates representing U k , and then apply the methods of the previous paragraph to evaluate σ k . One example of U k could be an operator applying the Hadamard gate H to each qubit corresponding to an X operator in σ k , and the sequence of gates HSZ to each qubit corresponding to a Y operator.

VIII. ERROR ANALYSIS AND MITIGATION
Quantum error correction, or more often error mitigation is essential for a reliable attainment of computations on a real QPU. [32][33][34] The quantum measurement, an integral element of any algorithm, is by itself probabilistic. In particular, expectation values of an operator are estimated over a large number (M) of discrete measurements. On a noiseless quantum computer, the variance in the expectation value of the Hamiltonian is limited by where E 2 is the average of the squared energy. It defines an uncertainty and can be resolved to an arbitrary level of precision by increasing the number of measurements. Importantly, the qubits may accumulate errors either due to the imprecise manipulation or interactions with environment. One of the major sources of errors that we have identified while collecting the data from the quantum processors is the readout error, emerging due to a certain probability that a qubit in a true |0 state is measured as a |1 or vice versa. Repeated measurements of prepared |0 or |1 states reveal transition rates w 01 and w 10 , defined as the probability that |0 is erroneously measured as |1 or |1 is measured as |0 , respectively. Moreover, the application of a particular circuit element may result in an imperfect transformation of the qubit state. The so-called gate noise is typically classified as a separate source of error but for the purpose of this study we have assumed it to be intrinsic to the readout error.
The procedure of error mitigation is based on the computation of the transition rates w 01 and w 10 and deriving an appropriate expression to correct the measured expectation values. In order to estimate these rates, we have explicitly prepared the state |0 (|1 ) 100,000 times and counted how many |1 s (|0 s) were measured, which determines the probability that a bit flip occurs on a readout for a given computational state of each qubit. The transition rates need to be measured and updated often to ensure that the correction scheme remains effective across the duration of the trials. In fact, they are calculated every time before the optimization step is reached to take into account changes in behavior of a specific qubit. Figure 6 reports the transition rates w 01 and w 10 evaluated for each qubit while computing band energies. The transition rates are sampled once per minute across the duration of a 50 minute run. The rates corresponding to a flip from |1 to |0 seem to oscillate with a period of roughly 18 minutes, suggesting that environmental effects indeed modulate the behavior of qubits.
FIG. 5. Transition rates estimated for a qubit on Rigetti's QPU. Blue circles denote the rates from state |1 to state |0 , while red circles correspond to the rates from state |0 to state |1 . The fitted trend in transitions suffering from less noise is marked with the gray line. We believe these transitions to arise due the environmental coupling.
The measured expectation value, on a single qubit, can be corrected using the following expression, derived in the Supplementary Material (SM): with p ± defined in terms of the transition probabilities for the single qubit, p ± = w 10 ± w 01 . The procedure can be easily generalized to any number of qubits measured in the computational basis, 35 as follows: where z i is the i th element of bitstring z, and z is among the set of bitstrings of length n (Z n 2 ). The fraction of measured bitstrings resulting in z is denoted as p(z). The correction have been successfully applied to the quantum computation of two-band electronic structure performed on Rigetti. Figure 2(b) shows a comparison between the corrected and uncorrected data points (closed and open circles, respectively), demonstrating that the errors have been significantly reduced.

IX. SUMMARY AND PERSPECTIVES
In summary, we have computed the band structure of silicon along high symmetry lines in the momentum space using quantum machine accessible via cloud. In order to perform quantum simulations beyond the tractability of modern supercomputers, we need to establish methods of translating a desired physical system to the language of qubits founded with quantum logic gates. The VQE algorithm adapted from quantum chemistry seems to be suitable for electronic structure computation and remarkably, is able to leverage even minimal quantum resources, as demonstrated by the results discussed in this work. In analogy to early quantum chemistry computation tackling the problems with known analytical solutions, we have selected the electronic structure of silicon which is considered trivial in materials science. The presented studies can be thus regarded as a first step towards scalable electronic structure quantum computation that would not be limited to a specific interaction or one particular quantum system. Even though the analyzed Hamiltonian was quite simple, we are convinced that adding interactions, field effects, or corrective terms will be possible in the nearest future.