Glassy Dynamics in a heavy ion irradiated NbSe2 crystal

Fascination with glassy states has persisted since Fisher introduced the vortex-glass as a new thermodynamic phase that is a true superconductor that lacks conventional long-range order. Though Fisher’s original model considered point disorder, it was later predicted that columnar defects (CDs) could also induce glassiness — specifically, a Bose-glass phase. In YBa2Cu3O7−x (YBCO), glassy states can cause distinct behavior in the temperature (T ) dependent rate of thermally activated vortex motion (S). The vortex-glass state produces a plateau in S(T ) whereas a Bose-glass can transition into a state hosting vortex excitations called double-kinks that can expand, creating a large peak in S(T ). Although glass phases have been well-studied in YBCO, few studies exist of other materials containing CDs that could contribute to distinguishing universal behavior. Here, we report on the effectiveness of CDs tilted ~30° from the c-axis in reducing S in a NbSe2 crystal. The magnetization is 5 times higher and S is minimized when the field is parallel to the defects versus aligned with the c-axis. We see signatures of glassiness in both field orientations, but do not observe a peak in S(T ) nor a plateau at values observed in YBCO. Finally, we discuss the possibility that competing disorder induces a field-orientation-driven transition from a Bose-glass to an anisotropic glass involving both point and columnar disorder.

associated with a vortex-glass state. Furthermore, when the field is misaligned with the CDs, various staircase structures 8 (see Fig. 1a) are known to form; a distinct signature of such structures has not yet been identified in creep measurements.
Many studies have characterized the effects of columnar defects on J c (θ H ) 5,[9][10][11][12][13][14][15][16][17] , where θ H is the angle of the applied field. Much less is known about the effect of field orientation on the creep rate (S) and, more generally, creep in materials besides YBCO that contain CDs. For example, it is unknown why the peak associated with rapid double-kink expansion in YBCO has not been observed in other materials [18][19][20][21][22] . Of particular interest is superconductors with low Ginzburg numbers (Gi), such as NbSe 2 , which can attain significantly lower creep rates 23 than superconductors with high Gi, such as YBCO (Gi ~ 10 −2 ). This evokes the question of whether glassy states in low Gi materials manifest as a plateau at such a high S ~ 0.02-0.04 and double-kink expansion creates a peak in S. More generally, it motivates garnering a better understanding of the dynamics of various vortex excitations and glassiness in materials with low Gi.
In this study, we characterize the effect of temperature, magnetic field and field orientation on vortex dynamics in a NbSe 2 crystal containing parallel CDs tilted ~30° from the c-axis. First, we observe the expected peak in J c (θ H ) when H is parallel to the CDs, and we find that this peak is indeed accompanied by a dip in S(θ H ). Second, we compare and characterize S(T) and S(H) when the field is parallel to the defects (H || CDs) versus the c-axis (H || c). Last, we find evidence of glassiness in both field orientations.

Sample Fabrication and Measurements
Our experiments are carried out on two undoped 2H-NbSe 2 crystals that were grown using iodine vapor transport 24 and have dimensions ~0.8 mm × 0.7 mm × 20 μm and ~1.5 mm × 0.3 mm × 8.5 μm (length L × width W × thickness δ). 2H-NbSe 2 is a layered transition metal dichalcogenide with an s-wave gap structure that has attracted intense interest 25 because it hosts a coexisting incommensurate charge density wave phase and superconductivity below T c ~ 7 K. Our primary motivation for studying NbSe 2 is that it is a clean system (few defects in as-grown crystals) that has a low Ginzburg number (Gi). Scanning tunneling microscopy studies have revealed a low density of Nb and Se vacancies and Nb interstitials in NbSe 2 crystals grown by iodine vapor transport [26][27][28] . One study found a defect density of ~0.4% 28 . Assuming a coherence length ξ ab ≈ 7.4 nm, penetration depth [29][30][31] λ ab ≈ 126 ± 3 nm, and upper critical field anisotropy 32  is the thermodynamic critical field.
One crystal (δ = 20 μm) was heavy-ion irradiated with 1.4 GeV 208 Pb 56+ ions at a dose of 1.45 × 10 11 ions/cm 2 corresponding to a matching field of 3 T (average distance between CDs ~ 26 nm) at the Argonne Tandem Linear Accelerator System (ATLAS) while mounted with the crystallographic c-axis ~30° from the incident beam. The sample underwent no additional processing steps post-irradiation. We chose to induce tracks at an angle of ~30° from rather than parallel to the c-axis to distinguish the effects of the CDs from those of mass anisotropy and intrinsic correlated defects (e.g., edge and screw dislocations) that are known to produce a peak in J c (θ H ) for H || c in YBCO 13 . Similarly, for tilted CDs, the mere existence of asymmetry between J c (θ H ) and J c (−θ H ) can provide evidence of correlated pinning.
Transmission Electron Microscopy (TEM) studies were performed on the irradiated crystal. The acquired image shown in Fig. 2a indicates that the columnar amorphous tracks are continuous and almost perfectly parallel to each other, consistent with previous studies 33 and with the small splay expected for 1.4 GeV Pb ions. Figure 2b is a higher magnification image showing an angle of ~29° between the radiation direction and the NbSe 2 [002] direction. From our TEM work, we measured an average CD diameter of about 4 to 6 nm. In addition to columnar tracks, heavy ion irradiation may induce secondary electrons that act inelastically with the material matrix, producing point defects in between the columnar tracks 34,35 . There is limited knowledge about the secondary damage produced by heavy ion irradiation. A recent scanning tunneling microscopy study of a heavy ion irradiated Fe(Se, Te) crystal showed that the superconducting order parameter was annihilated inside the columnar tracks and suppressed by the interlaying point defects 35 .
Magnetization (M) measurements were collected using a Quantum Design SQUID magnetometer with a rotating sample mount as well as transverse and longitudinal pick-up coils to measure each component of the magnetic moment. By measuring M versus T at 2 Oe, we find that the critical temperature of the irradiated crystal is T c ≈ 7 K, similar to that in pristine crystals 25

Results and Discussion
Magnetization in different field orientations. Figure 3 compares isothermal magnetic hysteresis loops, M(H), at T = 1.8 K for the pristine crystal for H || c (θ H = 0°), and the irradiated sample for both H || c and for the field aligned with the defects (H || CDs, θ H = θ CD = −31°). The pristine crystal demonstrates dramatically lower magnetization and irreversibility field than the irradiated crystal. This suggests a weak pinning landscape and that the columnar defects in the irradiated crystal are overwhelmingly the predominant source of pinning.
For the irradiated crystal, the magnetization is roughly 5 times higher when the field is aligned with the CDs than with the c-axis. A large enhancement was anticipated and had been observed in previous studies, though the magnitude was less 5 . This improvement could be attributed to the higher energy used during irradiation (1.4 GeV Pb 56+ versus 300 MeV Au 26+ in ref. 5 ), which might create straighter, more continuous tracks 40 .
The dip at low fields μ 0 H < 0.6 T is caused by the out-of-plane pinning anisotropy. That is, pinning by extended defects along the c-axis (or, in our case, tilted 30° from) should produce a weak dip in M(H) at zero field, while pinning along the crystallographic ab-plane is expected to produce a peak 41 . At fields below self-field H sf ≫ H, vortex lines over a large region of the sample peripheries are quite curved. As the applied field is increased (approaching self-field), this region decreases as vortices straighten over a wider portion of the sample center. Columnar defects are more effective at pinning straight vortices. Hence, the initial increase in M with increasing H is caused by a combination of the heightened effectiveness of individual CDs in pinning less curved vortices and growing portions of the sample containing straight vortices. Predicted theoretically 41 , the peak has been observed in irradiated YBCO 42 and Ba(Fe 0.93 Co 0.07 ) 2 As 2 crystals 43 .
Additional M(H) loops were collected at T = 4.5 K and at 20 different angles. Select curves are shown in Fig. 4a,b, capturing crossovers into different regimes. Note that the curves converge near zero field. This is because in the very dilute limit and for all field orientations, vortex lines will be oriented normal to the sample surface (aligned with the c-axis) to minimize their energy by shortening 5 .
As the field tilts away from alignment with the CDs (|θ H − θ CD | > ~6°), the low-field peak progressively shifts to lower fields and eventually disappears. In particular, at θ H = −24°, M(H) decreases nearly linearly with decreasing H. Further rotation of the field away from the CDs (θ H < −40°, θ H > −19°) changes the M(H) behavior. M initially abruptly decays with increasing H, showing similar shape to M when H||c (Fig. 3). As the field is increasingly tilted (θ H ≥ −2°), the M(H) curves exhibit a weak second magnetization peak (known as the fishtail effect) between 0.5 T and 1 T. This is most pronounced at θ H = 33°, as highlighted in Fig. 4b. The fishtail effect has been observed in a wide variety of materials, including low-temperature superconductors, cuprates, MgB 2 , and iron-based superconductors 44,45 and associated with an equally wide variety of effects, including elastic-to-plastic crossovers, vortex order-disorder phase transitions, and vortex lattice structural transitions 44 . In fact, a previous study 46 reported the appearance of a fishtail in a pristine NbSe 2 crystal when the applied field was tilted 30° from the c-axis and attributed it to a vortex order-disorder transition.  Extracted from the M(H) loops, the data is re-plotted as M(θ H ) at different fields in Fig. 5. The peak at θ H = θ CD is clear at all fields and M rapidly decays at the slightest field misalignment with the defects, corresponding to a large reduction in J c . It is important to note that this prominent peak indicates that pinning provided by the CDs is significantly greater than that from any point defects possibly introduced in between the CDs by secondary electrons during the irradiation process. If we compare critical currents when the field is aligned with the CDs versus the c-axis, we find that J c is ~240 kA/cm 2 compared to ~48 kA/cm 2 , respectively, at 0.6 T. Figure 6 shows such a comparison at 0.3 T over a broad temperature range, displaying an increase in J c by a factor of ~4 at 4.5 K and ~3 at 1.8 K. Note that the defects are effective even down to the lowest field of 0.2 T, where J c is only ~10% lower than at the maximum. This is consistent with all data in Fig. 5 being well above H sf ~ J c δ ≤ 550 Oe at this temperature. At most angles, lower fields produce higher M. However, for θ H > 0°, some low field curves cross, resulting in non-monotonic M(H) that is consistent with the regime in which the fishtail is observed (Fig. 4b).
Vortex creep when field is aligned with CDs versus c-axis. To analyze vortex excitations and the potential for glassy dynamics, we measured the dependence of the creep rate on temperature and field orientation. First, we consider two creep models: the Anderson-Kim model and collective creep theory. A defect (or collection of defects) can immobilize a vortex segment (or a bundle of vortex lines) by reducing the vortex line energy by the pinning energy U P (T, H), which is the energy barrier that must be overcome for vortices to move. The Lorentz force induced by the persistent current J then reduces U P to an activation barrier U act (T, H, J) and the vortex hopping rate is ~− e U k T / act B . The Anderson-Kim model 2 , which neglects vortex elasticity and therefore does not predict glassy behavior, often accurately describes creep at low temperatures T ≪ T c . It assumes U act (J) ∝ U P |1 − J/J c | for  Collective creep theory 2 predicts that the temperature dependence of the creep rate is where t 0 is the effective hopping attempt time and C ≡ ln(t/t 0 ) ~ 25-30. Here μ > 0 is the glassy exponent indicating the creep regime: μ = 1/7, 3/2 or 5/2, and 7/9 are predicted for creep of single vortices, small bundles (size less than the penetration depth λ ab ) and large bundles (size greater than λ ab ) of flux, respectively. At low temperatures T ≪ T c , U P ≫ μk B Tln(t/t 0 ) such that S(T) ≈ k B T/U P , coinciding with the Anderson-Kim prediction. We now compare creep data for the irradiated crystal in two different field orientations: H || CDs and H || c. Note that our measurements are restricted to low fields because at high temperatures and fields, the magnetic signal is quite small when H ||c. Figure 7a shows the measured creep rate versus field orientation at 4.5 K and 0.5 T. Creep is clearly minimized when the field is aligned with the defects; S is an order of magnitude smaller for H || CDs than for H || c. In fact, aligning the field with the defects suppresses creep at all fields and temperatures measured in our study, e.g., the comparison of S(H) in both field orientations at 1.8 K shown in Fig. 7b.
Comparing creep data for the irradiated sample to the pristine crystal can only be performed at very low fields because the measurement signal produced by the pristine crystal at higher fields is near the lower limit of our measurement sensitivity. The temperature dependence of the creep rate in the pristine crystal and the irradiated crystal at 0.02 T is shown in Fig. 7c. For both field orientations, S increases linearly with T up to 5.5 K, qualitatively adhering to the Anderson-Kim description. Despite the very low applied field, the CDs are effective at lowering creep when H || CDs, but not when H || c, seen from a comparison to the data from the pristine sample.
Considering collective creep theory, if U P ≪ Cμk B T, S(T) should plateau at S ~ 1/Cμ. Such a plateau is predicted in the case of glassiness, such that S ~ 0.02−0.04, equivalent to typical observations of plateaus in YBCO single crystals 7 and iron-based superconductors 39,[47][48][49][50][51][52] . For our NbSe 2 crystal, Fig. 7d shows S(T) at μ 0 H = 0.3-0.5 T for the two field orientations. In all cases, in Fig. 7c,d the creep rates are well below the usual collective creep plateau. The simplest interpretation is that U P is not negligible compared to Cμk B T (see eq. 1), which is in agreement with the pinning energy estimates described below. Although, consistent with this scenario, most of the  Fig. 7d also shows a broad temperature insensitive region in the 0.5 T data for H || c (S ~ 0.003) and a narrower one in the 0.3 T data for H || CDs (S ~ 0.002). However, interpretation of these data as indicative of a plateau at much lower than usual values would imply Cμ ~ 300-500, producing unphysically large values of either μ (10-17) or C (120-200); note that typical C and μ values 7,12 give Cμ < 75. Finally, quantum creep may be a significant component of our measured creep rates at these low temperatures, in which case, adding a temperature independent (and unfortunately unknown) contribution would imply an even smaller thermal creep contribution.
A plateau in S(T) is the most apparent manifestation of glassy vortex dynamics. In its absence, we need a different approach to assess the nature of the vortex depinning excitations. Analysis of the current dependence of the effective activation energy U * ≡ T/S can provide direct experimental access to μ without the need for assumptions regarding U P . According to collective creep theory 2 , the activation barrier depends on the current as where J c0 is the temperature-dependent critical current in the absence of flux creep. Considering the Arrhenius hopping rate ~− − t e U J k T 0 1 ()/ act B and equations (1) and (2), the effective pinning energy is where μ > 0 for glassy creep and μ → p < 0 for plastic creep 53 . Consequently, the exponent can easily be extracted from the slopes of U* vs 1/J on log-log plot. From Fig. 8, we see distinct elastic-to-plastic crossovers for all sets of data. At low T the dynamics is clearly glassy at both field orientations, with μ ~ 1. This is one of the main experimental findings of this study. As T increases the dynamics turns plastic, with p in agreement with the expectation for the motion of dislocations in the vortex lattice (p = −0.5) 54 .
For H || CDs, glassy dynamics with μ ~ 1 is expected for a Bose-glass state characterized by half-loop formation. However, glassiness was unforeseen for H || c. In this configuration, we expected to see evidence of staircase structures (see Fig. 1), which form when the field is tilted away from the CDs by an amount greater than the lock-in angle (θ L ), but less than the trapping angle (θ t ). Yet in the simplest scenario staircases should be non-glassy, as finite length kinks easily slide along CDs. So, several possibilities should now be considered: θ H = 0° is within the lock-in angle and half-loop excitations are responsible for μ ~ 1, the dynamics of the staircase vortices is glassy, or this orientation is beyond θ t and the CDs do not produce correlated pinning (so glassiness arises from standard random collective pinning).
A Bose-glass state formed when the field is aligned with CDs (and vortices are localized on these defects) will be robust to small changes in field orientation. That is, when the field tilted away from the CDs by an angle less than θ L , vortices will remain completely pinned by the CDs. This results in a plateau in M(θ H ) for |θ H − θ CD | < θ L that has been observed in cuprates [55][56][57][58] and Co-doped BaFe 2 As 2 14 . Though our data is too coarse to determine if there is a lock-in effect and identify θ L , we see from Fig. 5 that the magnetization is greatly reduced at θ H = 0° versus θ H = θ CD . So, θ H = 0° is clearly well beyond the lock-in angle. Consistently, θ L is expected to be very small in our NbSe 2 crystal (see estimate below). On the other hand, the asymmetry of M(θ H ) around θ H = 0°, which can only arise from the tilted CDs, suggests that staircases are present at this orientation 55 .
Having eliminated half-loops and random collective pinning as the cause of μ ~ 1 at H || c, we consider the possibility of a vortex-glass state or an anisotropic glass involving both columnar and point disorder, as predicted in ref. 4 . Segments of a single vortex line could be alternatingly pinned by adjacent CDs and interstitial point defects. As the current and thermal energy act on the vortex, the segments pinned by point defects might wander/ entangle (instead of sliding like kinks). Alternatively, interactions among weakly pinned kinks may create "kink bundles" that, by analogy with the 3D vortex bundles, should exhibit glassy collective creep with μ ~ 1. In either case, if the phase for H || CDs is indeed a Bose-glass then the system experiences a field-orientation-driven transition from a Bose-glass (H || CDs) to an anisotropic glass (H || c). As the expected exponent μ ~ 1 is identical for a vortex glass, Bose glass, and anisotropic glass, measurements of the exponent alone cannot distinguish between vortex configurations that lead to glassy dynamics. The real fingerprint of the Bose glass is the presence of a lock-in effect.
In light of this, we find it important to mention an alternate possible scenario: even for H nominally parallel to the CDs, a slight field misalignment θ H − θ CD > θ L could lead to staircase formation. Such a misalignment is challenging to avoid when θ ≈ Pinning energies. The effectiveness of CDs is typically assessed by evaluating the measured pinning energies, which can be calculated from the creep data. The scale of the pinning energy in a superconductor 59 is approximately the condensation energy . For NbSe 2 , we calculate that U P1 ~ 160-300 K within our measurement T range. From the Fig. 7a inset, we see that the effective activation energies U* extracted from our creep measurements plummets from being considerably greater than to comparable to U P1 as the field rotates from H || CDs to H ||c. This is because pinning energies larger than U P1 can be achieved through individual strong pinning by defects larger than V c , as is the case for our CDs.
Columnar defects are most effective at pinning vortices of smaller core size ξ ≤ R 2 ab (where R is the CD radius) [2][3][4]12,60 . This is not easily achieved in low-T c superconductors, which tend to have large coherence lengths. When ξ < R 2 ab (as is the case for our sample), under ideal pinning conditions ε r ≈ ε 0 (R/2ξ ab ) 2 . Considering an average R ~ 2.5 nm for the CDs in our crystal, at T = 1.8 K we obtain ε r /ε 0 ≈ 0.02, about twice our aforementioned experimental value determined simply from J c /J 0 . This demonstrates that the CDs in our crystal indeed behave as strong correlated disorder, producing about half of the ideal pinning. For comparison, analogous calculations predict that CDs in YBCO should ideally produce J c ~ J 0 , while experimental J c values fall short of that by a factor of ~3 to 4.
Ũ K (4 5 ) 1000 K h . This is somewhat smaller than our experimental U * (4.5 K, 0.5T) = T/S~3500 K, but consistent given the simplicity of the estimates. First, we note that the calculation of ε r (4.5 K) based on J c /J 0 is likely an underestimate, as J c may be reduced due to CDs discontinuities, vortex bending, and the possibility that some vortices may be occupying interstitial positions outside the CDs. Alternatively, if we use the estimate ε r ≈ ε 0 (R/2ξ ab ) 2 , for R = 2.5 nm we obtain ε r (4.5 K) ~ 10 K/nm; and . Ũ K (4 5 ) 1300 K h . We note that the calculation is highly sensitive to slight changes in the parameters, e.g., R ~ 2-3 nm yields In fact, the effective CD size may be larger because the irradiation induced tracks may depress the superconducting order parameter over a farther distance than the diameter measured by TEM due to, e.g. lattice strain. Second and perhaps more importantly, the above analysis neglects vortex-vortex interactions, which should be considered because the lateral dimension of the half-loops 4 ε ε . is not negligible compared to the vortex lattice parameter a 0 ~ 70 nm for μ 0 H = 0.5 T. Hence, repulsion of neighboring vortices produces a caging effect that increases the effective pinning energy, stiffening the lattice and reducing S.

Conclusions
In conclusion, we have studied the dependence of vortex dynamics on the orientation and magnitude of the applied magnetic field in a NbSe 2 crystal containing tilted columnar defects. As most studies of creep in samples containing columnar defects have been limited to heavy ion irradiated YBCO, studying NbSe 2 has allowed us to probe effects applicable to materials with lower Ginzburg numbers and larger vortex size to columnar track ratios. Specifically, we demonstrated that the critical current is maximized and creep is concomitantly minimized when the field is aligned with the defects (T = 4.5 K, μ 0 H = 0.5 T). This result was not necessarily intuitive, as the rapid expansion of double-kinks can promote fast creep when H || CDs in YBCO (at low temperatures and fields below the matching field). We also found that H || CDs preferentially produced lower creep rates than H || c over our entire measurement range, and that both field orientations resulted in glassy behavior. A Bose glass state is indeed expected when the field is aligned with the CDs. Yet the existence of glassiness when the field is misaligned is quite fascinating and suggestive that staircase structures might be able to entangle or localize in a way that leads to glassy behavior. Many open questions remain. First, it is unclear why a distinct, large peak in S(T) resulting from double-kink expansion has only been observed in YBCO. Second, do other materials containing CDs show glassiness when the field is oriented in a way that is favorable for staircase formation? In addition to testing this in other low Gi materials, it would be interesting to test in highly anisotropic samples in which pinning to the ab-plane is highly favorable over the c-axis. Third, is the potential anisotropic glass state enabled by secondary damage that appears in between the columnar tracks? These results motivate further studies of creep rates at various field orientations in other heavy ion irradiated materials.

Methods
TEM images. The TEM specimen of the irradiated NbSe 2 crystal was fabricated in a focused ion beam and the microstructure was characterized by using FEI Tecnai F30 Transmission electron microscopy (TEM, 300 kV).
Magnetization Measurements. Magnetization measurements were collected using a Quantum Design SQUID magnetometer with a rotating sample mount, and transverse and longitudinal pick-up coils to measure each component of the magnetic moment, m t and m l , respectively. The angle of the field was verified by calculating θ = − m m tan ( / ) H t l 1 , the total moment θ = m m / cos l H , and the magnetization M = m/δLW (where δ μm is the thickness, W is the width, and L mm is the length). Creep data were taken using standard methods 7 . Firstly, the field was swept high enough (ΔH > 4H * ) that the sample was fully penetrated with magnetic flux and in the critical state. Then, successive measurements of M were recorded every 15 s, capturing the decay in the magnetization (M ∝ J) over time (t). Last, the time was adjusted to account for the difference between the initial application of the field and the first measurement and S = |d ln M/d ln t| is calculated from the slope of a linear fit to ln M-ln t. T c was determined from the temperature-dependent magnetization at H = 2 Oe.