Multi-Dimensional Quantum Capacitance of the Two-Site Hubbard Model: The Role of Tunable Interdot Tunneling

Few-electron states confined in quantum-dot arrays are key objects in quantum computing. The discrimination between these states is essential for the readout of a (multi-)qubit state, and can be achieved through a measurement of the quantum capacitance within the gate-reflectometry approach. For a system controlled by several gates, the dependence of the measured capacitance on the direction of the oscillations in the voltage space is captured by the quantum capacitance matrix. Herein, we apply this tool to study a double quantum dot coupled to three gates, which enable the tuning of both the bias and the tunneling between the two dots. Analytical solutions for the two-electron case are derived within a Hubbard model, showing the overall dependence of the quantum capacitance matrix on the applied gate voltages. In particular, we investigate the role of the tunneling gate and reveal the possibility of exploiting interdot coherences in addition to charge displacements between the dots. Our results can be directly applied to double-dot experimental setups, and pave the way for further applications to larger arrays of quantum dots.


Introduction
Semiconductor quantum dots (QDs) populated by charge carriers (electron or holes) are among the most studied systems for quantum computing [1][2][3]. A set of QDs, such as a QD array [4], is coupled to a classical circuit through a set of voltage gates [5], by means of which the experimenter can manipulate and measure the quantum state of the confined particles. In particular, the qubit readout can be efficiently performed within the gate reflectometry scheme, exploiting the capacitive coupling between the dots and the neighboring metal gates. The readout is performed by applying small voltage oscillations to one or more gates, and by estimating-from the reflected signal-the state-dependent change in the impedance of the resonant circuit [6][7][8][9][10]. When the frequency of the oscillation is slower than the charge carriers' response to the varying potential (adiabatic regime), the impedance contributed by the QD(s) has a capacitive character [11]. In order to readout with high fidelity the state of particles in the dot, it is thus crucial to identify the regions in parameter space where such quantum capacitance (QC) presents large and state-dependent values. Typically, the largest values of the QC are achieved at the working points corresponding to transitions between different charge-stability regions, which correspond to voltage domains where the occupations of the QDs are relatively stable. At the boundaries between two such regions, the interplay between quantum tunneling on the one hand and QD potentials and Coulomb interactions on the other induces significant fluctuations in the QD occupations [12][13][14][15]. This can and has been exploited to implement the dispersive readout of spin qubits [16][17][18][19][20][21][22][23] and, more generally, for quantum-state discrimination [24].
In the adiabatic regime and in the presence of a single gate voltage V, the QC is related to the second derivative with respect to V of the relevant energy eigenvalue [25][26][27]. When M > 1 gates are coupled to the system, the energy depends on M distinct voltage variables. This opens novel possibilities for the gate-reflectometry readout. The measured quantity, numerically computed, and the general trends are discussed; in Section 5, we summarize our findings and mention possible future developments.

Model Quantum Capacitance Matrix
The QCM for a system in its k-th eigenstate is a matrix C k;i,j , whose elements are defined by the second-order derivatives of the energy eigenvalue E k with respect to the gate voltages V i [28].
Any measurable QC is a linear combination of the elements of the QCM, and the coefficients of such a combination are determined by the direction (chosen by the experimenter) of the voltage oscillation in the M-dimensional voltage space. In other words, the measured QC can always be expressed as a weighted average of the eigenvalues of the QCM. At a given working point V = (V 1 , . . . , V M ), the QC eigenvalue with the largest modulus represents the maximum measurable value of the QC, and the direction of the voltage oscillation that should be applied in order to obtain such a value is given by the corresponding QC eigenvector [28].
From a theoretical point of view, the eigenenergies E k and the corresponding eigenstates are expressed in terms of the parameters entering the Hamiltonian, which defines the adopted model. This suggests the possibility of defining a model quantum capacitance matrix (MQCM), whose elements are given by the mixed second derivatives of the relevant energy eigenvalue with respect to the M par model parameters p i : (1) The model parameters corresponding to a given physical system are controlled by the applied gate voltages V j . For a known functional dependence of the p i on the V j , one can express the QCM as a function of the model quantum capacitance matrix through the equation The matrix α is, in general, rectangular if the number of tunable parameters M par does not coincide with that of gate voltages M. The sector of α, which relates the plunger gates to the dots' chemical potentials is the lever-arm matrix [29]. In the following, we consider a case where α is a square matrix, i.e., M par = M. Additionally, we assume for simplicity that the matrix is diagonal, i.e., that each model parameter depends on a different gate voltage. A similar situation is experimentally implemented by identifying the virtual gates: these are suitable combinations of the gate voltages that modify only one model parameter, while leaving the others unaffected [5,29]. In general, one can derive the α matrix for each device either experimentally or by carrying out a multiscale simulation of the device [30].

Model Hamiltonian and Its Diagonalization
The considered system is represented by a DQD occupied by N = 2 electrons and controlled by three gate voltages (V 1 , V 2 , and V 3 ), which control the shape of the confinement potential. In a simple picture, these three gates modulate the energy levels of the two QDs and the height of the barrier between them. More specifically, the two QDs are represented as the sites of a Hubbard model, with interdot tunneling amplitude T (assumed to be real for the sake of simplicity) and onsite energies 1 , 2 . The Hamiltonian is given bŷ where U is the onsite interaction energy, andn i,σ =ĉ † i,σĉ i,σ is the number operator corresponding to the i-th dot and spin σ, whilen i = ∑ σni,σ .
With reference to the terminology introduced in Section 2, here, the tunable model parameters are 1 , 2 , and T; equivalently, we will use ≡ 2 + 1 and δ ≡ 2 − 1 in place of the couple ( 1 , 2 ). These parameters are modulated by the physical gates V 1 , V 2 , and V 3 applied by the experimenter; instead, we assume that U is unaffected. Therefore, The Hamiltonian (3) can be analytically diagonalized. In the case of N = 2, there are six eigenstates comprising three singlets and a triplet. It is a simple matter to show that the triplet energy is E tr ( , δ, T) = . The derivation of the three singlet eigenvalues is doable but more complicated, because it involves the diagonalization of a non-trivial 3 × 3 matrix. Using the resolving formula for a third-degree polynomial equation and some algebraic manipulations, we obtain where The subscripts in Equation (4) refer to the well-known values of the corresponding eigenvalues in the absence of on-site energies ( = δ = 0): where J is the exchange energy. To our knowledge, the general solutions in the presence of arbitrary on-site energies (Equation (4)) and the related quantities that we will derive in the following are novel results that can be directly applied to any DQD which can be represented by a Hubbard model. At all values of the model parameters, the singlet with energy E −J is the ground state, while the first excited multiplet is the triplet. The eigenvalues (reckoned from , in units of |T|) are plotted in Figure 1a, as functions of δ/|T| and for U = 10|T|. Relevant to the following discussion, the avoided crossings between the two lowest singlets (red and blue curves) at δ = ±U correspond to the classical charge transitions (0, 2) → (1, 1) and (1, 1) → (2, 0), respectively, at δ = −U and δ = U. The particular U/|T| ratio used in Figure 1 was chosen for visualization purposes, but it can be significantly different, even by orders of magnitude, as discussed in the following sections. In particular, T strongly depends on the interdot distance and on the height of the tunneling barrier, which, in the present case, can be modulated by a gate voltage.

Model Quantum Capacitance Matrix
In terms of the model parameters ( , δ, T), the MQCM for state k has three independent nonzero elements, namely In fact, since the eigenenergies are linear in , all second derivatives involving vanish. The analytical expressions for the quantities (7) related to the ground state (k = −J) are given by where while B p 1 ,p 2 is a sixth-degree polynomial and C p 1 ,p 2 is an eighth-degree polynomial in the energies δ, T and U. Namely, In Figure 1b, we plot the three inequivalent elements of the model quantum capacitance matrix (for U = 10|T|). Interestingly, where the modulus of the diagonal element related to the bias (purple curve) presents the two maxima, the modulus of the element related to tunneling (orange curve) displays two pronounced minima. The different physical origins of the two terms, that correspond to the quantum capacitance measured by applying the oscillating voltage either to the bias or to the tunneling voltage, is thus reflected in the qualitatively different behaviors. The presence of a significant off-diagonal element (brown curve) further suggests that a maximal quantum capacitance might result from the simultaneous oscillations of the two voltages (see below).

Quantum Capacitance Matrix and Its Diagonalization
In order to compute the QCM, we now need to fix a relation between the control and model parameters. We assume that the physical gates were virtualized in such a way that each of them controls only one model parameter, according to where α 1 , α 2 , and β are positive. The exponential dependence of the tunneling on the corresponding gate voltage is justified by the fact that |T| is exponentially enhanced (suppressed) by lowering (increasing) the barrier height, which can be assumed to be linear in the applied voltage V 3 . The accuracy of an exponential relation between V 3 and |T| has been experimentally demonstrated in Ref. [5]. With our convention, a large positive value of V 3 lowers the barrier and thus favors the interdot tunneling. Under these assumptions, and for identical values of the two lever arms (α 1 = α 2 ≡ α), one obtains the following QCM: where the matrix elements are given by The matrix C k has a null eigenvalue, corresponding to the eigenvector (1, 1, 0)/ √ 2 and related to the symmetric fluctuation of the on-site energies, assumed to linearly depend on the voltages V 1 and V 2 . The remaining eigenvalues are: From the above expressions, it follows that C k;− ≤ C k;+ . Moreover, for the system ground state, X −J < 0 and Z −J < 0 (see Figure 1b). Therefore, one always has that C −J;− ≤ C −J;+ ≤ 0 and that the QC eigenvalue with the largest modulus is C −J;− . At the working points where the mixed partial derivative vanishes (Y k = 0), the QCM is block-diagonal, and the non-zero QC eigenvalues become The normalized eigenvectors of C k corresponding to the non-zero eigenvalues can be written as and v k;

Singlet-Triplet Discrimination
The theoretical knowledge of the QC allows to discriminate between quantum states characterized by different QC matrices. A particular case of interest is the discrimination between ground states with different total spins. In the case at hand, as discussed above, the spin triplet has energy E tr ( , δ, T) = and, therefore, its QCM is a null matrix. It follows that any QC measurement, at any working point, yields zero if the state of the system is a spin triplet. This is not altered by the inclusion of a magnetic field, which can only induce the energy splitting of triplet states, but no mixing between the triplet and the singlets.
Conversely, QC measurements on a singlet state reveal a very different scenario, as we will show in Section 4. The ground-state QCM is strongly dependent on both the working point and the direction of oscillation of the applied potential. The marked difference between the singlet and triplet QC matrices allows to experimentally discriminate between the two cases.

Numerical Results
We now derive from the above equations the main features of the three-gate QCM for the two-electron ground state in the DQD. In order to make our analysis as general as possible, it is convenient to refer to dimensionless quantities. In this system, the natural units of energy, voltages and QCs are |T 0 |, |T 0 |/α, and α 2 /|T 0 |, respectively. In the following, the dimensionless counterpart of a quantity Q is denoted as Q. Thus, for example,

Typical Ranges of Parameters
In Ref. [5], Mills and coworkers investigate Si/SiGe QDs separated by gated barriers; the tunneling amplitude between neighbor dots is demonstrated to be in an exponential relation with the barrier voltage. We refer to that work to deduce the appropriate orders of magnitude for the parameters of our model. Although the functional form that we adopt here for the dependence of T on V 3 is slightly different, it reproduces reasonably well the measured range of |T| in Ref. [5] with |T 0 | ≈ 6 µeV, β ≈ 46 mV, and V 3 ranging from 0 to several tens of mV, corresponding to β ≈ 10 3 and V 3 up to some 10 3 .
The value of U cannot be estimated from Ref. [5], since that work considers transitions between single-particle states. Typical values, derived from previous calculations on Si and Ge DQDs [31], are of the order of some tens of meV. In particular, for an interdot distance of 16 nm, we find that U ≈ 24 meV, corresponding to U ≈ 4000. In the following, this is the value adopted for the onsite repulsion, while (V 2 − V 1 ) is varied in a range that allows to capture the U-dependent charge transitions at δ = ±U.

Directly Measurable QCs: The Diagonal Terms of the QCM
As discussed in Ref. [28], all the elements of the QCM can be directly or indirectly experimentally estimated by suitable procedures. In particular, each diagonal element of the matrix, C −J;i,i , is measured by applying an oscillation voltage only to gate i. Before discussing the QC eigenvalues, we show the dependence of such elements, and specifically C −J;1,1 = C −J;2,2 and C −J;3,3 , on δ and V 3 .
We start by investigating C −J;1,1 . We first observe that this quantity has a clear physical interpretation, which can be derived by the application of the Hellmann-Feynman theorem [28], namely, i.e., it quantifies the variation of the occupation of dot j under small variations of the applied (virtual) gate potential V j , which tunes the dot's chemical potential. The behavior of C −J;1,1 is reported in Figure 2. The contour lines corresponding to the highest absolute values of the QC are two drop-shaped domains which occur at δ ≈ ± U for a certain range of values of V 3 -approximately V 3 < 8000 for the considered values of the parameters β and U. When V 3 becomes large enough ( V 3 > 8000 in this case), the contour lines of the two drop-shaped domains merge into a single one. A further increase in V 3 results in a significant drop in the absolute value of the QC. This is qualitatively explained in terms of the semiclassical picture developed in Ref. [28]: the QC peaks occur in regions of the parameters domain where the classical charge transitions [in this case, (0, 2) → (1, 1) → (2, 0)] take place, and small variations in V j lead to large variations in the population of dot j. In these regions, the energies of the ground state and of the first excited singlet undergo an avoided crossing. The peaks' height and width are related to the parameters describing such an avoided crossing, in particular to the ratio |T|/U. When V 3 is increased from the reference value 0 to higher values (i.e., the working point moves vertically from the bottom to the top of Figure 2), the amplitude of the avoided crossing gap, which is proportional to |T|, increases. Correspondingly, the height of the QC peak decreases, because the charge transition becomes smoother, i.e., the region in the parameter space where the occupation numbers of the two QDs are not well defined becomes larger. It should be noted that, for the parameters used in Figure 2, the value of V 3 ≈ 8000, where the drop-shaped contour lines merge, corresponds to T ≈ e 8 : the tunneling amplitude has thus been increased by a factor ≈3000 with respect to its value at V 3 = 0, and has reached a value that is comparable to the onsite repulsion U. Physically, this corresponds to a regime where the validity of the physical assumptions underlying the Hubbard model-the description of the double QDs in terms of weakly coupled site-localized states-becomes questionable. The effect of tuning the tunneling parameter manifests itself in the energy spectrum, as shown in Figure 3. The increase in V 3 leads to an exponential increase in |T|, which stretches the avoided-crossing gap between the ground state and the first-excited singlet from being negligible for V 3 = 0 (because of the large ratio U/|T| = 4000) to being sizeable at V 3 = 7000 (corresponding to U/|T| ≈ 3.65) and V 3 = 8000 (corresponding to U/|T| ≈ 1.34). When V 3 is further increased, the three singlet energies tend to flatten, and the bias-induced charge transitions associated with the DQD structure are lost.
In order to visualize the charge configuration of the two dots, we compute the groundstate occupation numbers n i −J using the Hellmann-Feynman theorem, which yields [28] where the last equality is valid for the DQD. In Figure 4, we plot the occupation of dot 1 for the ground state, as a function of the bias δ, for selected values of V 3 . It can be seen that, for low values of V 3 (resulting in T ≈ T 0 U), the occupation number changes abruptly with δ, as expected within a semiclassical picture with U |T|. Instead, for values of V 3 yielding |T| ≈ U, the system undergoes a transition from a double QD to an effective single QD, as reflected in the plot of n i −J by the disappearance of the sharp charge transitions at δ = ± U, and by the flattening of the curve at the value 1, corresponding to the charge being evenly distributed between the two QDs.  Figure 5 shows the contour plot of the other diagonal element of the QCM, namely C −J;3,3 , for the same set of parameters considered in Figure 2. We notice a strikingly different behavior compared to C −J;1,1 , due to the fact that C −J;3,3 does not probe the charge transitions between the two QDs, as C −J;1,1 does instead (see Equation (19)). In fact, from the Hellmann-Feynman theorem [28], it follows that whereĥ is the hopping operator in the Hubbard model. Therefore, just like the an oscillating bias voltage probes the dot occupations (Equation (19)), this QC measurement probes a combination of the tunneling coherence ĥ T k and of its derivative with respect to T.
In the case considered herein, the highest peaks of C k;3,3 occur at values of V 3 for which the description of the double dot in terms of the Hubbard model becomes questionable. Nevertheless, these features are interesting, since they indicate the possibility of obtaining peaks in the quantum capacitance that are qualitatively different from the ones related to the bias, and that have a physically different origin.

Directly Measurable QCs: The Largest Eigenvalue of the QCM
If the oscillating voltage is applied along the direction in the voltage space defined by the eigenvectors of the QCM, then the measured QC is given by the corresponding eigenvalue. The most relevant outcome of the multidimensional quantum capacitance is represented by the eigenvalue with the largest modulus, corresponding to the largest measurable QC. This is obtained for the optimal direction in the voltage space, along which the signal should be applied. As above, we refer to the QC related to the Hamiltonian ground state, C −J;− . In Figure 6, C −J;− is plotted as a function of δ and V 3 , for β = 1000 and U = 4000 (the same values considered for the diagonal elements of the QCM). In Figure 7, we show the cross-sections of the contour plot of C −J;− at selected values of V 3 .
From the comparison between the eigenvalue and the diagonal elements (Figures 2 and 5), it is apparent that the peaks of |C −J;− | occur very closely to the peaks of |C −J;1,1 | and |C −J;3,3 |. This suggests that, in the regions where the peaks of |C −J,i,i | are located, the i-th component of the QC eigenvector v −J,− dominates over the other two components. In other words, the main features that show up in the plot of the QC eigenvalue do not seem to result from significant mixing between the bias and tunneling voltages, but rather from the sum of the two diagonal contributions. In order to quantitatively verify the above picture, we report the components of the QC eigenvector v −J,− as functions of δ and for selected values of V 3 ( Figure 8). As already shown (Equation (17)), the components 1 and 2 are always opposites. Additionally, the peaks of their absolute values are located close to the charge transition points, δ = ± U. There, in the parameter range corresponding to weak interdot coupling, the dominant components in the eigenvector are those related to the bias (1 and 2). The eigenvector is instead dominated by the tunneling-related component (3) away from the charge transitions due to the exponential relation between the virtual and real tunneling gates. When the hopping is sufficiently large that the two QDs merge into a single one (e.g., Figure 8c,d), then component 3 dominates at all values of δ, indicating that there are no charge transitions between the QDs. This is consistent with the scenario described by Figure 4.

Discussion and Conclusions
In conclusion, we applied the recently developed theory of multidimensional quantum capacitance to the case of a double QD, controlled by two bias gates and one tunneling gate. In particular, the use of the latter gate is shown to control the height and width of the peaks of the QC related to the change in charge configuration induced by the bias voltage. Besides, the tunneling voltage is responsible for the presence of additional peaks, with a qualitatively different physical origin. These peaks show up both in the relevant diagonal element of the QCM, and in its eigenvalue with largest modulus.
The analysis was performed within a biased Hubbard model, which allows the derivation of analytical solutions for all the relevant quantities, including the eigenvalues and eigenvectors of the QCM. On the other hand, the considered range of physical parameters partially exceeds the one where the Hubbard model represents an accurate description of the double-dot system. For these cases, more quantitative estimates will require the combination of the general approach for the derivation of the QCM with more accurate descriptions of the electron (or hole) eigenstates in the double quantum dot.
The procedure involves the determination of an effective Hamiltonian based on model parameters which enter the MQCM. The connection between the MQCM and the QCM also requires additional, device-specific information, determining the dependence of the model parameters on the real gates. This information might result from realistic multiscale device simulations, which quantitatively account for the coupling between the gate voltages on the one hand, and the tunneling and bias voltages on the other. In all cases, the possibility to quantitatively determine the optimal working points and directions of oscillation in the M-dimensional voltage space is expected to provide experimenters with more customizable tools for the manipulation and readout of quantum states.