Table Of ContentApproximate Ginzburg-Landau solution for the regular flux-line lattice.
Circular cell method
W. V. Pogosova, K. I. Kugelb, A. L. Rakhmanovb, and E. H. Brandtc
a Moscow Institute of Physics and Technology, 141700
Dolgoprudny, Moscow region, Russia
1
0 b Institute for Theoretical and Applied Electrodynamics,
0 Russian Academy of Sciences, 127412, Moscow, Russia
2 c Max-Planck-Institut fu¨r Metallforschung,
D-70506, Stuttgart, Germany
n
(February 1, 2008)
a
J
0
2 be constant. In this case, the flux-line lattice (FLL) can
A variational model is proposed to describe the magnetic betreatedasasetofindependentvortices: the magnetic
] propertiesoftype-IIsuperconductorsintheentirefieldrange
flux density is a linearsuperpositionofthe contributions
n
o between Hc1 and Hc2 for anyvaluesoftheGinzburg-Landau of individual vortices, and each contribution coincides
parameter κ > 1/√2. The hexagonal unit cell of the trian-
c with the field of an isolated vortex. The energy of the
gularflux-linelatticeis replaced byacircle ofthesamearea,
- system is the sum of the self-energies of the vortices and
r andtheperiodicsolutionstotheGinzburg-Landauequations
p their pairwise interaction [5]. At high fields, the Lon-
within this cell are approximated by rotationally symmetric
u donmodelloosesits applicabilitybecause the fractionof
solutions. The Ginzburg-Landau equations are solved by a
s the total volume of superconductor occupied by the vor-
. trialfunctionfortheorderparameter. Thecalculatedspatial
t texcoresis nolongersmall[5]. Differentapproximations
a distributions of the order parameter and the magnetic field
m are compared with the corresponding distributions obtained were proposed to include the vortex cores and to extend
- bynumericalsolutionoftheGinzburg-Landauequations. The the applicability limits of the London model [6].
d comparisonrevealsgoodagreementwithanaccuracyofafew The problem of solving the GL equations for the ide-
n percent for all κ values exceeding κ 1. The model can be ally periodic FLL can be simplified considerably by re-
o extended to anisotropic superconduc≈tors when the vortices placing the hexagonal unit cell of the vortex lattice by
c
are directed along one of the principal axes. The reversible a circle of equal area (Wigner-Seitz approximation). In
[
magnetization curve is calculated and an analytical formula this approach both the order parameter and the mag-
2 for the magnetization is proposed. At low fields, the theory netic flux density within the cell have axial symmetry.
v reduces to theLondon approach at κ 1, provided that the The presence of other vortices is taken into account by
7 exact value of Hc1 is used. At high fi≫elds, our model repro- the boundary condition: the supercurrentdensity equals
5
duces the main features of the well-known Abrikosov theory. zero at the cell boundary. This method was used in Ref.
0
The magnetic field dependences of the reversible magnetiza- [7],whereanexplicitexpressionforthe magnetizationin
1
1 tion found numerically and by our variational method prac- low fields was found, and the results are in good agree-
tically coincide. The model also refines the limits of some
0 ment with the predictions of the London model [5,8]. In
0 approximationswhichhavebeenwidelyused. Thecalculated Refs. [9–13], the GL equations and the equations of the
t/ magnetization curves are in a good agreement with experi- microscopic theory of superconductivity were solved nu-
a mental data on high-Tc superconductors. merically in the framework of circular cell approaches
m
anditwasshownthatthisapproximationnotonlyyields
PACS numbers: 74.20.De, 74.60.Ec, 74.25.Ha
d- good results at low induction but also at H ∼ Hc2. A
similarapproachwasused in Ref. [21],where the system
n
o of GL equations was reduced to a single equation which
c canbesolvednumerically. TheWigner-Seitzapproxima-
: tion can be extended even to the case of exotic pairing
v
i I. INTRODUCTION symmetries, when the order parameter has two complex
X
components [14].
r A numerical method to find the periodic solutions to
a The solution to the Ginzburg-Landau (GL) equations
the GL equations was developed in Refs. [15,16]. This
foundbyAbrikosov[1]wasusedwidely[2,3]tostudythe
exact method accounts for the actual symmetry of the
propertiesoftype-IIsuperconductorsatlowandhighap-
vortex lattice. It allows calculating the spatial distri-
plied magnetic fields H, i.e. close to the lowerandupper
butions of the magnetic field and the order parameter
criticalfieldsHc1 andHc2,respectively. AtH Hc1 the
∼ withintheunitcell,theelasticshearmodulusoftheFLL,
intervortex spacing is much larger than the vortex core
and the magnetization for any FLL symmetry and any
sizeiftheGLparameterislarge,κ 1. Thereforeinthe
≫ induction B and GL parameter κ with any desired ac-
London model, which is commonly used at low fields [4],
curacy. However, until recently there was no adequate
theorderparameterinthe superconductorisassumedto
1
approach allowing one to find the magnetization analyt- free energy could not be taken into account in this ap-
ically in the entire field range Hc1 < H < Hc2 and to proach. The Clem-Hao model was further extended to
obtain explicit formulas which may be used to analyze include anisotropy [23,26]. This approximation is now
experimental data. widely used for the analysis of the experimental data on
In Ref. [17] Clem proposed a model to solve the GL magnetization of type II superconductors [27–29].
equations using a trial function for the order parameter However, it has been recently shown in Ref. [30] that
(or GL function) ψ : theClem-Haomodelhassomedrawbacks. Itwasargued
| |
that the procedure of obtaining the local magnetic flux
r
ψ =f(r)= , (1) densitybyalinearsuperpositionofcontributionsofindi-
| | r2+ξv2 vidualvorticesintheformusedinRef.[23]isvalidonlyat
low fields. The applicationofthis approachto the entire
where ξ is a variational paprameter. This model yields
anapprovximateexplicitexpressionforthelocalmagnetic field range Hc1 < H < Hc2 leads to an appreciable dis-
agreementbetweentheClem-HaomodelandAbrikosov’s
field of an isolated vortex:
high-field result [30]. In the originalpapers of Hao et al.
1 [23,25,26], this disagreement was made up by the use of
h(r)= κξvK1(ξv)K0 r2+ξv2 , (2) a non-selfconsistent field dependence of the variational
(cid:16)p (cid:17) parameters. Then, in calculating the magnetic free en-
where K are modified Bessel functions. Here and be-
n ergy in Ref. [23] the lattice sum was approximated by
low the following dimensionless variables are used: dis-
an integral. As it was shown in Ref. [30], this procedure
tance r, magnetic flux density h, and order parameter f
leads to a noticeable error in the magnetization at low
are measuredin units of λ, Hc√2, α/β, respectively, fields. This error and the use of an inaccurate value of
−
whereλistheLondonpenetrationpdepth,Hc isthether- Hc1 have led the authors [23,25]to the conclusionabout
modynamic critical field, and α and β are the GL coef- thequantitativeincorrectnessoftheLondonapproachat
ficients. In this notation, we have Hc2 = κ, Φ0 = 2π/κ, low fields. Note that in Ref. [21] it was argued too that
where Φ0 is the magnetic flux quantum. For κ 1 the the Clem-Hao model overestimatesthe effect of suppres-
≫
minimization of the free energy gives ξ √2/κ and
v sion of the order parameter in the vortex cores at small
≈
fields.
1
Hc1 = (lnκ+ε), (3) In this paper, we propose a variational model for the
2κ
description of the regular flux-line lattice in a more con-
where ε =0.52. The exactvalue ε=0.50 was calculated sistent fashion as compared to Ref. [23]. Our variational
in Ref. [18] (see also Ref. [19]) from the numerical solu- procedure is based on Clem’s trial function (1). How-
tiontotheGLequationsforanisolatedvortex. Notethat ever, in contrast to the Clem-Hao model, we do not use
the lower critical field cannot be found self-consistently the superposition of vortex fields. Instead, we apply the
in the framework of the London model, because in this circularcellmethodandcalculatethemagneticfluxden-
approach the magnetic flux density diverges on the vor- sity directly from the second GL equation. The model
tex axis. Therefore, Hc1 should be regarded as a free enablesustofindanalyticalexpressionsforthelocalmag-
parameter in the London expression for the magnetiza- netic field and the order parameter. The results of our
tion [8]. However, the approximation for Hc1 found by variational procedure are compared with the results of
cutting offthe fieldofavortexata distanceequalto the the numerical solution of the GL equations. This com-
coherencelengthξ,isoftenusedintheLondonapproach parison reveals that the analytical formulas for the spa-
as well [20]. Note that recently Clem’s trial function (1) tialdistributionoftheorderparameterandthemagnetic
was applied to the study of the vortex core structure in fieldagreewiththe numericalresults toanaccuracyofa
superconductorswithmixed(d+s)two-componentorder fewpercentinawiderangeofκ>1andB. Byintroduc-
parameter [22]. ingtheeffective-masstensorthe∼theoryisextendedtoin-
Hao et al. [23] (see also Ref. [24]) extended the model cludeanisotropywhenthevorticesaredirectedalongone
[17]tolargermagneticfieldsuptoHc2throughthelinear oftheprincipalaxesofthecrystal. Theresultsforthelo-
superposition of the field profiles of individual vortices. calorderparameterandthe magneticfieldarethenused
Inthismodel,thetrialfunction(1)ismultipliedbyasec- to calculate the magnetization. The resulting expression
ondvariationalparameterf∞ toaccountforthesuppres- for the magnetization is in agreement with the London
sionoftheorderparameterduetotheoverlappingvortex model in small fields at κ 1, and with the Abrikosov
≫
cores. This model enabled the authors [23] to calculate approximation at H Hc2. The field dependences of
∼
the magnetization of type-II superconductors in the full themagnetizationfoundbythevariationalandnumerical
range Hc1 < H < Hc2. Their analytical formula is in a approaches practically coincide. At the same time, the
goodagreementwiththewell-knownAbrikosovhigh-field difference between the numerical result and the Clem-
result. For the case of low fields, Hao and Clem argued Hao approach is considerable. We also found the field
[25] that the London model is quantitatively incorrect dependenceofthemagnetizationintheWigner-Seitzap-
(it does not give the correct asymptotics at H Hc1) proximationclosetoHc2,wheretheGLequationscanbe
→
since the contribution of the vortex cores to the total linearized [5]. The calculated magnetization curves are
2
compared with the available experimental data on some system can be solved only numerically [9]. Nevertheless,
high-T superconductors. in high and small fields some results can be obtained
c
analytically. At small fields, the approximate solution
to the GL equations in the Wigner-Seitz approach was
II. THEORETICAL FORMALISM found in Ref. [7]. Here, the spatialvariationof the order
parameter at κ 1 can be neglected when calculating
≫
In the Wigner-Seitz approximationboth the order pa- the magnetic flux density. In this case, h can be found
rameter and the magnetic flux density within the cell analytically from the second GL equation (6) and the
have axial symmetry. In this case the order parameter boundary conditions (9) and (11). This yields the mag-
ψ can be presented as f(r)exp( iϕ), with radius vec- netization [7]:
| | −
tor r and phase angle ϕ. The free energy density of a
superconductor can be written as the sum of two contri- 1 K1(R) 1
butions: F =Fem+Fcore. Fem is relatedto the energies −4πM(B)=Hc1+ 2κ(cid:20) I1(R) + 2I12(R)(cid:21)−B. (12)
ofmagneticfieldandsupercurrent,andF to the sup-
core
pressionof ψ in the vortexcore. It is easy to show that When H Hc1, the radius of the cell is R 1, and
≫ ≪
| | Eq. (12) can be expanded in powers of R, yielding:
in the framework of the GL theory in the Wigner-Seitz
approach F and F are given by:
em core 1
4πM(B)=Hc1 [ln(2κ(H Hc1)+σ], (13)
2π R 1 2 − − 4κ −
2 2
F = f a +h rdr, (4)
em κAcell Z0 " (cid:18) − κr(cid:19) # with σ = 1.3456. A similar relationship was obtained in
the London limit in Ref. [8] for the regular FLL, with
σ =1.3431 [8] for the triangular lattice.
2π R 1 1 df 2 At high fields, the magnetization is given by [5]:
2 2
F = (1 f ) + rdr, (5)
core κAcell Z0 "2 − κ2 (cid:18)dr(cid:19) # H Hc2
M = − , (14)
4πβ (2κ2 1)
where a is the dimensionless vector potential, R and A −
A = πR2 are the cell radius and area, related to the
cell where β depends only on the symmetry of the FLL [5]:
A
magnetic induction B by A =2π/Bκ.
cell
The two GL equations can be written as: f4d2r
β = . (15)
1 d df 1 2 A Rf2d2r 2
3
r +f f +f a =0, (6)
− κ2rdr (cid:18) dr(cid:19) − (cid:18) − κr(cid:19) Heretheintegralsareta(cid:2)kRenover(cid:3)theareaoftheunitcell.
LetusfindthevalueofβA forthecircularcell. NearHc2,
themagneticfieldundergoesonlyslightspatialvariation.
dh 1
2
=f a . (7) The vector potential in this case is a Br/2. The order
dr − κr ≈
(cid:18) (cid:19) parameter is small at H Hc2, so the first GL equation
∼
(6) can be linearized. The resulting equation has the
The magnetic field and the vector potential are related
analytical solution [31]:
by
1d(ra) κrB B κ κrB
h= . (8) f(r)=sexp Φ − ,2, , (16)
r dr − 2 2 2
(cid:18) (cid:19) (cid:18) (cid:19)
These equationsmust be supplemented by the boundary where Φ is the Kummer function. Factor s depends on
conditions for the magnetic field and the order parame- the nonlinear term in Eq. (6) but it does not affect the
ter: value of β . Using Eqs. (15) and (16), we find β =
A A
1.1576 for the circular cell. This value is close to β =
A
h(R)=he, (9) 1.1596 calculated in Ref. [32] for the triangular lattice.
Thus,inbothlimits oflowandhighfields,the magne-
′ tization in the Wigner-Seitz approximation is in good
f(0)=0, f (R)=0, (10)
agreement with that for the regular triangular FLL,
whichhasthelowestenergy. Inthenextsection,wepro-
dh poseavariationalmodeltosolvetheGLequationsinthe
−2
rf (r)dr =−1/κ at r →0. (11) whole field rangebetweenHc1 andHc2 atany κ>1/√2
with good accuracy.
Condition (11) follows from Eqs. (7) and (8) [5]. The
system of Eqs. (6)-(11) is much simpler than the similar
equationsforthehexagonalunitcell. However,eventhis
3
III. VARIATIONAL PROCEDURE ξv(B,κ)=ξv0
×
6.3 1/2
B B
Using the trialfunction (1)multiplied bya variational 1 4.3 1.01
prefactor, f(r) = f∞r/ r2+ξv2, allows us to solve the × − (cid:18) − 1.05κ(cid:19) (cid:18)κ(cid:19)! ×
second GL equation (7)panalytically within the Wigner- B 0.9 1/2
Seitz cell: 1 0.56 , (22)
× − (cid:18)κ(cid:19) !
h(r)=uI0(f∞ r2+ξv2)+vK0(f∞ r2+ξv2), (17)
whereuandvcanpbefoundfromthebopundaryconditions B2 B 4 1/2
(9) and (11): f∞(B,κ)=(cid:18)1− 2.8κ2(cid:19)× 1−(cid:18)tκ(cid:19) !
u= f∞ K1(f∞ρ) , (18) 2 1/2
1.7B 1.4B
κξv K1(f∞ξv)I1(f∞ρ) I1(f∞ξv)K1(f∞ρ) 1+ 1 , (23)
−
× κ (cid:18) − κ (cid:19) !
v = f∞ I1(f∞ρ) , (19) where t = 0.985, ξv0 is the value of ξv at B = 0. The
κξv K1(f∞ξv)I1(f∞ρ)−I1(f∞ξv)K1(f∞ρ) latter can be calculated from the condition dF/dξv = 0
at B =0:
and we introduced the notation ρ = R2+ξ2. Note
tnhoattm,reigeotrtohueslcyosnpdeiatikoinng(,1t0h)esCinlecme ittrsiadlpefruinvacttiiovnev(c1a)ndnooets κξv0 =√2(cid:20)1− KK0122((ξξvv00))(cid:21). (24)
be equal to zero at the cell boundary. However, df/dr is
small at r = R and the comparison between the results Whenκ 1,Eq.(24)hasthesolutionξv0 =√2/κ(√2ξ
≫
ofvariationalandnumericalmethodsdemonstratesgood in dimensional variables). Eqs. (1), (17)-(19), (23), and
accuracy of the approach. (24)givethedistributionsoftheorderparameterandthe
The values of variational parameters should be found magnetic field within the Wigner-Seitz cell.
by minimization of the total free energy density F =
F + F . Using Eqs. (7) and (4), it is possible to
em core
obtain the following expression for the magnetic energy
density: F = Bh(0) [23]. Taking into account Eqs.
em
(17)-(19), we get: (cid:20)(cid:15)(cid:19)
(cid:18) (cid:3) (cid:3)(cid:19)(cid:17)(cid:19)(cid:24)
%+ (cid:21)
F
Bf∞
F =
em
κξ × (cid:18) (cid:3) (cid:3)(cid:19)(cid:17)(cid:24)
v %+ (cid:21)
F
K0(f∞ξv)I1(f∞ρ)+I0(f∞ξv)K1(f∞ρ). (20) I(cid:19)(cid:15)(cid:24)
×K1(f∞ξv)I1(f∞ρ)−I1(f∞ξv)K1(f∞ρ) %(cid:18)+ (cid:21)(cid:3) (cid:3)(cid:19)(cid:17)(cid:27)
F
The energy related to the spatial variation of the order
parameter in the vortex core is found from Eq. (5) by a
straightforwardintegration [23,26]: (cid:19)(cid:15)(cid:19)
(cid:19)(cid:15)(cid:19) (cid:19)(cid:15)(cid:24) (cid:20)(cid:15)(cid:19)
U(cid:18)G
1 2 2 FIG.1. Thespatialdistributionofthedimensionlessorder
F = (1 f ) +
core 2 − ∞ parameter in the unit cell of the flux-line lattice at different
1 2 magneticinductionsB=0.1Hc2,B =0.5Hc2,B=0.8Hc2for
2 2 2
+2Bκξvf∞(1−f∞)ln 1+ Bκξ2 + κ = 10. The solid lines correspond to the variational calcu-
(cid:20) v(cid:21) lations in the Wigner-Seitz approximation. The dashed lines
+f∞4 f∞4 + Bf∞2 1+Bκξv2 . (21) correspondtothenumericalsolutionoftheGinzburg-Landau
2 − 2+Bκξ2 κ(2+Bκξ2)2 equationsforthetriangularlattice(nearestneighborvortices
v (cid:0) v (cid:1)
are in the plane of the graph). The distance is measured in
unitsof theintervortex spacing d in the triangular lattice.
The field dependence of the variational parameters is
The upper critical field is determined as the field at
calculated numerically by minimization of the total free
which the order parameter in the superconductor be-
energy F(B,κ,ξv,f∞) with respect to ξv and f∞. This comes equal to zero. As can be seen from Eq. (23),
dependence can be approximated by the following ana-
one has f = 0 at H = 0.985κ. Thus, the difference be-
lytical expressions with an accuracy of about 0.5%:
tweentheexactHc2 =κanditscalculatedvalueisabout
4
1.5%. This result is quite natural, since variational pro-
ceduresingeneralgiveonlyapproximatesolutionstothe
GL equations. Similarly, the Clem value of Hc1 slightly (cid:20)(cid:15)(cid:19)(cid:19)
differs from the numerically calculated one [18,19].
Now we compare the obtained results for the order
(cid:19)(cid:15)(cid:26)(cid:24)
parameter and the magnetic field with the similar dis-
tributions computed for the triangularlattice by the nu-
merical method proposed in Ref. [16]. The dependence (cid:19)(cid:15)(cid:24)(cid:19)
w
ofthe orderparameteronthe distance fromthe cellcen-
ter is showninFig. 1 fordifferent values ofthe magnetic
(cid:19)(cid:15)(cid:21)(cid:24)
induction at κ = 10. The spatial distribution of the or-
der parameter in the triangular lattice along the nearest
(cid:19)(cid:15)(cid:19)(cid:19)
neighbor direction is also shown in Fig. 1. The results (cid:19) (cid:21)(cid:24) (cid:24)(cid:19) (cid:26)(cid:24) (cid:20)(cid:19)(cid:19)
of our variational calculations are close to the numerical +
onesatanyvaluesofthemagneticinduction. The differ- FIG. 3. Averaged order parameter squared ω vs B in di-
ence does not exceed several percent. Such an accuracy mensionless units for κ= 100. The solid line corresponds to
of our approach remains in a wide range of κ > 1. The thevariational calculations and thedotline tothenumerical
magnetic flux density as a function of the dista∼nce from solution.
the vortexaxisisshowninFig.2atκ=10,B =0.5Hc2.
Thereisgoodagreementbetweenthevariationalandnu- The anisotropyinthe GL functionalcanbe takeninto
merical results for h. In Fig. 3, the spatial average of account by introducing the phenomenological effective-
the order parameter squared, ω = ψ 2 , is plotted as mass tensor mj (j = 1,2,3), where mj are the effec-
a function of the magnetic inductiohn| a|tiκ = 100. The tive masses in the direction of the principal axes xj. It
comparisonwith the numericalresult shows that the de- was shown in Ref. [33] that the Ginzburg-Landau equa-
viation of variationally-calculated ω from the exact de- tions can be transformed to isotropic form by a simple
pendence does not exceed one percent in small and in- transformation if κ is replaced by κ∗ = κp1−1/2, where
termediate fields. Near Hc2 this deviation increases due p1 =m1/√m1m2m3,(thevorticesaredirectedalongthe
to a small difference between the calculated Hc2 and the x1 axis). Thus, in this case, the Wigner-Seitz approxi-
exact value. mation can be used. The Wigner-Seitz cell has an ellip-
tic shape and can be transformed to the circular cell by
a scaling transformation: the distance along the x axis
j
shouldbenormalizedbym1m2m3/√m1mj. Themagne-
tization of the anisotropic superconductor can be found
from the magnetization of the isotropic one by replacing
(cid:19)(cid:15)(cid:19)(cid:19)(cid:24)(cid:19)
∗
the GL parameter κ by κ .
(cid:19)(cid:15)(cid:19)(cid:19)(cid:21)(cid:24)
% IV. MAGNETIZATION
(cid:16)
(cid:12)
U
(cid:11)
K
(cid:19)(cid:15)(cid:19)(cid:19)(cid:19)(cid:19)
The magnetization is defined by the well-known rela-
tionship:
(cid:16)(cid:19)(cid:15)(cid:19)(cid:19)(cid:21)(cid:24)
(cid:19)(cid:15)(cid:19) (cid:19)(cid:15)(cid:23) (cid:19)(cid:15)(cid:27) 4πM =H B. (25)
U(cid:18)G − −
FIG.2. Thespatial distributionof thedimensionless mag- It can be calculated by two equivalent methods if the
netic field in the unit cell of flux-line lattice at the magnetic exact solution to the GL equations is known. First, the
induction B =0.5Hc2 for κ=10. The solid line corresponds magnetic field H may be calculated by minimization of
tothevariationalcalculationsintheWigner-Seitzapproxima- the Gibbs free energy G=F 2BH:
tion. The dashed line corresponds tothe numerical solution. −
1∂F
H = (26)
Thus,theresultsofourvariationalapproachagreewell 2∂B
with the exact numerical solution to the GL equations.
This derivative was calculated, e.g. in Refs. [8,23]. The
In the the next section we shall apply this approach to
second approach uses the virial theorem for the flux-line
the calculation of the magnetization curve.
lattice,whichwasproveninRef.[34],namelytheapplied
magnetic field H can be found from the local magnetic
field h and the order parameter f as
5
1
H = f2 f4+2h2 d2r, (27) magnetization (28) is shown in Fig. 4 for κ = 100. The
2BAcell Z − magnetization curves corresponding to the London and
(cid:0) (cid:1) Abrikosovapproximationsarealsoplotted. Atlowfields,
where the integral is taken over the area of the unit cell.
themagnetizationpracticallycoincideswiththeresultsof
Bothmethods areequivalentifthe exactsolutionsf and
the London approach, the difference between them does
h are used [34].
not exceed 0.5%. At H Hc2, the behavior of the mag-
The variational model gives the spatial distributions ∼
netizationisingoodagreementwiththeAbrikosovhigh-
ofthe orderparameterandthe magnetic fieldwithin the
fieldresult(14). AlthoughnearHc2 bothcurvesareclose
cell, which are close to the exact results. Let us find the
to each other, the error of our approximation is not so
magnetization by means of both methods. According to
small as in low fields because of the slight deviation of
Eq. (26) the magnetization M is:
the calculated Hc2 from the exact value. For example,
f∞ at H =0.8Hc2, κ=100 the error of our variational pro-
4πM = B+ cedure is about 5%. In order to improve the accuracy
− − κξ ×
v near Hc2, one may put the constant t in Eq. (23) to be
K0(f∞ξv)I1(f∞ρ)+I0(f∞ξv)K1(f∞ρ) equal to 1. As a result, the difference between the varia-
+
×K1(f∞ξv)I1(f∞ρ) I1(f∞ξv)K1(f∞ρ) tional and the Abrikosov M(H) curves decreases in the
−
1 vicinity of Hc2, whereas at low and intermediate fields
+2Bκ2ξ2{K1(f∞ξv)I1(f∞ρ)− the magnetization does not change.
v
In the inset of Fig. 4, we compare the calculated de-
f2 2+3Bκξ2
−I1(f∞ξv)K1(f∞ρ)}−2+ 2∞κ((cid:0)2+Bκξv2)v3(cid:1) + ptheendCelnecme-H−a4oπMap(pHro)awchithoftshuaptefropuonsditiionnRoeff.vo[3r0te],xwfiheeldres
κ2f2 ξ2 1 f2 2 [23]wasused. Itisclearlyseenthatthisapproachisvalid
+ 2∞ v{ −2 ∞ ln Bκξ2 +1 + only at small fields, and its use leads to an appreciable
(cid:20) v (cid:21) error in the magnetization even in the intermediate field
f2 2+3Bκξ2 f2 2+3Bκξ2
+ ∞ v ∞ v . (28) range. An additional error in M arises in small fields
2κ(2+Bκξ2)3 − 2κ(2+Bκξ2)3 } due to the approximate replacement of the lattice sums
(cid:0) v (cid:1) (cid:0) v (cid:1)
by integrals [23] in calculating the magnetic free energy;
Thus, in the former case the dependence of the magneti- for more details see Ref. [30].
zationonthe magneticinductionB isgivenbyEqs.(28)
and (22)-(24). The dependence of H on B is given by
Eq. (25). Thus, we find the implicit function M(H).
Within the second approach, the integral (27) can be
calculatedonlynumericallywhen f andh aredefined by
Eqs. (1) and (17). Our calculations show that not only
(cid:19)(cid:15)(cid:19)(cid:20)(cid:19)
thevaluesofH foundbybothmethodscoincide,butalso (cid:19)(cid:15)(cid:19)(cid:20)
the values of M, which is usually much smaller than H, 0
arepractically indistinguishable. The difference between p(cid:23)
(cid:16)
themis muchlessthanonepercentatanyinductionand 0
(cid:19)(cid:15)(cid:19)(cid:19)(cid:24)
κ > 1/√2. Note that the result for M in the second p(cid:23) (cid:19)(cid:15)(cid:19)(cid:19)
(cid:16) (cid:19) (cid:23)(cid:19) + (cid:27)(cid:19) (cid:20)(cid:21)(cid:19)
approach is very sensitive to the perturbations of f(r)
and h(r). For example, if one puts f∞ =1 near Hc1 and
minimizes the free energy only with respect to ξ this
v (cid:19)(cid:15)(cid:19)(cid:19)(cid:19)
does not change f(r) and h(r) considerably. However, (cid:19) (cid:21)(cid:19) (cid:23)(cid:19) (cid:25)(cid:19) (cid:27)(cid:19) (cid:20)(cid:19)(cid:19)
this procedure would lead to an appreciable change of +
M(H) when using Eq. (27), while according to Eq. (26)
FIG. 4. Calculated 4πM(H) using the variational
the magnetization practically remains the same. Below, method for κ = 100 (sol−id line). Also shown are the Lon-
we shall use Eq. (28) for the magnetization. don dependence 4πM(H) (13) with the exact value of Hc1
At low fields the variational parameters (22) and (23) (3)(dashedline)−andtheAbrikosovhigh-fieldresult(16)(dot
may be considered as constants independent of B when line). The inset refers to themagnetization found within the
κ 1. In this case, Eq. (28) can be expressed as a framework ofthevariationalmodel(solidline)andtheresult
≫
power series in terms of ξ . As a result, it is possible to of calculations [30] for the magnetization in the Clem-Hao
v
obtain Eq. (12) with Hc1 given by Eq. (3) at ε = 0.52. approximation [23] (dot line). The difference between these
Thus, in small fields the model reduces to the London curvesarises duetothesuperposition of vortexfieldsusedin
approximationprovidedthatthevariationally-calculated theClem-Hao model [23]. Dimensionless variables are used.
value of Hc1 is used, which is practically indistinguish-
able from the exact Hc1. Actually, the use of the exact Nowwecomparethemagnetizationfoundintheframe-
Hc1 in small fields is equivalent to taking into account work of the variational and numerical methods. At low
the effect of vortex cores. The field dependence of the and high fields the exact dependence M(H) coincides
6
with the results of the London (at κ 1) and the analyzedandreducedtothedimensionlessformbasedon
≫
Abrikosov approximations, respectively. As we found the Clem-Hao formulas with non-selfconsistent field de-
above, the results of the variational approach are in pendencesofthevariationalparameters[23]. Theresult-
agreement with these approximations. In the interme- ing magnetization curve is close to the Abrikosov high-
diate field range, where the London and the Abrikosov field result and to our variational dependence in the in-
approachesarenotapplicable,thedifferencebetweenthe termediate field range. The κ values obtained in Refs.
values of the magnetization calculated by numerical and [28,29] were κ = 70 for YBa2Cu4O8 and κ = 80 for
variationalmethodsisnotbiggerthan1%inawiderange Nd1.85Ce0.15CuO4−δ.
of κ 1 values. Thus, our results for the magnetization
≫
appear to be a good approximation to the exact numer-
ical solution of the GL equations at κ 1.
≫
Next we discuss the case of small κ values. In Fig. 5,
thefielddependencesofthemagnetizationareplottedfor
several small κ values. The solid and dotted lines corre- (cid:19)(cid:15)(cid:19)(cid:21)
(cid:19)(cid:15)(cid:19)(cid:20)(cid:19)
spond to the variational and numerical calculations, re-
0
spectively. The agreement between these results is good p(cid:23)(cid:19)(cid:15)(cid:19)(cid:20)
at κ > 1. At smaller κ values the variational and the (cid:16)
0
exact∼numericalresultsdiffernearthelowercriticalfield.
Inthiscase,theintervortexdistanceisoftheorderofthe p(cid:16)(cid:23)(cid:19)(cid:15)(cid:19)(cid:19)(cid:24) (cid:19)(cid:15)(cid:19)(cid:19) (cid:16)(cid:23) OQ(cid:19)(cid:3)+ (cid:23)
coherencelengthalmostintheentirefieldrange,andthe
variational approach based on an appropriate trial func-
tionforthe orderparameterinthe circularWigner-Seitz
(cid:19)(cid:15)(cid:19)(cid:19)(cid:19)
cell may lead to some deviation from the exact solution. (cid:19) (cid:21)(cid:19) (cid:23)(cid:19) (cid:25)(cid:19) (cid:27)(cid:19)
+
IntheinsetofFig.5wecomparethemagnetizationcal-
culatedbymeansofourvariationalprocedureandbythe FIG. 6. The field dependence of magnetization in dimen-
Clem-Hao model at κ = 1. It is seen that also for small sionlessunits. Thesolidanddottedlinesshowthetheoretical
(even as for large κ, see above) the calculation method variational dependences at κ = 70 and κ = 80. The circles
proposed in Ref. [23] leads to an inaccurate magnetiza- andtrianglesgivetheexperimentaldataforYBa2Cu4O8 [28]
tion. and Nd1.85Ce0.15CuO4−δ [29], respectively. The inset shows
the magnetization calculated by our variational model and
plottedversusthelogarithmoftheappliedfieldHatκ=100.
(cid:20)(cid:15)(cid:19) In Fig. 6 we compare these experimental curves with
(cid:19)(cid:15)(cid:25) our variationalresult. The circles and triangles in Fig. 6
k correspondto YBa2Cu4O8 and Nd1.85Ce0.15CuO4−δ, re-
(cid:19)(cid:15)(cid:22) (cid:3) (cid:3)(cid:19)(cid:17)(cid:26)(cid:24) spectively. The solid and dotted lines correspond to the
(cid:21)
+F (cid:19)(cid:17)(cid:27)(cid:24) theoretical curves at κ=70 and κ=80, respectively. It
0(cid:3)(cid:18)(cid:3)(cid:19)(cid:15)(cid:24) (cid:19)(cid:15)(cid:19)(cid:19)(cid:15)(cid:25) (cid:19)(cid:15)(cid:27) (cid:20)(cid:15)(cid:19) (cid:20)(cid:15)(cid:21) (cid:20) iosryclaenardlyexspeeenrimtheantt.gIonodthaegirneseemteonftFeigx.is6tsthbeettwheeoernettihcea-l
p(cid:23)
(cid:16) (cid:20)(cid:17)(cid:21) (variational) magnetization is plotted as a function of
lnH at κ = 100. This dependence is nearly linear in
(cid:20)(cid:17)(cid:24)
a wide range of intermediate fields. A similar behavior
(cid:21)
of the magnetization was observed in numerous experi-
(cid:19)(cid:15)(cid:19)
(cid:19)(cid:15)(cid:19) (cid:19)(cid:15)(cid:24) (cid:20)(cid:15)(cid:19) ments; see, for example, Ref. [20].
+(cid:3)(cid:18)(cid:3)+
F(cid:21)
FIG. 5. Magnetization curves calculated using the vari-
ational approach (solid lines) and exact numerical method V. CONCLUSIONS
(dotted lines) at various small κ values. The inset compares
the magnetization obtained by our variational model (solid We proposed an approximate method to solve the
line) with the Clem-Hao approximation [23] (dashed line) at Ginzburg-Landau equations for the regular flux-line lat-
κ=1. Dimensionless variables are used.
tice at any values of the magnetic induction and the
Ginzburg-Landau parameter κ > 1/√2. The Wigner-
Our formulas for the magnetization may be used for
Seitz approximationis used, and the hexagonal unit cell
the analysis of experimental data. In Fig. 6 the calcu-
ofthe vortexlattice is replacedbya circlewith the same
lated magnetization curves are compared with the mea-
area. Our model is based on Clem’s trial function for
suredmagnetizationofYBa2Cu4O8polycrystals[28]and the order parameter. The use of this function allows us
Nd1.85Ce0.15CuO4−δ singlecrystals[29]. Inthesepapers, to find the magnetic flux density self-consistently from
the magnetization curves at different temperatures were
the second Ginzburg-Landau equation. The comparison
7
between the variational results and the results of exact [11] J.Rammer,W.Pesch,andL.Kramer,Z.Phys.B68,49
numericalsolutionoftheGinzburg-Landauequationsre- (1987).
veals good accuracy of our approach: the difference be- [12] J. Rammer, J. Low Temp.Phys. 71, 323 (1988).
tween the spatial distributions of the order parameter [13] A. I. Larkin and Yu. N. Ovchinnikov, Phys. Rev. B 51,
and the magnetic field does not exceed several percent. 5965 (1995).
Such accuracy remains in a wide range of values of the [14] A. S. Mel’nikov, Zh. Eksp. Teor. Fiz. 101, 1979 (1992)
magneticinductionandκ 1. Themethodisappliedto [JETP 74 (6), (1992)].
the calculation of the field≫dependence of the reversible [15] E. H.Brandt, Phys. stat. solidi (b) 51, 345 (1972).
[16] E. H.Brandt, Phys. Rev.Lett. 78, 2208 (1997).
magnetization. An analytical expression for the mag-
[17] J. R.Clem, J. Low Temp. Phys. 18, 427 (1975).
netization is proposed. At low fields, the obtained de-
[18] C.-R. Hu,Phys. Rev.B 6, 1756 (1972).
pendence agrees with the predictions of London theory
at κ > 1. At high fields, it is in good agreement with [19] E. A.Shapoval, JETP Letters 69, 577 (1999).
[20] V. G. Kogan, M. M. Fang, and S. Mitra, Phys. Rev. B
the A∼brikosov result. It is shown that the values of the
38, 11958 (1988).
magnetization calculated within the framework of our
[21] A. E. Koshelev, Phys. Rev.B 50, 506 (1994).
variational model and of the numerical method of solu-
[22] A. S. Mel’nikov, I. M. Nefedov, D. A. Ryzhov, I. A.
tion of the Ginzburg-Landau equations are practically Shereshefskii, and P. P. Vysheslavtsev, Phys. Rev. B 62
indistinguishable, especially in small and intermediate
(2000).
magnetic fields at κ > 1. Our model yields the lim-
[23] Z. Hao, J. R. Clem, M. W. McElfresh, L. Civale, A. P.
its of the Clem-Hao m∼odel for the magnetization. The Malozemoff, and F. Holtzberg, Phys. Rev. B 43, 2844
presentedanalyticalformulas for the magnetizationmay (1991).
be used to analyze experimental data. As an illustra- [24] J. R. Clem, in Superconducting Electronics, edited by
tion we compared the experimental and calculated mag- H.Weinstock,andM.Nisenoff(Springer-Verlag,Berlin,
netization curves for different high-Tc superconductors 1989), p.1.
(YBa2Cu4O8 and Nd1.85Ce0.15CuO4−δ) and found good [25] Z.HaoandJ.R.Clem,Phys.Rev.Lett67,2371(1991).
agreement between theory and experiment. [26] Z. Hao and J. R.Clem, Phys. Rev.B 43, 7622 (1991).
[27] J.H.GohngandD.K.Finnemore,Phys.Rev.B46,398
(1992); Q. Li et al., Phys. Rev. B 46, 3195 (1992); D.
ACKNOWLEDGMENTS N. Zheng et al., Phys. Rev. B 48, 6519 (1993); V. C.
Kim et al., Phys. Rev. B 51, 11767 (1995); P. Pugnat
et al., Europhys. Lett. 29, 425 (1995); M.-S. Kim et al.,
TheauthorsacknowledgeusefuldiscussionswithL.G.
Phys. Rev. B 51, 3261 (1995); M. Xu et al., Phys. Rev.
Mamsurova, K. S. Pigalskiy, and N. G. Trusevich. This
B 53, 15313 (1996); Y. Zhuo et al., Phys. Rev. B 55,
work is supported by the Russian Foundation for Ba-
12719 (1997); M.-S. Kim et al., Phys. Rev. B 57, 6121
sic Research(RFBR), grants#00-02-18032and #00-15- (1998); M.Y. Cheon et al., Physica C 302, 215 (1998);
96570, by the joint INTAS-RFBR program, grant #IR- S.-I. Leeet al.,Physica C 341, 379 (2000); H.-J. Kim et
97-1394,and by the Russian State Program’Fundamen- al.,PhysicaC341,745(2000);A.Aburtoetal.,Physica
tal Problems in Condensed Matter Physics’. C 303,185 (1998); H.-J.Kimet al.,PhysicaC331, 292
(2000).
[28] W. Chen et al.,Phys. Rev.B 51, 6035 (1995).
[29] A. Nugroho et al.,Phys. Rev.B 60, 15384 (1999).
[30] W. V. Pogosov, A. L. Rakhmanov, and K. I. Kugel, Zh.
Eksp.Teor.Fiz.118, 908(2000) [JETP 91, 588(2000)].
[31] M. Abramowitz and F. A. Stegun, Handbook of Mathe-
[1] A. A. Abrikosov, Zh. Eksp. Teor. Fiz. 32, 1442 (1957) matical Functions (Dover,New York,1965).
[Sov.Phys. JETP 5, 1174 (1957)]. [32] W.M.Kleiner,L.M.Roth,andS.H.Autler,Phys.Rev.
[2] E.H.Brandt and U.Essman, Phys.stat.solidi (b)144, A 133, 1226 (1964).
13 (1987). [33] V. G. Kogan and J. R. Clem, Phys. Rev. B 24, 2497
[3] E. H.Brandt, Rep. Prog. Phys. 58, 1465 (1995). (1981).
[4] P.G.deGennes, Superconductivity of Metals and Alloys [34] M. M. Doria, J. E. Gubernatis, and D. Rainer, Phys.
(Benjamin, NewYork,1966). Rev. B 39, 9573 (1989).
[5] D. Saint-James, G. Sarma, and E. Thomas, Type-II Su-
perconductors (Oxford, Univ.Press, NewYork, 1969).
[6] I.G.deOliveiraandA.M.Thompson,Phys.Rev.B57,
7477 (1998).
[7] D.Ihle, Phys. stat. solidi (b) 47, 423 (1971).
[8] A.L. Fetter, Phys.Rev. 147, 153 (1966).
[9] D.Ihle, Phys. stat. solidi (b) 47, 429 (1971).
[10] R. J. Watts-Tobin, L. Kramer, and W. Pesch, J. Low
Temp.Phys. 17, 71 (1974).
8
