Effects of Lattice and Molecular Phonons on Photoinduced Neutral-to-Ionic Transition Dynamics in Tetrathiafulvalene--Chloranil\abst
For electronic states and photoinduced charge dynamics near the neutral-ionic transition in the mixed-stack charge-transfer complex tetrathiafulvalene--chloranil (TTF-CA), we review the effects of Peierls coupling to lattice phonons modulating transfer integrals and Holstein couplings to molecular vibrations modulating site energies. The former stabilizes the ionic phase and reduces discontinuities in the phase transition, while the latter stabilizes the neutral phase and enhances the discontinuities. To reproduce the experimentally observed ionicity, optical conductivity and photoinduced charge dynamics, both couplings are quantitatively important. In particular, strong Holstein couplings to form the highly-stabilized neutral phase are necessary for the ionic phase to be a Mott insulator with large ionicity. A comparison with the observed photoinduced charge dynamics indicates the presence of strings of lattice dimerization in the neutral phase above the transition temperature. \kword photoinduced phase transition, neutral-ionic transition, electron-lattice interaction, electron-molecular-vibration coupling
Nonequilibrium electronic states have drawn attention as a field where novel electronic functions and properties are being searched for. Correlated electron and electron-phonon systems with rich electronic phase diagrams are promising because synergy can be exploited, which is inherent to respective orders at low temperatures. Among nonequilibrium phenomena, photoinduced phase transitions occur synergistically since a low density of photons relative to the density of molecules or atoms can alter the electronic phase realized in such a system. A variety of photoinduced phase transitions are known to proceed very rapidly. [1, 2]
From the viewpoint of controlling functions in nonequilibrium states, it is important to learn a means by which a photoinduced phase transition proceeds. The most fundamental information is on relevant interactions that stabilize the ground state. It is generally true that, even if a few interaction parameters are sufficient to describe equilibrium properties, they are still insufficient for describing nonequilibrium dynamics. Additional interactions can be crucial in determining whether and how the phase transition is photoinduced. Among organic materials that have been regarded as strongly correlated electron systems, some have been realized to have substantial electron-molecular-vibration (EMV) couplings that contribute to the stabilization of the ground state and photoinduced dynamics.
For instance, the charge order in the quasi-two-dimensional organic salt -(BEDT-TTF)I [BEDT-TTF=bis(ethylenedithio)tetrathiafulvalene] is basically driven by a long-range Coulomb interaction. [3, 4] Relatively weak Peierls couplings to lattice phonons or molecular rotations are important for understanding the interrelation between the altered charge distribution and the structural deformation [5, 6] as well as the long-time behavior of photoinduced dynamics. [7, 8] From the early-stage dynamics after photoexcitation, however, substantially strong Holstein couplings to molecular vibrations are evident in the Fano destructive interference between vibrations containing C=C stretching and correlated electrons’ motion.  These molecular vibrations stabilize the charge order to a substantial extent.
Most organic compounds that show electronic phase transitions have C=C bonds in their constituent molecules. Their vibration frequencies are well known to depend on the charged state of the molecule. It is therefore expected that EMV couplings are generally strong even though the electronic state has strong electron correlations in terms of intermolecular spin and charge correlation functions. In this study, we focus on the mixed-stack charge-transfer complex TTF-CA, whose thermal and photoinduced neutral-ionic transitions have intensively been studied both experimentally [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24] and theoretically. [25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39] We review the origin of the discontinuities in ionicity and optical conductivity during the neutral-ionic transition. Their quantitative understanding will become important for controlling photoinduced phase transition dynamics.
2 Discontinuities during Neutral-Ionic Transition
The importance of the intersite Coulomb repulsion for discontinuities during the neutral-ionic transition has been pointed out in ref. \citennagaosa_jpsj86b. Without coupling to lattice phonons, i.e., without dimerization, the spin gap at the phase boundary is zero on the ionic (i.e., regular-Mott-insulator) side and on the neutral (i.e., band-insulator) side in the limit of vanishing transfer integrals. Here, the creation of a neutral-ionic domain wall requires of energy, which allows for the phase separation of neutral and ionic domains at the boundary. The transfer integral introduces quantum fluctuations. As long as is smaller than , however, the creation energy remains finite. The discontinuous neutral-ionic transition is demonstrated by quantum Monte Carlo simulations. 
The ionic phase is known to have dimerization. The strongest candidate for its mechanism is due to the spin-Peierls instability.  A regular ionic phase is a one-dimensional half-filled paramagnetic Mott insulator, so that it can be described by a spin-1/2 antiferromagnetic Heisenberg chain. If the superexchange interaction is modulated by lattice phonons, the system is stabilized by dimerization that produces a finite spin gap. The dimerization in the ionic phase is demonstrated by quantum Monte Carlo simulations.  As lattice phonons are strongly dimerized, i.e., as the electron-lattice coupling is strengthened, the phase boundary is shifted because only the ionic phase is stabilized, which decreases the ionicity on the ionic side of the phase boundary. The discontinuity in ionicity is strongly reduced by the electron-lattice coupling. As a consequence, the character of the ionic phase as a Mott insulator is weakened, and the transition becomes close to a continuous Peierls transition.
For the discontinuous ionicity change, Painelli and Girlando have pointed out that both EMV couplings and intersite Coulomb repulsion are important and that they collaborate.  For finite systems, they perform exact diagonalization in the limit of a large on-site repulsion and the second-order perturbation theory with respect to the electron-lattice (i.e., Peierls) and EMV (i.e., Holstein) couplings. The electron number configuration is roughly described as 2020 in the neutral phase and as 1111 in the ionic phase, so that the Holstein couplings stabilize the neutral phase, which are in contrast to the Peierls coupling that stabilizes the ionic phase through dimerization.
As noted above, both of the Peierls and Holstein couplings are suggested to be necessary for the reproduction of ionicity and optical conductivity near the neutral-ionic phase boundary. In spite of this, the Holstein couplings have often been ignored for the following reason. Without the Peierls coupling, the bond-charge densities cannot be dimerized, and the ionic phase cannot have the ferroelectric order. In contrast, the Holstein couplings only enlarge the site energy difference between the neighboring orbitals in the neutral phase, so that they can be absorbed into the renormalized site energy difference in the adiabatic limit. The situation is quite similar to a phase transition described within the mean field theory in the sense that the Holstein couplings stabilize the neutral phase in a synergistic manner and enhance the discontinuity during the transition. Namely, the degrees of renormalization are different between the two phases.
In this paper, we include both couplings in the model. The charge dynamics after photoexcitation of the neutral phase shows a large contribution from molecular vibrations. This is consistent with the ultrafast charge and molecular-vibration dynamics recently observed by the transient reflectivity measurements.  Numerical calculations show that the initial neutral state must be accompanied with a small dimerization to reproduce the experimentally observed, coherently oscillating dimerization immediately after photoexcitation. The presence of strings of lattice dimerization in the neutral phase has indeed been observed by the latest X-ray diffuse scattering. 
3 Model with Lattice and Molecular Phonons
We use the one-dimensional half-filled extended ionic Hubbard-Peierls-Holstein model
where creates an electron with spin at site , =, and =. The parameter denotes the transfer integral on a regular lattice (i.e., without lattice distortion), the site energy difference between neighboring orbitals when molecular distortions are absent, the on-site repulsion strength, and the nearest-neighbor repulsion strength. The lattice displacement at site modulates the transfer integral between the th and th orbitals and that between the th and th orbitals with the coefficient . The displacement in the th mode on the th molecule modulates the site energy with the coefficient . The quantities and are the time derivatives of and , respectively. The parameters and are their elastic coefficients, and and are their bare phonon energies, respectively.
For these model parameters, we take eV as the unit of energy and use =0.17,  =1.5, and =0.6; we vary around the phase boundary depending on the Peierls and Holstein coupling strengths. We define the strengths of these couplings as and . As long as the displacements are treated classically, the ground state is given by their static configuration, and the ground state is determined not by the distribution of , but by their sum, .  The displacements are scaled using ==1, so that we have and . For simplicity, we set .
As for phonons, we take one mode for the donor molecule and two modes for the acceptor molecule in addition to the lattice phonon mode, and use parameters that approximately reproduce the experimentally observed phonon energies: =0.013 (this energy is softened to be about 0.0066 when strongly coupled to electrons), =0.040, =0.055, and =0.12. Donor and acceptor molecules are specified by odd and even ’s, respectively. For simplicity, we set = because the conclusion does not depend on the detailed distribution of EMV-coupling strengths. The displacements , , and correspond to CA’s (320 cm), TTF’s (438 cm), and CA’s (957 cm) modes, respectively, in ref. \citenuemura_prl10, whereas we ignore TTF’s (740 cm) mode because its coupling to electrons is very weak.
Photoexcitation is introduced through the Peierls phase
The time-dependent vector potential for a pulse of an oscillating electric field is given by
where is the excitation energy, is the pulse width, and is the electric field amplitude.
The time-dependent Schrödinger equation for the exact many-electron wave function on the chain of =12 sites with periodic boundary condition is numerically solved by expanding the exponential evolution operator with a time slice =0.02 eV to the 15th order and by checking the conservation of the norm.  The initial state is set in the electronic ground state. The classical equations for the lattice and molecular displacements are solved by the leapfrog method, where the forces are derived from the Hellmann-Feynman theorem:
Unless otherwise stated, the initial displacements are at the minimum adiabatic potential for the electronic ground state, so that the initial velocities are zero.
4 Ground-State Properties
as a function of the site energy difference near the phase boundary with different combinations of and .
In the vicinity of the phase boundary, one phase is stable and the other is metastable. By comparing their energies, only the ionicity in the stable phase is plotted. It is clearly shown that, as increases, the neutral phase is stabilized, the discontinuity in ionicity is enlarged, the ionicity in the neutral phase on the large- side of the phase boundary becomes smaller, and the ionicity in the ionic phase on the small- side becomes larger (Fig. 1). As a consequence, in order for the ionic phase to be a typical Mott insulator, i.e., in order for the ionic phase to have nearly one electron per site, should be so large that the neutral phase is sufficiently stabilized. In particular, when is large, needs to be large [Fig. 1(b)]. Otherwise, the ionicity is too small even in the ionic phase.
Figure 2 has essentially the same information as Fig. 1. It is evident that, as increases, the ionic phase is stabilized, and the discontinuity at the transition is suppressed. When is large, the phase boundary does not shift monotonically as a function of [Fig. 2(b) and =0.33 (not shown)]. This is caused by the competition between the bond charge localized by and the site charge localized by . The coupling stabilizes only the ionic phase, while stabilizes both phases (i.e., the neutral phase strongly and the ionic phase weakly). When is small and the dimerization is weak (i.e., when lattice phonons only slightly lower the energy of the ionic phase), the increasing gains energy by localizing the bond charge, but it loses more energy by delocalizing the site charge. As a consequence, for large , the intermediate stabilizes the neutral phase more strongly than the ionic phase. This counterintuitive result cannot be obtained using the second-order perturbation theory.  Nevertheless, the shift of the phase boundary is small in this case.
Experimentally, it is known that the thermal transition from the neutral to ionic phases is accompanied by a large ionicity jump from 0.3 to 0.7.  This indicates a large . Of course, a finite is certain to induce dimerization and the three-dimensional ferroelectric order with a broken inversion symmetry.  As shown below, a large is also necessary for the reproduction of the optical conductivity alteration during the transition [Fig. 3(b)]. Therefore, both and are large in TTF-CA.
The optical conductivity in the ground state is calculated using
where is the current operator , is a peak-broadening parameter set at 0.05, and . Figure 3 shows in both phases close to the phase boundary.
When only is large, is largely altered during the transition [Fig. 3(a)]. The optical gap in the neutral phase strongly stabilized by is large, while that in the ionic phase is small in the present choice of for TTF-CA. Experimentally, such a large alteration is not observed at the transition. When only is large, the discontinuities in physical quantities at the transition are suppressed: the spectra in both phases are similar [Fig. 3(c)].
5 Photoinduced Dynamics from Ground State
In the case of a large and a large [Fig. 3(b)], the neutral phase near the phase boundary is photoexcited with an energy =0.65 just above the optical gap. The time evolution of the ionicity during and after photoexcitation is plotted in Fig. 4(a). For comparison, the displacement on the acceptor molecule (with a lower ) is shown in Fig. 4(b), the displacement on the donor molecule [ is multiplied so as to oscillate in the same phase with .] in Fig. 4(c), and the displacement on the acceptor molecule (with a higher ) in Fig. 4(d).
As increases, the electron density increases for the acceptor molecule and decreases for the donor molecule [eq. (6)], so that the displacement increases at the acceptor molecule and decreases at the donor molecule [eq. (5)]. Thus, the quantities shown in Figs. 4(b)-4(d) basically behave as times the cosine function. The ionicity receives a positive feedback from these molecular displacements and oscillates in the same phase with them. This cosine behavior is consistent with the experimental observation,  demonstrating that the neutral phase is stabilized by EMV couplings.
To investigate quantum effects, we treated the molecular vibration with the highest phonon energy quantum-mechanically as in ref. \citenkawakami_prl10. The difference between the wave profile when classically treated and that when quantum-mechanically treated was small (not shown). This is in contrast to the case reported in ref. \citenkawakami_prl10, where the photoinduced transition is from an insulator to a metal, so that low-energy electronic excitations exist and interfere quantum-mechanically with phonons. In such a case, the difference is large. In the present case, where the photoinduced transition is from an insulator to another insulator, the electronic excitations maintain a gap much larger than the phonon energies. The wavelet analysis does not show a trace of quantum interference.
In the present calculation, we did not adjust the distribution of so as to fit to the experimental data. Therefore, we can only roughly compare the experimental  and theoretical results for contribution of molecular vibrations to the photoinduced ionicity modulation. They are comparable in that the oscillation amplitude in the photoinduced ionicity change amounts to approximately one third in the present calculation and one fifth in the experimental result. We have calculated the photoinduced dynamics for different combinations of (, ), (0.167, 0.2), (0.167, 0.25), (0.167, 0.33), and (0.2, 0.33), all of which roughly reproduce the discontinuities in ionicity and optical conductivity at the transition. In these cases, the above contribution of molecular vibrations is roughly one third.
However, there is a qualitative difference between the experimentally observed and present numerical results, which originates from the fact that the initial state here is the ground state with the regular and static lattice configuration. It does not contain thermal fluctuations. In the experiment, the lattice phonons are ready to be dimerized.  The photoinduced ionicity dynamics has a large contribution from their slow oscillation, which fits to a damped oscillator. In the present case, the inversion symmetry of the neutral phase is broken by the oscillating electric field, so that lattice phonons begin to dimerize. However, its growth rate is quite low because it essentially corresponds to the spontaneous broken symmetry. The dimerization becomes significant only after t=600 in Fig. 4. If the initial state contains random numbers in , , , and , according to the Boltzmann distribution at a finite temperature of 0.01 eV, these fluctuations accelerate the growth of , but they do not oscillate it like times the cosine function as experimentally observed. It should be noted that this calculation cannot simulate thermally induced domains larger than the present system.
6 Dynamics with Precursory Dimerization
In this section, we introduce in the initial state a dimerization =0.01, which is much smaller than =0.058 of the ground state at =0.218 on the ionic side of the phase boundary. After obtaining the electronic ground state with fixed =0.01, we take it as the initial state and apply to it a pulse of an oscillating electric field with energy =0.65 again. The time evolution of the ionicity is plotted in Fig. 5(a). For comparison, the dimerization is shown in Fig. 5(b), and the displacement on the donor molecule in Fig. 5(c).
Now, rapidly increases with and oscillates like times the cosine function. Thus, receives a positive feedback from from an early stage. It is also clear that the oscillations of are temporally modulated by , as observed experimentally.  In our calculation, the molecular vibrations and lattice phonons are indirectly coupled through electrons.
In this calculation, fluctuations are not introduced into , , , or . This leads to the fact that an ultrafast charge transfer between neighboring donor and acceptor molecules continues to oscillate without dephasing. Although it appears as a very thick curve in Fig. 5(a), an ultrafast and large-amplitude oscillation is apparent if magnified on the time axis. In experiments performed at finite temperatures, such electronic motion is rapidly dephased. Indeed, if random numbers are introduced according to the Boltzmann distribution at a finite temperature of 0.01 eV in , , , and of the initial state, such an ultrafast and large-amplitude oscillation disappears but the slow oscillations caused by and survive. The resultant evolution of becomes closer to that observed in ref. \citenuemura_prl10. These fluctuations increase the ionicity on both sides of the phase boundary, which also becomes close to the experimentally observed behavior. As a consequence of the increased ionicity, the phase boundary is shifted to a larger . We show an example of transient ionicity in Fig. 6 in such a case of dephased charge transfer in a neutral state near the shifted boundary.
The ultrafast and large-amplitude oscillation [the very thick curve in Fig. 5(a)] is indeed suppressed here.
The above facts indicate the presence of strings of lattice dimerization before the photoexcitation of the neutral phase in the experiment performed at 90 K above the neutral-ionic transition temperature, i.e., 81 K.  They are also consistent with the X-ray diffuse scattering showing the characteristic length of 16 molecules for thermally induced ionic strings at 105 K.  At 90 K, the characteristic length would be larger and already beyond the size for which we can treat the exact many-electron wave function. Such one-dimensional precursors locally break the inversion symmetry and make the quick growth and the cosine-type oscillation of the dimerization possible. Considering the fact that the photoinduced neutral-to-ionic dynamics decays much faster than the photoinduced ionic-to-neutral dynamics,  the growth of such one-dimensional ionic domains is saturated before they form a three-dimensional ionic domain that is sufficiently metastable for longevity.
The picture of the photoinduced neutral-to-ionic transition obtained here is different from the previous one,  which suggested ultrafast growth of an ionic domain without dimerization, based on the static and regular neutral initial state without couplings to lattice phonons or to molecular vibrations. In reality, the neutral phase in equilibrium above the phase transition temperature contains strings of lattice dimerization, already before photoexcitation. Photoexcitation rapidly grows and oscillates the dimerization [Fig. 5(b)], which enhances the charge transfer between neighboring molecules [Fig. 5(a) compared with Fig. 4(a)]. In other words, precursors are transformed into ionic domains after photoexcitation.
The presence of a low-frequency component at 0.02 in Fig. 7(a) is due to the slowly growing dimerization contributing to after the period shown in Fig. 4(a). When the dimerization is introduced in the initial state, the dimerization oscillates from the beginning [Fig. 5(a)] and contributes more to the low-frequency component in Fig. 7(b). With the exception of such a low-frequency component, Figs. 7(a) and 7(b) are quite similar. All of the molecular vibrations evidently contribute to : a peak at =0.038 due to with =0.040, a peak at =0.050 due to with =0.055, and a peak at =0.11 due to with =0.12. All of them oscillate in the same phase as that in Figs. 4(b)-4(d).
Neutral-ionic transition and photoinduced dynamics have intensively been studied for the mixed-stack charge-transfer complex TTF-CA. The importance of the coupling to lattice phonons has been recognized from the dimerization and the consequent ferroelectric order in the ionic phase. The roles of the couplings to molecular vibrations have not been paid much attention to, except for a few theoretical studies. [27, 39] They are found to be important for the ionic phase to be a Mott insulator with large ionicity. The couplings to molecular vibrations stabilize the neutral phase, making the ionicity in the ionic (neutral) phase near the boundary large (small). The couplings to lattice phonons (molecular vibrations) reduce (enhance) the discontinuities in physical quantities. Both couplings are necessary for the reproduction of the experimentally observed ionicity, optical conductivity, and photoinduced charge dynamics.
The photoinduced ionicity dynamics also shows a large contribution from molecular vibrations. The comparison with the experimentally observed, photoinduced charge dynamics indicates the presence of strings of lattice dimerization as local and precursory symmetry breaking in the neutral phase above the transition temperature. Only when a small but finite dimerization is introduced in the initial state can we reproduce its cosine behavior and its large contribution to the photoinduced ionicity.
This work was supported by Grants-in-Aid for Scientific Research (C) (Grant Nos. 19540381 and 23540426), Scientific Research (B) (Grant No. 20340101) and Scientific Research (A) (Grant No. 23244062), and by “Grand Challenges in Next-Generation Integrated Nanoscience” from the Ministry of Education, Culture, Sports, Science and Technology of Japan, and the NINS program for cross-disciplinary study (NIFS10KEIN0160).
-  Y. Tokura: J. Phys. Soc. Jpn. 75 (2006) 011001.
-  K. Yonemitsu and K. Nasu: Phys. Rep. 465 (2008) 1.
-  R. H. McKenzie: Science 278 (1997) 820.
-  H. Seo: J. Phys. Soc. Jpn. 69 (2000) 805.
-  Y. Tanaka and K. Yonemitsu: J. Phys. Soc. Jpn. 77 (2008) 034708.
-  S. Miyashita and K. Yonemitsu: J. Phys. Soc. Jpn. 77 (2008) 094712.
-  Y. Tanaka and K. Yonemitsu: J. Phys. Soc. Jpn. 79 (2010) 024712.
-  S. Miyashita, Y. Tanaka, S. Iwai, and K. Yonemitsu: J. Phys. Soc. Jpn. 79 (2010) 034708.
-  Y. Kawakami, T. Fukatsu, Y. Sakurai, H. Unno, H. Itoh, S. Iwai, T. Sasaki, K. Yamamoto, K. Yakushi, and K. Yonemitsu: Phys. Rev. Lett. 105 (2010) 246402.
-  J. B. Torrance, J. E. Vazquez, J. J. Mayerle, and V. Y. Lee: Phys. Rev. Lett. 46 (1981) 253.
-  J. B. Torrance, A. Girlando, J. J. Mayerle, J. I. Crowley, V. Y. Lee, P. Batail, and S. J. LaPlaca: Phys. Rev. Lett. 47 (1981) 1747.
-  S. Koshihara, Y. Tokura, T. Mitani, G. Saito, and T. Koda: Phys. Rev. B 42 (1990) 6853(R).
-  M. Le Cointe, M.-H. Lemée-Cailleau, H. Cailleau, B. Toudic, L. Toupet, G. Heger, F. Moussa, P. Schweiss, K. H. Kraft, and N. Karl: Phys. Rev. B 51 (1995) 3374.
-  S. Koshihara, Y. Takahashi, H. Sakai, Y. Tokura, and T. Luty: J. Phys. Chem. B 103 (1999) 2592.
-  T. Suzuki, T. Sakamaki, K. Tanimura, S. Koshihara, and Y. Tokura: Phys. Rev. B 60 (1999) 6191.
-  S. Iwai, S. Tanaka, K. Fujinuma, H. Kishida, H. Okamoto, and Y. Tokura: Phys. Rev. Lett. 88 (2002) 057402.
-  E. Collet, M.-H. Lemée-Cailleau, M. Buron-Le Cointe, H. Cailleau, M. Wulff, T. Luty, S. Koshihara, M. Meyer, L. Toupet, P. Rabiller, and S. Techert: Science 300 (2003) 612.
-  L. Guérin, E. Collet, M.-H. Lemée-Cailleau, M. Buron-Le Cointe, H. Cailleau, A. Plech, M. Wulff, S. Koshihara, and T. Luty: Chem. Phys. 299 (2004) 163.
-  K. Tanimura: Phys. Rev. B 70 (2004) 144112.
-  H. Okamoto, Y. Ishige, S. Tanaka, H. Kishida, S. Iwai, and Y. Tokura: Phys. Rev. B 70 (2004) 165202.
-  H. Matsuzaki, H. Takamatsu, H. Kishida, and H. Okamoto: J. Phys. Soc. Jpn. 74 (2005) 2925.
-  G. D’Avino, A. Girlando, A. Painelli, M.-H. Lemée-Cailleau, and Z. G. Soos: Phys. Rev. Lett. 99 (2007) 156407.
-  L. Guérin, J. Hébert, M. Buron-Le Cointe, S. Adachi, S. Koshihara, H. Cailleau, and E. Collet: Phys. Rev. Lett. 105 (2010) 246101.
-  H. Uemura and H. Okamoto: Phys. Rev. Lett. 105 (2010) 258302.
-  N. Nagaosa and J. Takimoto: J. Phys. Soc. Jpn. 55 (1986) 2745.
-  N. Nagaosa: J. Phys. Soc. Jpn. 55 (1986) 2754.
-  A. Painelli and A. Girlando: Phys. Rev. B 37 (1988) 5748.
-  Y. Toyozawa: Solid State Commun. 84 (1992) 255.
-  T. Iizuka-Sakano and Y. Toyozawa: J. Phys. Soc. Jpn. 65 (1996) 671.
-  P. Huai, H. Zheng, and K. Nasu: J. Phys. Soc. Jpn. 69 (2000) 1788.
-  Y. Anusooya-Pati, Z. G. Soos, and A. Painelli: Phys. Rev. B 63 (2001) 205118.
-  T. Kawamoto, T. Iizuka-Sakano, Y. Shimoi, and S. Abe: Phys. Rev. B 64 (2001) 205107.
-  K. Yonemitsu: Phys. Rev. B 65 (2002) 205105.
-  L. Del Freo, A. Painelli, and Z. G. Soos: Phys. Rev. Lett. 89 (2002) 027402.
-  N. Miyashita, M. Kuwabara, and K. Yonemitsu: J. Phys. Soc. Jpn. 72 (2003) 2282.
-  J. Kishine, T. Luty, and K. Yonemitsu: Phys. Rev. B 69 (2004) 075115.
-  K. Yonemitsu: Phys. Rev. B 73 (2006) 155120.
-  K. Iwano: Phys. Rev. Lett. 97 (2006) 226404.
-  Z. G. Soos and A. Painelli: Phys. Rev. B 75 (2007) 155119.
-  K. Yonemitsu and N. Maeshima: Phys. Rev. B 79 (2009) 125118.