{'fraction_non_alphanumeric': 0.03074539226879791, 'fraction_numerical': 0.009664002609688468, 'mean_word_length': 4.814366998577525, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Discrete breathers in an array of self-excited oscillators: exact solutions and stability

Dynamics of array of coupled self-excited oscillators is considered. Model of Franklin bell is adopted as a mechanism for the self-excitation. The model allows derivation of exact analytic solutions for discrete breathers (DBs), and extensive exploration of their stability in the space of parameters. The DB solutions exist for all frequencies in the attenuation zone, but lose stability via Neimark-Sacker bifurcation near the boundary of propagation zone. Besides the well-known DBs with exponential localization, the considered system possesses additional and novel type of solutions - discrete breathers with main frequency in the propagation zone of the chain. The amplitude of oscillations in this solution is maximal at the localization site and then exponentially approaches constant value at infinity. We also derive these solutions in closed analytic form. They are stable in a narrow region of system parameters bounded by Neimark-Sacker and pitchfork bifurcations.


Introduction
Localized excitations in arrays of coupled oscillators are often referred to as discrete breathers (DBs) or intrinsic localized modes (ILM) [1], [2]. The appearance of the DB solution is triggered by the discreteness of the system and nonlinearity of the evolution equations [2]. Such dynamical excitations normally exhibit exponential localization; in models with purely nonlinear inter-particle interaction stronger hyperexponential localization is encountered [2]. These regimes have been observed and studied in various physical systems, and here we give several examples. Forced-damped coupled pendula and stability boundaries between different localized modes are studied in [3]. In [4], magnetic metamaterial breathers in energy conserving and dissipative systems are explored and coupling limit, stability and mobility properties are examined. Paper [5] demonstrates numerically and experimentally the excitation of the ILM as a result of subharmonic, spatially homogeneous driving. Several studies proved the existence of DBs in driven micromechanical silicon-nitride cantilever arrays [6][7][8]. In [9] it was examined how the activation of localized mode effects certain properties (thermal expansion, thermal and electrical conductivity) of metallic uranium and ionic sodium iodide crystals.
In many studies related to the DBs, conservative models are used [10]. In realistic systems the dissipation commonly cannot be neglected, and forcing of some kind should be applied to maintain the oscillations. When damping and forcing are not negligible, Hamiltonian approach is not applicable and stable localized excitation regimes become attractors of the dynamic flow. Mathematical analysis of such attractors in MDOF nonlinear systems is not a straightforward task. Consequently, most systems possessing the DBs are analysed by means of numerical simulations. Some recently studied models allow broader exploration of these dynamical regimes. In [11], [12] it was demonstrated that it is possible to obtain exact solutions for DBs in damped chains with vibro-impact nonlinearity in the case of a homogeneous external forcing. Moreover, it was found possible to determine a zone of existence of the DB in a space of parameters and to reveal its stability thresholds.
Instead of direct external forcing, one can counter-balance the dissipation by internal generation of energy inside the system. Such phenomenon, referred to as self-excitation, is well-known in physics [13]. This phenomenon is usually realized through an internal feedback mechanism. Few of the many examples of self-oscillatory systems are the Huyghens pendulum clock [14], flutter of suspension bridges (e.g. Tacoma Narrows Bridge in 1940 [15]), the human voice [16] and biochemical processes [17], [18].
Such systems often take advantage of energy transform in order to produce self-sustaining dynamics. Dynamics of various nanomechanical systems can be described in terms of the self-excitation. Optomechanical systems can be driven by bolometric [19] or radiation pressure feedback [20]. In nanoelectromechanical systems (NEMS), internal feedback has been achieved by field emission of vibrating nanowires subject to a DC voltage [21], the periodic charging of a nanowire [22], transport through a carbon nanotube quantum dot mediated by the backaction of tunnelling single electrons [23]. Recently, a thermodynamic feedback of a piezoresistive resonator was employed to achieve self-oscillation [24].
A famous classical example of generation of self-oscillatory motion driven by conversion of electricity to mechanical oscillation is Franklin bell [25]. The system consists of two bells, one connected to a lightning rod attached to a chimney, the other one -to the ground. A metal ball hung between the two bells was warning Franklin that storm was coming by swinging repeatedly from one bell to the other as result of repetitive charge and discharge periods.
A particularly significant study of the self-excitation in nano-mechanical systems realized the idea of Franklin bell on nano-level (Gorelik [26]). In this system, the Coloumb blockade effect in combination with soft links attaching a metallic grain to electrodes was found to result in significant oscillations. At sufficiently large value of the bias voltage, the oscillation amplitude increased until the balance between energy dissipation and energy absorption was achieved, and the system reached the stable self-excitation regime. A cyclic change in direction is caused by the repeated pick up of N electrons at the source electrode and their transfer to the drain electrode. As a result, the sign of the net grain charge alternated leading to the oscillatory motion of the grain. In this way, the "electron shuttle" mechanism for the charge transport was realized.
Electron-shuttles were further investigated in several studies and demonstrated rather ontersting and beneficial properties [27][28][29][30][31][32]. One of the most interesting opportunities is a suggestion of extremely accurate discrete current source [27], [28]. In [29] the regions of self-excitations were derived analytically and verified experimentally. Bifurcations driven by change of system parameters such as electromechanicalcoupling were demonstrated. Electron-shuttles were found to also operate at radio frequencies (RF), possess a set of resonance frequencies, and exhibit Coulomb blockade even at room temperature [30]. Spontaneous symmetry breaking has also been observed in these systems [31], revealing Arnold tongues in the frequency domain [32].
As a result of their diversity, electron-shuttles may possibly encourage progress in such fields as Metrology [33], Ferromagnetic island may encourage investigation of Kondo shuttling [34] and realization of spintronic devices (e.g. mechanical spin valves) [35], a superconducting shuttle may allow to target mechanically mediating phase coherence in a nanostructured device [36].
A self-excited system such as Franklin bell or the "electron-shuttle" may be mathematically represented as a single degree-of-freedom oscillator restricted by two inelastic constraints and subjected to constant force, which direction is reversed upon reaching the constraints. This work will consider dynamical behaviour of arrays of coupled self-excited oscillators of this type. Special attention will be paid to localization patterns and their stability. As the first problem to address, we consider the discrete breather, localized on a single excited oscillator. As a simple model, we adopt an infinite chain restricted by inelastic constraints, where the central particle is subjected to self-excitation of the specified type. The analysis presented in this paper focuses on derivation of analytic solutions for the discrete breathers, exploration of existence and stability boundaries in a parametric space and numerical verification of the results.

Description of the model and analytical solutions a. Model description
As it was mentioned before, the single self-excited model "Franklin bell" is a single particle which movement is restricted by two inelastic constraints. Interaction of the particle with the constraints is described by Newton impact model with restitution coefficient . These constraints are, in turn, positively and negatively charged electrodes. The particle is connected to a fixed location through linear elastic leafspring with stiffness . We suggest for simplicity that after each re-charging the particle carries the same absolute amount of electric charge -positive or negative, depending on which electrode was impacted last. If the applied voltage is constant, one can adopt that between the impacts with the electrodes the particle is subject to action of force with constant magnitude . Direction of the force is determined by the sign of the charge carried by the particle at the moment and is reversed after each impact. The sketch of the single-DOF basic model is shown in Figure 1. The distance between the electrodes is set to 2 withoutloss of generality. The stable self-excited regime of this single-DOF model is described by simple analytic solution; dependence of frequency of oscillations on the system parameters is expressed by Eq. (1).

Figure 1 -Sketch of the single particle system
The model of coupled Franklin bell oscillators addressed in this work is a chain of particles coupled by linear elastic springs with unit stiffness (without loss of generality). Similarly to the single-DOF model of the Franklin bell outlined above, the movement of each particle is restricted by two inelastic electrodes with restitution coefficient . To construct the simplest possible localized solution, we adopt that only the central particle ( = 0) impacts the constraints and, consequently, carries the charge. Similarly to the single-DOF system, this particle is subject to action of force with constant magnitude , which direction is reversed upon reaching the constraints. A sketch of the system is shown in Figure 2.

Figure 2 -Sketch of the model system
We suggest that the DB solution is periodic in time. In this case, the impacts with the electrodes can be expressed as periodic pulses, which act on the particle 0 n  . Consequently, equations of motion of the model can be cast in the following form:

2Δ
It is important to notice that equation ( 2 ) assumes symmetric solution and periodic response. Besides, it is adopted that only particle 0 n  reaches the constraints. Any solution that doesn't comply with these assumptions will require different formulation.
Without affecting the generality, in Eq.
( 2 ) the mass of each particle and clearance  (se Figure 2) are set to unity. One should supplement Eq. ( 2 ) by explicit formulation of the impact conditions. At the impact instance of the central particle t i,o , the following conditions are satisfied: In this model, the impact losses are the only source of damping in the system. The notations denote the velocities immediately after and before impacts, respectively. Naturally, the diapason of the restitution coefficient is The forcing in (  2  ) can be rewritten (in the sense of distributions) in terms of Fourier series with frequency

b. Localized breather in the attenuation zone of the chain
In this subsection we derive the analytical solution for the DB, similarly to the derivation in paper [11]. The periodic solution has the following form: is homogenous. The solution for n   should vanish, therefore one obtains: Therefore, by substituting Eq. ( 7 ) , (8) into Eq. ( 5 ) the general expression for the periodic localized solution is obtained: The unknowns of expression ( 9 ) are 0 , p  . These unknowns are found with the help of impact conditions . For the sake of simplicity, the first impact is set to occur at 0 t  : and ( 12 ) , one obtains expressions (one of them implicit) for the unknowns of expression ( 9 ) : it can be seen that the momentum transferred to the main particle is a function of two system parameters: , fk , and can be derived explicitly. An important property of the momentum is that it is directly related to the forcing, which results in the possibility of extremely weak impacts for very small forcing.
The oscillation frequency 0  , cannot be explicitly extracted from Eq.  A plot of 0 as a function of , is given in Figure 3. Once the periodic solution is found, it is natural to explore the extent of its existence. One can easily see that the maximal existence limit is the propagation zone (PZ) boundary 0 2   . In order to find the actual existence zone in the space of parameters and to ensure the self-consistency of the solution, it is necessary to satisfy the condition that the particles with 0 n  never reach the impact constraints. Thus, one should demand that for all particles 0 n  : is satisfied for any set parameters as long as Thus, this analysis leads to conclusion that the existence zone of the periodic solution is as broad as the AZ of the chain. This, however, by no means guarantees the stability of this solution. The latter is analyzed in Section 3.

c. Localized breather in the propagation zone of the chain
As will be further presented in the numerical section, at certain parameters, an additional localized solution in the PZ exists for this system, and it coexists with the DB solution derived above. To demonstrate this possibility, we construct below the solution, where only the first harmonics belongs to the PZ, while the rest infinite number of harmonics are situated in the AZ. The proposed solution is a superposition of the two essentially different forms.
First we present the system in an equivalent form which results from the symmetrical forces induced on the central particle by 11 , uu  .
The suggested solution consists of two parts -irradiated (travelling) and exponentially localized: it is easy to obtain the following expressions for wavenumber (it is, of course, common dispersion relation for the linear chain), and amplitude of the "propagating" harmonics : can be solved numerically, and one can directly obtain values of , by substitution into (  23  ) and ( 20 ) .
The existence of this solution depends on satisfying requirement ( 15 ) . It was verified numerically for a wide range of parameters that this solution is self-consistent.

d. Other solutions
So far, two, essentially different, exact solutions were derived analytically (and will be verified numerically in Section 4). This, naturally, raises the question about other stable attractors in the system. Theoretically, a number of additional solutions can be derived for the suggested model, e.g. the solutions with different number of harmonics in the propagation zone. We restrict ourselves by two solutions presented so far for two reasons: 1) These solutions are structurally different while additional solutions are expected to have structure similar to the solution in the PZ; 2) Numerical analysis implies that the solutions presented above are in a sense the most typical, while the other solutions are expected to be stable in considerably narrower regions of system parameters.

Stability analysis a. Localized AZ solution
In this section a linear stability analysis of the exponentially localized breather in the attenuation zone is carried out by means of Floquet theory. The theory requires computation of the monodromy matrix which describes the evolution of perturbation of the state vector ( 26 ) , for a single period. According to Floquet theory, the loss of stability occurs when at least one eigenvalue leaves the unit circle in the complex plane. In our case, each period of motion consists of 2 impacts and 2 periods of linear evolution between the impacts. This can be generally represented as product of 4 matrices ( 27 ) , where describes the evolution between impacts, and is the saltation matrix [37] which describes the evolution of the perturbed phase flow at the impact hypersurfaces: The evolution of the small perturbation between successive impacts is given by block construction in ( 28 ) , where is an 2 + 1 × 2 + 1 identity matrix, 0 is a 2 + 1 × 2 + 1 zero matrix, and Laplace adjacency matrix is given below in In current work both acceleration and velocity are discontinuous at impacts; this issue effects the structure of the saltation matrix. It is presented by block construction in ( 30 ) , where 0 is a 2 + 1 × 2 + 1 zero matrix, , are given in ( 31 ) . The only particle which is subjected to impacts is the central one, and for that reason the only non-zero term in matrix and the only non unity term in the diagonal of matrix is + 1, + 1 .                                                        To explore the stability of the DB solution, we compute numerically the eigenvalues of the monodromy matrix. It is impossible to consider the infinite number of particles; therefore, the monodromy matrix is truncated at certain number n. It was found that the stability thresholds converge as the number of particles increases. In all stability calculations for the localized breather, = 300 was found to be satisfactory. An example is given in Figure 4. The stability boundary in the − plane is presented in Figure 5. The loss of stability at the boundary is realized through Neimark-Saker bifurcation. It is important to notice that close to the stability boundary, this solution coexists with the solution in the PZ; stability of the latter will be discussed in the next subsection. Therefore, in order to numerically verify the boundary of stability, it is necessary to select initial conditions very close to the desired solution. The issue of coexistence will be further addressed in subsection 3c and in Section 4.  In the infinite chain, the propagating solution consists of a propagating wave and an exponentially localized component. Here, in order to perform the stability analysis, we substitute the initial infinite chain by a symmetric (with respect to the impacting particle) finite segment of the chain with oscillatory profile similar to the PZ DB (17)(18)(19)(20)(21)(22)(23)(24)(25). We also assume the boundary conditions at the edges of the selected segment that exactly reproduce the forces induced by the omitted parts of the chain. However, variation on these boundary conditions is equal to zero, so they have no effect on the monodromy matrix. Then, we can explore the stability of the PZ DB with respect to perturbations that are localized at the selected fragment of the infinite chain. Control of convergence of the stability boundary is achieved through exploration of longer fragments. The shape of the monodromy matrix is not different from the one derived for the AZ DB. Therefore, the only parameters which determine the Floquet multipliers are the DB frequency, and velocity and acceleration of particle 0 n  before and after the impact. All of these parameters are calculated from expressions which were given in section 2c.
The solution in PZ was found to be stable between two limit curves in the − plane. On the lower boundary, the solution loses stability through Neimark-Sacker bifurcation. On the upper boundary, the solution loses stability by pitchfork bifurcation for 0 < < 0.45 (example in Figure 7) and by Neimark-Sacker for 0.45 < < 1 ( Figure 6). Again, due to coexistence of the solutions, in order to verify results numerically, initial conditions were chosen in close vicinity of the analytical solution.
The value of for simulations was selected such that the results were not influenced by further increasing of (the process converged). Here, on the lower boundary value of = 300 was sufficient, while on the upper boundary it was required to use = 500.  In both cases, the loss of stability close to these boundaries occurs by Neimark-Sacker mechanism. The frequency band of the PZ DB is considerably narrower than of the AZ DB. It is reasonable to assume that other solutions in propagation zone will be stable in narrow regions and will also lose stability, when their harmonics will approach the PZ boundary, by the generic Neimark-Sacker scenario. It is also worthy to mention that there are zones of frequencies where no one of the explored solutions is stable.

Numerical Analysis
An example of the AZ DB obtained numerically, compared to the analytic solution, is presented in Figure 9. As one would expect from the exact solution, the analytical and numerical results match exactly. This solution is characterized by frequency in the AZ, and exponential decay of the amplitudes along the chain. In this example the amplitude of the nearest neighbour is 1 = 0.09.
The second attractor described was the PZ. As was shown in the stability analysis, this solution co-exists with other solutions and approach to the vicinity of this attractor depends on careful choice of the initial conditions. An example of the solution is shown in Figure 10. The PZ DB solution has longer period and nonzero asymptotic value of the neighbours amplitude (in this example -0.88). Another phenomenon illustrated by Figure 9 and Figure 10 is the coexistence in the system, as both solutions were obtained for the same system parameters with different initial conditions.  Another way to examine coexistence in the system is by presenting both solutions in a certain projection of the state space. In Figure 11 the two solutions are illustrated in a plane of 0 − 0 with the particular initial condition that led to each one of the solutions (localized solution -'o', propagating solution -'*'). The initial conditions for both solutions have 2 nonzero coordinates 0 , 0 (and can be represented in the 2D projection of Figure 11). Here, it is seen that stable attractors can be reached even if system is instantiated at a point distant from the analytic solution in the state space. It was shown that the localized solution loses its stability through Neimark-Sacker bifurcation. In order to explore the state of the system after the loss of stability, the forcing parameter was gradually decreased to the precise limit at which a modulated response, characteristic to Neimark-Sacker bifurcation, was achieved. The modulated response is most easily noticed in the response of the first neighbours. ( 1 , −1 ), an example is shown in Figure 12. The FFT of this response provides the frequency of modulation, which in the case of the parameters in Figure 12   In the stability analysis of the PZ DB it was found that the solution has a region where it loses stability through the pitchfork bifurcation. In Figure 13 a typical numerical response is presented at the pitchfork boundary.

Concluding remarks
Current work addressed an infinite symmetrical chain of coupled oscillators, subject to the selfexcitation. Analytical solutions for the DBs were derived for two qualitatively different dynamical regimes. The first solution was derived for the DB in the attenuation zone of the chain and is characterized by exponential decay of the amplitude to zero. The second solution has the main frequency in the PZ of the chain; the other harmonics belong to the AZ and are exponentially localized. In this regime, the energy input due to the self-excitation is balanced not only by the inelastic impacts, but also by the energy irradiation into the chain.
The loss of stability for both solutions in the vicinity of the propagation zone occurred through the Neimark-Sacker bifurcation; it is reasonable to assume that the Neimark-Sacker bifurcation is a generic scenario that is triggered by interaction with the boundary of the propagation zone. Besides, the PZ breather has another stability boundary, where both Neimark-Saker and pitchfork bifurcations are encountered. It seems reasonable to suggest that this stability boundary is not generic and appears due to particularities of the forcing mechanism.
Regime of the self-excitation in the single-DOF oscillator (Figure 1) is stable for all combinations of parameters. Expansion of the model to a system of coupled self-excited oscillators brings about the nontrivial stability patterns. Therefore one can control the charge transport in the system of coupled oscillators by modification of the parameters. This possibility outlines possible applications for the localized solutions revealed in the paper.
The authors are grateful to Israel Science Foundation (grant 838/13) for financial support.