Bond order wave (BOW) phase of the extended Hubbard model: Electronic solitons, paramagnetism, coupling to Peierls and Holstein phonons
aa r X i v : . [ c ond - m a t . s t r- e l ] J un Bond order wave (BOW) phase of the extended Hubbard model: Electronic solitons,paramagnetism, coupling to Peierls and Holstein phonons
Manoranjan Kumar and Zolt´an G. Soos
Department of Chemistry,Princeton University, Princeton NJ 08544 (Dated: November 11, 2018)The bond order wave (BOW) phase of the extended Hubbard model (EHM) in one dimension (1D)is characterized at intermediate correlation U = 4 t by exact treatment of N -site systems. Linearcoupling to lattice (Peierls) phonons and molecular (Holstein) vibrations are treated in the adiabaticapproximation. The molar magnetic susceptibility χ M ( T ) is obtained directly up to N = 10. Thegoal is to find the consequences of a doubly degenerate ground state (gs) and finite magnetic gap E m in a regular array. Degenerate gs with broken inversion symmetry are constructed for finite N for a range of V near the charge density wave (CDW) boundary at V ≈ . t where E m ≈ . t is large. The electronic amplitude B ( V ) of the BOW in the regular array is shown to mimic atight-binding band with small effective dimerization δ eff . Electronic spin and charge solitons areelementary excitations of the BOW phase and also resemble topological solitons with small δ eff .Strong infrared intensity of coupled molecular vibrations in dimerized 1D systems is shown to extendto the regular BOW phase, while its temperature dependence is related to spin solitons. The Peierlsinstability to dimerization has novel aspects for degenerate gs and substantial E m that suppressesthermal excitations. Finite E m implies exponentially small χ M ( T ) at low temperature followed by analmost linear increase with T . The EHM with U = 4 t is representative of intermediate correlationsin quasi-1D systems such as conjugated polymers or organic ion-radical and charge-transfer salts.The vibronic and thermal properties of correlated models with BOW phases are needed to identifypossible physical realizations.PACS numbers: 71.10.Fd, 73.22.Gk, 75.40.Cx, 78.30.-jEmail: [email protected] I. INTRODUCTION
Nakamura identified the bond order wave (BOW)phase [1] of the one-dimensional (1D) half-filled ex-tended Hubbard model [2] (EHM, Eq. 1). The keyfeatures are broken inversion symmetry, doubly degen-erate ground state (gs) and finite magnetic gap E m ina regular (equally spaced) array. Competition amongelectron delocalization t , on-site repulsion U >
V > V ≈ U/ t . The quantum transition to the chargedensity wave (CDW) phase at large V is first order[4, 5] for U > U ∗ ≈ t , continuous for U < U ∗ . TheBOW/CDW boundary is at V c ( U ) for U < U ∗ , whilethe boundary to the spin-fluid phase with E m = 0 isat V s ( U ) < V c ( U ). Other half-filled 1D Hubbard mod-els with spin-independent interactions have narrow BOWphases between the spin-fluid and CDW phases [10].Theoretical studies of the EHM have focused on thequantum phase diagram of an extended array withoutreference to possible physical realizations. Even for mod-els, however, 1D instabilities and finite temperature mustbe addressed. The principal goal of this paper is to char-acterize the BOW phase of the EHM at U = 4 t , a typ-ical choice for intermediate correlation, with special at-tention to the consequences of gs degeneracy and finite E m . Coupling to Peierls or Holstein phonons requires gsderivatives in the adiabatic (Born-Oppenheimer) approx-imation, and finite temperature properties are neededto assess physical realizations. It is clearly desirable tounderstand BOW phases prior to specific applications.Broadly similar properties are expected for other U orother 1D models with spin-independent interactions.Quasi-1D organic molecular crystals and conjugatedpolymers have strong electron-phonon (e-ph) couplingand electron-molecular-vibration (e-mv) coupling in ad-dition to intermediate correlation. The Su-Schrieffer-Heeger (SSH) model of polyacetylene has linear e-ph cou-pling and topological solitons as elementary excitations[11, 12]. The Peierls instability is driven by e-ph cou-pling and has spectacular e-mv consequences in infraredspectra when inversion symmetry is broken. Conjugatedpolymers and organic ion-radical salts have been a play-ground for 1D Hubbard models [12–17] with e-ph ande-mv coupling, variable electron or hole filling, degener-ate or nondegenerate gs, and either segregated or mixedstacks. The bandwidth is 4 t ≈
10 eV in polymers and4 t ≈ π -stacks, with comparable intermediate t/U .The Peierls instability of Hubbard-type models is a richseparate topic [18]. We are not aware of work on either e-ph or e-mv coupling in the BOW phase. As shown below,finite E m and degenerate gs lead to electronic solitons ina rigid regular array that nevertheless resemble SSH soli-tons. With suitable modification, extensive SSH analysis[12, 14] can be applied to models with a BOW phase.We recently proposed that a BOW phase is realizedin Rb-TCNQ(II), the second polymorph of a tetra-cyanoquinodimethane salt [19, 20]. The evidence is a100 K crystal structure ( P¯1 ) with regular TCNQ − stacksat inversion centers, negligible spin susceptibility below140 K that indicates a large E m , and infrared spectrathat demonstrates broken electronic inversion symmetry.Broken C i symmetry and finite E m in a regular 1Darray are precisely the signatures of a BOW phase [1].Large E m indicates proximity to the CDW boundary,and we will so choose V in the EHM at U = 4 t . We haveRb-TCNQ(II) and alkali-TCNQs in mind, but do notmodel them explicitly beyond invoking solitons for thetemperature dependence of infrared spectra. In additionto values of microscopic parameters, interactions be-tween chains must be addressed in actual models alongwith the Coulomb interactions and transfer integrals fordifferent stacking motifs. We consider electronic proper-ties of the EHM at U = 4 t with linear coupling to Peierlsand Holstein phonons. Vibrational degrees of freedomare introduced as needed in the adiabatic approximation.The EHM describes electronic degrees of freedom in aregular 1D array [2] H el = N X p =1 ,σ − t ( a † p,σ a p +1 ,σ + h.c )+ N X p =1 ( U n p ( n p − / V n p n p +1 ) (1)where t = 1 is the unit of energy, h.c. is the hermitianconjugate, a † p,σ ( a pσ ) creates (annihilates) an electronwith spin σ at site p , and n p is the number operator. H conserves total spin. The gs is a singlet ( S = 0) for U , V ≥ N . The half-filled band with N electrons and N sites has electron-hole (e-h) symmetry J = ± C i at sites that welabel as σ = ±
1. The correlated many-electron basisincreases as ≈ N at large N . Valence bond (VB)methods [21, 22] at present yield exact results up to N = 17 for low-energy states and the full spectrum upto N = 10. Nakamura identified [1] the BOW phaseusing field theory, symmetry arguments and numericalresults up to N = 12.To introduce the EHM phase diagram at 0 K, wedefine three threshold excitations from the singlet gs:the magnetic gap E m to the lowest triplet, the gap E J to the lowest singlet with opposite J , and the gap E σ to the lowest singlet with opposite C i . Increasing V atconstant U , t drives the system from a spin-fluid phaseto a CDW phase. Fig. 1 shows the evolution of E m , E J and E σ with V for U = 4 t , N = 12. The BOWphase spans V s < V < V c , where V s ( N ) and V c ( N ) aredefined by the excited-state crossovers E σ ( N ) = E m ( N )and E σ ( N ) = E J ( N ), respectively. Table I lists the V/t -0.200.20.40.60.8 E /t E J E m E σ E J E m V s E σ V V c EHM U=4t N=12 PBCBOWSpin Fluid CDW
FIG. 1: Excitation thresholds and crossovers of the 12-siteextended Hubbard model, Eq. 1, with periodic boundary con-ditions. E m , E σ and E J are energy gaps to the lowest tripletand the lowest singlets with opposite C i and e-h symmetry,respectively. V s is defined by E σ = E m , V by E σ = 0 and V c by E σ = E J .TABLE I: Excitation thresholds and crossovers of the ex-tended Hubbard model, Eq. 1, with N sites, U = 4, and t N = ± N = 4 n , 4 n + 2. N V s ( E σ = E m ) V ( E σ = 0) V c ( E σ = E J )8 1.8094 2.0597 2.159210 1.8190 2.0726 2.162412 1.8297 2.0840 2.164514 1.8311 2.0925 2.165116 1.8452 2.0981 2.1653 remarkably weak size dependence of V s ( N ) and V c ( N )for N = 4 n with periodic boundary conditions (PBC, t N = 1) and N = 4 n + 2 with antiperiodic boundaryconditions ( t N = − V s boundaryof a frustrated spin chain has been found [23] using E m ( N ) = E σ ( N ) up to N = 24. Multiple methods yieldthe boundary V c ( N ) of the CDW phase [7, 10]. Thepoints V = V ( N ) in Table I or Fig. 1 correspond to E σ = 0. They mark a gs degeneracy that is central toour discussion, where broken symmetry gs are readilyconstructed.The paper is organized as follows. Broken inversionsymmetry and elementary excitations are treated exactlyin Section II for finite N . Electronic solitons, both spinand charge, are found in regular chains with open bound-ary conditions (OBC, t N = 0) and compared to SSHsolitons. Section III deals with linear coupling to molecu-lar (Holstein) and lattice (Peierls) vibrations. The Berryphase formulation of polarization is applied to the in-frared activity of molecular vibrations when C i symme-try is broken. The Peierls instability of the BOW phaseis contrasted to the SSH model within the limitations ofan adiabatic approximation. The magnetic gap E m andspin susceptibility χ M of the BOW phase are obtainedin Section IV for large E m close to the CDW instabil-ity. Large E m reduces the thermal population of excitedstates and opens a new regime in which spin solitonsgovern χ M . The discussion in Section V summarizes theconsequences of broken symmetry and finite E m such asthe temperature dependence of the infrared intensity orof χ M . We briefly mention extensions to BOW phases ofrelated models. II. BROKEN SYMMETRY AND ELEMENTARYEXCITATIONS
We discuss the EHM, Eq. 1 with U = 4, t = 1, usingexact results for finite N with PBC( t N = 1) for N = 4 n and t N = − N = 4 n + 2. The gs kinetic energy isgiven by the bond orders p n of successive sites2 p n = h ψ | X σ ( a † n,σ a n +1 ,σ + h.c ) | ψ i . (2)Broken C i symmetry in the BOW phase leads to p n = p n − , while broken e-h symmetry in the CDW phaseleads to different electron count n p = n p − in the evenand odd sublattice. At constant U and t , the order pa-rameter B ( V ) of the BOW phase is B ( V ) = | p n ( V ) − p n − ( V ) | . (3) B ( V ) is large between sites 2 n , 2 n − n , 2 n + 1 in the other.Since both E m and E σ vanish rigorously in the spin-fluid phase with V < V s , finite gaps in Fig. 1 are dueto finite N . The dashed line E m ( V ) − E σ ( V ) is an ap-proximation for opening the magnetic gap. Similarly, thedashed line E J ( V ) − E σ ( V ) approximates the closing ofthe e-h gap at V c . Finite N limits E σ = 0 to points V ( N ) in Fig. 1 and Table I. At gs crossovers, we con-struct broken-symmetry states | ψ ± ( V ) i = ( | ψ σ =1 i ± | ψ σ = − i ) / √ N provide direct access to BOW systems that is not available for Monte Carlo, which yieldsthe gs energy, or density matrix renormalization group(DMRG) calculations. DMRG with OBC breaks C i symmetry for even N and returns a nondegenerate gs.We are not aware of DMRG with PBC that conserves C i and gives the energy and gs in both the σ = ± V ( N ) where E σ = 0 for finite N than to demonstrate a degenerate gsover the interval V s < V < V c in the extended system.On the other hand, we need the gs degeneracy at U =4 t beyond the single point V ( N ). To do so, we add thefollowing perturbation to Eq. 1 H ( J ) = J X p ~S p .~S p +2 . (5) H ( J ) acts on second neighbors with one electron each.The upper panel of Fig. 2 shows how V ( J ) scans E σ = 0across the BOW phase at U = 4 t , lowering V for antifer-romagnetic J > J <
0. The interval J = ± .
15 in Fig. 2 is sufficient to enforce E σ ( J ) = 0between V /t = 2 . N = 12and U/t = 4. Strictly degenerate gs in the σ = ± E σ V/t B ( V ) = | p + - p - | J =0.15 0.1 0.05 0 -0.05 -0.1 -0.15
8 12N=16N=12
FIG. 2: (top panel) Energy E σ of the lowest singlet withopposite C i symmetry in a 12-site EHM with U = 4 t and PBCas a function of V in Eq . 1 and the indicated J values in Eq.5. (bottom panel) Order parameter B ( V , J ) = | p + − p − | inEq. 3 at E σ = 0 for N = 12 (closed symbols) and B ( V ) for N = 8, 12, 16 for J = 0 (open symbols). The gs in Eq. 4 and the order parameter B in Eq.3 hold for V ( J ) where E σ = 0. The lower panel ofFig. 2 compares B ( V ( J )) at the indicated J valuesfor N = 12 with B ( V ) based on J = 0, when | ψ ± i in Eq. 4 are linear combinations of functions thatare not quite degenerate. B ( V, N ) depends weaklyon E σ and N in this interval. Of course, H ( J )also perturbs other properties, at least slightly, andbecomes a strong perturbation for V < t . We are inter-ested in V close to V c , where E σ is small and E m is large.Bond orders p ± ( V, N ) and B ( V, N ) are well approxi-mated near V ≈ V ( N ) without invoking J , and resultsin this Section are based on J = 0. For example, U = 4 t , V = 2 . t , N = 16 returns p + = 0 . p − = 0 . B = 0 . V c ,the average p = 0 .
574 is within 10% of 2 /π , the valuefor a tight binding (H¨uckel) band of free electrons with U = V = 0 in Eq. 1. The t term dominates near V c where U and V almost cancel.The SSH model [11, 12] is a tight-binding band withlinear e-ph coupling α = ( dt/du ) to the Peierls phonon,the optical mode of the 1D chain, and an adiabatic ap-proximation for a harmonic lattice. Its gs is dimerized,with t n = − (1 − δ ( − n ) along the stack. Bond ordersof the infinite chain are readily found analytically, with p + ( δ ) = 0 . p − ( δ ) = 0 .
452 and B ( δ ) = p + − p − =0 .
354 at δ = 0 .
10 that matches B ( V ) of the finite EHMabove. The origin of broken C i symmetry is quite dif-ferent: e-ph coupling in SSH, correlation U , V in EHM.Moreover, B ( V ) opens at V s ≈ . t and is almost con-stant for V > . t before vanishing abruptly at V c ≈ . t while B ( δ ) is monotonic in δ . The electronic gs are nev-ertheless similar and previous discussions of topologicalsolitons or domain walls between regions with opposite B can be applied to BOW phases [12, 14, 15, 26].We start with electronic domain walls and return laterto the Peierls instability. Domain walls are modeledas in SSH with odd N in Eq. 1 and t N = 0. Weretain PBC for V to minimize end effects. N e = N electrons correspond to a neutral soliton with S = 1 / N ± S = 0.Spin-charge relations are reversed just as in SSH. Welabel sites from r = 0 at the center to r = ± ( N − / N = 17exactly using a symmetry-adapted valence bond basisto compute p n in Eq. 2. As shown in Fig. 3 for N = 15 and 17, the bond orders of spin and chargesolitons are symmetric about the center. They areslightly larger at the chain ends than p − = 0 . p + = 0 . V = 2 . t .This end effect for t N = 0 also appears for even N and alternating bond orders with largest p + at either end.As context for these results, we note that Eq. 1 with -6 -4 -2 0 2 4 6 8 r B ond o r d e r , p r N e =N, N=17N e =N, N=15N e =N+1, N=17N e =N+1, N=15 U=4t, V=2.05t
FIG. 3: Bond orders p r in Eq. 2 of a spin soliton (closesymbols) or charge soliton (open symbols) in an N -site EHMwith t N = 0, U = 4 t and V = 2 . t in Eq.1. U = V = 0 and odd N can be read as a H¨uckel modelof an alternant hydrocarbon with a nonbonding orbitalat ǫ = 0 that is empty in the cation, singly occupied inthe radical and doubly occupied in the anion. Since the ǫ = 0 orbital has nodes at every other site for any N , itdoes not contribute to bond orders that are consequentlyequal for neutral and charge solitons in the SSH model.The correlated model has equal p r for charge solitonsthat differ slightly from the spin soliton. In every case,bond orders are p + at the ends and reverse smoothly inbetween.Next we compute spin densities ρ r = 2 h S zr i for theradical and charge densities q r = 1 − n r for ions.The upper panel of Fig. 4 has ρ r for N = 15 and17 while the lower panel has q r for the cation; theanion charges are − q r by e-h symmetry. The spinor charge density vanishes in alternant hydrocarbonswhere the nonbonding orbital has nodes, at odd r for N = 4 n + 1 and even r for N = 4 n −
1. As seen inFig. 4, correlation [26] generates small negative ρ r atthese sites and also small q r of opposite sign. The spinor charge density is large in the middle and decreases atthe ends. Electronic solitons in the BOW phase differin this respect from a regular tight-binding band, whichhas equal ρ r or q r at every other site for any odd N .To mimic the BOW results in Fig. 4 with SSH solitons,finite δ ≈ .
05 is needed for ρ r and δ ≈ .
08 for q r since large δ gives faster decrease. Electronic solitons inFigs. 3 and 4 connect regions with ± B in a regular chain.Domain walls or topological solitons are the elemen- -0.1500.150.3 < S r z > N=17N=15 -8 -6 -4 -2 0 2 4 6 8 r -0.200.2 q r = - n r U=4t, V=2.05tN e =N-1N e =N(a)(b) FIG. 4: (a) Soliton spin density 2 h S zr i and (b) soliton chargedensity q r = 1 − n r of an N -site EHM with t N = 0, U = 4 t and V = 2 . t in Eq. 1. tary excitations of the BOW phase. They resemble SSHsolitons in an appropriately dimerized lattice. The BOWphase has different energy for creating a pair of spin orcharge solitons. We take E ( N, N e ) as the gs energy ofEq. 1 for even N , N e electrons and t N = 0. The forma-tion energy of a pair of spin solitons ( S = 1 /
2) is2 W S ( N ) = E ( N + 1 , N + 1) + E ( N − , N − − E ( N, N ) . (6)In the limit of large N , parallel spins are the lowest tripletwith 2 W S = E m , the magnetic gap for PBC. The resultsfor N = 8, 12 and 16 in Table II are for V = V ( N ) inTable I, where E σ = 0, but the V dependence is weak.The formation energy for a pair of charge solitons is2 W C ( N ) = E ( N + 1 , N ) + E ( N − , N ) − E ( N, N ) . (7)The cation, anion and neutral system are singlets withdifferent charge. The natural comparison is to the chargegap,2 E C ( N ) = E ( N, N + 1) + E ( N, N − − E ( N, N ) . (8)For large N , E C ( N ) is the energy of separated ionradicals. As seen in Table II, 2 W C approaches E C frombelow, as expected since charge solitons are singletswhile the ions in E C are doublets. On the other hand,2 W C is larger than E J because V >
TABLE II: Representative EHM energies in units of t for U = 4 and V = V ( N ) in Eq. 1.Energy, t = 1 N = 16 N = 14 N = 122 W S ( N ) Eq . W C ( N ) Eq . E C ( N ) Eq . E m E J E a a Third lowest singlet sites with n = 0 (hole) and n = 2 (electron).Valence bond methods yield low-energy excitations inevery symmetry sector with fixed S , J and σ . Exceptfor total wave vector k = 0 or π , degeneracy in ± k isexpected and found. Finite-size effects increase with en-ergy and it becomes progressively more difficult to ex-tract more than 3-4 states for large N. The entries inTable II are from a much larger set. An “effective” δ eff ≈ . − .
10 is inferred from the order parameter B ( V ) or from spin or charge solitons in the EHM with U = 4 t and 2 . < V /t < V c . While BOW-phase resultsare for a regular stack, they are naturally related to SSHresults with δ < .
10 that is considerably smaller than δ = 0 .
18 based on the optical gap of polyacetylene[11].
III. COUPLING TO HOLSTEIN AND PEIERLSPHONONS
The SSH model invokes linear e-ph coupling, α = ( dt/du ) , to characterize the Peierls instabilityand elementary excitations of a half-filled tight-bindingband. Linear coupling to molecular vibrations is thebasis for interpreting polarized infrared spectra of π -radical stacks. The operators a † pσ , a p,σ create, anni-hilate electrons with spin σ in the lowest unoccupiedmolecular orbital of TCNQ, with equal energy ∆ = 0and n p = 1 for a TCNQ − stack. Charge fluctuations[27] modulate ∆ and illustrate linear Holstein coupling g n to the n th totally symmetric (ts) molecular vibration.When C i symmetry is broken, ts modes become stronglyIR allowed by borrowing intensity from the opticalcharge-transfer excitation and they are polarized alongthe chain. Accordingly, the appearance of ts modes inpolarized IR yields detailed information [28–30] aboute-mv coupling constants g n and has been widely usedto infer dimerization. Charge fluctuations break e-hsymmetry, and strict degeneracy at V ( N ) is criticalbecause ∆ modulation is a small energy.The Berry-phase formulation [31, 32] of polarizationmakes possible improved vibronic analysis of extendedsystems. It is directly applicable to quantum cell models[33] such as the EHM. The polarization P per unit chargeand unit length is a phase [32, 33] P N = 2 πN Im ( lnZ N ) (9) Z N = h ψ | exp (2 πiM/N ) | ψ i (10) M = N X p =1 p ( n p −
1) (11)where | ψ i is the exact gs of an N -site supercell and M is the conventional dipole operator for a regulararray with unit spacing. | Z N | = 0 is a metallic pointat which P is not defined [31]; it corresponds to V c for U < U ∗ and a continuous CDW transition [10]. TheEHM has real Z and P = 0 by either C i or e-h symmetry.The required generalization of Eq. 1 for Holstein cou-pling in the adiabatic approximation for molecular sitesis [33] H (∆) = H el + ∆ N X p =1 ( − p n p . (12)Finite ∆ breaks e-h symmetry, but not C i symmetry.The gs in the BOW phase is now | ψ ± ( V, ∆) i with J in Eq. 5 chosen to have E σ = 0 for ∆ = 0. Strictlydegenerate σ = ± E J in Fig. 1whose energy in the BOW phase decreases rapidly withincreasing V . We again vary V at U = 4, t = 1. Asidefrom a multiplicative constant, the IR intensity goes as[34] I IR ( V ) = (cid:0) ∂P ( V, ∆) ∂ ∆ (cid:1) . (13)Charge fluctuations give a finite derivative at ∆ = 0.The IR intensity in Eq. 13 is purely electronic. It canbe partitioned among ts modes as discussed [28–30] insystems where dimerization breaks C i symmetry.We compute the derivative in Eq. 13 in the BOWphase at V = V ( J ). The imaginary part of Z N isinitially proportional to ∆, as shown in the inset of Fig.5. I IR ( V ) increases rapidly with V up to V c , as expectedfor decreasing E J , and vanishes abruptly in the CDWphase where C i symmetry is restored and e-h symmetryis broken. The corresponding band result [33, 34] for( ∂P ( δ, ∆) /∂ ∆) increases as 1 /δ and diverges at themetallic point ∆ = δ = 0 where, however, P is notdefined. Since the real part of Z N ( V ) = 0 at V = V c ( N ),the finite system also has divergent ( ∂P ( V, ∆) /∂ ∆) at V c ( N ) and the BOW phase again resembles a bandwith small δ , with two major differences. First, I IR ( V ) V/t I I R ( V ) = ( ∂ P ( V , ∆ ) / ∂ ∆ ) -0.04 0 0.04 ∆ -0.100.1 P ( V , ∆ ) V=V U=4t, N=12V=V (J ) V c FIG. 5: Infrared intensity I IR ( V ), Eq. 13, of molecular vi-brations due to ∆ in Eq. 12 and broken electronic symmetryin the BOW phase of the EHM. The inset is the polarization P (∆) in Eq. 9. I IR vanishes by symmetry for V > V c or < V s . increases with V up to V c while B ( V ) is almost constant.Second, δ eff for IR intensity decreases with increasing V up to V c while δ eff for B ( V ) increases with V from V s . The SSH model has a single band gap of 4 δt insteadof the EHM’s multiple threshold excitations in Table Iand Fig. 1. For example, the metallic point V c has large E m ≈ . t .Domain walls introduce inversion centers that reduce I IR between regions with opposite B . More quantita-tively, we again consider systems with odd N , V = 2 . t N = 0. We break C i symmetry with ± ∆ at sites ± r from the center and evaluate ( ∂P/∂ ∆ r ) . Since ∆ inEq. 12 applies all sites in Fig. 5, we scale P ( V, ∆ r ) by N/ P ( V, ∆ r )at V = 2 . t for ∆ r = 0 .
05, with P = 0 at r = 0 andscaled P ≈ .
06 at the ends of a spin soliton. The samepattern is found for charge solitons (data not shown)with about 50% higher intensity. In either case, I IR ( V )at the ends is an order of magnitude less than the V ( J )values for degenerate gs in Fig. 5. This finite-size effectcan be traced to higher E J in the radical. We did notsolve N = 17 since finite ∆ doubles the dimensions ofthe many-electron basis. IR intensities are consistentwith a soliton and an inversion center between regionswith ± B ( V ).We consider next the Peierls instability. The SSH -6 -4 -2 0 2 4 6 r ( ∂ P ( V , ∆ r ) / ∂ ∆ r ) -6 -4 -2 0 2 4 6 r -0.0800.08 P ( V , ∆ r ) N=13U=4t, V=2.05t 15 FIG. 6: Infrared intensity I IR ( V ), Eq. 13, of molecular vi-brations due to ∆ in Eq. 12 in 13 and 15-site spin solitonswith V = 2 . t and t N = 0 in Eq. 1. The inset shows P ( V, ∆ r ) for ∆ r = 0.30 in Eq. 9 at site r , with r = 0 at thecenter and r = ± ( N − / model is the generic case for e-ph coupling α = ( dt/du ) in a harmonic 1D lattice with force constant k . Linearcoupling is retained in generalizations [18] that requirenumerical methods and include electron correlation,quantum fluctuations, nonadiabatic or 3D effects, spinchains and models with site energies. Quantum fluc-tuations are particularly important for small δ that iseasily reversed locally. Fluctuations reduce but do notwash out δ ≈ .
18 in the SSH model of polyacetylene [35].Now Eq. 1 reads H ( δ ) = H el + N X p =1 ,σ δ p t ( a † p,σ a p +1 ,σ + h.c )+ X p δ p / ǫ d . (14)with ǫ d = α /k and constant δ p = ( − p αu/k in the gs.The distribution δ p is subject to the constraint P p δ p = 0for a fixed chain length. The gs energy per site is ǫ T ( δ ) = ǫ ( δ ) + δ / ǫ d . (15)A minimum at δ = 0 implies a dimerized gs. A global adiabatic approximation has δ ( T ) > T < T P , thePeierls temperature, while a local adiabatic approxima-tion leads to domains walls between regions of opposite δ for T > ǫ T ( δ ) in Eq. 15. The gsof correlated models is unconditionally dimerized whenthe regular array has E m = 0 and χ d ( δ ) = − ( ∂ ǫ /∂δ )diverges at δ = 0. Examples [27] include the ionic phaseof organic charge transfer salts and the spin fluid phase ofHubbard models or Heisenberg spin chains. Dimerizationis conditional at ǫ d χ d (0) = 1 when χ d (0) is finite in theneutral phase of CT salts or the CDW phase of the EHM.Degenerate gs is different, as seen for N ǫ T ( δ ) in Fig. 7for V = V ( N ) in Table I and inverse stiffness ǫ d = 0 . ǫ ( δ ) goes as − B | δ | − χ d (0) δ /
2. The cusp B ( V ) isdue to degeneracy at δ = 0 while the energy gap E inthe singlet sector, listed in Table II, ensures finite χ d (0).Minimization of ǫ T yields δ eq = ± ǫ d B ( V ) / (1 − ǫ d χ d (0)) (16)with δ eq > δ > δ eq < δ <
0. Dimeriza-tion is unconditional and increases E m ( δ ) beyond E m (0)as shown in Fig. 7 for a vertical excitation. The gs cuspappears again because the lowest triplet is not degener-ate and hence evolves as δ . The slope of ( ∂E m /∂δ ) atthe origin is N B ( V , N ). It follows δ -0.2-0.100.10.20.30.40.5 E /t E m ( δ )-E m (0) E ( δ )-E (0)+N δ /2 ε d U=4t, V=V PBC N=88 FIG. 7: Ground state energy E in Eq. 15 and magnetic gap E m of an EHM with dimerization δ in Eq. 14 and ǫ d = 0 . that two domain walls that change the δ p patternwithout changing P p δ p lead to lower E m , and reduced E m are readily found for N = 12 or 16. Generalizationof 2 W S ( N ) in Eq. 6 is more suitable for spin solitonsbetween regions with opposite δ eq . Such simulations arebeyond the scope of the present study.An adiabatic approximation gives a dimerized gs with ± δ eq that has to be relaxed locally for solitons. X-raydetection requires that δ eq not be too small compared tothe zero point amplitude, h δ i / , of the Peierls phonon ~ ω P in Eq. 15, h δ i = ǫ d ~ ω P / t (17)Typical values [36] of ~ ω P ≈ cm − , ǫ d ≈ .
30 and t ≈ cm − in organic stacks return h δ i / ≈ . T = 0 that increases with T . Finite E m and small δ eq minimize contributions from thermal excitations, incontrast to SSH or correlated Peierls systems in whichsmall δ eq necessarily implies small E m that vanishes at δ = 0. An arbitrarily small perturbation selects one ofthe degenerate gs of a BOW phase at 0 K, but there isno long-range order in 1D at finite T . Solitons lower thefree energy and result in local adiabatic approximations.At low T , the mean separation R ( T ) between solitons islarge compared to their widths 2 ξ . Each soliton can thenbe centered on any of R ( T ) sites. The soliton density inan extended system with R > ξ and N → ∞ is ρ ( T ) = R ( T ) N = exp ( − Wk B T ) (18)where k B is the Boltzmann constant and W = W S forspin solitons or W C for charge solitons. The N = 17results in Fig. 4 are not sufficient for estimating 2 ξ , butindicate 2 ξ >
30 for spin solitons and 2 ξ ≈
30 for chargesolitons, consistent with the SSH estimate [11] of 2 ξ ≈ δ = 0 .
18. Since W S < W C , spin solitons arethermally accessible and the condition R ( T ) > ξ holdsup to ρ ( T ) ≈ IV. SPIN SUSCEPTIBILITY
The magnetic gap E m opens at the boundary V s ofthe spin-fluid and BOW phases, and it does so veryslowly at a Kosterlitz-Thouless transition [1]. DMRGwith PBC provides an independent calculation [10]showing that E m ( V ) opens at V s = 1 . t for the EHMwith U = 4 t . The gap is only 0 . t at V = 2 . . t at V = 2 .
10 and 0 . t at V = 2 . E m near the boundary V c of the BOW and CDWphases in Fig.1. As noted above, V c is a metallicpoint with bond orders p ( V ) in Eq. 2 close to 2 /π ,the band limit. The band limit returns E m = 0 for N = 4 n , when the degenerate orbitals at ǫ = 0 are halffilled, and E m = 4 tsin ( π/N ) for N = 4 n + 2. Openboundary conditions with t N = 0 have intermediate E m . The same pattern is seen in Fig. 8 for the EHMat V = 2 .
5, with smallest E m ( N ) for N = 4 n , t N = 1and N = 4 n + 2, t N = −
1. Quite unusually, E m increases with N at V = 2 .
5. DMRG calculations [10]of E m with PBC have minimum E m ( N ) at N ≈ V = 2 .
2. Exact E m ( N ) in Fig. 8 with t N = 0, ± V = 2 . N . Rapidly increasing E m at V c is also seen [10] for other potentials and hasimportant implications for modeling the spin suscepti-bility. In particular, larger N is not automatically better.The molar magnetic susceptibility χ M ( T ) allows quan-titative comparisons [13, 25, 37] for organic ion-radicalsolids with small spin-orbit coupling and g -factors closeto the free-electron value, g = 2 . χ M ( T )to be the spin susceptibility after standard corrections fordiamagnetism and impurities. The full spectrum of Eq.1 is now required, and charge degrees of freedom vastlyincrease the number of states. We extend exact resultsfor χ M ( T ) up to N = 10. The partition function of theEHM in Eq. 1 with even N is Q N ( T ) = N/ X S =0 X r =1 (2 S + 1) exp ( − E Sr ( N ) /k B T ) (19)The singlet gs is the zero of energy, E ( N ) = 0, and E Sr ( N ) refers to the state r with spin S . The molar spinsusceptibility is χ M ( T, N ) = N A g µ B tN Q N (cid:0) tk B T (cid:1) N/ X S =0 X r =1 S ( S + 1) × (2 S + 1) exp ( − E Sr ( N ) /k B T ) (20)where µ B is the Bohr magneton and N A is Avogadro’snumber. Finite E m = E leads to χ M (0) = 0, incontrast to finite χ M (0) for V < V s in the spin-fluidphase. We set E σ = 0 in the BOW phase and considerfinite N to be a coarse-grained approximation of anextended system with a dense spectrum for E ≥ E m .Figure 9 shows χ M ( T ) up to k B T = t for N = 8 withPBC ( t N = 1) and N = 10 with t N = − V = 2 . T , as has long been recognized in spin chains. Inaddition, the broad plateau around k B T ≈ t dependsweakly on V , which simply reflects the narrowness ofthe BOW phase. By contrast, the T ≈ E m ( V ) and is very sensitive to V at constant U , t . Although t ≈ K is a small electronic energy,experiment is limited to much lower T . The inset ofFig. 9 expands the relevant range. Finite E m suppresses χ M at low T after which χ M is almost linear in T . TheBOW phase of a frustrated spin chain [38] without chargedegrees of freedom also has an almost linear χ M ( T ) and E m /t V=2.2t
V=2.5t
V=2.2t DMRGV=2.0tV=2.1t t = 0 t = -1 1 1 -1-1t = 1t = 1 -1-1 -111 1-1t = 0 11 FIG. 8: Magnetic gap E m of the N -site EHM, Eq. 1, with U = 4 t and the indicated V’s. Boundary conditions with t N = 0 and ± E m around N ≈ V = 2 . t . a broad maximum.A BOW phase with E m > χ M ( T ) ≈ T ≈ T followed by a linear increase for T > T . Forexample, Rb-TCNQ(II) has [19] T ≈
140 K that roughlyfixes t ( V ) in the inset to Fig. 9 and completely specifies χ M ( T ) of the EHM with U = 4 t . Although χ M ( T )increases almost linearly to T = 300 K , the magnitudeat 300 K rules out the V > . t curves in the insetwhile the V = 2 .
10 curve with k B T /t ≈ .
15 at 300 Kfails for T . An EHM with U = 4 t is not quantitativefor Rb-TCNQ(II). Improved fits are possible for small δ eq ≈ .
01 but such modeling also entails variation ofCoulomb interactions.The EHM boundary V s between the spin fluid andBOW phases has been difficult to model and has beenreported [5] to be as high as V s ≈ . t at U = 4 t .Since χ M (0) is finite in the spin fluid phase and E m opens slowly for V > V s , very different χ M ( T ) are thencalculated in the BOW phase [5]. The evolution of E m ( V, N ) in Fig. 8 with N and DMRG results supportthe original estimate [1] of V s ≈ . t at U = 4 t based onexcited-state crossovers in Table I. Rapidly increasing E m ( V ) as V approaches the metallic point V c is foundin the BOW phase of related Hubbard models [10].Large E m is needed for the Rb-TCNQ(II) susceptibilityas well as for K and Na-TCNQ at high T where X-ray k B T/t χ t/ ( N A µ Β g ) k B T/t χ t/ ( N A µ Β g ) U=4tV=2.10tN=10 8 108V=2.15tV=2.20tN=8 FIG. 9: Temperature dependence of the molar spin suscepti-bility χ M ( T ) of N -site EHM with U = 4 t and V in Eq. 1 Thecurves in the inset are for the same U , V . structures [39] indicate regular TCNQ − stacks. V. DISCUSSION
Rice [28] recognized the possibility of measuring e-mvcoupling constants g n from polarized IR spectra when C i symmetry is broken on dimerization. Several groups[29, 30] extended the procedure to extracting transfer-able g n for selected π -donors and π -acceptors. K-TCNQwas a prime example [40] of a crystal with dimerizedTCNQ − stacks at 300 K. Polarized mid-IR spectrashow coupled ts modes that are shifted to the red fromthe corresponding Raman transitions [29, 40]. PowderRb-TCNQ(II) has virtually identical IR transitions [41]whose polarization along the stack has been confirmed insingle crystals [42], but strikingly different temperaturedependence. Raman spectra and g n of TCNQ − areexpected to be the same, since solid-state perturbationsare usually small [29, 30].The 100 and 295 K crystal structures of Rb-TCNQ(II)decisively indicate [19] a regular stack of TCNQ − at in-version centers and interplanar separation R = 3 . A .The 295 K structure is in excellent agreement withprevious data [43], and the triclinic space group P ¯1 isretained at 100 K . Low R factors and examination ofthermal ellipsoids at 100 K place a conservative limit ondimerization of R + − R − < . A [19]. Yet negligibly0small χ M ( T ) below 140 K implies a large E m and IRdata indicates broken electronic C i symmetry. Broken C i symmetry in a BOW phase accounts naturally forlarge E m near the CDW boundary and for IR intensityat 0 K. We take up the different temperature dependenceof the Rb salt.The IR intensity in Eq. 13 goes as ( ∂P ( V, ∆) /∂ ∆) [20, 27] and increases with V in the BOW phase as shownin Fig. 5. The intensities for a spin soliton in Fig. 6 aresmaller and vanish at the center. Each soliton introducesa C i center between regions with ± B ( V ) and reduces( ∂P ( V, ∆) /∂ ∆) over ≈ ξ sites. We approximate thetemperature dependence as I IR ( T ) I IR (0) = 1(1 + 2 ξρ S ( T )) . (21)Here ρ S ( T ) is the spin density given by χ M ( T ) /χ C ,where χ C = N A g µ B / k B T is the Curie susceptibility.The intensity is 50% lower at 2 ξρ S ( T ) = 1 where spinsolitons overlap. Eq. 21 relates two measured quantitiesthrough 2 ξ , as shown in Fig. 10 for 2 ξ = 60, whichis in the expected range. The χ M ( T ) data is for Rb-TCNQ(II) from ref. [19] and gives ρ S ( T ). The intensityratio of the 722 cm − mode is from Fig. 2 of ref. [44],with open and closed symbols on cooling and heating.Similar T dependence is seen for other mid-IR modesof crystals [42]. The fit supports a BOW phase inter-pretation and electronic solitons rather than a specificmicroscopic model or parameters. A microscopic modelmust account for χ M ( T ) in addition to the intensity ratio.We have examined the BOW phase of the EHM atintermediate correlation U = 4 t by direct solution of Eq.1 for finite N . We used degenerate gs at V = V ( N ) inTable I to break inversion symmetry in finite systemsand varied J in Eq. 5 to scan V ( N ) over the BOWphase. Exact degeneracy enforced by J turns out to beimportant for e-mv coupling to Holstein phonons but notfor the order parameter B ( V ) = p + − p − in Eq. 3. Wecompared BOW properties due to electronic correlationin a regular 1D chain to the SSH model of a dimerizedband. SSH results for topological solitons carry over formany aspects of spin and charge solitons in the BOWphase, albeit with different δ eff for different properties.We focused on the consequences of a degenerate gsand finite E m , especially large E m close to the metallicpoint V c . Broadly similar results are expected in BOWphases of other Hubbard-type models with intermediatecorrelation.In the adiabatic approximation for the lattice, lineare-ph coupling generates a dimerized gs in both theBOW phase of the EHM and the SSH model, but theyare different. The SSH model has a standard Peierlstransition at T P . Thermal population of excited statesstabilizes the regular array for T > T P , and low T P T I I R ( T ) / I I R ( ) CoolingHeatingEq.21 with 2 ξ =60 χ M ( e m u ) x 10 -4 FIG. 10: Intensity ratio I IR ( T ) /I IR (0) of the 722 cm − mode of Rb-TCNQ(II) on cooling (open symbols) and heat-ing (closed symbols) from ref [44]. The solid line is Eq. 21with 2 ξ = 60 and spin density ρ S ( T ) from χ M ( T ) in ref. [19]. necessarily implies weak coupling or a stiff lattice.The Peierls transition of polyacetylene is far above itsthermal stability, and we are not aware of evidence for δ eq ( T ) variations up to ≈
400 K. Spin-Peierls systemsillustrate decreasing δ eq ( T ) up to T SP <
20 K thatcan be modeled in the adiabatic approximation [45].The BOW phase of the EHM at intermediate U = 4 t samples a different sector of parameter space, one inwhich substantial E m up to ≈ . t suppresses thermalexcitations. Spin solitons in Eq. 18 with W S ( δ ) givea small χ M ( T ) in this range, and E m = 2 W S remainsfinite at δ = 0. The BOW phase has novel aspects thatneed further study. The principal theoretical issues arequantum fluctuations or nonadiabatic phonons that maysuppress a sharp Peierls transition when e-ph couplingis weak. The experimental problem is to detect smalldimerization against a background of zero-point motions.The present discussion is limited to the BOW phaseof the EHM at U = 4 t . We have developed theconsequences of coupling to lattice phonons and tomolecular vibrations in the adiabatic approximation.Similar results are expected [10] for other quantumcell models with electron-hole symmetry and will beneeded to model physical systems with BOW phases,starting with Rb-TCNQ(II). Since alkali-TCNQ salts aresemiconductors, they have Coulomb interactions ratherthan a Hubbard U and are stabilized close to the CDWboundary by the 3D electrostatic (Madelung) energy1[24]. The Na and K-TCNQ salts have dimerizationphase transitions [39] with some 3D character since thecations also dimerize. The regular structure at high T has small χ M ( T ) that increases with T in a manner thatsuggests a BOW phase. Both π − radical organic stacksand conjugated polymers are quasi-1D systems whoseinitial modeling is without interchain interactions.In summary, we have characterized the BOW phase ofthe EHM with intermediate U = 4 t by exact treatmentof finite systems with degenerate gs at V = V ( N ). Theelementary excitations are electronic solitons, both spinand charge, in a regular array. Solitons in the correlatedBOW phase resemble the familiar solitons of the SSHmodel for e-ph coupling in a tight-binding band. Severalmeasures indicate an “effective” dimerization δ < . δ = 0 .
18 for the SSHmodel of polyacetylene. Charge fluctuations are coupledto molecular (Holstein) phonons that become IR active in the BOW phase due to broken C i symmetry. The T dependence of IR modes is consistent with spin solitonswhose width is 2 ξ ≈
60 lattice constants. The BOWphase is dimerized at 0 K in the adiabatic approximation,but gs degeneracy and finite E m lead to novel aspectsfor a possible Peierls transition. Previous discussionsof SSH solitons greatly facilitate analysis of the BOWphase. So have previous treatments of e-ph and e-mvcoupling in 1D Hubbard models for conjugated polymersand organic ion-radical or charge-transfer crystals. TheBOW phase of Hubbard-type models with intermediatecorrelation has some unique aspects that invite furtherstudy as well as features that are common to such models. Acknowledgments:
ZGS thanks A. Girlando for ac-cess to unpublished IR spectra and A. Painelli for dis-cussions of polarization and BOW phases. We thank theNational Science Foundation for partial support of thiswork through the Princeton MRSEC (DMR-0819860). [1] M. Nakamura, Phys. Rev. B , 16377 (2000); J. Phys.Soc. Jpn. , 3123 (1999).[2] J.E. Hirsch, Phys. Rev. Lett. , 2327 (1984); Phys. Rev.B ,[3] P. Sengupta, A.W. Sandvik and D.K. Campbell, Phys.Rev. B , 155113 (2002).[4] Y.Z. Zhang, Phys. Rev. Lett. , 246404 (2004).[5] S. Glocke, A. Klumper and J. Sirker, Phys. Rev. B ,155121 (2007).[6] S. Ejima and S. Nishimoto, Phys. Rev. Lett. , 216403(2007).[7] G. P. Zhang, Phys. Rev. B , 153101 (2003).[8] K-M Tam, S-W Tsai and D.K. Campbell, Phys. Rev.Lett. , 036408 (2006).[9] A. Sandvik, L. Barents and D.K. Campbell, Phys. Rev.Lett. , 236401 (2004).[10] M. Kumar, S. Ramasesha and Z.G. Soos, Phys. Rev. B , 035102 (2009).[11] W.P. Su, J.R. Schrieffer and A.J. Heeger, Phys. Rev. B , 2099 (1980); Phys. Rev. Lett. , 1698 (1979)[12] A.J. Heeger, S. Kivelson, J.R. Schrieffer and W.P. Su,Rev. Mod. Phys. , 781 (1988).[13] Z.G. Soos, Annu. Rev. Phys. Chem. , 121 (1974); Z.G.Soos and D.J. Klein, in Treatise on Solid State Chemistry ,Vol. III (ed. N.B. Hannay, Plenum, New York, 1976) p.689.[14] T.J. Skotheim, Ed.
Handbook of Conducting Polymers ,Vol. 2 (Marcel Dekker, New York, 1986).[15] D. Jerome, Chem. Rev. , 5565 (2004).[16] N.S. Sariciftci (Ed.),
Primary Photoexcitations in Conju-gated Polymers: Molecular Exciton Versus Semiconduc-tor Band Model , World Scientific, Singapore, 1997.[17] Z.G. Soos, D. Mukhopadhyay, A. Painelli and A.Girlando, in
Handbook of Conducting Polymers , Sec.Edit. (Eds. T.A. Skotheim, R. Elsenbaumer and J. R.Reynolds, Marcel Dekker, New York, P997) p. 165.[18] P. Sengupta, A.W. Sandvik and D.K. Campbell, Phys.Rev. B , 245103 (2003).[19] T.M. McQueen, D.M. Ho, C. Jimenez Cahua, R.J. Cava, R.A. Pascal, Jr., and Z.G. Soos, Chem. Phys. Lett. 475, (2009).[20] Z.G. Soos, M. Kumar, S. Ramasesha and R.A. Pascal,Jr., Physica B , 5353 (2010).[21] S. Ramasesha and Z.G. Soos, Int. J. Quant. Chem. ,1003 (1984); Z.G. Soos and S. Ramasesha, Phys. Rev. B Valence Bond Theoryand Chemical Structure , D.J. Klein and N. Trinajstic,Eds. (Elsevier, Amsterdam,1989), p. 81; S. Ramaseshaand Z.G. Soos, in
Theoretical and Computational Chem-istry , Vol. 10, D.L. Cooper, Ed. (Elsevier, Amsterdam,2002) p. 635.[23] K. Okamoto and K. Namura, Phys. Lett. A , 433(1992).[24] M. Kumar, S. Ramasesha, R.A. Pascal, Jr. and Z.G.Soos, Europhys. Let. , 37001 (2008).[25] S.A. Bewick and Z.G. Soos, Chem. Phys. , 60 (2006);S.A. Bewick and Z.G. Soos, J. Phys Chem. B , 18748(2006).[26] Z.G. Soos and S. Ramashesha, Phys. Rev. Lett. , 2374(1983).[27] A. Girlando, A. Painelli, S.A. Bewick and Z.G. Soos,Synth. Met. , 129 (2004).[28] M.J. Rice, Solid State Commun. , 93 (1979); M.J. Rice,N.O. Lipari and S. Strassler, Phys. Rev. Lett. 21, 1359(1977).[29] R. Bozio and C. Pecile, in Spectroscopy of Advanced Ma-terials , Adv. Spectrosc. Vol. 19, R.J.H. Clark and R.E.Hester, Eds. (Wiley, New York, 1991), p. 1.[30] Painelli and A. Girlando, J. Chem. Phys. , 5655 (1986);C. Pecile, A. Painelli and A. Girlando, Mol. Cryst. Liq.Cryst. 171, (1989).[31] R. Resta, J. Phys.: Condens. Matter , R625 (2002).[32] R. Resta, Phys. Rev. Lett. , 1800 (1998).[33] A. Painelli and Z.G. Soos, Chem. Phys. , 48 (2006);Z.G. Soos, S.A Bewick, A. Peri and A. Painelli, J. Chem.Phys. , 6712 (2004).[34] L. Del Freo, A. Painelli and Z.G. Soos, Phys. Rev. Lett. , 27402 (2002).[35] W.P. Su, Solid State Commun. , 497 (1982).[36] M. Masino, A. Girlando and Z.G. Soos, Chem. Phys.Lett. , 428-33 (2003); G. D’Avino, A. Girlando, A.Painelli, M-H. Lemee-Cailleau and Z.G. Soos, Phys. Rev.Lett. , 156407 (2007).[37] J.S. Miller, Ed. Extended Linear Chain Compounds , Vol.3 (Plenum, New York 1983); H.T. Diep, Ed.
FrustratedSpin Systems (World Scientific, Hackensack, N.J. 2003).[38] M. Kumar, S. Ramasesha and Z.G. Soos, Phys. Rev. B , 054413 (2010).[39] M. Konno and Y. Saito, Acta Cryst. B , 2007 (1975);M. Konno, T. Ishii and Y. Saito, Acta Cryst. B , 763(1977). [40] D.B. Tanner, C.S. Jacobsen, A.A. Bright and A.J.Heeger, Phys. Rev. B , 3283 (1977).[41] R. Bozio, I. Zanon, A. Girlando and C. Pecile, J. Chem.Soc. Faraday Trans. , 235 (1978).[42] A. Girlando, private communication (2008).[43] H. Kobayashi, Bull. Chem. Soc. Jpn. , 3669 (1981); I.Shirotani and H. Kobayashi, Bull. Chem. Soc. Jpn. ,2595 (1973).[44] R. Bozio and C. Pecile, J. Chem. Phys. , 3864 (1977).[45] Z.G. Soos and S.A. Bewick, Chem. Phys. Lett. , 210(2006); I.S. Jacobs, J.W. Bray, H.R. Hart, L.V. Inter-rante, J.S. Kasper, G.D. Watkins, D.E. Prober and J.C.Bonner, Phys. Rev. B14