Frequency divider by underdamped Josephson transmission line

AbstractA parallel array of coupled short Josephson junctions
under boundary ac driving is theoretically investigated. The
driving response is monitored via numerical simulations on the
associated Frenkel-Kontorova model as an oscillating output
current at the other boundary junction. The theoretical analysis
is conducted assuming the phase locking of boundary driving with
exact stationary solution of sine-Gordon equation that can be
represented by a kink (fluxon) motion back and forth in the
restricted geometry of the array. As a consequence the results
suggest the conception of a device which could switch to the
conducting regime not according to intensity range but rather to
given (quantized) frequency range. Moreover our findings indicate
also a frequency converting scenario by choosing appropriately the
system length and injected fluxon number.



Introduction
A subject of the present paper is a theoretical investigation of boundary driven stationary regimes in parallel array of small Josephson junctions (often referred as Josephson transmission line (JTL)) which is successively described in terms Frenkel-Kontorova model or its continuous analogy (sine-Gordon equation) [1,2]. At the same time, these models have reappeared in a great number of physical contexts, including description of domain walls motion in ferro-and antiferromagnetic crystals, charge-density waves in solids, pseudospin model of quantum Hall bilayers and even nonlinear excitations in hydrogen-bonded molecules (see e.g. Refs. [3][4][5] and references therein). Thus the results derived in the present paper could be readily extended to those completely different physical systems. In particular, in the present paper we show the possibility of the division of the driving frequency by odd fractions, in other words, applying a fixed frequency (3Ω) driving at one boundary of JTL we shall demonstrate by numerical simulations the output current oscillations with frequency Ω at the other boundary. Such a regime is argued analytically as a phase locking of boundary excitation with an exact stationary solution of associated sine-Gordon equation describing fluxon motion forth and back inside the JTL.
The analysis of the fluxon motion in JTL mostly appears in literature in the context of ac or dc homogeneous driving by a bias current. One can quote the perturbation analysis [6] and modulation instability studies [7] of a e-mail: khomeriki@hotmail.com the process, theoretical [8,9] and experimental [10,11] investigations of the resonances between fluxons and linear plasma oscillations, fluxon dynamics studies in case of coupled JTL-s [12,13], fluxon induced directed transport under homogeneous ac driving with broken time symmetry [14] and many other interesting phenomena. At the same time, much less studies have been devoted to the dynamics of nonlinear excitations under the local (boundary) driving: In this connection one should mention the investigation of hysteretic regimes appearing due to the phase locking effect between the boundary driving and breather type localizations [15][16][17] (where the breathing frequency coincides with the driver frequency [18][19][20]) or, on the other hand, boundary driving can maintain stationary motion of the fluxon which creates a new frequency in the system being an odd fraction of the driving frequency [21][22][23]. However, in all of the previous studies the driving is applied from the both ends and in order to observe the effect one has to measure the averaged quantities or some emitted radiation. In the present paper I consider the scenario when the ac magnetic field is applied only at one boundary (input) while one can monitor the current at the other boundary (output) as displayed on the schematics of Figure 1. Then in stationary regime one can get the output signal with three times smaller frequency than at the input realizing thus frequency dividing device. For longer chain when the frequency division could not be observed because of the damping, one can think about a high precision frequency filter, since for given parameters of the junction the stationary fluxon motion could be induced by very narrow range of driving frequency. The aim of the present paper is to represent the analytical  . The response is detected at the right end as an output current Iout(t). Inset (the upper graph) shows the possibility of phase locking mechanism between different fractional frequencies. As it is clarified in the text below, only the locking with odd fractional frequencies could be realized.
description of such a stationary regime assuming the external driving locking with the exact stationary solution of sine-Gordon model and compare it with numerical simulations.

Analytical stationary solution for fluxon motion
Let us consider one dimensional array of N short Josephson junctions coupled through super-conducting wires as represented by Figure 1. It obeys the following set of equations [4,10] for the phase differences u n across n-th junctions: where n = 2, . . . , N − 1 and one has the input signal current I S (t) (induced by microwave radiation S(t)) at the left junction of the array and detects I out (t) at the right end. The time is normalized to the inverse plasma frequency ω p = 1/ √ L J C, C stands for the junction capacitance and L J = /(2eI c ) is the Josephson inductance and I c is a Josephson critical current in the single junction. The parameter λ J is defined by λ 2 J = L J /L S where L S is the inductance representing by the superconducting wires connecting the junctions. γ J = /(2eI c R 2 J C) is a damping parameter and R J the junction resistance. In typical experiments on Josephson junction parallel array [12], the parameters have the values: R J ≈ 100 Ω, C ≈ 300 fF, I c ≈ 10 μA, L S ≈10 pH and thus γ J ≈ 0.1 and λ J ≈ 3. However, by changing the critical current density and the temperature one can easily control the two parameters of the model γ J and λ J . For our numerical simulations we choose λ J = 3 and γ J = 0.02.
One can easily define the measured quantities from the phase variables u n in equations (1), particularly, longitudinal current I n (flowing from junction n − 1 to n) and the voltage across n-th junction could be expressed as: Now let us consider the ac driving source at the left in the form where the amplitude a is connected with a magnetic and/or electric coupling strength [24]. At the right, as far as I out (t) is unknown variable, in addition one has to use an ohmic law V N (t) = RI out (t), where R is the resistance of the external circuit (see Fig. 1). Then using the second equation from (2) one automatically gets: where γ = ω p L s /R is a damping parameter associated with external circuit. Numerical simulations are made on the model equations (1) putting there the relations (3) and (4). Making simulations for the different input parameters a and Ω and also in case of different length N of the system I monitor corresponding output current according to the relation (4).
For the analytical consideration we rewrite (1) in the form of Frenkel-Kontorova model [3]): where n = 1 . . . N and the following definitions for virtual junctions u 0 and u N +1 have been introduced according to the relations (3) and (4): Now let us work in the limit when the external resistance is small (R → 0), consequently γ → ∞ and it is easy to see from (6) that one has a vanishing boundary condition at the right endu N = u N = 0. Then taking into account that for strong coupling constants λ 2 J 1 one can proceed with continuous approximation of (5), i.e. introducing continuous n via u n ≡ u(n) and rewriting the combination u n+1 + u n−1 − 2u n ≈ ∂ 2 u(n)/∂n 2 , we get in case of negligible damping parameter γ J the equation with dirichlet boundary condition at the right end and Neumann condition at the left: where the boundary condition at the left automatically follows from the first relation of (6); time is scaled again in inverse plasma frequency units, space variable x is defined as x = n/λ J and L = N/λ J . Then the longitudinal current and voltage could be expressed according to (2) as follows: Our approach is based on the assumption that the boundary driving is locked to the exact stationary solution of the sine-Gordon equation (7). Using the analogical approach as in reference [4] one can find the following periodic stationary solution vanishing at x = L: where cn and sn are the standard Jacobi elliptic functions [25], K(μ) stands for a complete elliptic integral of the first kind of modulus μ and the solution is parameterized by the two free constants ω and ν ∈ [0, 1], the remaining parameters are given by Further let us assume that the frequency of the time dependent part T (t) of the stationary solutions (9) coincide with an odd integer fractions of the boundary driving frequency Ω. It should be especially mentioned that stationary motion in the system with some specific frequency could be excited only by excitations characterized with odd fractions of that frequency (and not even fractions). The schematics for justification of this statement is displayed in the inset of Figure 1. Indeed the positions of maximums and minimums of 3Ω and Ω coincide with each other, while this is not the case for the frequency oscillations 2Ω and Ω. This is why only the locking with only odd fractional frequency oscillations takes place in numerical simulations with (5).

Frequency division scenario
Recalling that the period of T (t) is 4K(ν)/ω we require thus ω = 2ΩK(ν)/(mπ), where m is an odd integer. For a given value of the parameter ν ∈ [0, 1], the above relation fixes the second parameter ω in terms of the driving frequency Ω and only one free parameter ν is left. Therefore fixing Ω (driver frequency) and varying ν one can plot the output current amplitude in terms of the input signal amplitude from the analytic expressions (9). Let us just note that signal and output current could be calculated according to (8) from (9)  use the value a = 0.2 which is depicted as vertical dashed line in main plot of Figure 2. However as it is seen from Figure 2 there are multiple possible solutions for ν from which (that is seen from numerical simulations) only two ones indicated by crosses are stable: the lower intersection point corresponds to the negligible output current amplitude (see the lower inset in the same Figure 2 for the time evolution of input-output current) and the upper intersection is responsible for the stationary fluxon motion regime (see the input-output current time evolution in the upper inset of Fig. 2). This latter solution (the space time dependence of which is displayed in the graph (a) of Fig. 3) is a subject of our interest and we will compare that analytical solution with the results of numerical simulations.
The numerical simulations have been made for the Frenkel-Kontorova chain (5) with the parameters γ J = 0.02, λ J = 3 and the conditions at the extremities (6) where I take the driving amplitude a = 0.2 and frequency Ω = 0.1843 (like in analytical calculations) and external damping γ = 50; the junction number is N = 18 such that the effective length L = N/λ J appears the same as in analytical calculations. Then I monitor longitudinal current according to the relation (2) and input and output currents are calculated using expressions (3) and (4), respectively.
As a result in the graph (b) of Figure 3 the space time dependence of the longitudinal current is displayed. As seen the graph is almost identical with the analytical solution presented in graph (a) of the same figure. The only difference is the method of driving: in numerical simulations I use harmonical driving while the analytic solution requires rather complicated time dependence of the input ∂ x u(0, t) (please compare dotted curves in the insets of the graphs (a) and (b)). However this difference does not play

Conclusions and comments
Concluding, it is suggested how the frequency of input source could be divided by odd fractions in underdamped parallel array of short Josephson junctions. As further numerical simulations show it is possible to change the range of frequencies which could be divided by changing the point of the connection of the external circuit, effectively decreasing by this the system length. For shorter chain of Josephson junctions one can achieve higher division numbers (e.g. the frequency could be divided by factors m = 1, 3, 5, 7...) realizing thus the frequency dividing device.
It should be mentioned here that one can detect the divided frequency measuring emitted radiation coming from the whole array of Josephson junctions. In the definitions of (5) the emitted power from the whole chain could be expressed as N −1 n=1 γ J (u n ) 2 , while the power release form the right extremity (see Fig. 1 The typical dependences of those powers versus time are displayed in Figure 4. As seen power released locally at the chain extremity exceeds one from a total chain. Note that part of the total power goes to the internal losses and detected radiation will be less than the value depicted as a dashed line in Figure 4. Besides that one can not increase γ J sufficiently in numerical simulations as far as this will cause destruction of the stationary motion of fluxon. Therefore it seems much more efficient to detect frequency division locally at the extremity (solid line in Fig. 4) than measuring the emitted radiation.
For longer chain (N > 24) the subharmonic locking can not provide the long living stationary regime and only fundamental frequency can support the motion of the kink forth and back. Thus for such a chain one can think about a frequency filter, i.e. only the frequencies from very narrow range will provide the transmission regime. It should be especially mentioned that in all the considered cases the driving frequency lies within the forbidden band gap, therefore one can switch the transmitting regime by injecting a fluxon from the right end of the chain and then support its motion by the driving from the left.
Obviously the analytical model (7) should work for long Josephson junction, although the numerical simulations require the proper discretization of the continuous model. Besides that, it seems more appropriate to work with discrete chain of junctions since with the latter model it will be easier to decouple input and output from the experimental point of view. Variation of the system parameters will affect the division scenarios, particularly, it is easy to see from (1) and its analytical approximation (7) that if one decreases parameter λ J this will lead to the rescaling of the effective length of the chain (it will increase) and for the same number of junctions in (1) we will have a "onger" system (7). Thus for the junctions with large critical currents I c one should take smaller number of junctions to observe frequency division scenario. We performed the numerical simulations also for this highly discrete case and although the effect is still there, one can not use the analytical model (7) since it works only in continuous approximation.
Considered effect is not restricted only to the sine-Gordon equation but generally speaking it can characterize any nonlinear system which permits the existence of moving localized solutions, particularly, the resent experimental observation of the similar stationary regime as a motion forth and back of solitary shallow water waves could be quoted as an example [26].