2 η$_{Q}$ = ( Q (Δ T = 30Γ) − | Q (Δ T = − 30Γ) | ) /Q (Δ T = 30Γ). We obtain η$_{Q}$ = 0 . 86 for E$_{B}$ = E$_{F}$ + 2Δ E/ 5 and 0.88 for E$_{B}$ = E$_{F}$ + 4Δ E/ 5. Fig. 2(c) shows DTC in units of Q$_{0}$k$_{B}$/ Γ. It is found that the rectification behavior is not very sensitive to the variation of E$_{B}$ . DTC is roughly linearly proportional to Δ T in the range − 20Γ < k$_{B}$ Δ T < 20Γ. In addition, we also find a small negative differential thermal conductance (NDTC) for E$_{B}$ = E$_{F}$ + 4Δ E/ 5. Similar behavior was reported in the phonon junction system. 10 Fig. 3 shows the heat current, differential thermal conductance and thermal power ( S = e Δ V/k$_{B}$ Δ T ) as functions of temperature difference Δ T for a three-QD case for various values of Γ$_{AR}$, while keeping Γ$_{B}$$_{(}$$_{C}$$_{)}$$_{,R}$ = Γ$_{B}$$_{(}$$_{C}$$_{)}$$_{,L}$ = Γ. Here, we adopt η$_{A}$ = | Γ$_{AL}$ − Γ$_{AR}$ | / (2Γ) instead of fixing η$_{A}$ at 0.3 to reflect the correlation of dot position with the asymmetric tunneling rates. We as- sume that the three QDs are roughly aligned with dot A in the middle. The energy levels of dots A, B and C are chosen to be E$_{A}$ = E$_{F}$ − Δ E/ 5, E$_{B}$ = E$_{F}$ + 2Δ E/ 5 and E$_{C}$ = E$_{F}$ + 3Δ E/ 5. U$_{AC}$ = U$_{BA}$ = 15 k$_{B}$T$_{0}$ , U$_{BC}$ = 8 k$_{B}$T$_{0}$ , U$_{C}$ = 30 k$_{B}$T$_{0}$ , and all other parameters are kept the same as in the two-dot case. The thermal rectification effect is most pronounced when Γ$_{AR}$ = 0 . as seen in Fig. 4(a). (Note that the heat current is not very sensitive to U$_{BC}$ ). In this case, we obtain a small heat current Q = 0 . 068 Q$_{0}$ at Δ T = − 30Γ, but a large heat current Q = 0 . 33 Q$_{0}$ at Δ T = 30Γ and the rectification efficiency η$_{Q}$ is 0.79. However, the heat cur- rent for Γ$_{AR}$ = 0 is small. For Γ$_{AR}$ = 0 . 1Γ, we obtain Q = 1 . 69 Q$_{0}$ at Δ T = − 30Γ, Q = 5 . 69 Q$_{0}$ at Δ T = 30Γ, and η$_{Q}$ = 0 . 69. We see that the heat current is sup- pressed for Δ T < 0 with decreasing Γ$_{AR}$. This im- plies that it is important to blockade the heat current through dot A to observe the rectification effect. Very clear NDTC is observed in Fig. 3(b) for the Γ$_{AR}$ = 0 . 1Γ case, while DTC is symmetric with respect to Δ T for the Γ$_{AR}$ = Γ$_{AL}$ case. From the experimental point of view, it is easier to measure the thermal power than the direction-dependent heat current. The thermal power as a function of Δ T is shown in Fig. 3(c). All curves except the dash-dotted line (which is for the symmetrical tunneling case) show highly asymmetrical behavior with respect to Δ T , yet it is not easy at all to judge the efficiency of the rectification effect from S for small | Δ T | ( k$_{B}$ | Δ T | / Γ < 10). Thus, it is not sufficient to determine whether a single QD can act as an efficient thermal rectifier based on results obtained in the linear response regime of Δ T/T$_{0}$ ≪ 1. 6 According to the thermal power values, the electrochemical potential e Δ V can be very large. Consequently, the shift of QD energy levels caused by Δ V is quite important. To illustrate the importance of this effect, we plot in Fig. 4 the heat cur- rent for various values of E$_{C}$ for the case with Γ$_{AR}$ = 0, U$_{BC}$ = 10 k$_{B}$T$_{0}$ and η$_{A}$ = 0 . 3. Other parameters are kept the same as those for Fig. 3. The solid (dashed) curves are obtained by including (excluding) the energy shift η$_{A}$ Δ V/ 2. It is seen that the shift of QD energy levels due Here N$_{A}$$_{(}$$_{B}$$_{)}$ is the average occupancy in dot A(B). There- fore, it is expected that the curve corresponding to E$_{B}$ = E$_{F}$ + 4Δ E/ 5 obtained with this delta function approximation is in good agreement with the full solu- tion, since E$_{B}$ is far away from the Fermi energy level. For cases when E$_{B}$ is close to E$_{F}$ , the approximation is not as good, but it still gives qualitatively correct behav- ior. Thus, it is convenient to use this simple expression to illustrate the thermal rectification behavior. The asym- metrical behavior of N$_{A}$ with respect to Δ T is mainly resulted from the condition Γ$_{AR}$ = 0 and Γ$_{AL}$ = 2Γ. The heat current is contributed from the resonant chan- nel with ϵ = E$_{B}$ , because the resonant channel with ϵ = E$_{B}$ + U$_{AB}$ is too high in energy compared with E$_{F}$ . The sign of Q is determined by f$_{LR}$ ( E$_{B}$ ), which indi- rectly depends on Coulomb interactions, tunneling rate ratio and QD energy levels. The rectification behavior of Q is dominated by the factor 1 − 2 N$_{A}$ , which explains why the energy level of dot-A should be chosen below E$_{F}$ and the presence of interdot Coulomb interactions is crucial. The negative sign of Q in the regime of Δ T < 0 indi- cates that the heat current is from the right electrode to the left electrode. We define the rectification efficiency as Q/γ$_{B}$ = π (1 − N$_{B}$ )[(1 − 2 N$_{A}$ )( E$_{B}$ − E$_{F}$ ) f$_{LR}$ ( E$_{B}$ )(3) + 2 N$_{A}$ ( E$_{B}$ + U$_{AB}$ − E$_{F}$ ) f$_{LR}$ ( E$_{B}$ + U$_{AB}$ )] , Fig. 2 shows the heat currents, occupation numbers, and differential thermal conductance (DTC) for the two- QD case, in which the energy levels of dot A and dot B are E$_{A}$ = E$_{F}$ − Δ E/ 5 and E$_{B}$ = E$_{F}$ + α$_{B}$ Δ E , where α$_{B}$ is tuned between 0 and 1. The heat currents are exporessed in units of Q$_{0}$ = Γ $^{2}$/ (2 h ) through out this article. The intradot and interdot Coulomb interactions used are U$_{ℓ}$ = 30 k$_{B}$T$_{0}$ and U$_{AB}$ = 15 k$_{B}$T$_{0}$ . The tunneling rates are Γ$_{AR}$ = 0, Γ$_{AL}$ = 2Γ, and Γ$_{BR}$ = Γ$_{BL}$ = Γ. k$_{B}$T$_{0}$ is chosen to be 25Γ throughout this article. Here, Γ = (Γ$_{AL}$ + Γ$_{AR}$) / 2 is the average tunneling rate in energy units, whose typical values of interest are between 0.1 and 0.5 meV. The dashed curves are obtained by using a simplified expression of Eq. (2) in which we set the average two particle occupation in dots A and B to zero (resulting from the large intradot Coulomb interactions) and taking the limit that Γ ≪ k$_{B}$T$_{0}$ so the Lorentzian function of resonant channels can be replaced by a delta function. We have QD junctions involving two QDs and three QDs for var- ious system parameters. We first determine Δ V by solv- ing Eq. (1) with J$_{e}$ = 0 (the open circuit condition) for a given Δ T , T$_{0}$ and an initial guess of the average one- particle and two-particle occupancy numbers, N$_{ℓ}$ and c$_{ℓ}$ for each QD. Those numbers are then updated accord- ing to Eqs. (5) and (6) in Ref. [7] until self-consistency is established. For the open circuit, the electrochemical potential will be formed due to charge transfer generated by the temperature gradient. This electrochemical po- tential is known as the Seebeck voltage (Seebeck effect). Once Δ V is solved, we then use Eq. (2) to compute the heat current.