- Email: [email protected]
Monte Carlo simulation of three-dimensional suspensions of superconducting granules
0038-1098/82/300487-03 $03.00/0 Pergamon Press Ltd.
Solid State Communicauons, Vol. 43, No. 6, pp. 487-489, 1982. Pnnted in Great Bntam.
MONTE CARLO SIMULATION OF THREE-DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES M. Hillen and D. Stauffer Institut fur Theoretische Physik, Universltat, 5000 Koln 41, West Germany
(Received 12 February 1982 by B. Miihlschlegel) We simulate on a computer the change from the superconducting to the normal state for a monodisperse suspension of granules, distributed randomly in a plate, as a function of magnetic field. In contrast to earlier monolayer approximations, our three-dimensional calculation gives an s-shaped curve, as a function of field, for the fraction of spheres which are still superconducting, as was observed experimentally. IF A LARGE NUMBER of superconducting spheres, suspended rigidly m a wax, is put into a magnetic field, they do not all become normal conducting at the same field. Instead, their diamagnetic interaction due to the Meissner effect produces a local field which differs for different spheres in the same external (applied) field He. Thus with increasing field He, first that sphere with the largest local field at its surface becomes normal, then that with the largest of all remaining spheres, etc. until all spheres have become normal in a sufficiently large applied field He. Of course, in this process the change of each sphere to the normal state changes the local field felt by all other spheres. Experimental studies [1,2] of this gradual transition of superconducting granules show an s-shaped curve if the fraction x, of spheres which are still superconducting is plotted vs the field H0. Prehmmary numerical studies [3] used only a monolayer of spheres and could not reproduce fully this s-shape. The present work uses more time on a better computer, the CDC Cyber 76, m order to take rote account fully the threedunensmnal nature of the expenment [1, 2]. Still, had we reed to simulate as many spheres on our computer as are present m the experimental sample, we would have needed more than ten years of computer t~me. Thus we have to extrapolate carefully from smaller systems to the desared system size. As m the monolayer work [3], we distribute sphers of umt radms randomly m our sample of size L x L × W, w~th a fillmg factor (volume fraction) of ten percent. If a sphere would overlap with another sphere, we Fmd another position for the new sphere by selecting a different set of three random numbers for its coordinates. (In the monolayer approximation, no significant changes were found if instead one distributes the spheres as in a hard sphere gas, as stmulated by Monte Carlo methods.) Every sphere m an external field HI adds to the total 487
field the contribution H(r) = ½r-a(l -- 3 cos20)H1
(1)
in z-direction, parallel to the field. Here r is the distance from the sphere center and 0 the angle of the distance vector with the applied field. (We neglect the field components perpendicular to the z-ax~s since they cancel in the average.) Summing up the field produced by all spheres, including the self-induced field H(r) = ~Ho and the external field He we obtain the field Hk at the center of an arbitrary sphere k in our system:
Hk = Ho+½Ho+ ~. ~r~(1--cos201h)Ht.
(2a)
i
In principle one would have to solve this set of equations, for example by iteration. In practice, we are interested m a rather small filling factor 0.1, where the sum in equation (2a) gives a contribution smaller by about one order of magnitude than the leading term 1.5//0 (see below). Therefore we neglect higher-order effects and replace Hi by the applied field H0 to get an exphclt expression:
Hk/I'lo = 1.5 + ~
~r~a(l--cos2Otk)
(2b)
l
which we use for all calculations (also those in [3]). For the same reason of using a rather small filling factor, we neglect in equation (2) the spatial variation of the field produced by one sphere across the volume of the other sphere; instead we worked with the field at the center posatlon of the given sphere k when calculating the influence of sphere i; only for the self-induced fie.ld i = k was the sphere geometry taken into account, and the field defined as that on the equatorial plane of the sphere [cos 0 = 0 in equation (1)]. Our spheres is assumed to be much larger in radius
488
THREE-DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES
than the London penetration depth. Each sphere is assumed to switch from the superconducting to the normal state as soon as the local field, equation (2b), on its equatorial plane becomes larger than a matenal dependent threshold value ~//,h. (For an isolated sphere this transitmn is reached when the applied field He reaches the superheating field H,h, which means that the local field on the equatorial plane Is ~H,h. In other words, due to demagnettzation the superheating field H,h is smaller by a factor 2/3 for spheres than it is for bulk matenal.) Spheres which are no longer supercon. ducting produce no field obeying the above formula and are omitted from the calculation in equation (2b). Figure 1 shows our numerical results m a L x L x W plate, together with the monolayer approximation of [3]. The calculation of the phase transitmn m the largest system, 259 x 259 x 15, containing more than 24,000 spheres, took 8052 sec computer time; in general the time increased as the square of the number of spheres involved. Plotting the characteristic points of the curve for each system vs the reciprocal length L at constant width W = 15, we obtain the extrapolated curve of Fig. 1 as our main result for the experimental system slze [1-5], 4050 x 4050 x 15. We see that, in contrast to the monolayer approximation [3], the curve is now rather symmetric. (In our calculation we keep the external field He fixed and look at the largest field observed at any of the sphere centers. That sphere is then onutted, the fields calculated agam, and the process repeated again and again. Thus theoretically one decreases H=hat fixed He. Experimentally, one keeps Hsh fixed and vanes H0. Of course, the results depend only on the ratmHsh/Ho,and the true approaches are thus equivalent apart from a vanable transformaUon.) Figure 2 compares our extrapolated curve and the monolayer curve voth expernnent. We normahze the fields such that for a fractmn x a = 30% of superconducting spheres (which corresponds t o r t e = Hmh m our calculation, i e. to the field where isolated spheres swRch) the two curves cross each other. This means that the superheating field in the expenment was 235 G (Hueber [2] gives 250G). The monolayer curves goes steeply to zero at He -- Hsh since the few remaining superconducting spheres for small concentrations x, are all Isolated. Instead, the experimental curve is smoother like an s-shaped function. And indeed our three-dimensmnal result reproduces this shape: For small x, m three dLmenstons, the few remaining spheres turned out to be grouped together into pasrs (few triplets) such that the d=stance vector between the two spheres is roughly parallel to the field. Then each sphere is sh]elding away some of the field for its neighbor (quantitatively, cos 0 _~ 1 m equation (1), Le. a negatwe contribution to
Vol. 4 3 , N o . 6
Xs
90 8O
50
30 2O
o 18
17
15
16
14 13 H,,mlHo
Fig. 1. Vanation of fraction x s of superconducting spheres w i t h s u p e r h e a t i n g f i e l d H . h . H m = = = ~ H . h is
the largest field observed at any superconducting sphere, and He ts the apphed field external field. At Hmex = 1.5He, isolated spheres switch. 100
"=
i
=
,
i
=
X$
90 S0 7O 6O 5O 40
--,°.°,
....
t~s
11
theory
. . . . . ma'~olc~er theory 30 20 10 I
0
025
I
I
05
075
10 125 Ho/Hsh
Fig. 2. Comparison of experimental curve with theoretical curves. the field. Therefore these pairs remain superconductmg even for fields larger than H,h, up to He - 1.09H, h. The average field ( H ) at the sphere centers was found to be 1.59tlo for large plates of thickness W = 15. We found this limiting value by first looking at systems of various lengths L at fixed ratio W/L and
Vol. 43, No. 6
THREE.DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES and
250
= -- 1.6(x. = 0.1)
489 (4)
for this 200 x 200 x 15 sample. (The width ( H 2) ( H ) 2 was extrapolated analogously to equatton (3) as 0.0019 + 0.0002 for the infinite sample andx, = 1.) Also the excess determined by the fourth moment changed very drastically.
2(30
150
( ( H - - ( H ) ) 4 ) ( ( H - - ( H ) ) 2 ) - 2 - - 3 = 0.18(x, = 1) 100
and
50
140
145
150
155
160
165
170
175 H/Ha
Fig. 3. Histogram for magnetic field H felt by 14,231 superconductmg spheres in 200 x 200 x 15 system for small external fields Ho. 100 8O 6O 4O 20
.... :.........~lllh
0 128 1'30 132
13l* 136
138
140
1/.2
11.4 146
HIHo
Fzg. 4. Histogram of the same system as in Fig. 3; but due to a strong field only ten percent of the spheres are still superconducting. extrapolating to L --, ae; then this extrapolated value was found to vary smoothly with W/L, roughly as (H)/Ho = 1.593 -- 0.210 W/L + O.141(W/L) 2.
(3)
Figure 3 shows the histogram of the fields observed by the spheres in the case where all the spheres are still superconducting (x, = 1;H (~ H,h), whereas Fig. 4 gives the histogram for the same system when only 10% of the spheres are still superconducting (x, = l ; H = 1.034H, h). We see that the fun distribution is roughly Gaussian whereas a strong asymmetry appears when most spheres became normal. Quantitatively, we found for the skewness of these two distributions
((H--(H))S)/((H--H)2) s/2
-- 0.15 (x,
=
l)
(5)
for the same sample size. (The absolute value of both excess and skewness increased with increasing L for both x, = O.1 and x, = 1. More details are gwen m [4].) The physical reason for this asymmetry ts obvious: The transition from the superconducting to the normal state hits mainly the spheres in large fields and thus gives a rather sharp cut-off as observed in Fig. 4. Szmilar to the monolayer calculation [3], the remaimng superconducting spheres avoid to be close to each other but showed no tendency to order periodically. Further details, mcludmg the computer program, are given elsewhere [4]. Future research along these lines could revolve higher concentrations instead of a filling factor of 0.1 only; different sphere radfi instead of unit radms for all; mhomogeneous distnbution of spheres instead of our pseudo-random procedure; and other improvements. Some of these effects could broaden the translUon regton towards the larger width observed experimentally (Fig. 2). In such future calculations, computer time could be saved by working wtth smaller systems than employed here. As our numerical evidence m Fig. 1 shows, already small systems like 100 x 100 x 15 could approximate well the desired 4050 x 4050 x 15 limit, and use 45 times less computer time than our largest system, 259 x 259 x 15. Acknowledgements - We thank G. Waysand for suggestmg this work, V. Wmkelmann form improvements m our computer program, and G. Deutscher, O. EntinWohlmann, Y. Kantor and A. Kapitulnik for discussions. REFERENCES I. 2. 3. 4.
=
= 3.83(x, = 0.1)
D. Hueber, C. Valette & G. Waysand, J. Phys. 41, L161 (1980). D. Hueber, Th~se, Universit~ de Paris Sud (1980). C. Valette, G. Waysand & D. Stauffer, Solid State Commun. 41,305 (1982). M. Hillen, Staatsexamensarbeit, Cologne University (1982).
Solid State Communicauons, Vol. 43, No. 6, pp. 487-489, 1982. Pnnted in Great Bntam.
MONTE CARLO SIMULATION OF THREE-DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES M. Hillen and D. Stauffer Institut fur Theoretische Physik, Universltat, 5000 Koln 41, West Germany
(Received 12 February 1982 by B. Miihlschlegel) We simulate on a computer the change from the superconducting to the normal state for a monodisperse suspension of granules, distributed randomly in a plate, as a function of magnetic field. In contrast to earlier monolayer approximations, our three-dimensional calculation gives an s-shaped curve, as a function of field, for the fraction of spheres which are still superconducting, as was observed experimentally. IF A LARGE NUMBER of superconducting spheres, suspended rigidly m a wax, is put into a magnetic field, they do not all become normal conducting at the same field. Instead, their diamagnetic interaction due to the Meissner effect produces a local field which differs for different spheres in the same external (applied) field He. Thus with increasing field He, first that sphere with the largest local field at its surface becomes normal, then that with the largest of all remaining spheres, etc. until all spheres have become normal in a sufficiently large applied field He. Of course, in this process the change of each sphere to the normal state changes the local field felt by all other spheres. Experimental studies [1,2] of this gradual transition of superconducting granules show an s-shaped curve if the fraction x, of spheres which are still superconducting is plotted vs the field H0. Prehmmary numerical studies [3] used only a monolayer of spheres and could not reproduce fully this s-shape. The present work uses more time on a better computer, the CDC Cyber 76, m order to take rote account fully the threedunensmnal nature of the expenment [1, 2]. Still, had we reed to simulate as many spheres on our computer as are present m the experimental sample, we would have needed more than ten years of computer t~me. Thus we have to extrapolate carefully from smaller systems to the desared system size. As m the monolayer work [3], we distribute sphers of umt radms randomly m our sample of size L x L × W, w~th a fillmg factor (volume fraction) of ten percent. If a sphere would overlap with another sphere, we Fmd another position for the new sphere by selecting a different set of three random numbers for its coordinates. (In the monolayer approximation, no significant changes were found if instead one distributes the spheres as in a hard sphere gas, as stmulated by Monte Carlo methods.) Every sphere m an external field HI adds to the total 487
field the contribution H(r) = ½r-a(l -- 3 cos20)H1
(1)
in z-direction, parallel to the field. Here r is the distance from the sphere center and 0 the angle of the distance vector with the applied field. (We neglect the field components perpendicular to the z-ax~s since they cancel in the average.) Summing up the field produced by all spheres, including the self-induced field H(r) = ~Ho and the external field He we obtain the field Hk at the center of an arbitrary sphere k in our system:
Hk = Ho+½Ho+ ~. ~r~(1--cos201h)Ht.
(2a)
i
In principle one would have to solve this set of equations, for example by iteration. In practice, we are interested m a rather small filling factor 0.1, where the sum in equation (2a) gives a contribution smaller by about one order of magnitude than the leading term 1.5//0 (see below). Therefore we neglect higher-order effects and replace Hi by the applied field H0 to get an exphclt expression:
Hk/I'lo = 1.5 + ~
~r~a(l--cos2Otk)
(2b)
l
which we use for all calculations (also those in [3]). For the same reason of using a rather small filling factor, we neglect in equation (2) the spatial variation of the field produced by one sphere across the volume of the other sphere; instead we worked with the field at the center posatlon of the given sphere k when calculating the influence of sphere i; only for the self-induced fie.ld i = k was the sphere geometry taken into account, and the field defined as that on the equatorial plane of the sphere [cos 0 = 0 in equation (1)]. Our spheres is assumed to be much larger in radius
488
THREE-DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES
than the London penetration depth. Each sphere is assumed to switch from the superconducting to the normal state as soon as the local field, equation (2b), on its equatorial plane becomes larger than a matenal dependent threshold value ~//,h. (For an isolated sphere this transitmn is reached when the applied field He reaches the superheating field H,h, which means that the local field on the equatorial plane Is ~H,h. In other words, due to demagnettzation the superheating field H,h is smaller by a factor 2/3 for spheres than it is for bulk matenal.) Spheres which are no longer supercon. ducting produce no field obeying the above formula and are omitted from the calculation in equation (2b). Figure 1 shows our numerical results m a L x L x W plate, together with the monolayer approximation of [3]. The calculation of the phase transitmn m the largest system, 259 x 259 x 15, containing more than 24,000 spheres, took 8052 sec computer time; in general the time increased as the square of the number of spheres involved. Plotting the characteristic points of the curve for each system vs the reciprocal length L at constant width W = 15, we obtain the extrapolated curve of Fig. 1 as our main result for the experimental system slze [1-5], 4050 x 4050 x 15. We see that, in contrast to the monolayer approximation [3], the curve is now rather symmetric. (In our calculation we keep the external field He fixed and look at the largest field observed at any of the sphere centers. That sphere is then onutted, the fields calculated agam, and the process repeated again and again. Thus theoretically one decreases H=hat fixed He. Experimentally, one keeps Hsh fixed and vanes H0. Of course, the results depend only on the ratmHsh/Ho,and the true approaches are thus equivalent apart from a vanable transformaUon.) Figure 2 compares our extrapolated curve and the monolayer curve voth expernnent. We normahze the fields such that for a fractmn x a = 30% of superconducting spheres (which corresponds t o r t e = Hmh m our calculation, i e. to the field where isolated spheres swRch) the two curves cross each other. This means that the superheating field in the expenment was 235 G (Hueber [2] gives 250G). The monolayer curves goes steeply to zero at He -- Hsh since the few remaining superconducting spheres for small concentrations x, are all Isolated. Instead, the experimental curve is smoother like an s-shaped function. And indeed our three-dimensmnal result reproduces this shape: For small x, m three dLmenstons, the few remaining spheres turned out to be grouped together into pasrs (few triplets) such that the d=stance vector between the two spheres is roughly parallel to the field. Then each sphere is sh]elding away some of the field for its neighbor (quantitatively, cos 0 _~ 1 m equation (1), Le. a negatwe contribution to
Vol. 4 3 , N o . 6
Xs
90 8O
50
30 2O
o 18
17
15
16
14 13 H,,mlHo
Fig. 1. Vanation of fraction x s of superconducting spheres w i t h s u p e r h e a t i n g f i e l d H . h . H m = = = ~ H . h is
the largest field observed at any superconducting sphere, and He ts the apphed field external field. At Hmex = 1.5He, isolated spheres switch. 100
"=
i
=
,
i
=
X$
90 S0 7O 6O 5O 40
--,°.°,
....
t~s
11
theory
. . . . . ma'~olc~er theory 30 20 10 I
0
025
I
I
05
075
10 125 Ho/Hsh
Fig. 2. Comparison of experimental curve with theoretical curves. the field. Therefore these pairs remain superconductmg even for fields larger than H,h, up to He - 1.09H, h. The average field ( H ) at the sphere centers was found to be 1.59tlo for large plates of thickness W = 15. We found this limiting value by first looking at systems of various lengths L at fixed ratio W/L and
Vol. 43, No. 6
THREE.DIMENSIONAL SUSPENSIONS OF SUPERCONDUCTING GRANULES and
250
= -- 1.6(x. = 0.1)
489 (4)
for this 200 x 200 x 15 sample. (The width ( H 2) ( H ) 2 was extrapolated analogously to equatton (3) as 0.0019 + 0.0002 for the infinite sample andx, = 1.) Also the excess determined by the fourth moment changed very drastically.
2(30
150
( ( H - - ( H ) ) 4 ) ( ( H - - ( H ) ) 2 ) - 2 - - 3 = 0.18(x, = 1) 100
and
50
140
145
150
155
160
165
170
175 H/Ha
Fig. 3. Histogram for magnetic field H felt by 14,231 superconductmg spheres in 200 x 200 x 15 system for small external fields Ho. 100 8O 6O 4O 20
.... :.........~lllh
0 128 1'30 132
13l* 136
138
140
1/.2
11.4 146
HIHo
Fzg. 4. Histogram of the same system as in Fig. 3; but due to a strong field only ten percent of the spheres are still superconducting. extrapolating to L --, ae; then this extrapolated value was found to vary smoothly with W/L, roughly as (H)/Ho = 1.593 -- 0.210 W/L + O.141(W/L) 2.
(3)
Figure 3 shows the histogram of the fields observed by the spheres in the case where all the spheres are still superconducting (x, = 1;H (~ H,h), whereas Fig. 4 gives the histogram for the same system when only 10% of the spheres are still superconducting (x, = l ; H = 1.034H, h). We see that the fun distribution is roughly Gaussian whereas a strong asymmetry appears when most spheres became normal. Quantitatively, we found for the skewness of these two distributions
((H--(H))S)/((H--H)2) s/2
-- 0.15 (x,
=
l)
(5)
for the same sample size. (The absolute value of both excess and skewness increased with increasing L for both x, = O.1 and x, = 1. More details are gwen m [4].) The physical reason for this asymmetry ts obvious: The transition from the superconducting to the normal state hits mainly the spheres in large fields and thus gives a rather sharp cut-off as observed in Fig. 4. Szmilar to the monolayer calculation [3], the remaimng superconducting spheres avoid to be close to each other but showed no tendency to order periodically. Further details, mcludmg the computer program, are given elsewhere [4]. Future research along these lines could revolve higher concentrations instead of a filling factor of 0.1 only; different sphere radfi instead of unit radms for all; mhomogeneous distnbution of spheres instead of our pseudo-random procedure; and other improvements. Some of these effects could broaden the translUon regton towards the larger width observed experimentally (Fig. 2). In such future calculations, computer time could be saved by working wtth smaller systems than employed here. As our numerical evidence m Fig. 1 shows, already small systems like 100 x 100 x 15 could approximate well the desired 4050 x 4050 x 15 limit, and use 45 times less computer time than our largest system, 259 x 259 x 15. Acknowledgements - We thank G. Waysand for suggestmg this work, V. Wmkelmann form improvements m our computer program, and G. Deutscher, O. EntinWohlmann, Y. Kantor and A. Kapitulnik for discussions. REFERENCES I. 2. 3. 4.
=
= 3.83(x, = 0.1)
D. Hueber, C. Valette & G. Waysand, J. Phys. 41, L161 (1980). D. Hueber, Th~se, Universit~ de Paris Sud (1980). C. Valette, G. Waysand & D. Stauffer, Solid State Commun. 41,305 (1982). M. Hillen, Staatsexamensarbeit, Cologne University (1982).