Superconducting Triplet Rim Currents in a Spin-Textured Ferromagnetic Disk

Since the discovery of the long-range superconducting proximity effect, the interaction between spin-triplet Cooper pairs and magnetic structures such as domain walls and vortices has been the subject of intense theoretical discussions, while the relevant experiments remain scarce. We have developed nanostructured Josephson junctions with highly controllable spin texture, based on a disk-shaped Nb/Co bilayer. Here, the vortex magnetization of Co and the Cooper pairs of Nb conspire to induce long-range triplet (LRT) superconductivity in the ferromagnet. Surprisingly, the LRT correlations emerge in highly localized (sub-80 nm) channels at the rim of the ferromagnet, despite its trivial band structure. We show that these robust rim currents arise from the magnetization texture acting as an effective spin–orbit coupling, which results in spin accumulation at the bilayer–vacuum boundary. Lastly, we demonstrate that by altering the spin texture of a single ferromagnet, both 0 and π channels can be realized in the same device.

of the main text. To establish a possible link between spin-diffusion length and the width of the rim current channels, several samples, including the one used for Figure   4 of the main text, were deposited with a copper (5 nm) buffer layer, inserted in between the Nb and Co (i.e., after FIB structuring, the layer configuration of the planar junction was Nb/Cu/Co/Cu/Nb). The Cu layer was introduced to enhance the spin-diffusion length at the S-F interface. However, the limited number of lobes in the interference patterns of such devices prevents the Fourier analysis from providing the necessary spatial resolution to quantify the variations in the size of the rim current channel. Qualitatively, no discernible changes could be attributed to the presence of the Cu buffer layer below the Nb electrodes.

S2. MICROMAGNETIC SIMULATIONS
A full description of the micromagnetic simulations of the cobalt layer is reported in a previous work. [1] Within the Object Oriented MicroMagnetic Framework (OOMMF) simulation software, we model the multilayer by dividing it into a three-dimensional mesh of 5 nm cubic cells. The exchange coefficient and saturation magnetization of Co are 3 × 10 −11 J/m and 1.4 × 10 −6 A/m respectively. The Gilbert damping constant α was set to 0.5 to allow for rapid convergence. In order to represent the polycrystalline nature of the sputtered films, we define the direction of anisotropy by a random vector field.

S3. CONTROL EXPERIMENTS
In this section we establish the vital role of spin texture for generating long range triplet (LRT) correlations. Furthermore, we show that junctions, where transport is dominated by S2 Supporting Figure S1. a. Resistance vs temperature obtained on five bar junction devices lacking any spin texture, using a measurement current of 10 µA. The weak link in devices A to D is fabricated following the same method as the triplet disk devices. The low T resistance of these devices equals the weak link resistance (0.2 Ω). As can be seen from the lower inset, only device E shows weak signs of proximity effect, with a heavily suppressed critical current. This device was fabricated using a lower Ga-ion dose, resulting in a Nb-residue in the weak link. The top inset shows a false-colored electron micrograph of a bar junction where the scale bar corresponds to 2 µm. b. shows the SQI pattern obtained on device E. The Fraunhofer pattern indicates a uniform current distribution throughout the weak link. The center lobe is shifted to -14 mT due to the magnetic stray fields of the ferromagnet. singlet correlations, are not affected by the spin texture and that the two-channel behavior is completely absent.
We fabricated control samples with rectangular configuration, where the Co layer is uniformly magnetized along its long axis (i.e., has no spin texture). The bar-shaped devices have equal width (1.05 µm) and layer thicknesses as the triplet disk devices (see top inset of Figure S1a). Due to the 5:1 aspect ratio, no magnetic vortex is stable in the cobalt layer. The same focused ion beam milling procedure was applied to form the weak link. Figure S1a shows the resistance vs temperature of five of such junctions. For devices A to D a similar Ga-ion dose is used for preparing the trench as the triplet disk devices (Nb is completely removed from the weak link). In device E, however, the weak link is fabricated by applying a lower dose, which results in the incomplete removal of the Nb layer in the trench. As shown in Figure S1a, only device E shows signs of short-range proximity effect, Supporting Figure S2. SQI pattern recorded on a disk-shaped sample containing a shallow trench (i.e., with a Nb-residue in the weak link) and therefore dominated by singlet transport.
The central peak is twice as wide as neighboring ones, and the peak height decreases in a 1/B fashion. Superconductivity in this sample is found to be unaffected by the removal of the spin texture of the ferromagnet by the application of an in-plane field.
with a heavily suppressed critical current. Figure S1b shows the superconducting quantum interference (SQI) pattern of device E. In contrast to the triplet disk junctions, we observe a distorted Fraunhofer pattern with a rapidly decaying I c (B), i.e., no two-channel behavior.
The pattern is shifted and distorted due to the magnetic dipole fields of the ferromagnet.
We additionally fabricated disk-shaped control samples where the Nb is not completely removed from the weak link to study the interaction of singlet supercurrents with the vortex spin texture. Figure S2 shows an SQI pattern of such a device. Contrary to the triplet devices, we observe a clear Fraunhofer SQI pattern with the central lobe being twice as wide as the subsequent ones and a 1/B decay of the peak height. The missing minimum on the negative field side of the pattern may be attributed to the presence of a nearby Abrikosov vortex in one of the superconducting leads. [2] This demonstrates that, despite the ferromagnetic vortex in the Co-layer, the singlet supercurrent is distributed uniformly across the junction. It should also be noted that conventional weak link devices maintain their critical current, even when the spin texture is removed by applying an in-plane field. Supporting Figure S3. In-plane I c (H y ) for a sample with triplet channels and a sample with a conventional weak link. a. The triplet sample loses its critical current above 180 mT in-plane field along the trench. At these fields, the spin texture of the disk has become uniform and therefore the generator for LRT correlations is absent. On the sample in a. we obtained the data shown in Figure 4 of the main text. b. The sample with a shallow trench retains a I c up to above 2 T, far above the fields required for stabilizing a uniform spin-texture.
In Figure S3 we compare the behavior of a device with triplet channels to that of a device with a shallow trench, i.e., a conventional weak link, under influence of an in-plane field H y along the trench. The triplet device shows variations in I c due to 0-π transitions and vice versa, and disappearance of the critical current around 180 mT. The conventional device containing the shallow trench ( Figure S2) retains an I c of around 10 µA up to 2 T.
Finally, in order to verify that the two-channel behavior is not a trivial consequence of the disk configuration, or a byproduct of the FIB milling used in forming the weak link, we examine the Josephson transport of disk-shaped junctions with a non-magnetic barrier. These control samples are structured from a MoGe (55 nm)/Ag (20 nm) bilayer.
The Ag weak link is formed through the same FIB treatment as the one applied to the Nb/Co/Nb junctions discussed in the main text. A false colored electron micrograph is shown in Figure S4a and Figure  experiments demonstrate that the rim channels only appear when transport is carried by the LRT correlations. In contrast to singlet or SRT transport, the LRT currents are also highly sensitive to the changes in cobalt spin texture, which are brought about by relatively small in-plane fields.

YSIS
In their 1971 paper, Dynes and Fulton [3] discuss a relation between the shape of the SQI pattern and the supercurrent distribution in a junction. They realized that, in the case of sinusoidal current phase relation,[4] the current distribution in the junction can be extracted by complex inverse Fourier transform of the field dependence of the critical current. Since the supercurrent can be found from integrating the critical current density along the width of the junction, we can write:

S6
Here w is the position along the junction, W the junction width, β = 2πL eff B/Φ 0 is the normalized field (L eff the effective junction length, B the applied field and Φ 0 = h/2e) and φ the phase difference at the center of the junction due to a voltage bias to the sample. In the last part, the integration bounds have been extended to infinity since outside the sample there are no supercurrents. The critical current is given by the absolute value of the complex expression: From this expression a Fourier transform can be recognized where position along the junctions width w and β form conjugate variables. The transform c is complex and therefore its real and imaginary parts encode for the even and odd components in I c (β) respectively.
Since the interference patterns measured in our experiments are relatively symmetric (i.e., an even function of the applied magnetic field), we can assume c to be dominantly real: Therefore the real part of c is an oscillating function that flips sign at each zero crossing.
The imaginary part is expected to be significantly smaller than the real part, except at the zero-crossing where the even part vanishes. Therefore, the imaginary part of c (I c,odd (β)) can be approximated by the critical current at the minima, flipping sign between each minimum and linearly interpolating between them. The inverse transform yielding J C (w) from I c (β) is then given by: We extract the critical current by a voltage cut-off in each subsequent IV -measurement (i.e., I c is the current for which we measure a certain voltage over the junction). We observe the minima of the SQI pattern to be non-zero, which might be originating from asymmetries in loop induction or junction capacitance between the two current channels. We vertically translate the extracted I c values such that the global minimum equals zero current (see figure S5a). This step in the data analysis prevents the overestimation of I c,odd (β), which would result in an overly anti-symmetric current density distribution. I c,even (β) is found by multiplying the translated I c by a flipping function that changes the sign of each subsequent lobe of the interference pattern, as can be observed in Figure S5b. We follow the above procedure for finding I c,odd (β), which is depicted in Figure S5c. The corresponding critical current density distribution is found by a numerical Fourier transform carried out in Python using the Numpy package, yielding the distribution in Figure S5d.
The original work of Dynes and Fulton describes an overlap junction between two superconductors that are relatively thick (i.e., the superconductors extend over multiple times λ S8 on each side of the junction. In this case, both superconductors can effectively shield out the magnetic field, yielding an effective junction length L = 2λ + d, where d is the barrier thickness. The disk devices, however, are planar junctions with a film thickness that is small in comparison to the London penetration depth (λ(0) = 140 nm). In this limit, the sine-Gordon equation is modified due to non-local electrodynamics and flux focusing effects, and the period of the SQI pattern is only determined by the geometry of the junction. [5-7] Generally, in this limit L eff = 0.56W is found. As expected, for our devices, the period of the SQI oscillations scales with the radius of the disk. We find for both radii that: where ∆B is the field spacing of the peaks in the SQI pattern. The effective length is then determined as: Resulting in a prefactor of L eff that is significantly lower for our disk junctions. However, if the shielding currents are also restricted in the lateral dimensions, it is shown that the prefactor can decrease. [8] Translated to an equivalent rectangular geometry described in [8], L eff = 0.26W corresponds to a junction with a total length of L equi = 0.7W . A comparison between such geometry and our disk devices is shown in Figure S6. Evidently, L equi is close to the diameter of the disk, giving a strong indication for the validity of the use of L eff = 0.26W in the Fourier analysis. The fact that L equi is slightly lower than the diameter can be explained by the circular shape of the sample boundaries, which causes a further restriction to the shielding currents.
It is important to note that the spatial resolution of current density distribution is deter- concerns a mere 6% of the total area used in the calculation. Therefore, we can assume the analysis remains valid up to a high degree for our devices.

A. Basic equations
To analyze a proximity effect in ferromagnet/superconductor system we use linearised Usadel equation for the anomalous function which takes into account triplet and singlet superconducting correlations. We presented anomalous function in the ferromagnetic region in the formf In this expansion the first term corresponds to the singlet component and the last three terms correspond to the triplet components of the anomalous function. The equations for S10 components where we include the spin-flip scattering are (for ω > 0): with the boundary conditions where γ is the transparency of the F/S interface, D is diffusion coefficient, h is the exchange field proportional to the magnetic moment and ω is the Matsubara frequency. We add the spin relaxation term with rate τ −1 sf . The sources in the right-hand side of Eq. 10 and 11 f s | S and f | S are the spin-singlet and spin-triplet components of the anomalous function in the superconducting electrode. In the simplest case it contains only the spin-singlet part f s | S =∆/ |∆| 2 + ω 2 and f | S = 0. However in general there are also the spin-triplet correlations generated due to the inverse proximity effect [9,10] or due to the magnetic proximity effect [11][12][13] through the spin-mixing mechanism.
If the distance between superconducting electrodes is large the current though F is me- (orange arrows), which adiabatically follows the local spin gauge field ∇θ v = (m×∇m) z . c. When the ESC is incident on a bilayer-vacuum boundary, the adiabatic approximation breaks down. As the ESC cannot propagate outside the bilayer, the total spin current at the boundary must remain zero. To compensate for the adiabatic spin current, a condensate of LRT pairs (in red and blue) with S x = ±1 emerges near the vacuum boundary to generate an opposing spin current.
expression for the elongated spatial derivative with the spin-dependent vector potential Z = ∇θ v . The termσ z Z plays the role of effective spin-orbital coupling.
The transformation (13) demonstrates the importance of edges in generating the spintriplet components which are not collinear with magnetization and therefore can penetrate into the ferromagnet. Indeed, far from the magnetic vortex core, the scale of changing S12 magnetization texture is large, that the emergent vector potential field can be considered as locally-homogeneous Z = const which corresponds to the magnetic texture changing in space with constant gradient. In this case, far from the edges only the correlations with f h exist in S and in F in the thin layer of the size ξ F = D/h near the S/F interface.
However even the constant spin-dependent gauge field Z generates the equilibrium spin current [14][15][16] At the superconductor/vacuum or ferromagnet/vacuum boundaries the spin current should vanish n j J j = 0. Using (14) this boundary condition can be written as The boundary condition (15) provides the generation of LRT both in the S and F electrodes. The generation in F electrode is described by the linearized Usadel Eqs. (8,9) and boundary conditions (10, 11) with r.h.s. f s | S = F 0 = ∆/ √ ∆ 2 + ω 2 and f | S = 0.
The generation in S electrode is in general determined by the non-linear Usadel equations with inhomogeneous Zeeman field h m which can be induced either by the inverse proximity effect or the magnetic proximity effect. In the case of inhomogeneous h we again obtain a non-zero LRT source f | S ⊥ h in the r.h.s. of the Eq. 11. In both scenarios, the generation of LRT in the presence of a magnetic vortex is described by the qualitatively similar equations derived below; a schematic of the LRT generation is shown in Figure S7. S13

Generation of LRT in S and F electrodes
In general the superconducting correlations both in S and in F are determined by the Usadel equation for the quasiclassical Green function (GF) Hereτ 1,2,3 are the Pauli matrices in Nambu space. The GF has the off-diagonal anomalous partτ 1f which describes the superconducting correlations. The field h can be either (i) the induced Zeeman field in S or (i) the real exchange field inside F. We consider these two cases separately.
(i) First, let us assume that h is the Zeeman field in S, homogeneous over the z coordinate.
The local approximation produces the collinear triplet f SRT m which in the simplest case is where ∇ ⊥ = (∇ x , ∇ y , 0), ξ −2 ω = 2(ω 2 + |∆| 2 )/D. The equation is simplified in case of the coplanar magnetization texture m z = 0. Then the transverse spin vector can be found in the form f LRT = (−m y , m x , 0)Ψ LRT . Solving Eqs. (17,18) we find the distribution of Ψ LRT in the superconductor. Then, using the boundary conditions (11) we obtain the distribution of LRT in the F from the equations Note that the Eqs.(19) is similar to (17) with the only difference in the decay scale determined by λ ω = 1/ ξ −2 ω + λ −2 sf where ξ ω = D/|ω| and λ sf is the normal metal spin diffusion length. S14 Supporting Figure S8. Sketch of the overlap S/F structure with magnetic vortex pattern in the lower magnetic layer. The notches can be seen on the side.
In case if the F thickness is smaller than the normal metal coherence length and spin diffusion length we can integrate Eqs.(19,20) by z to get the equation with the source Solving this equation we obtain the distributions of LRT show in Figure 5 of the main text.
(ii) The SRT can be generated by the exchange field in F determined by the linearized Usadel Eqs. (8,9). For large exchange field h T c they have solutions of two types which are the short-range and long-range modes with the scales ξ F = D/h and ξ ω = D/|ω| correspondingly. Hereafter we assume that ξ F is the smallest length of the problem such that the spatial dependencies of exchange field and geometrical factors are characterized by the scales ξ F .
The spin-triplet correlations can be written as the superposition of short and long-range The spin vector of short range modes is parallel to the exchange field f SRT = mΨ SRT . Under such conditions we search the short-range solutions of Eqs. (8,9) in the form where the scale is λ h = e iπ/4 /ξ F . In case if the thickness of F is smaller than the coherence and spin diffusion lengths d F the LRT distribution is determined by the equation similar to Eq. 17 where the z-averaged SRT is Ψ SRT z = (γ/d F )(T c /h)F bcs .

Example: LRT generated by the magnetic vortex
The distribution of magnetization in the magnetic vortex far from the vortex core is given by m = (−sinθ v , cos θ v , 0) where θ v is the polar angle with respect to the vortex center. In addition, we take into account the non-circular geometric shape, e.g., shown in Figure S8. We find the LRT distribution solving either the system (17,18,21) or (24, 25) using finite element numerical package FreeFem [17]. Both systems yield qualitatively similar LRT distributions shown in Figure 5 of the main text. S16 C. Generation of triplets at the trench Besides the boundary sources of LRTs, some sources come directly from the Usadel equation in the presence of non-magnetic inhomogeneities such as the order parameter profile. [18] Indeed, the inhomogeneous |f SRT | 2 yields the divergence of the second term in the spin current (14) ∇ j J j ∝ Z j ∇ j |f srt | 2 z. Such sources have to be compensated by the contribution from the LRT.
To support the above qualitative arguments we consider the generic 2D problem of a superconductor/ferromagnet bilayer shown in Figure S9. We assume that exchange field amplitude is h = 100T c0 and the magnetization pattern is M = M (a, x, 0)/ √ x 2 + a 2 with a = 2ξ 0 . The F-layer has a thickness of ξ 0 and is 8ξ 0 wide, where ξ 0 = D/T c . We consider two model overlap S/F/S junctions shown in Figure S9. The upper row corresponds to smooth edges of S electrodes parametrized by the coordinate-dependent transparency of the S/F interface. The lower row corresponds to the sharp edges of S electrode with γ(x) changing abruptly. By solving the linearised Usadel equations (8,9,10,11) with f s | S =∆/ |∆| 2 + ω 2 and f | S = 0, using the finite element numerical package FreeFem, [17] we can separate the SRT contribution f SRT = f s h/h and the LRT one f LRT = f − f SRT for several different configurations. This allows us to study formation of LRT correlations near the edges of F and S layers.
The fourth row showing the total LRT amplitude |f LRT |, indicates that they are indeed generated at the boundaries of F layer. Also, LRT correlations are generated near the edges of the overlap S layers which forms the trench. At the same time, this source of LRT generation is very sensitive to the presence of surface defects at the trench. The destructive influence of the defect can be seen in Figure S9 comparing the panels (c) and (f), with a flat surface and the rectangular well between the S electrodes. The suppression of LRT generation results in the strong suppression of the critical current j c shown in panel (d). The current is normalized to j 0 = σ n T c /e which is of the order of 10 2 /ξ 0 A/cm. Such current flowing through the cross section of F of the linear size R ∼ 1 µm yields j ∼ 10 −8 A, which is three orders of magnitude smaller than the experimentally observed value. S17