Imaging anisotropic vortex dynamics in FeSe

TL;DR

This paper investigates the anisotropic behavior of vortex pinning in FeSe using scanning SQUID susceptibility, revealing that vortices move easily along domain walls but are strongly pinned perpendicular to them, shedding light on vortex dynamics and nematicity.

Contribution

It introduces a novel experimental approach combining local magnetic field application with vortex imaging to characterize anisotropic vortex pinning in FeSe.

Findings

Vortices move easily along twin domain walls.

Vortices are strongly pinned perpendicular to domain walls.

Pinning anisotropy is linked to domain wall structures.

Abstract

Strong vortex pinning in FeSe could be useful for technological applications and could provide clues about the coexistence of superconductivity and nematicity. To characterize the pinning of individual, isolated vortices, we simultaneously apply a local magnetic field and image the vortex motion with scanning SQUID susceptibility. We find that the pinning is highly anisotropic: the vortices move easily along directions that are parallel to the orientations of twin domain walls and pin strongly in a perpendicular direction. These results are consistent with a scenario in which the anisotropy arises from vortex pinning on domain walls and quantify the dynamics of individual vortex pinning in FeSe.

Tables

Table 1. Table 1: Summary of measured k w subscript 𝑘 𝑤 k_{w} and k s subscript 𝑘 𝑠 k_{s} values.

Sample 1
Axis	$k_w$ (N/m)	$k_s$ (N/m)	$k_s/k_w$	T (K)
1	$4.1\pm 0.5\times {10}^{-9}$	$>8.5\times {10}^{-8}$	$>21$	7
1	$8.3\pm 0.9\times {10}^{-9}$	$>1.7\times {10}^{-7}$	$>22$	7
1	$2.7\pm 0.8\times {10}^{-8}$	$>1.2\times {10}^{-7}$	$>4.5$	7
2	$4.0\pm 0.9\times {10}^{-8}$	$>6.2\times {10}^{-7}$	$>16$	7.5
2	$8.6\pm 3.3\times {10}^{-8}$	$>3.6\times {10}^{-7}$	$>4.2$	6
2	$3.2\pm 0.7\times {10}^{-8}$	$>1.5\times {10}^{-7}$	$>4.7$	7
Sample 2
Axis	$k_w$ (N/m)	$k_s$ (N/m)	$k_s/k_w$	T (K)
1	$8.5\pm 1.9\times {10}^{-9}$	$>6.6\times {10}^{-8}$	$>7.8$	6
1	$1.7\pm 0.3\times {10}^{-8}$	$>6.4\times {10}^{-8}$	$>3.8$	6
1	$6.6\pm 1.4\times {10}^{-9}$	$>8.6\times {10}^{-8}$	$>13$	6
1	$8.0\pm 2.0\times {10}^{-9}$	$>3.6\times {10}^{-8}$	$>4.5$	4.5
1	$9.5\pm 2.4\times {10}^{-9}$	$>2.0\times {10}^{-7}$	$>21$	7
2	$1.2\pm 0.4\times {10}^{-8}$	$>7.4\times {10}^{-8}$	$>6.2$	4

Equations

$${\phi }_s({\vec{r}}_i,z_i)=-\sum _l\sum _j\frac{w{g}_{j,l}}{4\pi }\frac{z_i-z_{j,l}}{((z_i-z_{j,l}{)}^2+{\rho }_{i,j,l}^2{)}^{3/2}},$$

$${\phi }_1(\vec{r},z)=\frac{1}{(2\pi {)}^2}\int d^2\kappa ({\phi }_s(\vec{\kappa })e^{\kappa z}+{\phi }_r(\vec{\kappa })e^{-\kappa z})e^{i\vec{\kappa }\cdot \vec{r}}$$

$${\vec{H}}_2(\vec{r},z)=\frac{1}{(2\pi {)}^2}\int d^2\kappa {\vec{h}}_2(\vec{\kappa })e^{qz}{e}^{i\vec{\kappa }\cdot \vec{r}},$$

$$h_{2,z}(\vec{\kappa })=\frac{2{\kappa }^2}{\kappa +q}{\phi }_s(\vec{\kappa })$$

$${\vec{h}}_{2,||}(\vec{\kappa })=\frac{2i\vec{\kappa }q}{\kappa +q}{\phi }_s(\vec{\kappa })$$

$$\delta u(\vec{r})=\frac{1}{2}\int d^2{r}^{\prime }{\vec{H}}_2({\vec{r}}^{\prime })\cdot {\vec{B}}_V(\vec{r}-{\vec{r}}^{\prime })$$

$$\delta u(\vec{r})={\Phi }_0{H}_{2,z}(\vec{r})/2$$

$$\delta U(x,y)=\frac{{\Phi }_0}{2(2\pi {)}^2}{\int }_0^{-\infty }dz\int d^2\kappa h_{2,z}(\vec{\kappa })e^{qz}{e}^{i\vec{\kappa }\cdot \vec{r}}$$

$$\delta U(x,y)=\frac{{\Phi }_0}{2(2\pi {)}^2}\int d^2\kappa h_{2,z}(\vec{\kappa })e^{i\vec{\kappa }\cdot \vec{r}}/q$$

$${\vec{F}}_{\mathrm{SQUID}}(x,y)=\vec{\nabla }(\delta U(x,y))$$

$${\vec{F}}_{\mathrm{SQUID}}=k_{s}ds\hat{s}+k_{w}dw\hat{w}.$$

$${\Phi }_{ac}=\frac{d{\Phi }_{dc}}{dx}dx+\frac{d{\Phi }_{dc}}{dy}dy$$

Full text

I Introduction

FeSe is a particularly simple layered iron based superconductor (Fe SC), with an un-strained critical temperature of about 8 K [1], an increase in this critical temperature with pressure to 36.7 K at  8 GPa [2], and with a report of interface-induced high-temperature superconductivity above 50 K in single unit-cell films on SrTiO3 [3]. In general, iron-based superconductors exhibit complex interplay between superconducting, nematic, and magnetic orders, and undergo tetragonal-to-orthorhombic phase transitions close to magnetic ordering transitions [4]. Unlike other Fe SCs, FeSe does not magnetically order, providing an opportunity to study the superconductivity-nematicity relationship without the added complexity of an ordered magnetic state. The structural transition of FeSe occurs at 90 K and is understood to be driven by electronic nematic order [1, 5, 6, 7]. In the orthorhombic state, FeSe and other Fe SCs form domains separated by twin boundaries (TBs). The superconducting pairing mechanism has been discussed in terms of spin-fluctuation pairing [8, 9, 10] but the nodal character of the gap remains controversial. Some superfluid density [11], thermal conductivity [11], and tunneling spectroscopy measurements [12, 11] are consistent with line nodes in the orbital component of the superconducting order parameter; however, recent thermal conductivity [13, 14], STM [9], and London penetration depth [15] studies suggest that FeSe is fully gapped but with deep gap minima. The relationship between nematic order and superconductivity is likewise controversial, with NMR studies suggesting that nematic order competes with superconductivity [16, 17] and heat capacity and thermal expansion studies suggest that it enhances superconductivity [18].

Although vortex pinning in superconductors is of great practical as well as fundamental interest, its mechanism is still poorly understood [19]. Recently it has been possible to study this pinning directly by imaging individual vortices while manipulating them. Auslaender et al. [20] dragged vortices in a cuprate superconductor over distances of a few microns using a magnetic force microscope tip. They found an enhanced response of the vortex to pulling when the tip was oscillated transversely. They also found enhanced vortex pinning anisotropy, which they attributed to clustering of oxygen vacancies in their sample. Later work from Shapira et al. used MFM to drag vortices along twin boundaries in YBa2Cu3Ot-$\delta$ and demonstrated that the vortices moved in a series jumps, consistent with power-law behavior [21]. Kalisky et al. [22] showed that vortices, when dragged by a scanning SQUID microscope, avoided crossing twin boundaries in underdoped Ba(Fe1-xCox)2As2. Embon et al. [23] used a SQUID on a tip to image the movement over a few tens of nanometers of vortices driven by applied supercurrents in a thin Pb film. They were able to map out anisotropic and spatially inhomogeneous pinning forces on the vortices, which they attributed to multiple overlapping pinning sites.

The vortex pinning properties of FeSe have attracted interest [24, 25, 26, 27] and make it a potential competitor to high-$T_{c}$ cuprates for high field applications [26]. Critical current density studies of vortex pinning in FeSe have found that it is dominated by strong point-like pinning [24], while STM studies have shown that vortices preferentially pin on TBs in FeSe, where the superfluid density is reduced [28]. In this paper we present scanning SQUID magnetometry and susceptibility images of vortices trapped in single crystals of FeSe. The susceptibility images show structures that we attribute to motion of the pinned vortices driven by the magnetic fields applied by the field coil integrated into our SQUID susceptometers. Analysis of our data using a simple model is consistent with a quadratic dependence of the restoring force on displacement along the direction of the TBs, with a much larger restoring force in the orthogonal direction. The anisotropy in the in-plane restoring forces can be as large as a factor of 20.

II Materials and Methods

We used scanning SQUID microscopy to image magnetometry and susceptibility in bulk single crystal FeSe. Single crystal FeSe samples were grown by chemical vapor transport following the growth procedure outlined in Ref. 29 and exfoliated with Kapton tape and silver paint (Dupont 4929N) as in Ref. 30 to achieve a surface flat enough to be scanned using our susceptometers. The superconducting transition and structural transition temperatures of the batch of samples used were found to be 8.2 K and 88 K respectively. The bulk superconducting transition temperature was extracted from magnetic-susceptibility measurements taken using the vibrating sample mount option of the MPMS 3 from Quantum Design. The structural transition temperature was determined from resistivity measurements on a free standing crystal taken using the MPMS 3 paired with a Linear Research Model LR-700 AC resistance bridge. The susceptibility and resistivity data are shown in Fig. 1.

Our SQUID susceptometers contain two Nb pickup loop/field coil pairs arranged in a gradiometric layout [31]. The pickup loop and field coil are covered by Nb shielding so that flux passes only through the loop and not through the gaps between the leads. The inner radius of the pickup loop was 0.3 um, resulting in sub-micron spatial resolution. As the susceptometer scans across the surface of the sample, we record the magnetic flux passing through the pickup loop. The dc signal is recorded as magnetometry and is reported in units of the flux quantum $\Phi_{0}=h/2e$. We use an SR830 lock-in amplifier to pass an ac current through the field coil, creating a local magnetic field, and record the ac flux at that frequency to measure susceptibility. The gradiometric design cancels out the flux due to the field coil so that the pickup loop only measures the magnetic response of the sample. The susceptibility is normalized by the lock-in amplifier current and is reported in units of $\Phi_{0}$/A. This design allows us to image local magnetic fields and susceptibility at the surface of a sample simultaneously.

III Results

III.1 Imaging vortex motion

In four measured samples, we found two-lobed features (”butterflies”) in the susceptibility images, accompanying superconducting vortices. Not all vortices had these “butterflies” and their brightness varied from vortex to vortex. Increasing the sample temperature increases the brightness of the lobes up to Tc, after which both the vortex and accompanying butterflies disappear. The butterflies orient along one of two directions which are perpendicular to each other. Further, the relative brightness of the two lobes is consistent within butterflies of the same orientation, and the dimmer lobes are located in areas opposite to the location of the field coil shielding. The brightness of the butterflies was varied from vortex to vortex, ranging from around one tenth of a $\Phi_{0}/A$ to a few $\Phi_{0}/A$. Besides wiggling the vortex, the field coil could also push the vortex to a different pinning site. Not all vortices had butterflies, likely because they were pinned too strongly to show a signal in susceptibility. Due to unevennesss in the sample surface, it was difficult to control for the height and angle of the pickup loop/field coil pair relative to the surface area being scanned. We used x-ray diffraction to determine the tetragonal crystal axes in Samples 1 and 2 and found that a line cutting through the two lobes orients the butterflies either along or perpedicular to the tetragonal [100] direction. This orientation is consistent with the butterflies being aligned with TBs. Fig. 2 shows two representative susceptibility butterflies from Sample 1, along with the corresponding magnetometry image showing the vortex for the first butterfly. A background is subtracted from all susceptibility images by fitting a plane to the four corners of the scan.

In Fig. 3, we show a series of susceptibility scans for three vortices with butterflies at variouss temperatures, as well as a temperature series showing the development of stripes in diamagnetic susceptibility close to $T_{c}$. Near the superconducting transition, the vortex becomes easier to move and the features of the butterflies become sharper. In one sample we also observed striped variations in susceptibility run along the tetragonal [100] direction.

III.2 Modeling vortex motion

We model our susceptibility images by 1) calculating the magnetic fields inside the superconductor at the vortex position due to the applied currents through the field coil, 2) calculating the motion of the vortex in response to these fields using a simple model with an anisotropic pinning potential, and then 3) calculating the change in flux through the susceptometer due to the vortex motion. The SQUID susceptibility is given by the response flux $\Phi$ divided by the field coil current $I$.

III.2.1 Applied fields

Consider a geometry in which a scanning SQUID susceptometer, composed of niobium films with penetration depth $\lambda_{\rm Nb}$, with layout in the pickup loop/field coil region as illustrated if Fig. 2b, is assumed oriented parallel to the sample surface, has a spacing $z_{0}$ between the surface of the susceptometer and the surface of the sample, and is in the half-space $z>0$ (region 1). The superconducting sample, with penetration depth $\lambda$, is in the half-space $z<0$ (region 2). The fields generated by the susceptometer are calculated following Ref. 32 as described for our sensors in Ref. 31. For a thin superconducting film with thickness $t<\lambda_{Nb}$, one can define a stream function $g$ by ${\bf J}={\bf\hat{z}}\times\nabla g$, where ${\bf J}$ is the sheet current [32]. One interpretation of the stream function is that it defines a density of magnetization in the $\hat{z}$ direction, or equivalently a collection of small current loops in the plane. Once the stream functions $g_{j,l}$ are known for all of the grid points $j$ in all of the superconducting layers $l$ in the susceptometer, the source potential at any point $i$ in the half-space $z>0$ but outside of the susceptometer superconducting layers is given by [33]

[TABLE]

where $w$ is the area of the pixels used, and $\rho_{i,j,l}=\sqrt{(x_{i}-x_{j,l})^{2}+(y_{i}-y_{j,l})^{2}}$. For what follows we calculate the source potential $\varphi_{s}$ for all $\vec{r}_{i}$ and $z=0$.

Once the source fields are known, we match boundary conditions at $z=0$. This is done by expanding the scalar magnetic potential $\varphi_{1}(\vec{r},z)$ outside the sample ($z>0$) and the magnetic field $\vec{H_{2}}(\vec{r},z)$ inside the sample ($z<0$) in Fourier series that are constructed to satisfy Maxwell’s equations for $z>0$ and London’s equation for $z<0$ [33].

[TABLE]

[TABLE]

where $\vec{\kappa}=\kappa_{x}\hat{x}+\kappa_{y}\hat{y}$, $\kappa=\sqrt{\kappa_{x}^{2}+\kappa_{y}^{2}}$, $q^{2}=\kappa^{2}+1/\lambda^{2}$ and $\vec{H}_{1}=\vec{\nabla}\varphi_{1}$. Using the boundary conditions $\vec{B}\cdot\hat{z}$ and $\vec{H}\times\hat{z}$ continuous at $z=0$, as well as $\vec{\nabla}\cdot\vec{B_{2}}=0$ and $\vec{B}=\mu_{0}\vec{H}$ results in [33]

[TABLE]

[TABLE]

The magnetic field $\vec{H}_{2}(\vec{r},z)$ penetrating the superconductor is calculated taking the Fourier transform of $\vec{h}_{2}(\vec{\kappa})$ (Eq. 3).

III.2.2 Vortex motion

The additional magnetic energy $\delta u$ of the vortex per unit length due to the susceptometer applied field is given by

[TABLE]

where $\vec{B}_{V}(\vec{r}^{\prime}-\vec{r})$ is the field of the vortex centered at $x,y$. If we assume the vortex is oriented in the $\hat{z}$ direction and neglect the finite extent of the vortex fields, this becomes simply

[TABLE]

The total additional energy $\delta U$ of the vortex due to the field penetration into the superconductor is then

[TABLE]

The integration over $z$ can be done analytically, resulting in

[TABLE]

The force on the vortex due to the susceptometer fields is the gradient of this extra potential:

[TABLE]

We find that the force is reduced when the field coil shield is above all or part of the vortex, meaning that when comparing scans to the corresponding susceptometer layout, the areas opposite to the field coil shield will show reduced signal.

Because we do not know the exact form or mechanism of the vortex pinning potential in FeSe, we use a simple quadratic model with spring constants $k_{w}$ and $k_{s}$ ($k_{s}\geq k_{w}$) associated with orthogonal axes $\hat{w}$ and $\hat{s}$ rotated by an angle $\theta$ in the $ab$ plane relative to the scan axes $\hat{x}$ and $\hat{y}$ . When the susceptometer scans in the $xy$ plane relative to the vortex position, the susceptometer applied fields pull the vortex towards (or away) from the pickup loop/field coil. Since we modulate the current through the field coil at about 1 kHz, we assume that the vortex response is much faster than the applied force, and therefore that the displacement of the vortex ${\vec{dr}}=dw\,{\hat{w}}+ds\,{\hat{s}}$ from its equilibrium position can be calculated from the balances of forces condition

[TABLE]

III.2.3 Response flux due to vortex motion

The expected ac flux from the vortex motion is calculated using the gradient of the dc flux through the susceptometer pickup loop due to the vortex:

[TABLE]

The flux due to a vortex is calculated following the procedure in Ref. 31; the derivatives in Eq. 12 are taken numerically, and the susceptibility is calculated by dividing $\Phi_{ac}$ by the applied field coil current.

III.2.4 Results of the model

Example predictions for the ac flux due to vortex motion with the susceptometer geometry in Fig. 2 is shown in Fig. 4 for various assumed values of $k$ along the x- and y-axes of the images, with $z_{0}$=2 $\mu$m. The resulting shapes for the isotropic cases ($k_{s}=k_{w}$) are similar to an incomplete torus, while the shapes for the anisotropic cases ($k_{s}\neq k_{w}$) are distinctly lobed. The weaker axis ($k_{w}$) spring constant determines the intensity of the susceptibility signal, while the stronger axis spring constant ($k_{s}$) largely determines the shape. As the strength of the strong axis spring constant is increased, the lobes become more axially symmetric, and the dark region in the center fades. Consistent with the force profile, the lobes located across from the field coil shielding are dimmer than their partners. Crucially, we find that the apparent axis of the anisotropy is not qualitatively changed by the asymmetric susceptometer geometry.

Using reduced chi-square fitting, we calculated the optimal spring constant in the weak axis (“$k_{w}$”) direction and a lower bound for the spring constant in the strong axis direction (“$k_{s}$”) for 12 butterflies from two samples. We find that in all cases the signal along the hard axis is so weak compared to the noise level that the optimal $k_{s}$ value approaches infinity. In Fig. 5 we show the fitting process for the butterflies shown in Fig. 2. We fit $k_{w}$ and $k_{s}$ separately by taking cuts along the weak and strong axes, respectively. Full images for the data, model, and difference between data and model are shown in Fig. 6 The difference in brightness between model and experiment in the dimmer lobe can be explained by the limited spatial extent of the field coil shielding in the model SQUID geometry.

Taking the optimal $k_{w}$ value and the lower limit for $k_{s}$ results in a lower limit for the ratio $k_{s}/k_{w}$, which we use to characterize the anisotropy of vortex motion. The range of $k_{s}/k_{w}$ varies from around 4 to over 20, and three vortices show ratios above 20. The results are summarized in Table 1, along with the temperature at which the data were taken.

Another way to infer the vortex dynamics in these samples is to directly extract the vortex displacement from scans. This approach avoids the assumptions of a toy model potential entirely. Providing that the direction of the vortex motion is known, we can use Equation 12 and the gradient from the magnetometry image to obtain the maximum vortex displacement. In this analysis we exclude the values close to the vortex center and far from the vortex, as the magnetometry gradient is very small in these areas and amplifies any noise in the susceptibility. This is the case along the lines where $F_{s}=0$ and $F_{w}=0$, respectively. Furthermore, since the motion in the $\hat{s}$ direction is severely limited, we can assume the vortex moves mostly in the $\hat{w}$ direction. The motion along lines through $F_{s}=0$ and with constant but small $F_{w}$ is therefore mostly radial. In Fig. 7 we show line cuts of the experimental vortex displacement and calculated susceptometer force along the weak and strong axes of a butterfly from Sample 1. Because the signal along the strong axis is dominated by the background, we subtracted a linear fit from susceptibility line cuts along this axis. The absolute magnitude of the vortex displacement is at least an order of magnitude smaller in the strong-direction compared to the weak-direction. Furthermore, while the displacement along the weak axis tracks with susceptometer force, the displacement along the strong axis appears to be largely independent of the force applied by the field coil. A linear fit to the weak axis data is plotted in Fig 1(b). The effective $k_{w}$ in this case is $1.0\pm 0.06\times 10^{-8}$ N/m, about 2.5 times the optimal $k_{w}$ found by fitting the data to the quadratic model. These results are consistent with a much larger spring constant along the strong axis than along the weak axis, and also consistent with a linear force-displacement relation in the weak direction (Fig. 7(b)).

IV Discussion

The ac current through the field coil generates an oscillating local magnetic field and superconducting screening currents that pull and push the vortex towards and away from the center of the pickup loop/field coil. As a result, as long as the vortex response to the applied force is faster than the time variation of the field coil current, the flux through the pickup loop is in-phase with the field coil current, as reflected in the real part of the susceptibility image. This scenario is rather general to scanning SQUID measurements with simultaneous magnetometry and susceptibility imaging and does not rely on the particulars of the FeSe samples, except that in FeSe the motion is both large and highly anisotropic.

We can first rule out that the motion is isotropic but that the observed lobed shape is simply due to the SQUID susceptometer geometry based on the simulations in Fig 2. The expected ac signal in the isotropic case is proportional to the susceptometer force, which is mostly radially symmetric except for shielded regions. The force is small when the pickup loop/field coil center is above the vortex, increases to a maximum at a distance equal to the field coil radius and then decreases again as the field coil moves away from the vortex. This force profile results in the “doughnut” shape discussed earlier. The effect of increasing the spring potential along one axis is to weaken the effect of the force along that axis so that as k becomes very large the signal along that axis becomes indiscernible, producing lobed features.

As noted in the introduction, previous MFM and scanning SQUID microscopy measurements have measured vortex dynamics in superconductors[20, 22, 23]. In contrast to the Auslaender et al. and Kalisky et al. studies, in which the vortices were dragged to new locations on the sample, in our case and for small displacements in the Embon et al. [23] study the vortex returns to its original location and is instead oscillated in place. In the Embon et al. study the displacements are typically a few tens of nanometers, the forces are pico-Newtons, and the spring constants are $\sim 10^{-5}$ N/m, whereas in the present study the displacements are typically one micron, the forces are typically femto-newtons, and the weak spring constants are typically $10^{-8}$ N/m, 3 orders of magnitude weaker in FeSe than in Pb.

A 2012 STM study [28] found that in thin film FeSe, superconductivity is suppressed along twin domain boundaries and that vortices tend to pin on these boundaries. The cause of this suppression is thought to be due to the increased height of the Se atoms along the boundary. We did not observe chains of vortices along any direction; however, it is possible our vortex density was not high enough to observe such patterns. In three samples for which we were able to determine the direction of the crystalline axes, the vortex motion was highly anisotropic along the TB directions, suggesting that it is easier to pull vortices along the TBs than acorss them. Further evidence of this can be found in the striped variations in susceptibility signal close to $T_{c}$ which run along the same direction as the weak axis of the butterflies. Given that the only other known symmetry-breaking through lattice orthorhombicity, for which the crystal axes are oriented at 45 degrees to the butterfly and stripe direction, the only reasonable conclusion is that the vortices are pinned on TBs and the variations susceptibility reflect the suppressed superconductivity on the boundaries.

V Conclusion

In conclusion, scanning SQUID susceptibility can be used to image the motion of vortices under the influence of locally applied magnetic fields. We have applied this method in the instance of vortices in FeSe. Detailed calculations of the magnetic fields generated by the susceptometer, combined with a simple model for the pinning forces on the vortex, show that these pinning forces can be highly anisotropic and is consistent with a quadratic dependence of the restoring force on displacement. We calculate an effective spring constant for the weak and strong axes and show that this is consistent with vortex pinning that is strong across TBs and weaker along them.

Acknowledgements.

We thank J. Straquadine and E. Rosenberg, and C. Watson for experimental assistance and P. Massat for feedback on the manuscript. This work was supported by the Department of Energy, Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division, under Contract No. DE-AC02-76SF00515