Giant vortex states in type I superconductors simulated by Ginzburg-Landau equations
The quantization of magnetic flux in superconductors is usually seen as vortices penetrating the sample. While vortices are unstable in bulk type I superconductors, restricting the superconductor causes a variety of vortex structures to appear. We present a systematic study of giant vortex states in type I superconductors obtained by numerically solving the Ginzburg-Landau equations. The size of the vortices is seen to increase with decreasing film thickness. In type I superconductors, giant vortices appear at intermediate thicknesses but they do not form a well-defined vortex lattice. In the thinnest type I films, singly quantized vortices seem to be stabilized by the geometry of the sample instead of an increase in the effective Ginzburg-Landau parameter.
pacs:74.78.-w, 74.25.Wx, 74.25.Ha
A well-known property of the superconductors is the quantization of the magnetic flux which is an essential property for many of the applications of superconductivity. In the typical case of a type II superconductor in the mixed state, the magnetic field penetrates the sample as vortices that form an Abrikosov lattice each vortex carrying one quantum of flux. With the recent advances in imaging techniques vortices can be resolved in greater detail compared to the traditional Bitter decoration. Perhaps the most impressive is the work of Hasegawa et al.  in which the internal magnetic field structure of vortices was measured. They report vortices in thin Pb films with single and multiple flux quanta. These exotic vortex structures with more than one flux quantum have also more recently been reported in confined superconducting systems and in thin Pb films that had triangular or quasiperiodic artificial pinning sites [2, 3, 4, 5]. Multiquantum vortices are seen in the experiments despite the fact that a vortex with flux quanta are generally thought to be energetically unfavourable. The magnetic energy of a giant vortex scales as while it scales as if the flux quanta are separated into singly quantized vortices . In this paper, we use the convention that a giant vortex is a cylindrically symmetric vortex with more than one flux quantum while the non-symmetric case is called a multivortex.
Tinkham was the first to propose that vortices could be the reason behind the increase in the critical field of thin type I superconductors . Later on Maki showed that the Abrikosov lattice is a stable solution to the Ginzburg-Landau equations for arbitrary provided that the film thickness is sufficiently small . The work of Maki was continued by Lasher who predicted that, between thick type I films that are in the intermediate state and thin films that are in the mixed state, there is a transitional area where giant vortices form honeycomb lattices and possibly more complex structures . The theoretical work on these giant vortex structures is mostly based on linearized Ginzburg-Landau equations with the exception of a few recent numerical papers [10, 11, 12, 13]. The current methods in numerical modelling allow us to solve the complete Ginzburg-Landau equations to address the problem of the complex giant vortex structures in type-I superconductors.
In this paper we present a systematic study of vortex states in low- superconducting thin films in low magnetic fields. We have numerically determined the local energy minima of the Ginzburg-Landau functional for a series of thin films of different thicknesses. The results show stable giant vortex configurations at certain thicknesses and a change in the size of vortices with film thickness.
2 Computational model
The static Ginzburg-Landau equations were solved by finding a (local) minimum of the associated energy functional. It is convenient to write the energy in a dimensionless form leaving only the overall energy scale dimensional which does not affect the solutions. The dimensionless energy is
where is the dimensionless Ginzburg-Landau parameter and and . Here is the order parameter, is the vector potential, is the permeability of free space, is the mass of an electron and is the coefficient of the term as usual in the context of Ginzburg-Landau theory. The penetration depth is absorbed in the overall energy scale and therefore . The method used to minimize (2), the used boundary conditions and the details of restricting to model pinning sites are described in our earlier work . In short, the boundaries are so that direction is periodic and the and directions end in vacuum with magnetic field parallel to the .
The Ginzburg-Landau functional was minimized for and in thin film geometry (27 m x 37 m x ) with thicknesses and . The magnetic field values used were 0.2, 0.3, 0.4, 0.5 and 0.6 mT. Any lengths or magnetic field values that are given in dimensional form in this paper assume that 150 nm. The simulations were chained so that for each thickness the result with the lower magnetic field was used as the starting point for the higher field value simulation. The Meissner state was used as the initial configuration in the 0.2 mT simulations. In all the simulations there were two pinning sites at the same positions to break the symmetry caused by periodic boundaries. Additionally, the sample was enclosed in vacuum extending 4.5 in the non-periodic directions which allows to include the effects of magnetic stray field energy in the simulations.
3 Results and discussion
The vortex configurations obtained by the simulations can be seen from the absolute value of the order parameter that is shown in figure 1 for both and in two cases, the thinnest and the thickest film at mT. The singly quantized vortices in the thinnest films (figure 1a and b) are visibly larger in size than in the thick films (figure 1c and d) for both values of . In the thick film with the larger vortices are giant vortices with vorticities of all values up to six. All the vortex lattices have a strongly disordered triangular lattice which is partly due to the low value of . As we have shown in our earlier work (ref. ), the disorder is also increased by the pinning sites that can be seen as the two larger vortices in figure 1d. The two pinning sites are present at the same position in all the simulations. The giant vortex in the pinning site of which an enlarged version is also shown in figure 1c is actually a multivortex since the cylindrical symmetry of the vortex is broken. A closer look at the phase information reveals the individual vortex cores that are marked with white squares. The distinction between a cylindrically symmetric giant vortex with the vortex cores all exactly centered and a multivortex with the vortex cores very closely packed together is quite subtle. However, all multiply quantized vortices in pinning sites seem to be multivortices while the free vortices are giant vortices. This could possibly explain the two types of vortices seen in the experiments of Hasegawa et al. .
The phases of the order parameters corresponding to the absolute values in figure 1 are shown in figure 2. It is clearly visible that the phase contours leading to a giant vortex are clustered together all the way to the sample border. In contrast, the phase contours of the singly quantized vortices are evenly spaced to fill the whole length of the sample. The clustering of the phase contours reflects the nature of the vortex-vortex interaction which is attractive in the case of . The clustering is prominent in the thick film (figure 2c) with while the thin film (figure 2a) with the same value shows a similar spread in phase contours as the ones with repulsive vortex-vortex interaction.
To visualize the order of the vortex lattices we define the two-dimensional case of the radial distribution function as
where is the average density of vortices in the sample, is the area of the ring from to and is the number of vortices inside it . The radial distribution functions were calculated from the positions of the vortices in cross-sections similar to the images in figure 1 and are shown in figure 3. The vortices in simulations form a better defined vortex lattice than the vortices in the case. In the case of , the nearest neighbour peak broadens and moves to larger distances with increasing film thickness. The shift in the peak position is caused by the giant vortices that have larger separation than singly quantized vortices. The broadening of the peak means that the giant vortices do not form a lattice with a well defined nearest neighbour distance. In the case of , the thinnest film has stronger peaks than the thickest which indicates more robust vortex lattice in the thin sample.
The change of the vortex size with film thickness was quantified by considering an effective coherence length and an effective penetration depth . The effective coherence length was defined by looking at how the order parameter changes near the sample surface, i.e. the absolute value of the order parameter should approach zero as and the magnetic field should decrease into the superconductor from the surface as . These were fitted to the average order parameter and field profiles of the samples where the averaging was done to the middlemost splice perpendicular to the magnetic field averaged over the periodic direction. Examples of the fits to the order parameter profile and magnetic field profile are shown in figure 4. The fits were limited to the region where the profile is linear in semi-log plot near the sample boundary. The accuracy of the fits is better in the thicker samples where the vortices are not so near the boundary as in the thinner samples. For consistency, these fits were also checked against fitting the single vortex profiles which yielded similar results but proved to be difficult to do systematically.
No systematic magnetic field dependence of the fitted characteristic lengths was found. Thus, the fit results for a fixed sample thickness with different magnetic field values were averaged and are shown in figure 5. The standard deviations of the fit results for the same thickness are indicated with the error bars. The effective lengths and , given by the fits, decrease with increasing film thickness. The effective lengths for both -values fall on the same curve when scaled with the bulk values and . The inset shows how the effective Ginzburg-Landau parameter is independent of the sample thickness. Therefore, the effective lengths have the same thickness dependence which is shown as the black lines that have the form
where is the only free parameter obtained by fitting. For equation 2 can also be written as
where . Equation 3 is quite similar to the given by Pearl  for the interaction length of two Pearl vortices in a thin superconducting film. Both effective lengths diverge at and at thicknesses less than the simulation gives normal state as the minimum energy. The divergence of the effective lengths is obviously unphysical since, for example, Pb thin films are superconducting down to 5 monolayer thicknesses . However, the general trend that both effective lengths increase as the thickness decreases is in agreement with experiments that show how the critical field decreases () as the thickness decreases . The decrease of the critical field with decreasing thickness continues up to the point when the thickness of the sample starts to restrict the electron mean free path after which the critical field starts to increase [18, 19]. Since the electron mean free path related effects are not included in our model it is natural that the simulation results start to deviate from the experiments when the mean free path becomes important for the real sample.
It is the stray field outside the sample that is responsible for the change in effective lengths and vortex-vortex interaction with film thickness. This is in good agreement with the work of Brandt in ref.  that shows when taking the stray fields into account the vortex lattice behaves like it had an effective that tends to infinity as thickness goes to zero. The effective in Brandt’s work is defined through the behaviour of the bulk shear modulus of the vortex lattice which reflects the change in the vortex-vortex interaction to repulsive. In this work stays constant as it is defined by the ratio of the effective lengths which differs from Brandt’s definition. The result that stays constant is interesting, because the vortices seen in thin type I films are explained with electron mean free path increasing which drives the material into type II behaviour. Our simulations predict that singly quantized vortices could be stabilized just by the geometry of the sample. An alternative viewpoint to effective is to consider how the expansion of vortices near the surface strengthens the repulsive interaction of vortices which leads to enhanced repulsion in thinner films .
Naturally, the giant vortices grow monotonically in size with the increasing number of flux quanta. Examples of the radial profiles of the vortices with different quantum numbers are shown in figure 6a. A practical estimate for the radial shape of a vortex, known as the Clem model , is given by
where is the radial coordinate and is the radius of the vortex. The Clem model gives a rather good estimate of the single vortex shape near the core. The Clem model can be extended to give estimates of multiple quantum vortices with the help of Lasher’s equation  which gives the order parameter of the giant vortex as
which is also plotted as the coloured dashed lines in figure 6a for . While (6) catches the general flattening of the vortex core with increasing number of quanta, it deviates quite a lot from the simulation results outside the vortex core. The deviation is not surprising given that the analytical expression is based on linearized Ginzburg-Landau equations and our simulations are quite far from the situation when the linearization is valid. A more practical function that could be used in fitting experimental data would be
which is similar to (6) except that and the scaling is changed from to . These are shown as black lines in figure 6a. The size of the giant vortex increases almost linearly with increasing number of quanta which is shown in figure 6b. The size of the vortex was characterized by the cross-sectional area of the vortex where . The linear behaviour is to be expected because the magnetic flux density is approximately constant in the vortex core which means that each flux quantum takes up the same amount of space.
In our work the giant vortices do not form a honeycomb lattice of vortices of the same vorticity as was predicted by Lasher . Instead, giant vortices of different vorticities are mixed into the same lattice quite thoroughly, the smaller vortices filling in the gaps between the larger ones. The distributions of the giant vortices are shown as histograms in figure 7 for the case. The thickest film with is in the Meissner state at the lowest magnetic field. Just above the first critical field ( mT) there are giant vortices with and a few with . The one vortex with is sitting in a pinning site. Increasing the magnetic field brings in the higher vorticities up to . Figure 7b shows how the singly quantized vortices change into giant vortices when increasing the film thickness at mT. In the two thinnest films there are only singly quantized vortices present. The total flux decreases with increasing film thickness as expected from the change in the demagnetization factor with film thickness. The giant vortices seen in the simulations are a natural extension to the flux tube structure seen in the experiments on lead films in the intermediate state . The magnetic flux in a normal state spot in the intermediate state is quantized and a giant vortex is simply a very small flux tube in the intermediate state. Interestingly, the intermediate state continuously shrinks into a singly quantized vortex structure that is similar to the mixed state when the film thickness is decreased. The giant vortex structure serves as a transitional phase in between them.
We have simulated low- superconducting thin films of different thicknesses. The results show that giant vortex structures are a stable solution to the Ginzburg-Landau equations even in low magnetic fields provided that the superconductor is of type I and of intermediate thickness. No complex geometries or pinning sites were needed for these giant vortices to form. On the contrary, the pinning sites broke the cylindrical symmetry of the giant vortices which turned them into multivortices. In the case of , no giant vortices were seen. Based on the distribution of the phase contours one can speculate that the reason behind the giant vortex state could be a short range attraction and long range repulsion between the flux lines. Smaller vortices were seen to coexist with giant vortices, the smaller ones filling the gaps. In addition, the size of the vortices was seen to increase with decreasing film thickness in a manner bearing some similarity with Pearl vortices. The results predict that singly quantized vortices could be stabilized in type I superconductors by the thin film geometry alone which was unexpected because the model does not take electron mean free path effects into account.
-  Hasegawa S, Matsuda T, Endo J, Osakabe N, Igarashi M, Kobayashi T, Naito M, Tonomura A and Aoki R 1991 Phys. Rev. B 43 7631
-  Cren T, Serrier-Garcia L, Debontridder F and Roditchev D 2011 Phys. Rev. Lett. 107 097202
-  Engbarth M A, Bending S J and Milošević M V 2011 Phys. Rev. B 83 224504
-  Silhanek A V, Gutierrez J, Kramer R B G, Ataklti G W, Van de Vondel J, Moshchalkov V V and Sanchez A 2011 Phys. Rev. B 83 024509
-  Kramer R B G, Silhanek A V, Van de Vondel J, Raes B and Moshchalkov V V 2010 Physica C 470 758
-  Poole Jr. C P, Farach H A, Creswick R J and Prozorov R 2007 Superconductivity, Second Edition. Academic Press
-  Tinkham M 1963 Phys. Rev. 129 2413
-  Maki K 1965 Ann. Phys. (N.Y.) 34 363
-  Lasher G 1967 Phys. Rev. 154 345
-  Schweigert V A, Peeters F M and Deo P S 1998 Phys. Rev. Lett. 81 2783
-  Berdiyorov G R, Milošević M V and Peeters F M 2006 Physica C 437-438 25
-  Shi L M, Zhang L F, Meng H, Zhao H W, Zha G Q and Zhou S P 2009 Phys. Rev. B 79 184518
-  Sweeney M C and Gelfand M P 2010 Phys. Rev. B 82 214508
-  Palonen H, Jäykkä J and Paturi P 2012 Phys. Rev. B 85 024510
-  Duan F and Guojun J 2005 Introduction to Condensed Matter Physics, Volume 1. World Scientific Publishing, Singapore
-  Pearl J 1964 Appl. Phys. Lett. 5 65
-  Özer M M, Thompson J R and Weitering H H 2006 Phys. Rev. B 74 235427
-  Cody G D and Miller R E 1966 Phys. Rev. Lett. 16 697
-  Lischke B and Rodewald W 1974 Phys. Status Solidi B 63 97
-  Brandt E H 2005 Phys. Rev. B 71 014521
-  Carneiro G and Brandt E H 2000 Phys. Rev. B 61 6370
-  Clem J R 1975 J. Low Temp. Phys. 18 427
-  Lasher G 1965 Phys. Rev. 140 A523
-  Prozorov R, Giannetta R W, Polyanskii A A and Perkins G K 2005 Phys. Rev. B 72 212508