Devil's staircase structures in $\varphi_{0}$ junction

TL;DR

This paper investigates the emergence of fractional subharmonic steps in the IV-characteristics of a $\

Contribution

It reveals the origin of subharmonic steps in $\

Findings

Subharmonic steps follow a continued fraction pattern.

Steps are caused by magnetic moment precession locking.

External radiation induces fractional steps in the junction.

Abstract

The superconductor-ferromagnet-superconductor $\varphi_{0}$ junction provides a direct coupling between Josephson phase and magnetic moment of ferromagnetic barrier. We demonstrate an appearance of additional fractional subharmonic steps in the IV-characteristics of $\varphi_{0}$ junction under external electromagnetic radiation due to spin-orbit coupling. An origin of subharmonic steps is related to the locking of magnetic moment precession to the Josephson oscillations. We prove that the positions of those steps follow a continued fraction algorithm.

Equations

$$I+A\sin \Omega t$$

$$\displaystyle\frac{d\textbf{m}}{dt}=-\frac{\Omega_{F}}{(1+\alpha^{2})}\bigg{(}\textbf{m}\times\textbf{}\textbf{h}_{eff}+\alpha\left[\textbf{m}\times(\textbf{m}\times\textbf{h}_{eff})\right]\bigg{)},$$

$$E_s$$

$$E_M$$

$${\text{h}}_e=[Gr\sin (\phi -r{m}_y)\widehat{{\text{e}}_y}+m_z\widehat{{\text{e}}_z}],$$

$$V=\left(N\pm \frac{1}{n\pm \frac{1}{m\pm \frac{1}{p\pm ...}}}\right)\Omega$$

$$m_y(t)$$

$$\phi (t)=\sum _n{ϵ}^n{\phi }_n(t),I=\sum _{n=0}^{\infty }{ϵ}^n{I}_n.$$

$$\dot{\phi }(t)$$

$${\dot{\phi }}_n(t)=I_n+f_n(t).$$

$$f_0(t)$$

$$f_1(t)$$

$${\dot{\phi }}_0(t)$$

$$I_s^{(0)}$$

$${\phi }_0(t)$$

$$I_s^{(0)}$$

$$\frac{d\mathbf{M}}{dt}=-\gamma \mathbf{M}\times {\mathbf{H}}_e+\frac{\alpha }{\text{Arrowvert}\mathbf{M}\text{Arrowvert}}\left(\mathbf{M}\times \frac{d\mathbf{M}}{dt}\right),$$

$${\mathbf{H}}_e={\mathbf{H}}_0+\tilde{\mathbf{H}},\mathbf{M}={\mathbf{M}}_s+\tilde{\mathbf{M}},$$

$$\frac{d\tilde{\mathbf{m}}}{dt}+{\Omega }_F\tilde{\mathbf{m}}\times {\mathbf{h}}_0+\alpha \left(\frac{d\tilde{\mathbf{m}}}{dt}\times {\mathbf{m}}_s\right)=-{\Omega }_F{\mathbf{m}}_s\times \tilde{\mathbf{h}},$$

$$i{\Omega }_p\text{m}+{\Omega }_F\text{m}\times {\mathbf{h}}_0+i{\Omega }_p\alpha \text{m}\times {\mathbf{m}}_s=-{\Omega }_F{\mathbf{m}}_s\times \text{h}.$$

$$Re\{m_y(t)\}$$

$$m_y(t)=\frac{{\tilde{\gamma }}_2\sin \phi (t)}{D+r{\tilde{\gamma }}_2\cos \phi (t)},$$

$$m_y(t)$$

$$\Omega >>1,{\beta }_c\Omega >>1,A>>1$$

$$\dot{\phi }(t)$$

$$\phi (t)=\sum _n{ϵ}^n{\phi }_n(t),I=\sum _{n=0}^{\infty }{ϵ}^n{I}_n.$$

$${\dot{\phi }}_n(t)=I_n+f_n(t),$$

$$\displaystyle\bigg{[}\left(\frac{r\tilde{\gamma}_{2}}{2D}\right)^{2}-1\bigg{]}\sin\varphi(t)+\frac{r\tilde{\gamma}_{2}}{2D}\bigg{(}\sin 2\varphi(t)-\frac{r\tilde{\gamma}_{2}}{2D}\sin 3\varphi(t)\bigg{)}\approx$$

$$f_0(t)$$

$$f_1(t)$$

$$f_2(t)$$

$${\dot{\phi }}_0(t)$$

Full text

I Introduction

Josephson junction (JJ) with ferromagnetic barrier opens an additional gate for investigation of the interplay between ferromagnetism and superconductivity. Its current-phase relation is very sensitive to the characteristic of the ferromagnet Braude2007 ; buzdin2008direct ; Linder2015 ; Silaev2017 ; Bobkova2017 ; Rabinovich2018 . Particularly, the superconductor-ferromagnet-superconductor (SFS) Josephson junction oscillates between [math]- and $\pi$ states with an increase in ferromagnetic layer thickness Ryazanov2001 ; Oboznov2006 ; Robinson2006 ; Iver2008 . SFS JJ holds a great possibility for the realization of the qubit for quantum computation Yamashita2005 ; Massarotti2015 ; Yamashita2017 and can be used as electronic and spintronic components, such as phase shifters and superconducting spin valves Golod2010 ; Leksin2011 ; Gingrich2016 ; Soloviev ; Zhu2017 ; Niedzielski2018 ; Golod2019 .

In case of non-centrosymmetric magnetic barrier with broken inversion symmetry like $MnSi$ or $FeGe$, the spin-orbit coupling leads to a general current-phase relation $J=J_{c}\sin(\varphi-\varphi_{0})$, where $\varphi_{0}$ is proportional to the magnetic moment perpendicular to the gradient of the asymmetric spin-orbit potential buzdin2008direct ; buzdin2009 . In this case the spin-orbit coupling in $\varphi_{0}$-junction provides a direct coupling between the suppercurrent and the magnetic moment. This coupling leads to the supercurrent induced magnetization dynamics buzdin2009 ; Linder2014 ; Chudnovsky2016 . Recently, in Ref.Assouline2019 the authors report the observation of anomalous phase shift $\varphi_{0}$ in $Bi_{2}Se_{3}$ Josephson junctions (JJ) and provide a direct measurement of the spin-orbit coupling strength. Also, in $\varphi_{0}$ junction a full magnetization reversal is demonstrated by using electric current pulse shukrinov2017magnetization .

An increase of the Shapiro steps amplitude due to spin-orbit coupling near the ferromagnetic resonance, the appearance of the half-integer Shapiro steps, and precession of the magnetization vector with radiation frequency buzdin2009 are predicted, if the $\varphi_{0}$ Josephson junction is exposed to microwave radiation. An external electromagnetic field can control qualitative features of the magnetic moment dynamics in a current interval which corresponds to the Shapiro step. The radiation can also produce topological transformations of precession trajectories cond-mat18 , which depend on the amplitude of the applied electromagnetic radiation. Such results might be used for developing novel experimental resonance methods for determination of the spin-orbit interaction in the non-centrosymmetric materials.

Generally, the IV-characteristic of a conventional JJ in the underdamped case demonstrates harmonic and subharmonic steps in the presence of external electromagnetic radiationshukrinov2013devil ; azbel , while it shows harmonic steps only in the overdamped junctions RENNE ; WALDRAM . The steps for the underdamped JJ can form the so-called devil’s staircase (DS) structure as a consequence of the interplay between Josephson and applied frequencies ben-jacob ; shukrinov2013devil ; shukrinov2014structured ; sokolovic2017devil . The DS structure is a universal phenomenon and appears in a wide variety of different systems, including infinite spin chains with long-range interactionsnebendahl13 , frustrated quasi-two-dimensional spin-dimer systems in magnetic fieldstakigawa13 , and even in the fractional quantum Hall effecthriscu13 . In Ref.chen2017 a series of fractional integer size steps is observed experimentally in the Kondo lattice CeSbSe.

Coupling between the superconducting current and magnetic moment in the SFS Josephson junction is one of the most interesting topics nowadays Linder2015 . The connection between the staircase structure and current-phase relation provides an intriguing point in this field of science. Especially, the manifestation of the staircase structure in the IV-characteristics of the junction with the corresponding information on current-phase relationGolubov ; Sellier2004 ; Maiti2015 ; Pic2017 ; Shukrinov2018 , thus, serves as a novel method for its determination. The appearance of the DS structure and its connection with the current-phase relation in experimental situations still in need for detailed investigations.

In Ref.Shukrinov2018 we have demonstrated an appearance of subharmonic steps in the IV-characteristic for overdamped SFS junction with spin wave excitations (magnons). It was found that the width of the current steps at $V=2\Omega$ (where $\Omega$ is the frequency of the applied magnetic field) are larger than the width of the odd and fractional steps. Therefore, magnetization dynamics can be manifested in the IV-characteristic since this behavior is different from that of the conventional Shapiro stepShukrinov2018 ; hikino2011ferromagnetic . The origin of the even steps is related to the interaction of Cooper pairs with even number of magnons hikino2011ferromagnetic . In Ref.Mori2014 the authors show that the breathing of the domain wall leads to the appearance of staircase structure in the IV-characteristic of SFS junction, but they did not investigate its relation with DS structure. The IV-characteristics of superconductor-quantum spin Hall insulator-superconductor system Hang2018 in the presence of microwave radiation exhibit a structure where odd steps are completely suppressed, implying a fractional Josephson effect.

Analytic treatment of the coupled JJ-nanomagnet system driven by a time-dependent magnetic field both without and with an external ac field is studied in Ref.Roopayan2017 . The authors have shown the existence of Shapiro-type steps in the IV-characteristics of the JJ subjected to a voltage bias for a constant or periodically varying magnetic field and explore the effect of rotation of the magnetic field and the presence of an external ac drive on these stepsRoopayan2017 . Here we clarify the possibility for the appearance of fractional steps in the overdamped $\varphi_{0}$ junction using a modified RSJ model.

In this paper, we focus on the overdamped case ($\beta_{c}=0$, $\beta_{c}$ is the McCumber parameter) to reflect clearly the appearance of the subharmonic steps due to the presence of spin-orbit coupling in $\varphi_{0}$ junction. It is found that due to the coupling between Josephson phase and magnetic moment through the spin-orbit coupling, the additional fractional subharmonic steps appear in the IV-characteristic. The appearance and the positions of those steps depend directly on the spin-orbit coupling and the ratio of Josephson and magnetic energies. An analytical consideration of the linearized model justifies the appearance of the fractional steps in IV-characteristics, in agreement with our numerical results.

The plan of the rest part of the paper is as follows. In Sec.I, we describe the self-consistent modified RSJ and LLG dynamical equations for $\varphi_{0}$ junction. This is followed by a discussion of the IV-characteristics of such systems in Sec II. We demonstrate the effect of spin-orbit coupling in the IV-characteristic. Then, we confirm that the appearance of the subharmonic steps under external radiation are due to the presence of spin-orbit coupling. In addition to this we justify analytically the conditions of the frequency locking and discuss the experimental realization of the found effects. Finally, we conclude in Sec. III. Some details of our calculations are specified in the supplemental material.

II Model and Methods

The geometry of the considered $\varphi_{0}$ junction is shown in figure 1. The ferromagnet easy-axis and the gradient of the spin-orbit potential ($n$) are directed along the z-axis. In this case, $\varphi_{0}=rm_{y}$ and $r$ characterizes the relative strength of the spin-orbit interactionbuzdin2009 .

The electric current through JJs in the dimensionless form is determined by RSJ equation stewart1968current (we consider overdamped case with $\beta_{c}=0$ )

[TABLE]

where $I$ is normalized to the critical current $I_{c}$, $t$ is normalized to $\omega^{-1}_{c}$, $\omega_{c}$=$2eRI_{c}/\hbar$, $m_{y}=M_{y}/M_{s}$ satisfies the constraint $\sum_{i=x,y,z}m_{i}^{2}(t)=1$, $M_{s}=\|\textbf{M}\|$, and the voltage $V$ is normalized to $I_{c}R$. The external electromagnetic radiation is characterized by the amplitude $I_{ac}$ and frequency $\omega$, $A=I_{ac}/I_{c}$. The dynamics of the magnetic moment $m_{y}$ is determined by the LLG ounadjela2003spin :

[TABLE]

where $\Omega_{F}=\omega_{F}/\omega_{c}$, $\omega_{F}$ is the ferromagnetic resonance frequency, $\alpha$ is the Gilbert damping ( here we consider $\alpha=0.1$) and $\textbf{h}_{eff}$ is the total effective field of the system which is given by $-\bm{\nabla}_{M}E/v$, $v$ is the volume of ferromagnet. The total energy of the system for the current biased junction is $E=E_{s}+E_{M}$ with

[TABLE]

where $E_{J}=\Phi_{0}I_{c}/2\pi$ is the Josephson energy, $E_{M}$ is the anisotropy energy and $K_{an}=\omega_{F}M_{s}/\gamma$ is the anisotropy constant buzdin2009 ; shukrinov2017magnetization . So, the effective field in normalized form is given by

[TABLE]

with $G=E_{J}/K_{an}v$. The magnetic moment and the phase dynamics of the considered S/F/S Josephson junction is determined by Eqs. (1) and (2).

To compute the IV-characteristics, we study the temporal dependence of $V(t)=\hbar\dot{\varphi}(t)/(2e)$ obtained by numerical solution of equation (1). The dc bias current $I$ is normalized to $I_{c}$, and the voltage is normalized to $\hbar\omega_{c}/(2e)$. We assume a constant bias current and calculate the average voltage. We employ a fourth-order Runge-Kutta integration scheme. As a result, we find the temporal dependence of the voltage in the junction at a fixed value of bias current $I$. Then, the current value is increased or decreased by a small amount $\delta I$ (the bias current step), to calculate the voltage at the next point of the IV-characteristics. We use the final voltage achieved at the previous point of the IV-characteristics as the initial condition for the next current point. The average of the voltage $V_{av}$ is given by $V_{av}=\frac{1}{T_{f}-T_{i}}\int^{T_{f}}_{T_{i}}V(t){\rm d}t$, where $T_{i}$ and $T_{f}$ determine the interval for the temporal averaging. The initial conditions for the magnetization components are assumed to be $m_{x}=0$, $m_{y}=0$ and $m_{z}=1$, while for the voltage and phase we take zeros.

To analyze the positions of subharmonic steps in the IV-characteristics and reflect the DS structure, we use the continued fractions. They consist of number of levels. The first-level is denoted by N and the second-level gives two groups of subharmonics (N-1+1/n) and (N-(1/n)). The algorithm of continued fractions is presented in figure 2 where the red circles represent the Shapiro step number (first level). The green rectangles represent the second level subharmonic steps and the arrows represent the approaching direction, for example, the group $(N-1)+(1/n)$ approaches $(N-1)^{th}$ Shapiro step and $N-(1/n)$ approaches the $N^{th}$ Shapiro step. The third level subharmonics is represented by the blue diamond. This level is determined by fixing $n$, $n+1$ and changing $m$, for example, the second level group $(N-1)+(1/n)$ gives rise to two third level groups between $n=1$ and $n=2$, also two other third level group are a rise from the second group $N-(1/n)$ with $n=1$ and $n=2$, etc shukrinov2013devil ; shukrinov2014structured ; Shukrinov2018 . The positions of the current steps follow continued fraction formula shukrinov2013devil ,

[TABLE]

where N,n,m,p,… are positive integers. Terms with only N form harmonics, while other terms describe subharmonics or fractional steps.

In what follows, we first obtain an approximate analytical solution of equation (1) in Sec.B.1, which demonstrates the existence of subharmonic Shapiro steps for $\varphi_{0}$ junction. Then, in Sec.B.2, we carry out a detailed numerical study of equation (1) where the analytical results are verified and the devil’s staircase structure of the Shapiro steps is studied.

III Results and discussion

III.1 Perturbative analytical solution

First, we find the expression for $m_{y}$. If the deviation of the magnetic moment from the equilibrium point due to Josephson energy is small (i.e., $G<1$), we can linearize the LLG equation. In this case supplement , the magnetic moment $m_{y}(t)$ reads

[TABLE]

where $\tilde{\gamma}_{2}=Gr\gamma_{2}$, $\gamma_{2}=\left(1-(1-\alpha^{2})\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}\right)$, and $D=\left(1-(1+\alpha^{2})\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}\right)^{2}+4\alpha^{2}\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}$.

Then we demonstrate the results of the perturbative analysis of equation (1) for overdamped $\varphi_{0}$ junctions in high-frequency limit Kornev . We analyze this equation for $\Omega$ and $A>>1$ in which we can do the following expansions Shukrinov2011 ; Likharev

[TABLE]

where $I_{0}$ is the biased current, $\epsilon<<1$ and $I_{n}$ for $n>0$ are determined self-consistently from the condition of the absence of additional dc voltage: $\lim_{T\longrightarrow\infty}$ $\int_{0}^{T}$ $\dot{\varphi}_{n}dt=0$Likharev . Using equation (6), the RSJ up to terms $\sim O(rm_{y})$ has the form

[TABLE]

Using (7), the equations for $\dot{\varphi}_{n}$ can be obtained by equating terms in the same order of $\epsilon$. Then, the expression of $\dot{\varphi}_{n}$ is given by

[TABLE]

For $n=0$ and $n=1$, we have

[TABLE]

The $0^{th}$ order with $n=0$ represents the autonomous IV-characteristic of the junction. In this case we have

[TABLE]

while the supercurrent is given by

[TABLE]

After integrating equation (11) with respect to time, we have

[TABLE]

Inserting equation (13) into equation (12) and after some algebra supplement , we come to:

[TABLE]

Shapiro step within this order appears when the ac component of the supercurrent vanishes. This means we have three possible cases; $I_{0}=n\Omega$, $I_{0}=n\Omega/2$, and $I_{0}=n\Omega/3$. Furthermore, in the supplement we show that the $1^{st}$ order leads to fractional subharmonic steps which occur at different set of integers. These locking conditions appear only in the presence of spin-orbit coupling.

So, we find the complete frequency locking conditions for subharmonic steps, which are shown in figure 3. In the next section, we demonstrate these conditions for some steps.

III.2 DS structure in the IV-characteristics of $\varphi_{0}$ junction

Here we show the DS structure in the IV-characteristics of $\varphi_{0}$ junction under external electromagnetic radiation with $A=0.5$ and $\Omega=0.5$ and prove their correspondence to the continued fraction formula (5). Figure 4(a) demonstrates two IV-characteristics without and with spin-orbit coupling. At $r=0$ the IV-characteristic shows only harmonic Shapiro steps at $V=n\Omega$ with $n$ integer. However, the additional fractional subharmonics appear at $r=0.5$ between the harmonic steps as a result of spin-orbit coupling. Figures  4(b) and 4(c) demonstrate the enlarged parts of the IV-characteristic shown in figure 4(a). We see in figure 4(b) the fractional steps between $V=0$ and $V=0.5$ which can be described by the continued fractions of second level $(N-1)+1/n$ and $N-1/n$ with $N=1$ in both cases. In figure 4(c) we see the manifestation of second level continued fractions $N-1/n$ and $(N-1)+1/n$ with $N=2$ between voltage steps $V=0.5$ and $V=1$.

III.2.1 Effect of $r$ on the width of the subharmonic steps

Let us now inspect the effect of spin-orbit coupling on the presence of subharmonic steps in details. Figure 5(a) compares the IV-characteristics at different value of $r$ and we see that the enhanced subharmonic steps appear at large value of $r$.

In figure 5(b) we compare the harmonic Shapiro steps at $V=\Omega$ and $V=2\Omega$ at different value of $r$. We find qualitatively different behavior: the width of the first harmonic Shapiro step at $V=\Omega$ decreases slightly with increasing $r$, while the Shapiro step at $V=2\Omega$ demonstrates some horizontal shift with almost the same width for different values of $r$.

The corresponding width dependence of these steps as a function of amplitude of external radiation $A$ is shown figure 5(c). For the given simulation parameters, we observe a significant change of the Bessel dependence of the Shapiro step width with increase in $r$ in the range $A\sim 0.5-2$.

However, as we notice above in figure 5(a), the width of fractional Shapiro steps increases with $r$. The $r$-dependence of the harmonic step width at $V=0.5$ and subharmonic step at $V=0.75$ is shown in figure 5(d). We see the qualitatively different behavior: harmonic step width is decreased, while a reversal dependence occurs for subharmonic step.

The perterbative analysis shows that the width of all steps, harmonic and subharmonic, is not Bessel function of $r$ supplement , moreover, we have found that the width of subharmonic steps is proportional to $B_{x,y}$=$J_{x}(kA/\Omega)J_{y}(A/\Omega)$, where $k=1,2,3$.

III.2.2 Effect of $G$ on the appearance of the DS structure

Another important parameter which can control the appearance of the subharmonic steps is a ratio of Josephson and magnetic anisotropy energies $G=E_{J}/Kv$. In Ref.Shukrinovepls2018 it is shown that a reorientation of easy axis from $m_{z}=1$ to $m_{y}=1$ occurs in the regime of small Josephson frequency and large $G>20$. Here we show that the subharmonic steps are enhanced in the regime before a complete reorientation occurs.

In figure 6(a) and (b) we show the enlarged parts of the IV-characteristics at $G=0.2$ and $G=5\pi$. At $G=0.2$ the subharmonic steps appear mostly between $V=0$ and $V=0.5$ [see figure 6(a)], while much less number of steps with smaller width appear between $V=0.5$ and $V=1$ [see figure 6(b)]. All these fractional steps disappear at $G=5\pi$ when a reorientation of easy axis occurs to $m_{y}(t)=1$ (see insets).

Figure 6(c) shows the temporal dependence of $V(t)$ and $m_{y}(t)$ at $I=1.16$ (step at $V=\Omega$). The temporal dependence for both $V(t)$ and $m_{y}(t)$ is regular. Results of the corresponding Fast Fourier Transform (FFT) analysis for $V(t)$ and $m_{y}(t)$ are shown in figure 6(d). From the FFT analysis it is clear that the oscillation frequency of $m_{y}(t)$ is locked to the external frequency $\Omega=0.5$ as well as to $\Omega_{J}=0.5$. In addition to this, the FFT analysis for $V(t)$ show harmonics of $n\Omega$ with integer $n$.

Next, we have recorded the temporal dependence of $V(t)$ and $m_{y}(t)$ at $I=1.284$ for fractional step at $V=0.75$, presented in figure 6(e) and the corresponding result of the FFT analysis is shown in figure 6(f). Here the temporal dependence for both $V(t)$ and $m_{y}(t)$ is more complex. The oscillation amplitude for $m_{y}(t)$ is smaller than in case $V=0.5$ [see figure 6(c)]. The FFT analysis show that a locking between the precession frequency for $m_{y}$ and Josephson frequency occurs at fractions of $\Omega$.

Based on the presented results, we can come to the following conclusions. The first is locking of the magnetic moment precession at $I=1.16$ and $G=0.5$ to the external radiation with frequency $\Omega=0.5$. At $I=1.284$ we also observe the locking of magnetic moment precession to the subharmonic frequency $3\Omega/2$ which correspond to the step at $V=0.75$. The appearance of Shapiro steps at $V=n\Omega$ is not related to the magnetic moment precession. Shapiro step at $V=0.5$ appears at $G=5\pi$ with $m_{y}(t)=1$ when the precession frequency is zero.

Finally, we present some estimations for the experimental realization based on data in Ref’s.buzdin2009 ; Assouline2019 ; exper_data1 ; exper_data2 ; exper_data3 ; exper_data4 . The main parameter which controls the appearance of the current steps is $G=E_{J}/(Kv)$ which represents the ratio of the Josephson and the magnetic anisotropy energy. With flux quantum $\Phi_{0}=2.067833\times 10^{-15}(J/A)$, critical current $I_{c}=4\mu A$, volume of the ferromagnet $v=d\times l\times w\sim 150nm\times 2\mu m\times 2\mu m$, and anisotropy constant $K=500$ ($J/m^{3}$), we obtain $G=\Phi_{0}I_{c}/(Kv)=0.4$. Actually, for stronger anisotropy, parameter G is smaller. Play with the value of anisotropy parameter, we may have G as much less or much large than one.

IV Conclusion

In this work we solve the dynamical equations for $\varphi_{0}$ junction which describe the coupling between Josephson phase and magnetic moment through spin-orbit coupling. Using the high-frequency perterbative method we predict the appearance of subharmonic steps in the presence of spin-orbit coupling. Furthermore we confirm these prediction using the exact numerical scheme. The IV-characteristic of the $\varphi_{0}$ junction demonstrates devil’s staircase structures of Shapiro steps. The origin of the found steps related to the synchronization between Josephson and magnetic moment oscillations. The position of the found steps follow continued fraction formula. The structure and width of these steps depends on the ratio of Josephson and magnetic energies and spin-orbit coupling. Our calculations predict that spin-orbit coupling manifests itself through the appearance of subharmonic steps in the IV-characteristics, thus providing new insight into the precise nature of the current-phase relation and new opportunities for potential applications.

V Acknowledgments

The authors thank I. Rahmonov and A. Mazanik for fruitful discussions. The reported study was partially funded by the RFBR research projects 18-02-00318, 18-52-45011-IND. Numerical calculations have been made in the framework of the RSF project 18-71-10095.

VI References

Supplemental Material to “Devil’s staircase structure in $\varphi_{0}$ junction”

I Linearized Landau-Lifshitz-Gilbert equation

If the deviation of the magnetic moment from the equilibrium point due to Josephson energy is small( $G<1$), we can linearized the Landau-Lifshitz-Gilbert (LLG) equation. In its general form the LLG equation reads

[TABLE]

where $\alpha$ is the Gilbert damping and $\gamma$ is the gyromagnetic ratio. We assume that the effective magnetic field and magnetization can be written as sums of constant and alternating parts

[TABLE]

where the components of $\bm{H}_{0}$ are $(0,0,H_{0})$ with $H_{0}=K_{an}/M_{s}=\omega_{F}/\gamma$, $K_{an}$ is the magnetic anisotropy constant, $M_{s}=\Arrowvert\textbf{M}\Arrowvert$ is considered as a constant and equal to the saturation value of the magnetization, $\omega_{F}$ is the ferromagnetic resonance frequency and $\gamma$ is the gyromagnetic ratio. The components of $\bm{\tilde{H}}$ are $(\tilde{H}_{x},\tilde{H}_{y},0)$, the components of $\bm{M}_{s}$ are $(0,0,M_{z})$ and those of $\bm{\tilde{M}}$ are $(\tilde{M}_{x},\tilde{M}_{y},0)$. The magnitude of alternating parts are considered smaller than the steady parts, i.e. $\tilde{H}<<H_{a}$, $\tilde{M}<<M_{z}$. In what follow we use the dimensionless formula. We normalize $H_{e}$ to $H_{0}$, $M_{i}$ ($i=x,y$) is normalized to $M_{s}$, time is normalized to $\omega_{c}^{-1}$ ($t\longrightarrow t\omega_{c}^{-1}$) where $\omega_{c}$ is the characteristic frequency of JJ, and $\omega_{F}$ is normalized to $\omega_{c}$. The linearization of (1) can be found by inserting (2) into (1) and neglecting the products of the alternating parts. In the dimensionless form the linearized LLG is then reads as

[TABLE]

where the components of $\bm{h}_{0}$ are $(0,0,1)$, $\bm{\tilde{h}}$ are $(h_{x},h_{y},0)$, $\bm{m}_{s}$ are $(0,0,1)$, and $\bm{\tilde{m}}$ are $(m_{x},m_{y},0)$. We assume a harmonic time dependence for the effective field $\bm{\tilde{h}}$ and $\bm{\tilde{m}}$ in the form of $\bm{\tilde{m}}=\textit{{m}}\ e^{i\Omega_{p}t}$ and $\bm{\tilde{h}}=\textit{{h}}\ e^{i\Omega_{p}t}$, where $\Omega_{p}$ is the precession frequency normalized to $\omega_{c}$ ($\Omega_{p}=\omega_{p}/\omega_{c}$). According to this, the linearized LLG reads as

[TABLE]

Projecting Eq. (4) on the axes of Cartesian coordinate system we find the real part for $m_{y}(t)$ (we will take $\tilde{m}_{y}\equiv m_{y}(t)$, $\tilde{h}_{i}\equiv h_{i}(t)$ and $i=x,y$.)

[TABLE]

where $\gamma_{1}=2\alpha\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}$ and $\gamma_{2}=\left(1-(1-\alpha^{2})\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}\right)$. To find $m_{y}(t)$, we need the x and y-component of the effective field. Here we have $h_{x}=0$ and $h_{y}=Gr\sin(\varphi-rm_{y})$. Since $m_{y}<<1$, we can write $\sin(\varphi(t)-rm_{y}(t))\approx\sin\varphi(t)-rm_{y}(t)\cos\varphi(t)$. Using Eq.(5), the expression of $m_{y}(t)$ is given by

[TABLE]

where $\tilde{\gamma}_{2}=Gr\gamma_{2}$, and $D=\left(1-(1+\alpha^{2})\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}\right)^{2}+4\alpha^{2}\frac{\Omega_{p}^{2}}{\Omega_{F}^{2}}$. After using the expansion $1/(1+x)\approx 1-x+..$, we come to

[TABLE]

II Origin of subharmonic step in $\varphi_{0}$ junction

To find the origin of subharmonic step, we use high-frequency limit Kornev :

[TABLE]

where $\Omega$ is the frequency of the external electromagnetic radiation normalized to $\omega_{c}$ (characteristic frequency of JJ), $\beta_{c}$ is the McCumber parameter, and $A$ is the amplitude of the external electromagnetic radiation normalized to $I_{ac}$.

Using Eq.(7), the RSJ equation reads as

[TABLE]

Expend $I$ and $\varphi$ as Shukrinov2011

[TABLE]

where $I_{0}$ is the bias current, $\epsilon<<1$ and $I_{n}$ for $n>0$ are determined from the condition of the absence of additional dc voltage: $\lim_{T\longrightarrow\infty}$ $\int_{0}^{T}$ $\dot{\varphi}_{n}dt=0$Likharev . Using (10), the equations for $\dot{\varphi}_{n}$ can be obtained by equating terms in the same order of $\epsilon$. Then, the expression of $\dot{\varphi}_{n}$ is given by

[TABLE]

The second and third terms in Eq.(9) after using $\varphi(t)=\varphi_{0}(t)+\epsilon\varphi_{1}(t)$ (to the first order) read as

[TABLE]

Next, we consider

[TABLE]

The zeroth order with $n=0$ represents the autonomous IV-characteristic of the junction. In this case

[TABLE]

and the suppercurrent is given by $I^{(0)}_{s}=-f_{1}(t)$ (see Eq.(11))

[TABLE]

After integrating (14) with respect to time, we come to

[TABLE]

Next, we insert (16) into (15) and use $\sin(-x)=-\sin x$. The suppercurrent is given by

[TABLE]

and using $e^{iz\cos\theta}=\sum_{n=-\infty}^{\infty}i^{n}J_{n}(z)e^{in\theta}$, $J_{n}(z)$ is the Bessel function of first kind, and n is integer, we come to

[TABLE]

We see that Shapiro step appears when the ac component of the supercurrent vanishes. This means we have three possible cases: $I_{0}=n\Omega$, $I_{0}=n\Omega/2$, and $I_{0}=n\Omega/3$.

Now, for the first order we have

[TABLE]

which has the solution $<\dot{\varphi}_{1}>=0$ and $I_{1}=0$, so $\dot{\varphi}_{1}(t)=I_{s}$. Next, we get the imaginary part of (18) by separating odd and even parts using

[TABLE]

with $i^{2n}=(-1)^{n}$, $Im\{i^{2n}e^{i\theta}\}=(-1)^{n}\sin\theta$, $Re\{i^{2n}e^{i\theta}\}=(-1)^{n}\cos\theta$, $Im\{i^{2n+1}e^{i\theta}\}=(-1)^{n}\cos\theta$, and $Re\{i^{2n+1}e^{i\theta}\}=(-1)^{n+1}\sin\theta$. The imaginary part reads

[TABLE]

where $\Omega_{1}=2n\Omega-I_{0}$, $\Omega_{2}=2n\Omega-2I_{0}$, $\Omega_{3}=2n\Omega-3I_{0}$, $\tilde{\Omega}_{1}=(2n+1)\Omega-I_{0}$, $\tilde{\Omega}_{2}=(2n+1)\Omega-2I_{0}$, and $\tilde{\Omega}_{3}=(2n+1)\Omega-3I_{0}$. By integrating Eq.(21) with respect to time, we get

[TABLE]

For the $1^{st}$ order, the supercurrent is given by ($I^{(1)}_{s}=-f_{2}(t)$)

[TABLE]

Now, we need to insert Eq.(22) into Eq.(23) using

[TABLE]

where $q=1,2,3$, $\Omega_{m}=2m\Omega-I_{0}$,$\tilde{\Omega}_{m}=(2m+1)\Omega-I_{0}$ , and m is integer. The final expression for $I_{s}$ is given by

[TABLE]

with

[TABLE]

where $k$ takes the value of 1,2,3, $B^{1}_{x,y}$=$J_{x}(A/\Omega)J_{y}(A/\Omega)$,$B^{2}_{x,y}$=$J_{x}(2A/\Omega)J_{y}(A/\Omega)$ and $B^{3}_{x,y}=J_{x}(3A/\Omega)J_{y}(A/\Omega)$. The ac component of the supercurrent vanishes at $\Omega_{1}=\pm k\Omega_{m}$, $\Omega_{2}=\pm k\Omega_{m}$, $\Omega_{3}=\pm k\Omega_{m}$,$\tilde{\Omega}_{1}=\pm k\Omega_{m}$, $\tilde{\Omega}_{2}=\pm k\Omega_{m}$, $\tilde{\Omega}_{3}=\pm k\Omega_{m}$, $\Omega_{1}=\pm k\tilde{\Omega}_{m}$, $\Omega_{2}=\pm k\tilde{\Omega}_{m}$, $\Omega_{3}=\pm k\tilde{\Omega}_{m}$,$\tilde{\Omega}_{1}=\pm k\tilde{\Omega}_{m}$, $\tilde{\Omega}_{2}=\pm k\tilde{\Omega}_{m}$, $\tilde{\Omega}_{3}=\pm k\tilde{\Omega}_{m}$.

Finally, we confirm the appearance of subharmonic steps which are predicted from the perturbative analysis. This perturbative analysis can be applied when at least one of these condition is satisfied $A>>1$, or $\Omega>>1$, or $\beta_{c}\Omega^{2}>>1$ Kornev .

In figures 1(a) and (b) we show an analogous plot for the IV-characteristics of JJ without and with spin-orbit coupling ($r=0.3$) at $A=5$. As it is predicted by the perturbative analysis, a set of subharmonic steps appear when $r\neq 0$ (see figures 1(a) and (b)).