A semiclassical model for charge transfer along ion chains in silicates

Journal of Physics: Conference Series PAPER • OPEN ACCESS A semiclassical model for charge transfer along ion chains in silicates To cite this article: Juan F R Archilla et al 2024 J. Phys.: Conf. Ser. 2769 012015 View the article online for updates and enhancements. You may also like Nonlinear interactions of two-kink-breather solution in Yu-Toda-Sasa-Fukuyama equation by modulated phase shift Mi Chen and Zhen Wang - Revisiting multi-breathers in the discrete Klein–Gordon equation: a spatial dynamics approach Ross Parker, Jesús Cuevas-Maraver, P G Kevrekidis et al. - Quantum quench in the attractive regime of the sine-Gordon model Axel Cortés Cubero and Dirk Schuricht - This content was downloaded from IP address 150.214.182.234 on 17/12/2024 at 07:09 Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. Published under licence by IOP Publishing Ltd XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 1 A semiclassical model for charge transfer along ion chains in silicates Juan F R Archilla1,∗, J¯anis Baj¯ars2, Yusuke Doi3and Masayuki Kimura4 1Group of Nonlinear Physics, Universidad de Sevilla, ETSII, Avda Reina Mercedes s/n, 41012-Sevilla, Spain 2Faculty of Physics, Mathematics and Optometry, University of Latvia, Jelgavas Street 3, Riga, LV-1004, Latvia 3Division of Mechanical Engineering, Graduate School of Engineering, Osaka University, 2-1 Yamadaoka, Suita, Osaka 565-0871, Japan 4Department of Electrical and Electronic Engineering, Faculty of Science and Engineering, Setsunan University, 17-8 Ikeda-Nakamachi, Neyagawa, Osaka 572-8508, Japan E-mail: ∗Corresponding author: [email protected] Abstract. It has been observed in fossil tracks and experiments in the layered silicate mica muscovite the transport of charge through the cation layers sandwiched between the layers of tetrahedra-octahedra-tetrahedra. A classical model for the propagation of anharmonic vibrations along the cation chains has been proposed based on first principles and empirical functions. In that model, several propagating entities have been found as kinks or crowdions and breathers, both with or without wings, the latter for specific velocities and energies. Crowdions are equivalent to moving interstitials and transport electric charge if the moving particle is an ion, but they also imply the movement of mass, which was not observed in the experiments. Breathers, being just vibrational entities, do not transport charge. In this work, we present a semiclassical model obtained by adding a quantum particle, electron or hole to the previous model. We present the construction of the model based on the physics of the system. In particular, the strongly nonlinear vibronic interaction between the nuclei and the extra electron or hole is essential to explain the localized charge transport, which is not compatible with the adiabatic approximation. The formation of vibrational localized charge carriers breaks the lattice symmetry group in a similar fashion to the Jahn-Teller Effect, providing a new stable dynamical state. We study the properties and the coherence of the model through numerical simulations from initial conditions obtained by tail analysis and other means. We observe that although the charge spreads from an initial localization in a lattice at equilibrium, it can be confined basically to a single particle when coupled to a chaotic quasiperiodic breather. This is coherent with the observation that experiments imply that a population of charge is formed due to the decay of potassium unstable isotopes. 1. Introduction Tracks in mica muscovite were observed due to the special properties of the material. It can be exfoliated into very thin sheets, which are also semitransparent and therefore allow for direct observation. Some of those dark tracks were due to positively charged particles such as positrons, protons, or antipions [1–3]. Many other tracks were along the potassium hexagonal lattice layer and were therefore attributed to quasi-one-dimensional lattice excitations called quodons [4, 5]. XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 2 They were experimentally observed in an experiment where alpha particles were sent onto a monocrystal side, and subsequently, it was observed the ejection of an atom from the other side along the lattice close-packed directions [6]. Fossil tracks were most probably produced by the recoil of potassium ions after beta decay, which is 99% electron emission, leaving behind a positive charge. Therefore, it was deduced that quodons could transport electric charge [7, 8], which opened the possibility of the experimental measurement of electric current. This was achieved by bombarding a mica monocrystal with alpha particles which were expected to produce a large number of nonlinear excitations and measuring the current in the absence of an electric field, a phenomenon called hyperconductivity. There was an initial peak of current that after some seconds would diminish to the current transported by the flux of alpha particles [9]. It was interpreted as the nonlinear excitations trapping the accumulated electric charge left behind from electron beta emission. When this reservoir was exhausted, the only charge available was the one provided by the alpha particles. More experiments were done with other silicates [10] and other materials as it was developed a test to distinguish hyperconductivity and to separate quodon currents from Ohmic currents in semiconductors or conductors [11, 12]. The authors recommend a recent review on the subject [13]. A classical model for a cation chain in muscovite was developed based on first principles and empirical potentials. Interestingly, it was found the existence of crowdions or lattice kinks with energies of 26 eV, which is below the energy provided by the nucleus recoil after beta emission and above the energy to eject an atom, therefore coherent with the mica ejection experiment [14, 15]. There were also found kinks with other energies but with wings, also called nanopterons [16, 17]. Interestingly, a crowdion is a moving interstitial, which within an ionic crystal implies the movement of electric charge, making it a candidate for hyperconductivity. There are, however, observations of primary fossil tracks, which scatter in many other fainter ones, which should have much smaller energies, implying that less energetic nonlinear excitations as breathers were also of interest. Breathers are nonlinear localized solutions with also a vibration [18]. They are well described mathematically [19], and there are methods to construct numerically exact ones [20]. They can also appear in systems with long-range interaction [21] as alpha-helix proteins [22]. Breathers, also known as intrinsic localized modes (ILMs), were found in 3D in Si [23]. The theory of ILMs in 3D molecular dynamics was developed in Ref. [24], where it was found that they appear after X-ray recoil in molecular dynamics of ionic crystals but also in metals, such as Ni, Nb, and Fe. It was demonstrated that ILM frequencies can be within gaps in the phonon band or above it. The ILMs were highly mobile with energy of the order of 1 eV. ILMs were later found in other bcc metals, such as V and W [25], in fcc crystals, such as Cu [26, 27], and in hcp Be [28]. ILMs were also found in covalent crystals with the diamond structure, such as the insulator C, and the semiconductors Si and Ge [26, 27]. They were also constructed in graphene both with classical [29], and ab initio molecular dynamics [30]. To summarize, it was demonstrated that breathers or ILMs can appear in many materials with different electric behavior and with many different crystal structures. The spectral theory of breathers in the moving frame was developed in Ref. [31]. In the same work, it was applied to the muscovite model and traveling exact breathers with small energy were found. A related phenomenological model for muscovite was also used with the property that it was very easy to obtain traveling breathers in two dimensions [32]. The theory of exact breathers in the moving frame was extended to two and more dimensions using that model as an example [33]. The theory was also extended to a variation of the latter model with the addition of a quantum particle. Numerical methods were developed to deal with the differences in the time scale of the charge and the atoms that were able to conserve the charge probability [34]. Their spectral properties were deduced and described in Ref. [35]. In this paper, we construct a semiclassical model for the specific model based on potassium XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 3 chains in mica muscovite constructed from first principles and empirical potentials in Ref. [15], by adding a quantum particle. In many aspects, the model is similar to the phenomenological model used in Ref. [34, 35], but it is more complicated as the diagonal terms of the quantum charge Hamiltonian are not constant but correspond to the interaction with the other electric charges in the crystal. The vibronic interaction between the nuclei and the extra electron or hole is strongly nonlinear with the consequence that the tunneling or an electron between nearest–neighbor ions is only probable when the nuclei are close enough and the transition probability increases strongly as they approach. The consequence is the formation of quodons, dynamical states that transport electric charge in a localized manner, breaking the lattice discrete translation invariance as happens with the spontaneous symmetry breaking in the Jahn-Teller or pseudo Jahn-Teller effect[36, 37]. Certainly, the space of quodons has to keep the lattice translation invariance but given the low probability for the formation of quodons, their population at any given time should also break the translational invariance. The paper is structured as follows: after the introduction, the model is described in Section 2, including the hole or electron potential; the Hamiltonian and dynamical equations are then obtained in Section 3, which are linearized in Section 4. Section 5 presents the results of numerical integration for different initial variables, such as localized traveling and stationary trial solutions, and localized charge in a system at equilibrium; also it is observed a breather rebounding in a charge and a kink with an extra charge, and finally, a charge trapped by a chaotic quasiperiodic breather. The main part of the article finishes with the conclusions. There are also two appendixes: Appendix A describes the transformation of the semiclassical system into a real canonical Hamiltonian one, and Appendix B describes methods for numerical integration that preserve the charge probability at each integration step. 2. Model We propose a tight-binding model for a positive charge, a hole, or an electron, that we will call a charge with Q= 1 for the hole and Q=−1 for the electron. The model is an extension of the model already used for muscovite by Archilla et al. [14, 15, 17, 38]. The ket |n⟩represents the state in which the extra charge is located at the lattice position n, being ⟨n|the corresponding bra or adjoint operator, with ⟨m|n⟩=δm,n. An extra charge state is given by |ϕ(t)⟩=Pncn(t)|n⟩, where cn(t) represent the time-dependent possibility that the extra charge is located at site n. We will omit in what follows the explicit time dependence of |ϕ⟩and cnwhen convenient. The adjoint operator to |ϕ⟩is ⟨ϕ|=Pnc∗ n⟨n|, where c∗ nis the complex conjugate of cnand |c2 n|=c∗ ncnis the probability of finding the extra charge in the site nin the state |ϕ⟩. The total probability is one as there is an electron or hole in the system, i.e.: N X n=1 |cn|2= 1.(1) The cations K+are subjected both to an on-site potential U(un) which represents the interaction of K+with the ions of the surrounding lattice except for the K+ions in the same row, which are explicitly included as the interatomic interaction V(un−un−1) = VC(un−un−1) + VZ(un−un−1), where VCis the Coulomb repulsion; and VZ, a Yukawa type potential is a simplification of the ZBL potential [39], and corresponds to the repulsion between nuclei screened by the electron cloud. The variables unrepresent the n-th ion separation from the equilibrium position. We will use scaled units convenient for the modeled system: uL=a= 5.19 ˚ A is the equilibrium distance between potassium ions; uEis the Coulomb energy corresponding to two units of charge eat distance uL, i.e., uE=kee2/a ≃2.77 eV; the unit of mass is the mass of a potassium atom uM=mK= 39.1 amu; and the unit of time becomes a derived quantity: XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 4 uT= (uMu2 L/uE)1/2= (mKa3/kee2)1/2≃0.2 ps. The unit of angular frequency is equivalent to 5 Trad/s or uf= 1/(2πuT)≃0.8 THz≃26.7 cm−1. In the scaled units, the interaction potential energies become: V=VC+VZ(2) VC(un−un−1) = 1 1 + un−un−1 (3) VZ(un−un−1) = B 1 + un−un−1 exp(−1 + un−un−1 ρ) = C 1 + un−un−1 exp(−β(un−un−1)) ,(4) with B= 184.1 and ρ= 0.05690 β= 1/ρ = 17.58 and C=Bexp(−β) = 4.285 ×10−6. The on-site potential is obtained with the use of empirical potentials and electrostatic potentials [40], considering the interaction with the ions in the layers above, the first one composed of O−2, and the second of a mixed species between Si and Al with a positive charge +3.1. The resulting potential is a Fourier series that can be truncated at the fourth term: U= 4 X m=0 vmcos(2πx),with {vm}={2.4474,−3.3490,1.0997,−0.2302,0.0321}.(5) The on-site potential Uis represented among other potentials in Fig. 1. It has a potential barrier of 20 eV, as obtained with molecular dynamics [41], and a small amplitude frequency of 110 cm−1, as observed experimentally [42]. It is soft for x≲0.3a, and hard for larger distances as corresponds to a bounded system. Both the interaction and the on-site potentials are described in detail in Refs. [15, 31]. In the following, we consider for the first time some of the consequences of adding an extra unit of charge as a hole or an electron. 2.1. Hole or electron potential Localized energy transport in muscovite has been observed experimentally [6], and it has been deduced that dark tracks are produced by positive charge [7, 8]. Subsequently, charge transport has been observed experimentally in muscovite [9] and other silicates [10–13]. Different tracks suggest different types of carriers [43]. In this article, we attempt to model the transport of positive or negative charge attached to the K+ions and coupled with the lattice. An extra charge will appear due to β−decay of 40K [44], for which the nucleus transforms in 40Ca with an extra proton. The ion becomes Ca++, with the extra hole migrating to the neighboring ion. The far less frequent β−decay will transform 40K into 40Ar, and the ion will become Ar0, with an electron that can migrate to neighboring K+transforming them into neutral K [8]. An extra charge in site nwill experience the electrostatic interaction with the surrounding ions, given by UQ(un), and also with the other K+in the same row, given also by VC(un−un−1) in (3). UQis different from Uin (5), because the short-range interaction is already taken into account, and therefore only the extra electrostatic interaction of the extra charge has to be considered [40]. Considering the interaction with the oxygen ions at the immediate layers above and below and with the Si-Al mixed species in the following layer with charge +3.1 leads to a potential described also by a Fourier series truncated at the fourth harmonic. We denote it as XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 5 0 0.2 0.4 0 2 4 6 8 x/a Hole U (scaled) 0 0.2 0.4 0 2 4 6 8 x/a U (scaled) Electron UK Uh UK+Uh UK Ue UK+Ue Figure 1. (Left) Substrate potentials UKexperienced by a potassium ion K+;Uh, the electric potential experienced by a positive hole h+; and Uk+Uh, is the potential experienced by the double cation K++.UKis a sum of electrical and Buckingham terms [40]. (Right) Also, UK is represented with Ue, the potential experienced by an extra electron e−, and UK+Ue, the potential experienced by the neutral atom K0. The latter is therefore due only to the Buckingham potentials in Ref. [40]. The unit of energy is uE≃2.77 eV. UQ, being Uh=UQfor a hole and Ue=−Uhfor an extra electron. They are given by: Uh= 4 X m=0 hmcos(2πx),with {hm}={−0.6160,0.6941,−0.0930,0.0167,−0.0018},(6) UQ=QUh. Figure 1-Left shows the different potentials. Note that Uhhas a maximum at the equilibrium position for a hole (Q= +1) because it is energetically favorable to become closer to one of the negative oxygen ions. On the contrary, for an electron Q=−1, the potential energy Uehas a minimum. The potential Uhhas a maximum at x= 0 and a minimum between at x= 0.5 with a potential well of ≃ −3.74 eV, which diminishes the height of the potential barrier. The potential Ue=−Uhhas the opposite properties. The potential VC(un−un−1) is equal to (3), as is the Coulomb repulsion between an extra charge in site nand the nearest K+ions. In this case, the interaction with the extra charge with the K+ions in the same row does not include a Yukawa or ZBL potential because there is not an extra nuclei. Therefore the Hamiltonian operator for the extra charge will be ˆ HQ=X n En|n⟩⟨n|−Jn,n+1|n⟩⟨n+ 1|−Jn,n−1|n⟩⟨n−1|.(7) The transfer integrals Jn,n−1=Jn−1,n are related to the probability of a transition from the state |n−1⟩to the state |n⟩and vice versa. The term Encan be obtained easily as the expected XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 6 value of the charge Hamiltonian ⟨n|ˆ HQ|n⟩when the nondiagonal terms Jn,n−1are zero. It will be composed of the classical energy of the extra charge in the site n, that is, the electrostatic interaction with the lattice, plus the electrostatic interaction with the nearest K+ions. Then: En=QUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0.(8) We have subtracted the hole electrostatic energy at equilibrium −2Qand added a reference value E0, so as En=E0at the equilibrium distance. The value of E0has no physical consequences, and it will be generally taken as zero, however, some other values may be more convenient than others for numerical integration and to obtain periodic solutions [35]. The terms Jn−1,n have to be negligible for the equilibrium distance and become large only at the distance where the two electronic clouds of the K+ions interact. Therefore, a reasonable assumption is that they are exponentials with a decay rate similar to the ZBL repulsion between nuclei. That is, Jn−1,n =Jn,n−1=J0exp(−α(1 + un−un−1)) = I0exp(−α(un−un−1)) .(9) A first guess of αis α=βas both terms are the consequence of the overlapping of the electron shells. In principle, I0is the same for a hole or an electron, but this assumption might be revised. Note that I0is the only parameter for which we do not have an approximate value at the moment. We expect to deduce it from the band structure and experimental mobilities in muscovite, but at this stage, it will be taken as an adjustable parameter. Also, the value of αinitially equal to βmight have to be reconsidered. Figure 2 shows the interchange of probability when the particles approach as a consequence of the functional form of the transfer integrals. 2 2.2 2.4 2.6 2.8 3 14.5 15 15.5 16 16.5 time un+n (blue), |cn|2+n (red) Figure 2. Displacements (blue lines) and charge probability (red lines) for a system after an initial compression of particles 15 and 16 to a distance of 0.4. It brings about the oscillations of the particles in the antiphase for some time. The interchange of probability when the ions approach can be seen, and also the high frequency of the charge transfer. XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 7 3. Hamiltonian and dynamical equations We obtain the dynamical equations for the charge amplitudes from the Schr¨odinger equation iℏ∂/∂t|ϕ⟩=ˆ HQ|ϕ⟩collecting together the coefficients of the basis states |n⟩. The Planck constant in scaled units will be denoted τ=ℏ/uE/uT, where uEand uTare the scaled units of energy and time, therefore τ= 0.0011968. Then: iτ˙cn=QUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0cn−[Jn,n−1cn−1+Jn,n+1cn+1]. (10) The equations of movement for the variables unare obtained from the Hamiltonian as ˙pn=−∂Htot/∂unand ˙un=∂Htot/∂pn=pn. Where Htot is the Hamiltonian obtained as Htot =Hlat +⟨ϕ|HQ|ϕ⟩.(11) The first component of the Hamiltonian is the lattice classical Hamiltonian: Hlat =X n 1 2p2 n+U(un) + V(un−un−1),(12) with Vand Ugiven in (2) and (5). The second component of the Hamiltonian is the expected value of the charge Hamiltonian in a generic state |ϕ⟩=Pncn|n⟩, which is given by: HQ=⟨ϕ|ˆ HQ|ϕ⟩=X n Enc∗ ncn−[Jn,n+1c∗ ncn+1 +Jn,n−1c∗ ncn−1],(13) with Enfrom (8) and Jn,n+1 from (9). 3.1. Final equations Calculating the derivatives of Enand Jn,n+1, and collecting together the different terms, we obtain the dynamical equations for un: ¨un=−U′(un)−QU′ h(un)|cn|2 −1 (1 + un+1 −un)21 + Cexp(−β(un+1 −un)) + Q|cn|2+Q|cn+1|2 +1 (1 + un−un−1)21 + Cexp(−β(un−un−1)) + Q|cn|2+Q|cn−1|2 +Cβ exp(−β(un−un−1)) 1 + un−un−1−Cβ exp(−β(un+1 −un)) 1 + un+1 −un +αI0exp(−α(un+1 −un))(c∗ n+1cn+c∗ ncn+1) −αI0exp(−α(un−un−1))(c∗ ncn−1+c∗ n−1cn).(14) Note that the last two lines are real. The dynamical equations for the extra charge are given by: iτ˙cn=QUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0cn −I0exp(−α(un−un−1))cn−1−I0exp(−α(un+1 −un))cn+1 .(15) These equations can be written in real form and also as canonical Hamiltonian equations as explained in Appendix A. XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 8 4. Linearization We expand the terms in the dynamical equations (14)-(15), using 1/(1 + x)≃1−x+x2and 1/(1 + x)2≃1−2x+ 3x3, and we can neglect the ZBL potential because of small displacements un,C≃10−6and βC ≃10−4. We obtain the linearized dynamical equations: ¨un=−ω2 0un+ 2(un+1 +un−1−2un),(16) iτ˙cn=E0cn−I0[cn+1 +cn−1].(17) The frequency of the lattice homogeneous oscillations ω0is given by ω0= (− 4 X m=1 (2πm)2vm)1/2≃4.4800 .(18) A value of E0= 0 implies only a shift in the cnfrequencies of E0/τ, which can be convenient for integration purposes but has no physical consequences as always the products cnc∗ mthat appear in the dynamical equations are invariant with respect to a global frequency shift [35]. Note that the variables unand cnbecome decoupled at the linear limit. The dispersion relations are independent and are given by ω2=ω2 0+ 4c2 ssin2(q/2) ,for the variables un; (19) ω=E0 τ−2I0 τcos(q),for the variables cn.(20) The constant cs=√2 in (19) is the sound velocity in the lattice system without onsite potential and it is written as a symbol for comparison with other scalings. The second equation (20) multiplied by the scaled Planck constant τprovides the charge energy [45], i.e.: HQ=τω =E0−2I0cos(q).(21) 5. Simulation tests In this section, we test some physically interesting initial solutions and observe the result of the integration of the full system and some of its properties to check the model proposed. We limit ourselves to simulations for an extra hole, that is, Q= +1 in the previous sections. The preferred numerical methods are those that preserve the physical properties of the system at each integration step, in particular, charge probability conservation. They are described in Appendix B. 5.1. Extended solutions Linear solutions of the linearized equations are extended ones cn=1 √Nexp(i[qn −ωt]). Substitution in (16) and (17) leads to ¨un= 0 and τω =−2I0cos(q) and the charge Hamiltonian and frequency are Hh=τω =E0−2I0cos(q) and ω=E0/τ −2I0 τcos(q) as seen above. The lattice Hamiltonian is zero because un= 0 and pn= 0. Note that the unit energy in scaled variables is exactly uE=kce2/a, the electrostatic potential energy of a unit charge at the lattice unit distance. So an extra charge would provide twice that amount if the charge is localized in a single site, but it is diminished for the extended solution. The velocity of the waves in cnshould be the phase velocity as there is, in principle, a single plane wave, that is, Vteo =ω q=−2 q I0 τcos(q). The physical reason for the lattice to remain frozen is that the charge density is constant because |cn|2= 1/N so each charge is subjected to opposite repulsive forces from each neighbor with the same modulus that cancel themselves, i.e.: ¨u=−Q|c2 n+1|+Q|c2 n−1|=−Q/N +Q/N = 0. These solutions are somewhat irrelevant since nothing happens. However, they are very useful for coherence. XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 15 6. Conclusions In this work, we have presented a model for charge transport along K+chains in silicates mediated by nonlinear excitations. The motivation for this work is the experimental observation of hyperconductivity, the phenomenon of charge transport in the absence of an electric field when a side of the crystal is bombarded with alpha particles. The model relies heavily upon the work of previous publications, but for the first time, an electric charge is considered through a quantum Hamiltonian. The vibronic interaction between ions and an extra electron or hole, described by the transfer integral, is strongly nonlinear, increasing as neighboring ions become closer enhancing the probability of charge transmission. This allows for the existence of localized charge states that break the discrete translational invariance of the lattice. These states, when mobile, are the proposed charge carriers in hyperconductivity experiments. Most of the parameters are obtained through physical deduction, although some are yet not well known, particularly the transfer integral. Extensive work on that subject is being done, and it will be published elsewhere. In this work, we analyze the coherence of the model, its behavior, and spectra, for different initial conditions based on different ans¨atze, obtained from the tail analysis of extended and exponentially localized profiles, isolated charges, and other means. We have found interesting phenomena, such as the self-localization of some stationary solutions and the trap of a charge by a chaotic breather. The obtention of exact traveling solutions is the object of the present research, as well as the estimation of the missing parameters. Although in this work we have developed the model both for holes and electrons, we have limited the simulations to holes as they are the best candidate for the phenomenon of hyperconductivity because 99% of the 40K decay leaves a positive charge behind. New development may require a modification or refinement of our model, but, at present, it seems physically sound. At present, the propagation of charge has not been achieved due to the large number of frequencies that appear due to the different amplitudes of the particle oscillations. We expect to solve this problem both numerically and conceptually for the physical system. That work will be reported in future publications. Acknowledgments JFRA thanks projects MICINN PID2022-138321NB-C22 and VII PPIT-US 2024. He also acknowledges the Universities of Osaka and Latvia for hospitality. JB acknowledges financial support by the Faculty of Physics, Mathematics and Optometry of the University of Latvia. YD acknowledges JSPS Kakenhi (C) No. 19K03654. MK acknowledges support from JSPS Kakenhi (C) No. 24K07393 and 21K03935. References [1] Price P B and Walker R M 1962 Nature 196 732–734 [2] Russell F M 1967 Phys. Lett. 25B 298–300 [3] Russell F M 1967 Nature 216 907–909 [4] Russell F M and Collins D R 1995 Rad. Meas. 25 67–70 [5] Russell F M 2022 Nature Rev. Phys. 3147 [6] Russell F M and Eilbeck J C 2007 EPL 78 10004 [7] Russell F M 2015 Springer Ser. Mater. Sci. 221 3–33 [8] Archilla J F R and Russell F M 2016 Lett. Mater. 63–8 [9] Russell F M, Archilla J F R, Frutos F and Medina-Carrasco S 2017 EPL 120 46001 [10] Russell F M, Russell A W and Archilla J F R 2019 EPL 127 16001 [11] Russell F R and Archilla J F R 2022 Phys. Status Solidi RRL 19 2100420 [12] Russell F R and Archilla J F R 2022 Low Temp. Phys. 48 1009–1014 [13] Russell F M, Archilla J F R and Medina-Carrasco S 2021 13th Chaotic Modeling and Simulation International Conference Springer Proc. Complex. ed Skiadas C H and Dimotikalis Y (Cham: Springer) pp 721–734 [14] Archilla J F R, Kosevich Yu A, Jim´enez N, S´anchez-Morcillo V J and Garc´ıa-Raffi L M 2015 Springer Ser. Mater. Sci. 221 69–96 [15] Archilla J F R, Kosevich Yu A, Jim´enez N, S´anchez-Morcillo V J and Garc´ıa-Raffi L M 2015 Phys. Rev. E 91(2) 022912 XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 16 [16] Remoissenet M 1999 Mechanical Solitons (Berlin, Heidelberg: Springer) pp 143–189 [17] Archilla J F R, Zolotaryuk Y, Yu A Kosevich and Doi Y 2018 Chaos 28 083119 [18] Flach S and Gorbach A V 2008 Phys. Rep. 467 1–116 [19] MacKay R S and Aubry S 1994 Nonlinearity 71623 [20] Mar´ın J L and Aubry S 1996 Nonlinearity 91501 [21] Archilla J F R, Christiansen P L, Mingaleev S and Gaididei Y B 2001 J. Phys. A 34 6363–6373 [22] Archilla J F R, Gaididei Y B, Christiansen P L and Cuevas J 2022 J. Phys. A 35 8885–8902 [23] Voulgarakis N K, Hadjisavvas G, Kelires P C and Tsironis G P 2004 Phys. Rev. B 69(11) 113201 [24] Hizhnyakov V, Haas M, Shelkan A and Klopov M 2014 Physica Scripta 89 [25] Murzaev R T, Kistanov A A, Dubinko V I, Terentyev D A and Dmitriev S V 2015 Comp. Mater. Sci. 98 88–92 [26] Hizhnyakov V, Haas M, Shelkan A and Klopov M 2015 Springer Ser. Mater. Sci. 221 229–245 [27] Hizhnyakov V, Haas M, Klopov M and Shelkan A 2016 Lett. Mater. 661–72 [28] Bachurina O V, Murzaev R T, Semenov A S, Korznikova E A and Dmitriev S V 2018 Phys. Solid State 60 989–994 [29] Hizhnyakov V, Klopov M and Shelkan A 2016 Phys. Lett. A 380 1075–1081 [30] Lobzenko I P, Chechin G M, G S Bezuglova Y A B, Korznikova E A and Dmitriev S V 2016 Phys. Solid State 58 633–639 [31] Archilla J F R, Doi Y and Kimura M 2019 Phys. Rev. E 100 022206 [32] Bajars J, Eilbeck J C and Leimkuhler B 2015 Springer Ser. Mater. Sci 221 35–67 [33] Baj¯ars J and Archilla J F R 2022 Physica D 441 133497 [34] Baj¯ars J and Archilla J F R 2022 Mathematics 10 3460 [35] Archilla J F R and Baj¯ars J 2023 Axioms 12(5) 437 [36] Bersuker I B 2006 The Jahn–Teller effect (Cambridge: Cambridge University Press) [37] Bersuker I B 2017 J. Phys. Conf. Ser. 833 012001 [38] Archilla J F R, Doi Y and Kimura M 2019 Phys. Rev. E 100 022206 [39] Biersack J, Ziegler J P and Ziegler M D 2008 SRIM - The Stopping and Range of Ions in Matter (Chester, Maryland: Published by J.P. Ziegler) [40] Gedeon O, Machacek J and Liska M 2002 Phys. Chem. Glass. 43 241–246 [41] Collins D R and Catlow C R A 1992 Am. Mineral. 77 1172–1181 [42] Diaz M, Farmer V C and Prost R 2000 Clays Clay Miner. 48 433–438 [43] Russell F M 2018 Nonlinear Systems, Vol. 2: Nonlinear Phenomena in Biology, Optics and Condensed Matter ed Archilla J F R et al. (Cham: Springer) pp 241–260 [44] Mougeot X and Helmer R G 2012 LNE-LNHB/CEA–Table de Radionucl´eides, K-40 tables http://www.lnhb.fr/nuclides/K-40 tables.pdf [45] Griffiths D J and Schroeter D F 2018 Introduction to Quantum Mechanics 3rd ed (Cambridge: Cambridge University Press) [46] Ikeda K, Doi Y, Feng B F and Kawahara T 2007 Physica D 225 184–196 [47] Hairer E, Lubich C and Wanner G 2006 Geometrical Numerical Integration Springer Ser. Comp. Math. (Berlin Heidelberg New York: Springer) ISBN 978-3-540-30663-4 Appendix A. Real canonical Hamiltonian equations The equations in Section 3.1 can be converted into canonical Hamiltonian equations for the non-separable Hamiltonian H=Hlat +Hh, with Hlat and Hh, given by Eqs. (12)–(13). Let us denote anand bnas the real and imaginary part of cn, and let us define the scaled variables An=√2τan,Bn=√2τbn, then Bnis the conjugate momentum of An. The coordinate variables become z= [u1, . . . , uN, A1, . . . , AN] and the momenta variables Π = [p1, . . . , pN, B1, . . . , BN]. With this notation: ˙un=∂H ∂pn ;˙ An=∂H ∂Bn ; ˙pn=−∂H ∂un ;˙ Bn=−∂H ∂An .(A.1) Or ˙zn=∂H ∂Πn ;˙ Πn=−∂H ∂zn .(A.2) This system is real and in the form of canonical Hamiltonian equations, and it is convenient for numerical integration. XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 17 With this scaling, the new evolution equations for the displacements are changed, while the equations for the charge show no change: ˙pn=−U′(un)−QU′ h(un)(A2 n+B2 n) 2τ −1 (1 + un+1 −un)21 + Cexp(−β(un+1 −un)) + Q(A2 n+B2 n+A2 n+1 +B2 n+1) 2τ +1 (1 + un−un−1)21 + Cexp(−β(un−un−1)) + Q(A2 n+B2 n+A2 n−1+B2 n−1) 2τ +Cβ exp(−β(un−un−1)) 1 + un−un−1−Cβ exp(−β(un+1 −un)) 1 + un+1 −un +2αI0exp(−α(un+1 −un))(An+1An+Bn+1Bn) 2τ −2αI0exp(−α(un−un−1))(AnAn−1+BnBn−1) 2τ,(A.3) and: ˙ An=1 τQUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0Bn(A.4) −I0 τexp(−α(un−un−1))Bn−1−I0 τexp(−α(un+1 −un))Bn+1 ; ˙ Bn=−1 τQUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0An +I0 τexp(−α[un−un−1])An−1+I0 τexp(−α[un+1 −un])An+1 .(A.5) The lattice Hamiltonian remains unchanged, and the charge Hamiltonian changes in the obvious way. We reproduce it here for completeness: Hlat =X n 1 2p2 n+U(un) + V(un−un−1),with (A.6) V(un−un−1) = 1 1 + un−un−1+C 1 + un−un−1exp(−un−un−1 ρ).(A.7) For energy calculations, it will be convenient to subtract (1 + C) = V(0) from V, i.e., so as the energy is zero at equilibrium, but we do not include it here for clarity and because the dynamical equations do not depend on a constant term. The expected value of the charge Hamiltonian in a generic state |ϕ⟩=Pncn|n⟩is given by: HQ=X n En (A2 n+B2 n) 2τ−2Jn,n+1τAnAn+1 +BnBn+1 2τ,(A.8) or explicitly: HQ=X nQUh(un) + Q 1 + un−un−1 +Q 1 + un+1 −un−2Q+E0(A2 n+B2 n) 2τ −2I0exp(−α(un+1 −un))(AnAn+1 +BnBn+1) 2τ.(A.9) XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 18 With this change of variables, we have obtained canonical Hamiltonian equations, but also, the system is a bit more amenable as the scaling has moved from the time to the probability amplitudes (A, B). As seen in Sect. 4, the parameter ΩQ= 2I0/τ is the maximum linear frequency of the charge of the order of magnitude of the lattice minimal frequency. The eigenfrequency En/τ should be relatively small for breathers because Enis zero at equilibrium and, therefore, small for the relatively small displacements corresponding to breathers. Appendix B. Numerical integration In this section, we use it as the main reference Ref. [47]. Our system of ordinary differential equations (ODE) is canonical Hamiltonian, and most important, it conserves the charge probability, which is a property of the Schr¨odinger equation. Therefore, it is suitable to construct a numerical integration scheme that holds the same important physical properties. Another interesting properties of integrators is whether they are explicit or implicit. Explicit methods are much more efficient, but they are possible only when the Hamiltonian is separable, i.e., it can be written as a sum of parts with only generalized coordinates, but our Hamiltonian is not separable: H=H(Π) + H(z). With these limitations, let us write the equations of motion in condensed form to clarify their structure: ˙u=∂H ∂p =p(B.1) ˙ A=∂H ∂B =M(u)B(B.2) ˙p=−∂H ∂u =f(u, A, B) (B.3) ˙ B=−∂H ∂A =−M(u)A , (B.4) where the first two equations correspond to the generalized coordinates and the other two to the generalized momenta of the Hamiltonian system. The matrix M(u) and the function for forces fcan be obtained trivially from Eqs. A.3–A.5. This structure allows the use of a splitting method [47], considering the two systems of ODEs: System I        ˙u=p ˙ A= 0 ˙p=f(u, A, B) ˙ B= 0 System II        ˙u= 0 ˙ A=M(u)B ˙p= 0 ˙ B=−M(u)A (B.5) System I corresponds to an ODE for the lattice variables while the charge variables are kept constant, and system II corresponds to an ODE for the charge variable while the lattice variables are kept constant. System I consists of a separable Hamiltonian H=H(p) + H(q), and can be written as ¨u=f(u, A, B), which allows it to be integrated with the explicit Verlet method of second-order, also known as St¨ormer or leap-frog method. It is given by: un+1/2=un+h 2pn(B.6) pn+1 =pn+hf(un+1/2, An, Bn) (B.7) un+1 =un+1/2+h 2pn+1 (B.8) An+1 =An(B.9) Bn+1 =Bn,(B.10) XXV-JAHN-TELLER IOPConf.Series:JournalofPhysics:Conf.Series 2769(2024)012015 IOPPublishing doi:10.1088/1742-6596/2769/1/012015 19 where his the step in time and un+1/2is an intermediate point to calculate. We denote its flow by Φh V. As the charge variables do not change, the charge probability is trivially conserved. System II is also a separable Hamiltonian for the charge, i.e., Hq=H(B) + H(A). It is also linear and can be solved using the implicit midpoint rule given by: un+1 =un(B.11) pn+1 =pn(B.12) An+1 =An+hMun+1 +un 2Bn+1 +Bn 2(B.13) Bn+1 =Bn−hMun+1 +un 2An+1 +An 2.(B.14) The values An+1 and Bn+1 are obtained by solving the second pair of linear equations above. The lattice Hamiltonian is conserved in System II because (u, p) do not change, so we only have to consider the part corresponding to the charge. The charge Hamiltonian and probability are given by quadratic expressions Hq=1 2ATMA +1 2BTMB and |c|2=ATIA +BTIB, where Trepresents the transpose and Ithe identity matrix. Both Mand I(trivially) are symmetric. The implicit midpoint method conserves quadratic functions and, therefore, the Hamiltonian and the charge probability. Therefore, the charge probability is exactly conserved in both steps I and II and the Hamiltonian is exactly conserved in step II and in step I up to the second order in the integration step h. The Strang or Marchus splitting [47] obtains the global map in a symmetric way as: Φh= Φh/2 V◦Φh IM ◦Φh/2 V.(B.15) The numerical map is time-reversible, conserves exactly the charge probability and the Hamiltonian to second order in the integration step. The global map conserves exactly the charge probability and the Hamiltonian to second order in the integration step. Note that, if convenient, the splitting can be done by composing a larger number of flows. For example, it is possible to divide the IM step into as many IM steps as convenient. The same can be done with the V steps. There are other splitting methods that could be used. For example, the splitting into five different steps is developed and applied to a related system in Refs. [33, 35].