Phase slip lines in superconducting few-layer NbSe$_2$ crystals

We show the results of two-terminal and four-terminal transport measurements on few-layer NbSe$_2$ devices at large current bias. In all the samples measured, transport characteristics at high bias are dominated by a series of resistance jumps due to nucleation of phase slip lines, the two dimensional analogue of phase slip centers. In point contact devices the relatively simple and homogeneous geometry enables a quantitative comparison with the model of Skocpol, Beasley and Tinkham. In extended crystals the nucleation of a single phase slip line can be induced by mechanical stress of a region whose width is comparable to the charge imbalance equilibration length.

At sufficiently large current bias any superconductor (SC) shows dissipation. The ultimate limit is the pairbreaking current described by the Ginzburg-Landau theory. 1,2 Before that limit, dissipation may emerge as a result of the motion of Abrikosov vortices, which can be introduced by an external field or by the sample current itself. However, in the late '60s experiments on tin whiskers 3 showed that, even in absence of vortices, dissipation is observed before the pair-breaking current density is reached owing to the nucleation of phase slip centers (PSCs).
PSCs were first observed in 1D superconducting systems, where the transverse size is much smaller than both the coherence length ξ GL and the magnetic penetration depth λ. Accordingly, the models used to describe them were all inherently one dimensional. [4][5][6] It was soon found, however, that 2D systems can support an analogue of PSCs as well, the so called phase slip lines (PSLs). [7][8][9] The physical mechanism of the PSLs is markedly different from that of PSCs. It was found from numerical simulations based on time-dependent Ginzburg Landau (TDGL) equations 10,11 that PSLs are a sort of rivers of fast vortices and antivortices (called kinematic vortices 11 ) that annihilate in the middle of the sample. 12,13 Despite the physical difference, PSLs produce similar current-voltage characteristics (IVC) as PSCs, namely, discontinuous voltage jumps separating linear portions whose extrapolation on the current axis leads to a finite excess current. 9,[14][15][16][17] In ordinary polycrystalline or amorphous superconducting films, it is not easy to observe signatures of PSLs. The main experimental difficulty 15 is to prevent the formation of Abrikosov vortices by keeping a sufficiently high energetic barrier at the edges. 18 Moreover, it is crucial to ensure a good heat dissipation since at the large current densities close to the PSL nucleation Joule heating can easily smear the IVC. 9 All these limitations are absent in van der Waals SCs. 19,20 The nowadays standard exfoliation techniques make it possible to obtain clean monocrystalline devices consisting of one or few atomic layers, thus with a large Pearl screening length λ ⊥ = λ 2 /d. The edge roughness is on the nanometer scale, i.e. much sharper than the typical roughness of devices patterned by electron beam lithography (EBL). Encapsulation in hBN, besides preventing oxidation and contamination, provides an efficient thermal sink. 21,22 Therefore, one would expect that PSLs should dominate the IVC of van der Waals SC for current values well below the pair-breaking value.
The most studied van der Waals SC is the 2H polytype of NbSe 2 (in the following simply indicated as NbSe 2 ). [23][24][25][26][27][28][29][30][31][32][33][34][35][36][37] This material can be easily exfoliated in few-layer-thick flakes and then encapsulated in hBN. 25 Typical devices for transport measurements have a size of several micrometers. Several key results have been obtained in the recent years on this material. It was shown 25,28 that both the superconducting and the charge-density-wave phase survive even in monolayer crystals and that Ising superconductivity is observed owing to the large spin-orbit pseudo-magnetic field. 27,31 However, such experiments have focused on the low bias regime. To date, very little is known about the high current bias regime. We expect that the 2D character of exfoliated NbSe 2 crystals allows for the observation of PSLs possibly even below the critical temperature.
In this Letter we present the results of finite-bias measurements on devices based on few-layer (bi-or trilayer) NbSe 2 crystals. We shall first discuss two-terminal measurements on a point contact between Au electrode and trilayer NbSe 2 . Such measurements allow us to quantitatively understand the high bias transport features in terms of sequential nucleation of PSLs. In the second part, we discuss four-terminal measurements on extended flakes, where clear evidences of PSL nucleation have been observed in all devices measured so far. Finally, we show that it is possible to create an artificial nucleation site for PSLs by mechanically stressing a NbSe 2 flake. The first sample here described (sample A) was fabricated starting from a doped Si substrate capped with a 285 nm-thick thermal SiO 2 layer. Six ohmic contacts were defined by EBL followed by thermal evaporation of a 17 nm-thick Ti/Au film (2 nm Ti as sticking layer followed by 15 nm Au as top layer). Residuals of resist were removed by a plasma oxygen treatment. A trilayer NbSe 2 flake (thickness d = 1.8 nm) was then transferred onto the contacts using a dry-transfer tech-  nique. 38 NbSe 2 commercial crystals (HQgraphene) were exfoliated onto scotch tape and then transferred onto a polydimethylsiloxane (PDMS) film. The latter was then stamped on the prepatterned electrodes. Accurate control of the alignment during stamping was obtained using a micromanipulator. Finally, a several layer-thick flake of hBN was stamped on the NbSe 2 flake using the same technique. The hBN layer is needed to protect the NbSe 2 crystal underneath from oxidation and to improve the thermal coupling with the substrate.
An optical image of sample A is shown in Fig. 1(a). The NbSe 2 flake consists of three areas: the thinnest (three layer-thick) area overlaps three contacts, namely the two bottom contacts (5 and 6 ) and the one on the right (4 ); the medium-thick (10 layers) area overlaps the two contacts on top (2 and 3 ); finally, a large and bulky region adheres on the left contact (1 ) over a macroscopic area. Owing to its low resistance, this latter contact is used as a source in all two-terminal measurements. As a drain contact, we used one of the other contacts. All but the bottom left contact show a very low contact resistance in the few hundred ohm range at low temperatures. In contrast, contact 6, possibly due to resist residuals on the electrode surface, behaves like an opaque tunnel junction whose finite-bias differential conductance reflects the NbSe 2 density of states (see Supporting In- formation). We focused on two-terminal measurements at finite bias between contacts 1 and 5. The overlap area between NbSe 2 and the Au electrode of contact 5 is only about 0.1-0.2 µm 2 . Being so small makes this Au-NbSe 2 interface much more homogeneous than the others. Hence, transport data are quantitatively easier to interpret, although qualitatively similar to what observed for the other contacts of this sample and for other samples as well (details on these additional measurements are reported in the Supporting Information).
Transport measurements have been performed in two different cryostats with base temperature of 2 K and 300 mK, respectively. The latter one is equipped with a superconducting coil for measurements in finite magnetic field. Two-terminal differential conductance is measured by lock-in technique. A small AC excitation δV AC = 20 µV is superimposed to the swept DC voltage bias V DC . Differential (AC) and direct current signal are then simultaneously recorded.
The color plot in Fig. 1(b) shows the differential conductance measured between contacts 1 -5 as a function of finite voltage bias and temperature. Fig. 1(c) displays the differential conductance as a function of current for positive bias and for T = 3.3 K. This corresponds to the orange-white dashed segment in panel (b). We can identify three main features: a slow parabolic increase of the conductance at relatively low bias is followed by a sharp dip and then by a series of downward steps. The first step is much larger than the others. The spacing between the edges of the small steps increases with the bias while their amplitude slightly decreases. It is important to notice that all the features here observed occur at voltage bias V DC much larger than the superconducting gap of NbSe 2 , which is of the order of 1.764k B T c ≈ 0.76 meV for T c = 5 K. For V DC of the order of the gap, the conductance trace is featureless. The latter observation indicates that the observed features are not related to any quasi-particle tunneling process.
A clear clue on the origin of the steps comes from the observation of the temperature dependence of the step edges. We have highlighted one of them in Fig. 1(d) (dashed line), which shows the differential conductance as a function of current and temperature for positive bias. All step edges follow the same dependence as a function of temperature. It is instructive to compare the current at the step edges to the overall critical current for the onset of dissipation (the latter defined as the current at which the voltage drop exceeds 1 µV, see Supporting Information for further details) in few-layer NbSe 2 devices. Figure 1(f) shows (empty symbols) the value of the critical current for several few-layer NbSe 2 flakes measured in four-terminal configuration. For each sample, the current and temperature values have been rescaled to keep into account the individual values of T c and the different critical current values at low temperature. All these curves have the same temperature dependence, which is similar to that of the step edges in sample A (black dots). Interestingly, except deviations near T c , the temperature dependence of all the curves seems to reproduce the temperature dependence of the BCS gap (see discussion in the Supporting Information). This fact is non-trivial since, e.g., for dirty SCs the depairing critical current density far from T c is expected to scale as 2,39 What is important for the present discussion is the fact that the temperature dependence of the critical current in plain NbSe 2 devices is the same as that of the step edges in Fig. 1(d). This fact, together with the peculiar staircase pattern, indicates that conductance steps are a signature of the sequential PSL nucleation near the contact interface.
A further demonstration of the PSL nature of the steps is given by conductance measurement at finite perpendicular magnetic field. As mentioned above, the introduction of Abrikosov vortices rapidly disrupts the PSLinduced steps in the IVC. Numerical simulations based on TDGL equations 12 confirm that, already at fields of the order of 10% of the critical value, the step structure in the IVC is almost completely washed out. Figure 1(e) shows a color plot of the differential conductance as a function of current and perpendicular magnetic field measured at T = 300 mK. We notice that, except for the first (and largest) step, all the other steps are washed out already at 300 mT, which is only 10% of the critical field. On the one hand, this confirms that all the steps except the first are signatures of PSLs. On the other hand, the persistence of the first step at high magnetic fields indicates that its nature is different.
The above evidences allow us to construct a quantitative description of the high bias behavior of the point contact. A sketch of the region around contact 5 is shown in Fig. 2(a). The low tunnel resistance of the point contact implies that the overlap region between Au and NbSe 2 (dashed rectangle in the sketch) must be at least partially proximitized, since a low transparency contact would lead to a much lower subgap conductance. Therefore, the critical current in this region is lower than that in the rest of the crystal. By gradually increasing the current, this region will then turn normal first. This corresponds to the first step in the conductance plot. The critical current I p corresponding to this step is suppressed for perpendicular magnetic fields of the order of the critical field H c2 ≈ 3 T. 25 The first step is therefore only weakly affected by moderate magnetic fields, as observed in Fig. 1(e).
When the overlap region turns normal, the timeaveraged supercurrent density J s depends nonmonotonically on the distance x from the normalsuperconducting (NS) interface, as sketched in Fig. 2(b). At the NS interface (x = 0) the supercurrent is zero: all of the current is carried by unpaired electrons. The conversion of normal current into supercurrent at a NS interface is mediated by inelastic scattering processes that equilibrate the local charge imbalance. 4,40,41 As a result, the normal current density decreases nearly exponentially as a function of x with a decay length scale Λ q , the charge imbalance equilibration length. Correspondingly, J s increases and within a distance of the order of Λ q almost all the current is carried as supercurrent. At larger distances the total current (and thus the supercurrent) density decreases because the cross section of the SC increases. Let us consider a line connecting points at the same distance x from the NS interface (red line in Fig. 2(a)). The cross section of the SC corresponding to this line is the SC thickness d times the length of the line w+2h+πx, where w ≈ 1 µm and h ≈ 0.2 µm are the width and height of the proximitized region, respectively (see Fig. 2(a)). Hence, J s as a function of x scales approximately as where the numerator describes the normal-tosupercurrent conversion near x = 0 and the denominator accounts for the increase of the cross section at large x. This behavior is sketched in Fig. 2(b): J s as a function of x looks like a distorted bell with a maximum near x = Λ q and a tail that scales as (w + 2h + πx) −1 . The width of the bell is approximately L ≡ (w + 2h)/π ≈ 445 nm, where L is defined as the value of x such that the length (w + 2h + πx) is twice that for x = 0. By increasing the current bias, a PSL will nucleate when J s (x) locally exceeds the critical current density (see Supporting Information for a more detailed illustration of the process). By further increasing the current, more PSLs will sequentially nucleate. They will mostly concentrate within a region of size L, since J s (x) drops very quickly for large x. The convergence of the current lines near the narrow NS interface tends to confine the PSLs in a region whose size is comparable to the contact size. This allows us to adapt the analytical model of Skocpol, Beasley and Tinkham 4 (SBT) for the nucleation of PSCs in a homogeneous wire of length L.
The differential resistance as a function of current for a system of n interacting PSLs confined in a length L is given by 5 where R 0 = ρL/A, with ρ the low temperature limit of the resistivity in the normal state and A the transverse section of the crystal, which is approximately A ≈ d(w + 2h) = dπL. Therefore R 0 = ρ/(dπ) ≡ R /π, where R is the sheet resistance of trilayer NbSe 2 , which is expected to be 26,27,42 approximately 66 Ω. The nucleation of the n-th PSL occurs when the current equals (see Eq. 3 in Ref. 5) where b ≈ 0.5 is ratio between time-averaged and maximum supercurrent, 5,41 and I c is the critical current of the SC. Finally, we need to take into account a series resistance R b due to the finite interface transparency plus the resistance of the cryostat cables. Also, for current higher than the critical value I p , we add the normal resistance of the proximitized region R p . The differential conductance as a function of current and temperature is thus given by where the current dependence of R (n) P SL is given by Eq. 3 and the temperature dependence of the critical current values I c (T ) and I p (T ) is simply proportional to that of the BCS gap ∆(T ), as a result of the empirical observation discussed above. The resulting two-terminal differential conductance is plotted in Fig. 2(c). The input parameters R b = 431 Ω, R p = 80.2 Ω, I c (0) = 16 µA, I p (0) = 11 µA and T c = 4.95 K can be immediately obtained from the experimental data. The last two parameters, i.e. 2Λ q /L = 0.355 and R = 84 Ω are chosen in order to fit the distribution and the amplitude of the PSL steps. In Fig. 2(d) we highlight the conductance step edges as found in the model calculation (dashed lines) on top of the experimental data. We observe that, despite the drastic approximations behind our simple model, the agreement is satisfactory, especially for low number of PSLs. As the number of PSLs increases, the assumption that PSLs are rigidly confined in a length L becomes more and more inaccurate. This produces the deviations observed in Fig. 2(d) between the experimental and the calculated position of the conductance step edges for high number of PSLs.
The ratio 2Λ q /L obtained from the fit allows us to estimate the charge imbalance equilibration length Λ q ≈ 0.355L/2 = 80 nm. Assuming a mean free path 26,35 of about 30 nm, a Fermi velocity of about 43 10 5 m/s and a diffusion length 35 D = v F /3 = 10 −3 m 2 /s one can deduce a charge imbalance relaxation time τ q = Λ 2 q /D ≈ 6 ps.
Having found that transport in the vicinity of point contacts is dominated by the nucleation of PSLs, it is natural to expect that PSLs should also occur in fourterminal measurements on extended flakes. Figure 3(a) sketches the geometry of a typical four-terminal measurement. The contact electrodes can be either prepatterned Au strips or many-layer graphite flakes stamped on top of NbSe 2 . The graphs in Fig. 3(b,c) show IVC of two selected samples for different temperatures (more data is shown in the Supporting Information). In all the samples we measured so far, we have observed signatures of PSL nucleation: above a threshold current the IVC shows voltage jumps with distinct linear segments. The extrapolation point on the current axis of such linear segments reveals a finite and large excess current. While the qualitative behavior of all the measured samples is similar, the details of the IVC are sample-dependent, as e.g. the position and amplitude of the steps. Such variance in the IVC is closely related to the diverse geometry of the NbSe 2 flakes obtained by exfoliation. Different geometries clearly lead to different PSL distributions and thus to different step patterns in the IVC. In addition, the effective geometry of the supercurrent flow also depends on the position of the contacts. From the point contact measurement discussed above, we learned that, for an ultra-thin SC, normal contacts suppress the critical current density owing to the proximity effect. In devices based on exfoliated NbSe 2 the contacts occupy a substantial fraction of the flake area, 25-27 rather than being point-like. Therefore, voltage probes distort the supercurrent path when, at high current bias, the NbSe 2 region underneath is forced normal conducting. For example, in the sketch of Fig. 3(a) if the current bias is large enough to turn the contact regions normal, then the supercurrent will concentrate in the region between the black arrows. This will clearly be a preferred nucleation site for PSLs. Another factor relevant for the PSL arrangement is the different distribution of defects in each exfoliated flake. In clean van der Waals materials such defects might be due to bubbles of air, water or hydrocarbon contaminants 44 trapped underneath the flake during the stamping process. Such bubbles might stress the crystal and affect the arrangement of PSLs. For the reasons just exposed, it is clear that in plain NbSe 2 crystals the position where PSLs nucleate is sample-dependent and cannot be easily predicted. On the other hand, the mechanical flexibility of van der Waals SCs allows one to engineer mechanical stress in a narrow region across the whole sample width. Such region can act as a preferred nucleation site for a PSL. We chose to stamp the flake on a substrate with an etched trench, see Fig. 4(c). Van der Waals interaction with the substrate tends to pull the suspended part, which will result in tensile stress. If the width of the stressed region is smaller or comparable to Λ q , we expect that only one PSL will nucleate there, see Fig. 4(d). New PSLs will then nucleate far from the former one, and at considerably higher current densities. Figure 4 shows the results of four-terminal measurements on a trilayer NbSe 2 flake stamped on a substrate with prepatterned contacts and a 100 nm-wide etched trench in middle. Panel (a) shows the IVC for different temperatures, while panel (b) shows the color plot of the differential resistance as a function of current and temperature. For temperatures not too close to T c , only one step in the IVC is visible. After this step the IVC is linear up to the experimentally accessible bias. Such linear regime extrapolates to a finite and large excess current. No other features are visible up to I max = 300 µA, the maximum current bias we applied. Only close to the critical temperature (for T ≥ 4.9 K in Fig. 4(a)) it is possible to observe a second step within I max , which indicates the nucleation of a second PSL. After this second step, many PSLs precipitously nucleate and drive the sample into the normal state. The single isolated step is clearly evident in Fig. 4(b). The linear regime after the step is visible as a light blue region in the differential resistance color plot. This behavior is markedly different from that observed in plain devices, where several steps are visible immediately after the first (see Fig. 3). Interestingly, the isolated PSL nucleating at the stressed line behaves as a dynamic Josephson junction 14 within the same crystal lattice. This scheme can therefore be useful to study the impact of the lattice orientation on the Josephson supercurrent, e.g. by controlling the angle between the trench and the crystal axes. Moreover, it represents a fundamental building block for the implementation of more complex coherent devices within a single crystal as, e.g., a SQUID.
In conclusion, we have investigated transport in fewlayer NbSe 2 crystals under large current bias. Both twoand four-terminal measurements reveal that in this material phase slip lines are ubiquitous, as opposed to conventional metal films. Indeed, a direct current-driven transition to the normal state is never observed. In a submicron point contact the finite-bias conductance as a function of temperature can be successfully described with an adaptation of the Skocpol-Beasley-Tinkham model. In extended flakes current-voltage characteristics also show strong signatures of phase slip line nucleation, whose details are, however, sample dependent. Finally, we demonstrate that mechanical strain can induce an artificial nucleation site for an isolated phase slip line.

SUPPORTING INFORMATION AVAILABLE
PSL nucleation near the normal contact. Tunneling through the opaque contact 6. Parabolic behavior of the differential conductance at low bias. Supplementary data: four-terminal and two-terminal measure-ments. Temperature dependence of the critical current.

ACKNOWLEDGMENTS
The work was funded by the Deutsche Forschungsgemeinschaft within Grants DFG SFB1277 (B04) and GRK1570. Bulk NbSe 2 and hBN was purchased from HQ Graphene.