Critical Forces at Fractional Matching Fields in Superconducting Thin Films with Triangular Pinning Lattice

We numerically study the commensurability effects in the critical forces in a type II superconducting thin film with vortices under the influence of a triangular pinning lattice at sub matching fields and zero temperature. The analysis is made applying a transport force in two mutually perpendicular directions and magnetic field perpendicular to the film surface. The results show critical force peaks at fractional matching fields, as a direct consequence of the commensurability effects. For the first time, we simulate the sequence of force peaks detected in recent experiments. Anisotropic effects in the critical forces are also reported, in analogy with that found at higher fields.


Introduction
Commensurability effects in type II superconducting films under the influence of periodic pinning arrays, has been investigated during the last years  . As nnotechnology advances, pinning center modeling has become viable due to the control of size, shape and composition of the defects.Several types of periodical pinning lattices have been studied both experimentally and theoretically, such as triangular, square, honeycomb, Kagomé and others.In these works, a wide variety of commensurability effects have been found in dynamical properties of vortices, such as critical forces, dynamic phases, vortex-flow resistance, and magnetization  .When the vortex lattice matches the pinning lattice we say that the system is commensurate, in this case the pinning potential felt by the vortices is maximum.This occurs when the rate between the number of vortices and pinning centers (N v /N p ) or analogously the rate between the applied magnetic field and the first matching field (B/B ϕ ), is an integer or a rational fraction.The matching effects have been extensively analyzed, in most of cases for B=nB ϕ (n=1,2,3...), finding pronounced peaks of critical current vs (B/B ϕ ) at these values of field.For the particular case of triangular pinning arrays, numerical and experimental works have been conducted concerning the commensurability effects in the critical forces 1,2,11,12,15 .However, only few works 2,11,12 have analyzed the matching effects at sub matching fields.The results of these works show a disagreement between the sequence of critical force peaks predicted by theoretical simulations and experimental results, which needs a better understanding.S. Ooi et.al. 11 studied experimentally the critical currents and vortex behavior in 2D triangular pinning potentials in single-crystals Bi 2 Sr 2 CaCu 2 O 8+y films, finding fractional matching peaks at B/B ϕ =1/4, 1/3, 1/2, 2/3, 3/4, 6/7 and 1.Additionally, recent experiment conducted by D. Bothner et al. 12 investigated the commensurability effects in the critical currents for superconducting Nb thin films with periodic and quasi-periodic antidot arrays.In this work, critical current peaks were found for the triangular pinning lattice at B/B ϕ =1/4, 1/3, 2/3, 3/4 and 1. Theoretical simulations, conducted by C. Reichhardt and N. Grønbech-Jensen 2 , have also analyzed films with triangular pinning lattice as a function of pinning strength.They found current peaks at B/B ϕ =1/7, 1/4, 1/3, 2/3, 6/7, and 1 for high values of pinning strength, and additional force peaks at B/B ϕ =1/9, and 1/2 for lower values of pinning strength.However, the current peak at B/B ϕ = 3/4, which was found in recent experimental results 11,12 , was not predicted.
In this work, a numerical study of the critical forces in a two-dimensional superconducting film with a triangular pinning lattice is performed, focusing on fractional matching fields and its respective vortex states.The results show a sequence of critical force peaks which has some differences from previous simulations 2 , but are in good agreement with experimental data 11,12 , including the force peak at B/B ϕ =3/4.Two ingredients may be important for the appearance or not of critical force peaks; the pinning strength and the pinning radius.Concerning the pinning strength, numerical analysis 2 showed that small values of pinning strength may lead to additional force peaks, however, even for small values, they did not find the force peak at B/B ϕ =3/4.On the other hand, V. V. Moschalkov et.al. 20 has experimentally shown that a decrease in the pinning radius may provoke the attenuation a Programa de Pós-Graduação em Ciência e Tecnologia de Materiais, Faculdade de Ciências, Universidade Estadual Paulista -UNESP, Bauru, SP, Brazil b Department of Chemistry, University of Southampton, Highfield, Southampton SO17 1BJ, United Kingdom c Departamento de Física, Faculdade de Ciências, Universidade Estadual Paulista -UNESP, Bauru, SP, Brazil or even disappearance of some current peaks.It is not clear which pinning radius was used by C. Reichhardt and N. Grønbech-Jensen 2 , but we think that the force peak at B/ B ϕ =3/4 is evident in our results because we use a pinning radius close to the experimental value, and probably larger than in that work.Moreover, our calculations show that different pinning radius may lead to different ground state vortex configurations.In some cases, small pinning radius may lead to a disordered ground state, while larger pinning radius to ordered configurations.Then, changes in the size of pinning centers may dramatically change the commensurability effects, explaining the appearance of additional force peaks in our simulations.

Model
In this simulation, we consider an infinite superconducting thin film placed in the xy plane with an external magnetic field B applied perpendicularly to the film surface and zero temperature.
A set of Langevin equations describe the dynamical properties of the system 22 : In our case, we minimize the energy of the vortex system until finding the optimal vortex positions.Then, the set of Langevin equations are integrated numerically by Molecular Dynamics technique using the optimal vortex configuration as initial boundary value.After that, the driving force is increased slowly, in steps of ΔF=0.01 up to values just above the critical force F c .For each value of force, we use 30000 time steps for evaluation of the time averages.To obtain the critical force values, it is necessary the calculation of the time average vortex velocity , which is proportional to the macroscopically measured voltage-current V(I) curve, for x direction and analogously for y direction 22 .To determine these critical forces, we use the criterion that velocities , .Vx y 0 001 < corresponds to vortex fluctuations around the equilibrium position.

Results and Discussion
As the commensurability effects strongly depend if a given vortex is sited on top of a pinning site or not, it is essential to obtain the vortex positions of the ground state at zero transport current.Finding the right initial vortex positions is a very delicate work, because each vortex arrangement has slightly different energy value.For that task, we use the GSA method described before.The equilibrium configuration of vortices at zero transport current for B/B ϕ =1/9, 1/4, 1/3, 1/2, 2/3, 3/4, 8/9 and 1, are illustrated in Figure 1.It is important to note that the pronounced force peaks in Figure 2 always correspond to regular vortex lattices, including that for B/ B ϕ =3/4, which was not predicted before 2 .Other minor peaks were found at B/B ϕ =1/9 and 8/9, where all vortices are pinned but do not form a regular lattice.Another configuration that is necessary to highlight is for B/B ϕ =1/2.The vortex configuration in Figure 1 (d) is ordered, with a lowest energy and no interstitial vortices, different to that proposed 2 .
The results of the simulation for the transport forces applied in two mutually perpendicular directions are plotted in Figure 2. In that figure, pronounced peaks appear at B/B ϕ =1/4, 1/3, 1/2, 2/3, 3/4 and 1, when the vortex and pinning lattice are commensurate, in good agreement with experimental results 11,12 .
Furthermore, it is possible to see well defined anisotropic behavior in the critical forces, mainly for the peaks at B/B ϕ =1/3, 1/2, and 2/3.Anisotropic behavior was previously observed in experimental and theoretical for infinite films, but for higher magnetic fields and different system parameters [13][14][15][16]19,21 . Our smulations also show that, in general, F c y > F c x for fractional matching fields, with the exception of B/B ϕ =1/3.This behavior is in agreement with previous works 14,15,21 .
As pointed out before, previous theoretical 2 and experimental 11,12 works show discrepancies between the sequence of force peaks.As we suggested before, this may be attributed to the pinning size.To that end, we performed simulations for B/B ϕ =1/2 with smaller pinning size to corroborate The term on the left of Eq. ( 1) is the velocity of the vortex.In the right side, the first term is the vortex-vortex interaction, described as U vv (r ij )=C v K 0 (r ij /Λ), where r ij =|r j -r i | is the distance between the vortex i and j, C v is the vortexvortex strength defined as C v = Φ 0 2 /8πΛ and K 0 is the modified Bessel function.The second term in the right side is the vortex-pinning interaction which is described as a Gaussian function 22,23 , where r ip =|r p -r i | is the distance between the vortex i and the pinning center p, a p =ξ is the pinning center range and C p =0.2C v is the strength of the pinning potential.The last term is the driving force due to the transport current .In our simulations, the length scales are normalized by 4ξ, time scale by τ=16ηξ 2 /C v and energy scales by C v .Where ξ is the coherence length, η the Bardeen-Stephen friction and ϕ 0 =h/2e is the elementary flux quantum.The results can be compared to real thick samples through Λ, the effective penetration depth defined as Λ=2λ L 2 /d, where d is the sample thickness and λ L the London penetration depth 22,23 .
We simulate the vortex system under the influence of a triangular pinning array.The transport forces are applied in two mutually perpendicular directions, x and y.We consider a fixed pinning density of n p =N p /L x L y =0.114/ξ 2 .The initial vortex ground states are prepared using Generalized Simulated Annealing (GSA) minimization method 24 .In analogy with the annealing used in Metallurgy, where a metal is cooled slowly until it reaches its lowest energy state, the Simulated Annealing technique introduces one or more artificial temperatures that are gradually cooled to reach the extreme of a function, which corresponds to the global minimum.our argument.The static vortex configuration is shown in Figure 3, using a p =ξ/2 and same pinning strength, C p =0.2C v .As a result, vortices did not organize themselves in a regular lattice, moreover, part of the vortices stabilize in interstitial positions.This arrangement is no more commensurate, resulting in a much lower critical force.

Conclusions
In this work, we numerically studied the behavior of the critical forces for an infinite type II superconducting film under the influence of a triangular pinning array.The sequence of force peaks found in the simulation agree with  recent experimental data 11,12 at sub matching fields.In the absence of a driving force, we found ordered vortex states for some fractional matching fields, including for B/B ϕ =3/4, which was not predicted in previous simulations 2 .These ordered states lead to a higher energy barrier that avoid vortex movement and lead to peaks in the critical force.According to our simulations, the agreement between our results and experimental data may be attributed to the use of larger pinning size when compared with previous simulations 2 .
As the transport force was applied in two different directions, we were also able to see anisotropic effects in the critical forces.In general, our system shows that when the transport force is applied along the short diagonal length, vortices move easily than when it is applied along the longer diagonal direction.This anisotropic behavior agrees with other experimental and theoretical works performed for higher fields 14,15,21 .

Figure 2 .
Figure 2. Critical force, F c y and F c x , vs B/B ϕ for the triangular pinning array with C p =0.2C v and a p =ξ.

Figure 3 .
Figure 3. Static vortex configuration for smaller pinning radius, where the black dots represent the vortex lattice and the black circles the triangular pinning lattice at B/B ϕ =1/2 obtained from the minimization, with C p =0.2C v and a p =ξ/2.