CHAPTER 2 The relaxation method investigation of the static state of the Abrikosov lat-
2.2 Discrete TDGL on the hexagonal grid with Wilson link variables
2.2.3 Discretized TDGL equation
Substitute Eq.(2.24) into Eq.(2.1), the space part for TDGL equation is discretized, the formula is
t . In order to simplify the equation, the t + t terms put in left hand side, and t terms put in right hand side. Finally, The TDGL equation or the equation of motion of order parameters as following
n1;n2(t + t) = n1;n2(t) + tF n1;n2; Un1;n2 ; (2.26) where
F n1;n2;Un1;n2 = (2.27)
I assume that the superconducting sample is large enough. In this case there is a well developed periodic vortex structure and the sample’s boundaries very weakly in‡uence the bulk. Conse-quently the periodic boundary condition (PBC), see Fig.2-2, is suitable for simulation of such a system. The PBC are a set of boundary conditions that are often used to simulate a large system by modelling a small part of it which is located far from its edge. The present system is more complicated due to local gauge invariance so that the magnetic translation group [34]
should be considered. Periodicity is only up to a phase factor (AB phase) which is di¤erent for di¤erent locations and directions. The relation between the order parameters in the boundary as following
In order to simplify the boundary conditions, we let
Figure 2-2: Periodic boundary conditions. The magnetic translations between the opposite along the boundaries on the grid are indicated. In actual computation four additional lines of
"images" are used points. These are outside of the sample which is shown as a gray area.
L = 2N s2;
where N is a integer number, therefore, the boundary conditions become
0;n2 = L;n2; L+1;n2 = 1;n2; 0;0 = L;L: (2.29)
Discrete TDGL equation and simulation result in rectangular grid
2.2.5 Rectangular grid and boundary conduction
The points on the sample are also described by two integers(in unit of ) n = (n1; n2), where n1; n2 = 1~nmax.The sample is following direction
rn= a4
s (n1 1+ n2 2) ; (2.30)
with unit vector 1 = (1; 0) and 2 = 0;
p3
2 . The unit vectors show in Fig.2-3 as following and the shape of grid is rectangle.
The boundary condition for rectangle is periodic boundary condition with magnetic
trans-Figure 2-3: The unit vectors and vortex distance. a is the distance between vortex, where de…ne in Eq.(2.14). 1; 2 are and unit vectors, where de…ne in Eq.(2.13)
formation group which is similar as hexagonal grid. The relation between order parameters on the boundary is
0;n2 = exph
i s2L i
L;n2; L+1;n2 =h
i s2L i
1;n2 (2.31)
n1;0 = exph
i s2L i
n1;L; n1;0=h
i s2L i
n1;L+1: Chose sample L = 2N s2, the boundary condition became following
0;n2 = L;n2; L+1;n2 = 1;n2 (2.32)
n1;0 = n1;L; n1;0= n1;L+1:
In rectangle gird, link variables only have two directions, one for x direction another for y direction. The formula for link variable U = exp i n1;n2 , where = 1; 2. The must di¤erent of link variables between hexagonal grid and rectangle grid is n1;n2 ,for the link variable
represent the phase di¤erent between order parameters, in the other word, link variables are di¤erent in di¤erent shape of grid. The formula of n1;n2 as following(detail see Appendix )
1
Although the formulas seem the same, but there are di¤erent link directions. Substituting Eq.(2.14) into Eq.(2.33), the …nal formulas are
1
2.2.6 Free energy and TDGL equation in rectangular grid
Similar as hexagonal grid, varied the di¤erential terms with link variable dxd (x) !
NP
The factor 4=3 cause by the lattice distance for x direction(d) and y direction p23d are di¤erent.
And p
3=2 come form volume integration.
Replace continuum free energy by discretized free energy in Eq.(4.3), the formula is
d Next, discrete time part
d
dt n1;n2 = n1;n2(t+ t) n1;n2
t : (2.37)
Similar as discretized TDGL equation in hexagonal grid, the TDGL equation can be written as equation of motion of order parameter …nally. The formula is
n1;n2(t + t) = n1;n2(t) + tF n1;n2; Un1;n2 ; (2.38)
2.3 Simulation result and comparison of two grids.
2.3.1 Rectangular grid
The theory of 2LLL was de…ned by D. Li, A. M. Malkin, and B. Rosensten[15], the de…nition as following
2
LLL = ah
A
; (2.40)
0.2 0.4 0.6 0.8 1.0 b 0.1
0.2 0.3 0.4
2
Figure 2-4: Average super‡uid density D j j2E
as function of magnetic …eld(rectangular grid).
The red line is values for the analytic expression Eq.(2.40).The blue points are simulation values for rectangular grid.
where A= h
4i
hj 2ji2 = 1:16 for hexagonal structure[32].When the magnetic …eld is high enough(ah 1), the order parameters belong in LLL. Fig2-4 show relation between 2 and b in t = 0,
The red line is theory value for Eq.(2.40), and the blue points are simulation values. When ah become smaller, the
D j j2
E
for theory and simulation are closer. Superconductivity would totally break down when ah < 0, since there is no condensation potential. The vortex structures in di¤erent magnetic …elds are shown in following (Fig.2-5~Fig.2-7).
It’s clearly that the distribution of vortex structure is hexagonal(vortex solid). The vortex density become larger when magnetic …eld become larger. Moreover, the values for 2 become smaller. The Cooper pairs are broken by the magnetic …eld.
2.3.2 Hexagonal grid
Hexagonal gird is similar as Arbikosove lattice, the relation between average order parameters
2 and magnetic …eld show in Fig.2-8
Figure 2-5: Vortex structure( b = 0:1). The picture (a) is vortex lattice in real space, the vortex structure have the hexagonal symmetry. The picture (b) is the vortex spectrum in qusi-monentum space, the peaks are sharp and the peak distribution has hexagonal symmetry, therefore, this tell that the vortices in the solid phase.
Figure 2-6: Vortex structure( b = 0:3). Similar as Fig.2-5, the vortex lattice in real space(picture(a)) has hexagonal symmetry, and the peak distribution for the vortex spec-trum(picture(b)) still has hexagonal symmetry. The order parameters(j j2) in large magnetic
…le(b = 0:3) are smaller than small magnetic …eld(b = 0:7), because of the Cooper pairs are broke by the magnetic …eld
Figure 2-7: Vortex structure( b = 0:7). Similar as Fig.2-5, the vortex lattice in real space(picture(a)) has hexagonal symmetry, and the peak distribution for the vortex spec-trum(picture(b)) still has hexagonal symmetry. The order parameters(j j2) in large magnetic
…le(b = 0:7) are samller than small magnetic …eld(b = 0:3).
Figure 2-8: Average super‡uid density D j j2E
as function of magnetic …eld(hexagonal grid).
The red line is values for the analytic expression Eq.(2.40). The purper points are simulation values in hexagonal grid, and almost the same as values for rectangular grid.
The red line is theory value where come form Eq.(2.40), the purple points are simulation value. The value for 2 rectangular grid and for hexagonal grid are almost the same. The vortex structure in di¤erent magnetic …eld as shown in the following, see Fig.2-8~Fig.2-10.
Similar as vortex structure on rectangular grid, the vortex structure are also has hexagonal symmetry. The hexagonal grid is more symmetric than rectangular grid, however, the vortex lattice structure on the rectangular grid is more beautiful than on the hexagonal grid, since the relaxation time for hexagonal grid is larger than rectangular grid.
Figure 2-9: Vortex structure( b = 0:1). The vortex lattice has hexagonal symmetry.
Figure 2-10: Vortex structure( b = 0:3). The vortex lattice also has hexagonal symmetry.
Figure 2-11: Vortex structure( b = 0:7). The vortex lattice also has hexagonal symmetry, and the values of order parameters become samller, the Cooper pairs are broken by the magnetic
…eld.
Chapter 3
Model and simulation method for
the vortex dynamics in clean system
3.1 Time Dependent Ginzburg-Landau theory in continuum
3.1.1 Electric …eld in a mixed state superconductor and the ‡ux ‡ow
In previous chapter I discussed the static vortex system. The system would eventually relax to a lowest free energy state. The relaxation dynamics was not realistic in a sense that e¤ects of the electric …eld were neglected (which is not essential for the static properties of the system since the …nal state is the same). In this chapter I start to consider dynamics of a type II superconductor in external magnetic and electric …elds. When the electric …eld applied to the system, it acts on vortices as an external force. Consequently the TDGL equation should be modi…ed as following
2
~2 2m
@
@t +ie
~ = F [ ; A]
; (3.1)
where the free energy was de…ned in Eq.(chapter II). Here, is an electric potential and is the inverse di¤usion constant. The electric …eld is E = r 1c@A@t, while the magnetic …eld is given by B = r A. The TDGL equation is invariant under the gauge transformation
= exp ie
where is an arbitrary function of space and time. Choosing the zero scalar potential gauge is convenient for the simulations which follow:
= 0 (3.3)
A = 1
2 Bxbi + 1
2By cEt bj:
Here both magnetic and electric …elds are assumed to be constant. Dynamics of electromagnetic
…eld should in principle also taken into account by the Maxwell quietens. However as explained in chapter II, for strongly type II superconductors magnetization is of order 1= 2 and hence negligible. Vortices overlap and their magnetic …elds become homogeneous. Electric …elds in the ‡ux ‡ow state are also homogeneous (except at very low values in the presence of pinning, see chapter IV) and their dynamics can be neglected as well except for the normal conductivity
J = nE: (3.4)
Substituting the free energy F [ ; A] of Eq.(chapter II) into equation, the TDGL equation becomes The formula is formally identical to the static TDGL equation, however the vector potential now contains a time dependent part. The time derivative term is purely dissipative (relaxational) and carries an implicit the information about the vortex velocity. When vortices move faster the frequency of phase change becomes larger, and …nally, the vortex structure would break down and the order parameter goes to zero. Because of the energy dissipation, resistance of
type ll superconductors is not zero. The formula for supercurrent density is the same as in the static case.
Js = ie ~
2m ( r r ) + e 2
mscj j2A: (3.6)
Note that the supercurrent density is also invariant under gauge transformation.
3.1.2 Dimensionless electric …eld and conductivity
Using the same units as in static case, ; tGL, the derived units of electric …eld is EGL =
ctGLHc2= c2 Hc2. In Y BCO this …eld is very large. Therefore the dimensionless electric …eld will be de…ned as
= E=EGL: (3.7)
Consequently the dimensionless vector potential is
a= 1
2 by;1
2bx t ; (3.8)
here x and y in unit of , and t in unit of tGL. The velocity of vortices, which in physical units is V = cE=B, in dimensionless units reads:
v = VtGL The dimensionless TDGL equation takes a form
@
@t = bH ah + ( ) ; (3.9)
where bH = 12D2 2b,
ah = 1 t0 b
2 ; (3.10)
and the dimensionless supercurrent density
js= i
2( r r ) + j j2a; (3.11)
are the same as in the static case, with the only di¤erence being the electric …eld term in vector potential, Eq.(3.8).
3.2 Discrete TDGL equation on rectangular grid
3.2.1 Rectangular grid
The discretized version of the dynamic TDGL equation is similar to the statistic one, the only di¤erence being the link variables for the gauge …eld. The link variable
Un1;n2 = exp n1;n2 ; (3.12) For the nonequlibrum system, the link variables change with time, the link variables in 1 and
2 directions are
where the distance between vortex aM =q
p4
The discretized dynamic TDGL equation is similar as static case, the only di¤erence is the link variables. The formula show as following
d 3.2.2 Discretized supercurrent density
The the de…nition of supercurrent density
j = fgrad
a ; (3.17)
where = 1; 2, f is the dimensionless free energy. The supercurrent density in x and y direction as following
Next, I expand discretized supercurrent density formula to get continuum formula. Use Taylor expansion to expand link variables Un1;n2 and order parameters /n
1;n2,
where = 1; 2, both 1n1;n2 and 2n1;n2 were de…ned in Eq.(3.14). The discretized supercurrent density in x direction can be write as
jx = i Similarly, The discretized supercurrent density in y direction
jy = i
Both the continuum limit for discretized supercurrent density in x and y direction are the some as continuum formulas Eq.(3.11).
3.3 Simulation result and discussion
3.3.1 Parameters for Y BCO
To the following parameters of Y BCO sample as used. Tc = 100K, Hc2 = 140T , m = 2me. The coherence length and the penetration depth:
= ~
The Lattice spacing is magnetic …eld b = 0:1 in this simulation a ' 13nm. The unit of time is:
n= 8c2 2 = 5 103ohm 1cm 1= 1:15 103 1012sec 1 = 4:545 1015
! = 8c22 n= 8 (1:53140)2(4:545 1015)
(3 1010)2 = 1:06 s=cm2
tGL= 2=2 = 1:24 10 14s; (3.22)
and unit of electric …eld
EGL= Hc2
3.3.2 Super‡uid density at ‡ux ‡ow
In contrast to the static case, under the in‡uence of the Lorentz force the vortices move.
At low electric …elds, the vortex structure approximately remains hexagonal for an isotropic material (assumed in the present study) in directions perpendicular to magnetic …eld. See Fig.
3.1(a) compared to perfect hexagonal lattice in statics at the same …eld b = 0:4. The moving lattice is not exactly hexagonal since electric …eld breaks explicitly the hexagonal symmetry,
Figure 3-1: Super‡uid density with small electric …eld. The picture(a) show the vortex lattice without electric …eld. As show in previous chapter, the vortex lattice has hexagonal symmetry.
The picture(b) show the vortex lattice with small electric …eld. The electric …eld apply in y direction, and the vortices move in x direction. The moving vortex lattice is not exactly hexagonal.
see theoretical symmetry considerations in[7]. In the large electric …elds, the vortex matter
‡ows in my simulation like a liquid(In fact, the theoretical results for moving lattice still has hexagonal). For example, when I consider the vortex ‡ow for ah = 0:45 and electric …eld applied in the y direction, the vortex structures in di¤erent electric …eld show in Fig.3-1.
3.3.3 De…nition of nonlinear conductivity and comparison of simulation with the analytic results
Using expansion in small parameter
ah(v) = 1
2 1 t b v2 = 1
2 1 t b
2
b2 : (3.27)
D. Li, B. Rosenstein and M. Malkin [15] calculated the nonlinear conductivity for clean type-ll superconductors. In the ‡ux ‡ow regime, in addition to the normal state conductivity, there is a large contribution form the Cooper pairs represented by order parameter …eld. The nonlinear
Figure 3-2: Moving lattices with pinning. Two snapshots at times t=2000 (a) and t=2100 (b).
The vortices move in the x direction
conductivity will be de…ned as follows
s = Js
E = i~e
2mE r + ie
~cA + cc: (3.28)
Using the unit 0 = 4c2 2, the dimensionless nonlinear conductivity in the LLL approximation is proportional to the super‡uid density,
LLL= i
2vh LLL@y LLL LLL@y LLLi =D
j LLLj2E
= ah(v)
A(v)ev2: (3.29) Numerical simulation allows to check the analytic theory and extend it beyond perturbation theory applicability range.
I simulated the case b = 0:9, t = 0, for which the parameter ah(v) = 0:05 is small, the order parameters are within the LLL approximation applicability range. For velocity is low (corresponding to low electric …eld), so that the Abrikosov beta does not change from its static value, A(v) ' 1:16. The comparison of the analytic theory and the simulation is shown in Fig.3-3. For A(v) 6= 1:16 at large electric …eld , the value for the theoretical line have to be modi…ed. As expected, the di¤erence between the theory and the simulation becomes smaller when ah(v) become smaller. That is, when number of Cooper pairs in becomes smaller, most of them stay in LLL.
0.05 0.10 0.15 0.20 0.25 0.30 0.01
0.02 0.03 0.04 0.05
Figure 3-3: Comparison of simulated ‡ux ‡ow nonlinear conductivity with the analytic results.
The blue points are simulation values for j= . The red points are simulation value for 2 =
LLL which provides the …rst approximation and the green line is the full analytic expression Eq.(3.29).
0.05 0.10 0.15 0.001
0.002 0.003 0.004 0.005 j
Figure 3-4: Comparison of simulated nonlinear j curve with analytic result. The green line is theory value Eq.(3.29). The red points are simulation values for LLL .The blue points are simulation for supercurrent density values.
Next, the j characteristic is shown in Fig.3.4. The blue points are simulation result for supercurrent density, and the red points are simulation result for super‡uid density multiplied by electric …eld, 2 . Note that 2 is the LLL conductivity according to Eq.(3.29). The green line …t well with red points when ah become small.
The theory is inapplicable when ah(v) is large, so I consider two cases with large ah(v). In the …rst the magnetic …eld b is …xed at b = 0:4, while temperature t0 has two values 0:0 and 0:3 (just change of ah, without thermal noise on the mesoscopic level). Electric …eld ranges from 0 to 0:2. The parameters ah(0) for this case are therefore, 0:3 and 0:15 respectively. See Fig.3-5. The red triangles are values for ah = 0:3, while the blue squares are values for ah = 0:15.
The supercurrent current density at higher temperature is smaller since the temperature breaks Cooper pairs.
In the case II, the magnetic …eld b has two values 0:4 and 0:7, while the temperature t0 is
…x, t0 = 0:0. Electric …eld from 0 to 0:2. The parameters ah(0) for case Il are, therefore, 0:3 and 0:15 respectively. See Fig3-6. The supercurrent density in large magnetic …eld(b = 0:7) is
Figure 3-5: J-E curve for di¤ erent ah (…x b). The red triangles are values for ah(0) = 0:3.
The blue squares are values for ah(0) = 0:15.
Figure 3-6: J-E curve for di¤ erent ah (…x t). The red triangles are values for ah(0) = 0:3. The blue squares are values for ah(0) = 0:15:
smaller the small magnetic …eld(b = 0:3) under the same electric …eld, for Cooper pairs density are less in large magnetic. In addition, the j curve in b = 0:7 is more linear then in b = 0:3, the reason is that when vortex move faster, the Cooper pairs are more di¢ cult to transport.
Chapter 4
Vortex statics and dynamics in superconductor with periodic pinning.
4.1 Time dependent Ginzburg-Landau theory with periodic pin-ning
The di¤erence between the clean superconductor’s GL free energy and the energy in the presence of the arti…cial pinning centers located ra is the pinning potential:
V (r) = Tc
X
a
U (r ra) : (4.1)
The TDGL equation should be modi…ed as follows
2
~2 2m
@
@t = ~
2m r ie
~cA
2
+ Tc 1 t0+ U (r ra) b0j j2 : (4.2)
Using the same unitsas in the clean case, the dimensionless TDGL equation becomes
@
@t = 1
2(r ia)2+ 1
2 1 t0+ u (r ra) ( ) ; (4.3)
where u(r ra) is the dimensionless pinning potential. The dimensionless current density as de…ned before j = jn+ js, where
jn = GL (4.4)
js = i
2( r r ) + j j2a:
As discussed in previous, the dissipation mainly comes from the vortex ‡ow. The electric resistance comes also from normal electrons, when the electric …eld is "allowed" to enter the sample due to the ‡ux ‡ow. Pinning slows and eventually stops vortices from moving. This is the reason why there exists the critical current.
4.2 Pinning distributions
Two kinds of pinning array are discussed in this thesis. The …rst one assumes that small parts of the sample are normal or even insulating. An example is the columnar defect or simply a hole. Obviously the order parameter inside such a center (neglecting the proximity e¤ects) is always zero. The weaker Tc pinning is just a local modi…cation of the critical temperature at the center’s site. It in‡uences the condensation potential, Eq.(4.6). In other words, the vortices will not be always con…ned to the pinning centers. The simulation results are shown in following subsections.
The shape of pinning sites I simulated are rectangles. The pining distribution is hexagonal, as shown in Fig.4-1. I assume that vortices outnumber the pinning sites, namely, the magnetic
…eld is above the so called matching …eld. Therefore the …lling factor parameter de…ned as,
f = pinning sites0 number
vortices0 number ; (4.5)
. For example the f = 1=2 case is shown in Fig.4-1
Figure 4-1: Pinning distribution. The black rectangles are pinning sites. The distribution of their centers is hexagonal, commensurate with the Abrikosov lattice.
4.3 Strong pinning array
The magnetic …eld b = 0:4 and the temperature t = 0 in the following simulation.
4.3.1 Super‡uid density
The vortex distribution is in‡uenced by pinning. In this simulation all the pinning distributions are hexagonal, however the size of the pinning sites varies (pinning size 4; 9; 16). The pinning factor f = 0:5. The results are shown in Fig.4-2~Fig.4-4. In Fig.4-2 for the size of 4, the vortex structure still has hexagonal symmetry, while for larger sited in Fig.4-3 and Fig.4-4, the vortex structure does not remain hexagonal.
4.3.2 Dynamics in the presence of pinning.
a. Interstitial vortices
The vortices in the superconductors with arti…cial pinning array can be separated into two di¤erent sets: the pinned vortices which are trapped on the pinning centers, and the interstitial vortices which can be considered as "free vortices". If the number of vortices is smaller than the number of pinning sites, all of them are likely to be pinned. If magnetic …eld is equal to the matching …eld, there are no interstitial vortices, however when vortices outnumber the
Figure 4-2: Vortex structure with periodic pinning (pinning size = 4 ). The vortex lattice still has hexagonal symmetry, while the vortex core in the pinning centers are larger then that of the intersticial vortices.
Figure 4-3: Vortex structure with periodic pinning(pinning size = 9 ). The vortex structure is di¤erent compared to the clean system. Pinning leads to expansion of the vortex cores.
Figure 4-4: Vortex Structure with periodic pinning(pinning size = 16 ). Similar to Fig.4-3, the pinning vortices lead the vortex distortiontion increasing the unit cell
sites, some of them will be "liberated" due to repulsion from the vortices already pinned. On Fig.4-5 the dynamics of the interstitial vortices is shown. The electric …eld is applied in the y direction, thus the vortices move along the x direction. The parameters for the simulation are:
the electric …eld = 0:0001. Pinning centers size is 4. Two snapshots at times t = 2000 (a) and t = 2100 (b) for the case when for each rectangular pinning center there are two vortices.
Two di¤erent kinds of vortices are seen: the interstitial vortices have smaller cores and move, and the pinned vortices which have larger cores (like in the static case) and don’t move.
For the pinning force is larger that the driving force, the pinned vortices are not moving,
For the pinning force is larger that the driving force, the pinned vortices are not moving,