Molecular wires of the acene-family can be viewed as a physical realization of a two-rung ladder Hamiltonian. For acene-ladders, closed-shell ab initio calculations and elementary zone-folding arguments predict incommensurate gap oscillations as a function of the number of repetitive ring units, , exhibiting a period of about ten rings. Results employing open-shell calculations and a mean-field treatment of interactions suggest anti-ferromagnetic correlations that could potentially open a large gap and wash out the gap oscillations. Within the framework of a Hubbard model with repulsive on-site interaction, U, we employ a Hartree-Fock analysis and the density matrix renormalization group to investigate the interplay of gap oscillations and interactions. We confirm the persistence of incommensurate oscillations in acene-type ladder systems for a significant fraction of parameter space spanned by U and .
I. INTRODUCTION
The question how properties of a macroscopic system emerge when more and more atoms or molecules accumulate exhibits many different facets and, for that reason, reappears every once in a while in different contexts. On the simplest level of “emergence,” one could consider the appearance of band-structures in crystal growth. A first qualitative description of this phenomenon can already be given within a picture of non-interacting quasi-particles. A conceptually analogous situation arises in mesoscopic physics: when a quantum dot couples to an electrode, the dot states hybridize with contact states, thereby acquiring the level broadenings, Γ, that signal finite lifetime effects. Qualitatively new physics emerges upon accounting for interaction effects. For instance, in quantum-dots, the Kondo phenomenon results from the interplay of hybridization and the Coulomb-energy. Another example displaying the appearance of strong correlation effects with growing system size is superconductivity.1,2 There, one is facing a crossover scenario, where the super-conducting gap, , competes with the single particle level spacing, .
In this paper, we investigate an incarnation of a similar motif as it appears in molecular electronics: we ask how the electronic spectral properties of a molecular wire built out of repetitive units (“rings”) evolve with increasing . Motivated by earlier experiments,3–5 our specific example will be the acene-family such as benzene, naphthalene, anthracene, and others. Due to their peculiar electronic structure, oligo-acenes have received a considerable amount of attention in organic chemistry over the past 50 years.6–8 A recent technological interest is based on suggestions to use them as light harvesters in organic photo-voltaic devices.9 Besides being a simple realization of nano-wires, acenes can be viewed as the narrowest graphene nano-ribbons with zig-zag terminated edges. Due to the recent progress in synthesis methodology, oligoacenes up to length (nonacene) have already been produced.10,11
Previously, we have found that linear (oligo)acenes exhibit oscillations in the (or length) dependence of optical excitation gaps.12 Oscillations of electronic properties with the geometric size of a nano-system are per se frequently encountered. An established example is given by the threefold gap-oscillation of zigzag carbon nanotubes (armchair-ribbons) as a function of tube diameter (width).13–15 In this context, the three-fold periodicity reflects the symmetry of the hexagonal Brillouin zone, with the Dirac points located at the corners. A remarkable aspect of the acene oscillations is that they are incommensurate in the sense that the oscillation period is not dictated by lattice symmetry; it can reach periods of ten times the length of the unit cell. The oscillations come with an important consequence, namely, a close-to closing gap for molecules of certain (periodically repeating) finite length.
We have demonstrated that similar to the case of the nanotubes, the acene-oscillations can likewise be described (nearly) quantitatively by a selection rule (zone-folding) for the wave-vectors applied to the band structure of the underlying extended tight-binding model (Fig. 1). The procedure is analogous to the threefold periodicity in band-gaps of nanotubes13 and nanoribbons;14,15 in the acene case, however, incommensurate periodicity is reached because of the non-universal position of labelling the longitudinal momentum at the acene-analogue of the Dirac point. Based on the computational analysis of the acene-like Hubbard model, we find that the zone-folding argument is not invalidated by interactions beyond the independent quasi-particle picture; the oscillations persist given the interaction strength U does not exceed a threshold of the order of the bandwidth—at least as long as interactions are screened, i.e., short-ranged.
(a) Parts of the oligoacene backbone in the tight-binding representation with a nearest neighbor hopping (t, yellow), hoppings to the second neighbor (, dashed blue) and third neighbor (t″, dotted black). On the left, we draw these hoppings for the “outer” carbon site, on the right for the “inner” carbon. denotes the second neighbor hooping of the outer sites, of the inner sites. In the present study, we adopt and, in contrast to Ref. 12, . (b) The tight binding band-structure corresponding to the long-wire limit for t = 1, t′ = 0, and t″ = 0.3. The lines are the analytical result according to (2) with . The crosses are the results from a tight binding system consisting of NR = 31 rings and periodic boundary conditions. The valence- and conduction-band are crossing at the Γ-point kD.
(a) Parts of the oligoacene backbone in the tight-binding representation with a nearest neighbor hopping (t, yellow), hoppings to the second neighbor (, dashed blue) and third neighbor (t″, dotted black). On the left, we draw these hoppings for the “outer” carbon site, on the right for the “inner” carbon. denotes the second neighbor hooping of the outer sites, of the inner sites. In the present study, we adopt and, in contrast to Ref. 12, . (b) The tight binding band-structure corresponding to the long-wire limit for t = 1, t′ = 0, and t″ = 0.3. The lines are the analytical result according to (2) with . The crosses are the results from a tight binding system consisting of NR = 31 rings and periodic boundary conditions. The valence- and conduction-band are crossing at the Γ-point kD.
Building upon our previous work, we intend to give a comprehensive study of interactions and incommensurate oscillations (IO) of acenes by contrasting Hartree-Fock (HF) mean field treatments against microscopic numerical studies. We model acene-type molecular wires by finite length Hubbard ladders. A principal physical incentive is given by the potentially rich phase diagram of interacting acene wires.16 In this context, one could hope to manipulate acenes via chemical synthesis, or engineering of the environment so as to explore some part of this phase diagram. (Admittedly before this is experimentally feasible, major technological difficulties have to be overcome that are related to chemical stability. Momentarily, long acene wires need to be kept within a noble-gas matrix to prevent further reactions.11)
Further methodological incentive is provided by recent open-shell -calculations for the acene-series;17 they indicated that in short wires the IO might disappear in the presence of anti-ferromagnetic order in the ground state. At first sight, this result may look very plausible. However, effects of quantum-fluctuations are neglected in its derivation and in quasi-one dimensional systems their effect typically is very strong. Specifically, in the long-wire limit, quantum fluctuations tend to destabilize mean-field ordered phases that break a continuous symmetry.18 In this (quantum-disordered) phase, the order-parameter exhibits significant short-range correlations that leave a characteristic signature in local probe measurements.
In this work, we clarify the fate of the anti-ferromagnetic correlations in acene-type ladders in the presence of quantum fluctuations, i.e., by providing a joint analysis from the Hartree-Fock (HF) mean field and density matrix renormalization group (DMRG). Within HF mean field analysis, we reproduce the observation made in the ab initio study, Ref. 17, for acene-molecules with Coulomb interactions: acene-ladders with onsite Hubbard interactions exhibit anti-ferromagnetic order on the mean field level along with a large charge gap, thereby diminishing IO. By microscopic DMRG analysis, however, we show that also in short ladders quantum fluctuations destroy the long-range order. In particular, they significantly reduce the HF excitation gap, restoring the IO signature and up to intermediate interaction strengths of the order of the bandwidth.
II. MODEL AND METHODS
A. Model definition
We study the acene Hubbard-ladder
at half filling, where the first term describes the tight-binding dynamics, and the second term in Eq. (1) the onsite Hubbard repulsion with strength U. In the case of periodic boundary conditions (PBCs), we have a linear length of L = 2NR sites, where NR gives the number of hexagon (benzene) rings. For hard wall boundary conditions (HWBCs), we complete the last ring leading to L = NR + 1. The total number of sites is 2L. The nearest neighbor hopping is characterized by an amplitude , and the longer-range hoppings accordingly by t′, , and (Fig. 1). If not stated otherwise, t will be taken as the unit of energy and =. Note that and give rise to a non-avoided crossing of conduction- and valence band at non-commensurate wave-vectors predicted for the case U = 0.19 For , the system features a quadratic band touching.16 For U = 0, the energy spectrum of (1) reduces to
where denotes a band and parity index, respectively, yielding four bands resulting from the four-site unit cell.
In the following, we exclude next-nearest neighbor hopping, , so there is no hopping between sites sharing the same sublattice; the model becomes bipartite and a chiral symmetry is imposed connected to the particle-hole symmetry of the spectrum (Fig. 1). This choice can be made without loss of generality regarding IO, as the main effect of breaking the chiral symmetry in more general models with non-vanishing would be to shift the position of the Dirac point . The corresponding analytical calculations are given in Appendix A.
General comments on Hubbard-models. We briefly comment about the applicability of Hubbard models to real systems. Indeed, models along the lines of (1) have been investigated intensively over the last fifty years in physics and in chemistry. They provide substantial simplifications by (a) ignoring the coupling to molecular vibrations, (b) including only one orbital per atom, and (c) including all other degrees of freedom in terms of effective parameters, only, such as a short range interaction that mimics screening. In fact, information about the molecular geometry (e.g., the relative position of the atoms) is absorbed in hopping parameters and interaction matrix elements. These parameters span a multi-dimensional space which has a hyperplane that features those particular combination of parameters that can be reached in realistic situations.
In return, because Hubbard models are in a certain sense minimal, they allow for a numerically exact treatment of the many-body problem, with a size of the computational Hilbert space that is much larger as compared to most other many-body treatments that attempt to be more “first principles” on the Hamiltonian level.
Keeping this in mind, Hubbard models like (1) are by no means as oversimplistic as they may appear at first sight. First, if properly parametrized such models have a good chance to give a realistic, sometimes even a quantitative description of low-energy excitations also in those situations where quantum-fluctuations are dominant, so that most competing methods fail. Second, one often is interested in the system behavior under varying external conditions, i.e., in the phase diagram of the model; one would like to figure out what phases exist, what are their properties and, ideally, also where are the phase-boundaries located at. Since the qualitative aspects of phase diagrams tend to be very robust under a change of the Hamiltonian, one often adopts convenient choices for the model parametrization, even though they are sometimes far from realistic.
Once the phase-diagram has been mapped out, it may be possible to predict the qualitative behavior for a given experiment, even though the actual parameters representative for this situation may not be known with high accuracy, at least in those situations where very stable phases have been identified.
B. Computational details
In this spirit, we continue our analysis and treat the model Hamiltonian (1) numerically within DMRG.20,21 Similar to previous authors,22,23 we also observe that the ground state is a spin singlet independent of the wire length. Thereby, we track the many-body Hamiltonian at a fixed particle number (equal to the number of lattice sites) and in the spin sector. By targeting not only the ground state but also excited states, we obtain the optical gap by taking the energy difference between the first excited and the ground state. For non-interacting systems, the optical gap corresponds to the gap between the highest occupied molecular orbital (HoMO) and the lowest unoccupied molecular orbital (LuMO) as it is given by the lowest particle-hole excitation. Due to spin-rotational invariance, the excited states of the model in (1) come in spin-multiplets according to irreducible representations of SU(2). By performing reference calculations in the spin-sectors , we have ascertained that the lowest excitation is a triplet state, and that the ground state is indeed a singlet state. Details of the DMRG implementation are given in Appendix C. To investigate the effect of quantum fluctuations and to link with previous mean field treatments, it is useful to compare the DMRG results to a HF description. We have also implemented a self-consistent HF cycle to obtain mean field estimates for the energy gaps, ground state, and magnetization ( Appendix C). In the HF theory, we look for the lowest lying (charge-neutral) excitation energy to determine the energy gap, for which we compute the difference in energy between the lowest unoccupied and the highest occupied eigenvalue of the Fock operator.
III. RESULTS AND DISCUSSION
A. Hartree-Fock analysis
Before investigating finite size features of the spectrum, we analyze the infinite wire limit. In Fig. 2(a), we display the evolution of the HF gap for repulsive interaction strength U. We extract the bulk gap from finite-size calculations by extrapolating the finite size gap vs. the inverse linear length of the system ( Appendix B). While the non-interacting wire exhibits a band crossing at the Fermi-energy and therefore is metallic, the repulsive interaction U opens up a correlation gap.
(a): Evolution of the excitation gaps (+) and (×) in units of the interaction U vs. the increasing strength of the repulsive Hubbard interaction U. For large U, we recover the expected asymptotic behavior. For , the numerical extrapolation is impaired by the IO; the asymptotic behavior is indicated by the (yellow) dashed line; see text and for comparison Fig. 9. (b) Corresponding staggered magnetization (○) over the band-gap in units of U displaying a linear dependence in the complete regime.
(a): Evolution of the excitation gaps (+) and (×) in units of the interaction U vs. the increasing strength of the repulsive Hubbard interaction U. For large U, we recover the expected asymptotic behavior. For , the numerical extrapolation is impaired by the IO; the asymptotic behavior is indicated by the (yellow) dashed line; see text and for comparison Fig. 9. (b) Corresponding staggered magnetization (○) over the band-gap in units of U displaying a linear dependence in the complete regime.
Within the HF theory, the ground state is found to be unstable against antiferromagnetic order; the corresponding staggered magnetization is shown in Fig. 2(b). As we will see in our DMRG-calculations below, this order will be destroyed by quantum fluctuations. We here witness the established fact that mean-field theories overestimate the tendency for symmetry breaking, at least if the interactions are short ranged. In passing, we mention that broken symmetries refer to the long-distance properties of a quantum state, which are not well described on the HF-level. In contrast, the total energy is predominantly dictated by the pair-correlation function at short-distances. The latter is largely insensitive to the global symmetries. Therefore, with respect to total energy calculations, HF can often give useful answers even though the long-range physics is described qualitatively incorrect (Löwdin’s “symmetry breaking dilemma” of HF24).
In the perturbative limit of weak repulsion, , we can restrict ourselves to states with anti-ferromagnetic (AF) ordering on top of the non-interacting ground state. In this situation, the gap opens in long (but finite-length) wires as
where is the finite size gap of the non-interacting wire, ; MHF denotes the staggered magnetization amplitude. We recall that at the response of energies to the staggered field is linear in U, implying a quadratic gap opening in the thermodynamic limit (where ).
At intermediate interactions, the HF values for the gap opening allow for a phenomenological two-parameter fit,
where D is of the order of the non-interacting bandwidth; we obtain and , provided U is not too small. As expected, in the large U limit, the HF particle-hole gap is set by U. Its evolution follows closely one of the anti-ferromagnetic order parameters , Fig. 2(b).
We now turn to the HF analysis of acene wires with finite length. In this case, one expects that the opening of the HF gap decreases the ground-state energy only if exceeds the non-interacting (“native”) level-spacing at the Fermi energy significantly. Hence, for shorter wires, the gap opening is most effective at particular values of NR that have an anomalously small single-particle gap—and ineffective at the others. In other words, only at NR-values situated close to a minimum of the IO an AF-order develops and at the others the non-interacting ground-state prevails. In this situation, IO remain visible even in .
For longer wires, however, the IO fade away because the native band-gap is suppressed as and eventually no longer overcomes . As soon as exceeds the maxima of the native band gaps, the IO disappear from the HF spectra. Indeed, the IO in are seen to become weaker with increasing NR (and U) (Fig. 3). This explains why large system sizes are needed in order to obtain the HF gap in the small U regime: the finite size induced gap has to fall below the HF gap.
Evolution of the HF-magnetization MHF and HF-gap with increasing number of rings NR at U = 0.8 and HWBC. As is seen, the HF magnetization MHF of the acene-wire (×) is non-vanishing whenever the HF-gap (○) exceeds the non-interacting gap (+). For this reason, even the mean-field magnetization exhibits the IO, in principle. They fade away only for large system sizes, when the non-interacting finite site gap falls below the interaction induced HF gap.
Evolution of the HF-magnetization MHF and HF-gap with increasing number of rings NR at U = 0.8 and HWBC. As is seen, the HF magnetization MHF of the acene-wire (×) is non-vanishing whenever the HF-gap (○) exceeds the non-interacting gap (+). For this reason, even the mean-field magnetization exhibits the IO, in principle. They fade away only for large system sizes, when the non-interacting finite site gap falls below the interaction induced HF gap.
As the magnetization closely follows , it emerges preferably in acene wires that exhibit a relatively small (native) band-gap (Fig. 3). Vice versa, the HF-ground state remains non-magnetized at wire-lengths, such as , where the native gap is particularly large; there mean-field magnetism only appears at larger interaction values. In contrast, for large interaction values, the HF gap surpasses the non-interacting gap already for a small system size (Fig. 4), and the IO disappear in that regime.
Evolution of the HF-magnetization (×) and HF-gap () with increasing number of rings NR, L = 2NR + 1, at U = 1.0, 1.5, 2.0. As in Fig. 3, the IO vanish if the HF gap is larger than the non-interacting gap. Accordingly, at these instances the staggered magnetization MHF stays finite.
Evolution of the HF-magnetization (×) and HF-gap () with increasing number of rings NR, L = 2NR + 1, at U = 1.0, 1.5, 2.0. As in Fig. 3, the IO vanish if the HF gap is larger than the non-interacting gap. Accordingly, at these instances the staggered magnetization MHF stays finite.
B. DMRG calculations
We reignite the discussion from the viewpoint of DMRG by the infinite wire limit. Within our DMRG-implementation, this limit can only be achieved through finite size extrapolation: for every value of U, the extrapolation is performed to yield one data point of the blue line for the excitation gap in Fig. 2(a), where it is compared against the HF result. We mention that at small values of U numerical calculations are more challenging to converge in DMRG, so that the low U asymptotics is not fully resolved in Fig. 2(a). The difficulty is that due to the metallic character of the wire, the entanglement content is spread over a large manifold of modes that have non-negligible weight in the reduced density matrix.
By comparing to the HF result, one directly observes the effect of strong quantum fluctuations in the DMRG trace shown in Fig. 2(a). Fluctuation effects are visible by a reduction of the excitation gap but also by a modified gap evolution with U as compared to the HF trajectory. The asymptotic power-law at large U seen in Fig. 2 is readily understood in terms of the mapping to the Heisenberg spin-chain. It proceeds via the effective coupling constant J(U) = 4t2/U that sets the scale for in this limit.25–28 As is seen in Fig. 2(a), the Heisenberg limit is reached at values of U exceeding the bandwidth, D, by about a factor of two.33
Fig. 5 displays the most important result of our work. It confronts the mean-field data for the excitation gaps with the results obtained from the DMRG-calculations for finite-length wires. Specifically, the evolution of their DMRG-excitation gaps, , is shown for the increasing length of the wire NR at different U. It is clearly seen there, that the DMRG traces exhibit pronounced IO even for those values of U at which the residual oscillations in have almost faded away, already. Very much in the spirit of the long-wire limit, we understand this result as a manifestation of quantum-fluctuations in finite-size wires. They destroy the tendency of mean-field theories to open up excitation gaps, thereby restoring the IO as a hallmark of the symmetry-unbroken phase.
Evolution of the lowest lying excitation gap for acene-ladders with increasing number of rings NR, L = 2NR + 1, for varying interaction strengths U (). For comparison also the HF-result () is shown. In contrast to the strongly reduced IO on HF-level, they survive for not too large U in the DMRG calculations.
Evolution of the lowest lying excitation gap for acene-ladders with increasing number of rings NR, L = 2NR + 1, for varying interaction strengths U (). For comparison also the HF-result () is shown. In contrast to the strongly reduced IO on HF-level, they survive for not too large U in the DMRG calculations.
IV. CONCLUSIONS
We have studied the acene Hubbard model for a large range of on-site repulsion U. The gap oscillations, that are the characteristics of the weakly correlated phase, survive as a periodic modulation also in the strongly correlated phase. The latter still persist for a sizable interaction strength, as quantum fluctuations significantly reduce the formation of the correlation gap.
With an eye on acene-molecules, we expect that our results have the following implications. We have performed model-calculations employing short-range interactions. They may apply to acenes in the presence of a screening environment, such as substrates or electrochemical solutions. With screened Coulomb-interactions, the effective value of U should be significantly smaller than the molecular bandwidth, implying that our weakly correlated scenario should apply. In contrast, for gas-phase molecules screening is very weak and the applicability of our model-calculation is not immediately guaranteed. This is the situation investigated in Refs. 7, 8, and 17; we conclude with a brief overview.
Concerning Ref. 17, it remains to be seen if quantum-fluctuations can restore the spin-rotational symmetry also in the presence of a long-range Coulomb interaction. Very recent quantum-chemistry calculations on the level of density-functional multi-reference theory (DFT/MRCI) for oligoacenes (up to NR = 9) and pp-RPA (up to NR = 12) may give an indication that this indeed could be the case.7,8
The authors of Ref. 7 confirm earlier findings29 that best agreement between DFT/MRCI results and experimental IR-excitation energies is achieved for closed-shell reference states (with B3LYP functional) even though spin-unrestricted calculations may produce a lower total energy. It is thus suggested that quantum-fluctuations are indeed strong in long acene-wires so that single-determinant approximations of the ground-state are not representative. In agreement with this interpretation, these authors observe that the weight of the reference state in the DFT/MRCI ground state is progressively decreasing with growing wire length, so that more and more Slater determinants have to be included to accurately represent this state.
The authors of Ref. 8 come to a qualitatively similar conclusion. The pp-RPA method is an interesting variant of RPA-type approaches, since it incorporates vertex-corrections that are missed, e.g., in standard RPA bubble-summations. Thus, it allows to include a certain amount of quantum fluctuations when treating unscreened wires (i.e., with Coulomb-interactions) that would otherwise be too long for more accurate approaches. A certain drawback of the method is that it is uncontrolled: while it can signalize qualitatively that quantum fluctuations are likely to be important, there is no way of telling with confidence what they actually do in that case, quantitatively.
On the qualitative level, the conclusions of Refs. 7 and 8 are consistent. In particular, also the authors of Ref. 8 emphasize that the multi-determinant character of the ground state increases rapidly with growing wire length. It hinders the formation of long-range magnetic order and implies that mean-field results (i.e., unrestricted open-shell calculations) are to be interpreted with a grain of salt.
It remains to be seen, however, if restoring the spin-symmetry in the ground state is enough to also restore the IO. Namely, even assuming that long-range and short-range models of molecular wires share the same phases, it is still unclear whether the effective value of U in the gas-phase calculations is small or large as compared to the (non-interacting) bandwidth.
Interestingly, the authors of Ref. 8 propose an exponential decay of the lowest lying spin-triplet gap with increasing NR saturating near NR = 11 at roughly 0.13 eV without a trace of IO. If indeed correct, this finding could still be interpreted within the short-range model. It could be taken as a hint that the effective interaction U as it appears in the study, Ref. 8, should be considered large. However, the authors’ data underlying their exponential fit (Table 1 in Ref. 8) appear to be inconsistent with an exponential asymptotics. In fact, the data for the spin-restricted cases exhibit even non-monotonic behaviour (upturn at large wire length), which is qualitatively consistent with the onset of IO as predicted by us.
ACKNOWLEDGMENTS
We express our gratitude to Holger Bettinger, Reinhold Fink, Steve Kivelson, and Jan Wilhelm for inspiring discussions. Also, we thank Jan Wilhelm for sharing a manuscript17 prior to publication. P.S. and R.T. are supported by DFG-SFB 1170 and the ERC starting grant TOPOLECTRICS (No. ERC-StG-336012). F.E. acknowledges support of the DFG under Project Nos. EV30/8-1 and EV30/11-1. Most of the calculations were performed on the compute cluster of the YIG group of Peter Orth.
APPENDIX A: NON-INTERACTING BAND-STRUCTURE AND AN EFFECTIVE MODEL
We consider a non-interacting (U = 0) polyacene wire with the lattice and single particle Hamiltonian of (1) explained in Fig. 1. For this model, we derive analytic expressions for the wave number of the nodal point kD, Fermi velocity υF, and the oscillation period of gaps of finite chains. The ladder-like lattice features a symmetry with a line group D∞h. An important member of the group is a mirror symmetry with respect to the interchange of stringers which we call parity. The third-neighbor hoppings can be divided into two groups: a hopping that connects atoms related by the parity and the two hoppings that connect different rungs (Fig. 1). We denote the former by and the latter by . We label the 4 Bloch bands by the wave-number in units of a−1 (the inverse lattice spacing), a parity and a remaining band index . The explicit expression for the band dispersions is given in (2). The bands are particle-hole symmetric only if t′ = 0, whence .
The two bands with b = −1 give rise to a linear dispersion bυF(k −kD) at the Fermi level if the condition is met, as we will show below. We derive analytic expressions for kD and υF under the condition t′ = 0. The nodal (“Dirac”) point kD for t′ = 0 is the solution of the equation
In the limit of small, the nodal (Dirac) point and the Fermi velocity are expressed by
For oligoacenes within the tight-binding approximation parameterized as above, the band crossing implies gap oscillations with the unit cell period of
In Equations (A2) and (A3), the symbol denotes all terms that are of second order in and . We see that the band crossing sets in with , provided that the signs of t and are equal.
The hopping t′ preserves the band crossing and weakly renormalizes kD and υF. We remark that the low-energy behavior of the non-interacting polyacene model used in the present paper can be mapped onto a simpler effective model with and t and non-zero, i.e., a ladder model with alternating rungs. This model has been introduced originally by Kivelson and Chapman.19
APPENDIX B: EXTRACTING THE BULK GAP FROM FINITE-SIZE WIRES
Assuming that the gap opening only affects the band structure close to the Dirac point, one would expect that the gaps are given by a minimal model of a type
that follows the zone-folding argument of Ref. 12. The square root term in Eq. (B1) accounts for the level crossing with parameters , , and ϕ; and ϕ are fit-parameters that accommodate additional finite-size effects. Somewhat unexpectedly, our calculations indicate that the band-structure underlying Eq. (B1) does not provide a faithful description of the true bandstructure of the infinite wire with the range of U-values that we have been able to investigate. At values , the band-structure deviates significantly from the massive Dirac-shape and Eq. (B1) no longer applies. (For instance, we find a stronger damping of the oscillations as expected from Eq. (B1) indicating an interaction induced flattening of the band.)
In the limit of large interaction U, the finite size oscillations are completely washed out and the thermodynamic limit of the interaction induced gap can be extracted from a linear fit in the inverse system size 1/L, see Fig. 6. However note that for U = 2.5 an oscillatory part is already visible for not too large system sizes. Therefore we take a heuristic approach modeling the decay of the oscillations with a power law decay
as shown in Fig. 7. Here, a/L describes the decay of the non-oscillatory part, while the cosine term describes a decay of the oscillation faster than 1/L, i.e., . b gives the amplitude of the oscillations, kD gives the periodicity, and η allows for shifting the phase of the oscillation. α is a heuristic exponent modelling the damping of the oscillations. In our fits, we always find and α approaches zero for .
Extraction of the DMRG gap for large interaction values U (□: 2.5, ○: 5.0, ▽: 10.0, ⋄: 50.0) via a linear fit of the finite site excitation gap vs. the inverse system length 1/L, L = 2NR + 1.
Extraction of the DMRG gap for large interaction values U (□: 2.5, ○: 5.0, ▽: 10.0, ⋄: 50.0) via a linear fit of the finite site excitation gap vs. the inverse system length 1/L, L = 2NR + 1.
Extraction of the DMRG gap for not too large interaction values U (×: 1.0, +: 1.5, ▽: 2.0). The extrapolation is performed using Eq. (B3) for vs. the inverse length 1/L, L = 2NR + 1, and displayed by lines, while the symbols correspond to the DMRG results.
Extraction of the DMRG gap for not too large interaction values U (×: 1.0, +: 1.5, ▽: 2.0). The extrapolation is performed using Eq. (B3) for vs. the inverse length 1/L, L = 2NR + 1, and displayed by lines, while the symbols correspond to the DMRG results.
The gaps shown in Fig. 2 are obtained from the linear fits for , see Fig. 6, and from fits of Eq. (B3) for . We mention that the fit for U = 1.0 is still sensitive to the fit range, and the result should be taken with care. For this reason, we refrained from including results for . In addition, in that regime, it is not sufficient to fit a single cosine, as the results are getting close to the saw tooth like behavior as in the non-interacting case. In order to fit the results for small U, we introduce a (truncated) saw tooth function
and use it instead of the cosine in Eq. (B1). Here we imposed a limit of n = 17. The results for U = 0.2, 0.5 are displayed in Fig. 8. While the fits do capture the oscillation quiet well, the finite size oscillations are orders of magnitudes larger than the extracted gaps. In return, we do not trust the extracted values and did not include them in Fig. 2.
Extraction of the DMRG gap for small interaction values from vs. the inverse length 1/L, L = 2NR + 1. For not too large U, ○: 0.2, ▢: 0.5. The fits correspond to Eq. (B1) where the cosine is replaced by a saw tooth equation (B4).
APPENDIX C: METHODOLOGICAL DETAILS
1. Numerical Hartree-Fock
Due to SU(2) symmetry combined with a local on-site interaction, we can restrict our HF equations to the Hartree terms
where denotes the non-interacting system according to Eq. (1), σ denotes up and down spins, the opposite spin to σ, and the expectation value of the local density operator at site x with spin σ. We have checked that by breaking the remaining conservation; therefore allowing for Fock terms of the form , we only find solutions where the spin quantization axis is rotated with respect to the Hartree solution. Note that in this case the solution depends on the initialization of the self consistency loop.
Since up and down spins do not mix in (C1), we can represent the up and down spin sectors by independent matrices. We start the HF calculations by first diagonalizing the non-interacting system in order to obtain the ground state energy E0 as a reference. Note that in this case we obtain a homogeneous solution . We then explicitly introduce a staggered magnetic order by setting , where we take “+” on one sublattice of the bipartite system and “−” for the other sublattice (Fig. 1). With these initial values, we diagonalize and calculate the locate densities which we then iteratively insert into (C1). We perform at least five iterations and continue until the ground state energy EHFȂ of Hamiltonian changes by less then 10−10t. In order to avoid getting stuck in an oscillation between two solutions, we damp the HF self consistency loop by taking an average of 0.7 times the new density and 0.3 times the densities of the preceding iteration.
Finally, we check whether the ground state energy EHF including the double counting corrections
is smaller than the ground state energy E0 of the non-interacting system. Indeed we always find . The spin-averaged local density of our HF solution is homogenous and the total staggered magnetisation MHF is given by the difference between the magnetization of the two sub lattices
2. Hartree-Fock: Commensurate case
The main problem of obtaining the HF solutions in Fig. 2 roots in the finite size oscillations being incommensurable. Therefore we need very large system sizes in order to start with an U = 0 finite size gap being smaller than the interaction induced gap in the HF calculations. This extraction of the HF gap gets simplified if we switch to the case of commensurate oscillations by choosing a fine-tuned t″ from our solution Eq. (A1). In Fig. 9, we display the results in the case of a periodicity of twenty rings. We obtained the corresponding t″ = 0.13 by solving Eq. (A1) numerically.34Since we now always hit a minimum of the gap oscillations precisely, i.e., a zero in the non-interacting case, we can extract the HF gap without difficulties. As expected, the HF gap indeed opens . In addition, we see the same exponential behavior for intermediate to large U as displayed in Fig. 2.
The HF gap as shown in Fig.2, left, in the case of commensurate oscillations t″ = 0.13 where the periodicity is given by 20 rings. The gap is obtained by HF calculations for .
The HF gap as shown in Fig.2, left, in the case of commensurate oscillations t″ = 0.13 where the periodicity is given by 20 rings. The gap is obtained by HF calculations for .
3. DMRG
Within our DMRG calculations, we keep up to 8000 states per block; the dimension of the target space grows up to , the discarded entropies are considerably below 10−4 for , while for the discarded entropies reach for systems with in each DMRG step, typically significantly smaller. We employ at least 7 finite lattice sweeps with an optimized infinite lattice warm up for the smallest system sizes. In order to reach large system sizes, we iteratively restart a NR simulation in order to investigate a NR + 1 ring. Note that for this we only have to deactivate the wave function prediction in the initial step of a restart. In order to initialize the simulations, we apply a sliding environment block approach that we already employed successfully to fractional quantum hall systems.30–32 Instead of using a reflection or a growing right block during the DMRG infinite lattice warm up sweeps, we employ a right block of up to 5 sites, which we can always build exactly in such a way that we always have 4NR sites for PBCs and 4NR + 2 sites for HWBCs. With this initial warm up, we avoid running into excited states in the initial sweep.
REFERENCES
For reference, the precise value used in the numerics is t″ = 0.130 064 666 287 555.