- Email: [email protected]
Numerical investigations on influence of B-dependent flow resistivity on third harmonics of generated magnetic field
Physica C 412–414 (2004) 723–728 www.elsevier.com/locate/physc
Numerical investigations on influence of B-dependent flow resistivity on third harmonics of generated magnetic field Atsushi Kamitani a
a,*
, Ayumu Saitoh a, Soichiro Ikuno
b
Department of Informatics, Faculty of Engineering, Yamagata University, 4-3-16, Johnan, Yonezawa, Yamagata 992-8510, Japan b School of Computer Science, Tokyo University of Technology, 1404-1, Katakura, Hachioji, Tokyo 192-0982, Japan Received 29 October 2003; accepted 5 January 2004 Available online 6 May 2004
Abstract The numerical code for simulating the time evolution of the shielding current density in the high-temperature superconductor has been developed on the basis of the element-free Galerkin method. The magnetic flux density generated by the shielding current density is calculated by use of the code and its spectral analysis is performed. The results of computations show that an increase in the amplitude of the applied ac magnetic field will cause the appearance of the third harmonics of the magnetic flux density. Furthermore, it is found that the rapid growth of the third harmonics arises not from the B-dependence of the critical current density but from that of the flow resistivity. 2004 Elsevier B.V. All rights reserved. PACS: 02.70.D; 07.05.T; 74.60.J; 74.72.B Keywords: Galerkin method; Computer-modeling and simulation; Critical currents; YBCO
1. Introduction Claassen et al. have developed the new method for measuring the critical current density of the high-temperature superconductor (HTS) without making contact to the sample or modifying it in any way [1]. After placing a multi-turn coil near the HTS surface and driving it with an audio-frequency ac current, the shielding current density flows in the HTS sample. By observing the harmonics of the generated magnetic field, they found that the third harmonics develop rapidly with an increasing coil current. On the basis of the results,
they have developed the novel method for measuring the critical current density. The purpose of the present study is to develop the numerical code for analyzing the time evolution of the shielding current density on the basis of the element-free Galerkin (EFG) method [2–4] and to numerically reproduce Claassen’s experimental results by use of the code. In particular, the influence of the B-dependence of the critical current density and the flow resistivity on the third harmonics is investigated in detail. 2. Governing equations
*
Corresponding author. Tel./fax: +81-238-26-3331. E-mail address: [email protected] (A. Kamitani).
We assume for simplicity that the HTS has the same cross section over the thickness and that it is
0921-4534/$ - see front matter 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.physc.2004.01.097
724
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
placed in the homogeneous magnetic flux density parallel to its thickness direction. By taking z-axis along the thickness direction, the applied magnetic flux density is written as B 0 ¼ B0 sin 2pft ez . Throughout the present paper, the unit vectors along the x-, y-, and z-directions are denoted by ex ; ey , and ez , respectively. In order to reflect the crystallographic anisotropy [5] in the critical current density jC , we further assume the multiplethin-layer approximation [6,7]: the HTS is composed of M pieces of thin layers of thickness 2e and the shielding current density never flows across the interface between every two adjacent layers. In the following, X and oX are the cross section of the HTS and its boundary, respectively. In the pth layer of the HTS, the electric field E p and the shielding current density j p are closely related by the J –E constitutive relation, E p ¼ Eðjj p j; BÞ j p =jj p j
ð1Þ
where B denotes a magnitude of the applied magnetic flux density. As the relation, we adopt the flux-flow and flux-creep model. For the model, the explicit form of Eðj; BÞ is given in Ref. [8]. The B-dependence of the critical current density jC and that of the flow resistivity qf are also assumed as follows: jC ðBÞ ¼ jC0 ½BK =ðB þ BK Þ
ð2Þ
qf ðBÞ ¼ qf0 ðB=Bq Þd
ð3Þ
where qf0 , jC0 , Bq , BK , and d are all constants. Eq. (3) becomes the Bardeen–Stephan model for d ¼ 1, whereas qf is a constant for d ¼ 0. Thus, the Bdependence of the flow resistivity is controlled by the parameter d. On the other hand, the Bdependence of the critical current density is characterized by the parameter B0 =BK . Under the above assumptions, the state of the shielding current density can be expressed by the M-dimensional vector-valued function: ~ Sðx; tÞ ½S1 ðx; tÞ; . . . ; SM ðx; tÞ
T
B and ~ ER , are dewhere M-dimensional vectors, ~ ~ fined by B ½ez B 0 ; . . . ; ez B 0 T and ~ ER ½ez ðr E 1 Þ; . . . ; ez ðr E M ÞT . In addition, bI b is defined by denotes an identity operator and Q Z b~ Qðjx x0 jÞ~ Sðx0 ; tÞ dx0 : ð6Þ Q S X
Here, QðrÞ is the M M matrix whose (i; j) element qij ðrÞ is given in Ref. [7]. The initial and the boundary conditions to Eq. (5) are assumed as follows: ~ S ¼~ 0 at t ¼ 0 and ~ S ¼~ 0 on oX. By solving Eq. (5) together with the conditions, we can investigate the time evolution of the shielding current density.
3. Numerical methods By applying the backward Euler method to Eq. (5), we get b þ 1 bI ~ l0 Q S n þ Dt~ EnR e b þ 1 bI ~ S n 1 ~ ¼ l0 Q Bn ~ Bn 1 ð7Þ e where Dt is the time step size and the superscript n indicates the value at the nth time step t ¼ tn ð nDtÞ. Apparently, the boundary condition to Eq. (7) is written as ~ Sn ¼ ~ 0 on oX. Before the derivation of the weak form equivalent to the boundary-value problem of Eq. (7), let us first define two functional spaces as H ðXÞ f½S1 ; . . . ; SM T jSp 2 L2 ðXÞð1 6 p 6 MÞg T
where Sp ðx; tÞ=e is the z-component of the current vector potential defined by j p ¼ r ðSp ez =eÞ:
Throughout the present study, both x and x0 denote the position vectors of points on the xy plane. The time evolution of ~ Sðx; tÞ is governed by the following integro-differential equation: o b 1b ~ o ~ ~ l0 Q þ I S þ B þ ER ¼ ~ 0 ð5Þ ot e ot
ð4Þ
H ðoXÞ f½S1 ; . . . ; SM jSp 2 L2 ðoXÞð1 6 p 6 MÞg where L2 ðXÞ and L2 ðoXÞ are Hilbert spaces. In addition, the inner products in H ðXÞ and H ðoXÞ are defined by
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
h~ S; ~ S0i
M Z X p¼1
s~ T;~ T 0t
Sp ðxÞSp0 ðxÞ dx
ð8Þ
Tp ðsÞTp0 ðsÞ ds
ð9Þ
X
M I X p¼1
oX
T;~ T 0 2 H ðoXÞ. for 8~ S; ~ S 0 2 H ðXÞ and 8~ After the straightforward calculations, we can derive the weak form that is equivalent to the boundary-value problem of Eq. (7). The resulting weak form can be written in the form
1 b ~n ~ b l0 dS; Q þ I S þ sd~ kt S; ~ kt þ s~ S n ; d~ e
o ~ ~n o ~ ~n dS; Ey þ dS; Ex þ Dt ox oy
D E b þ 1 bI ~ S n 1 d~ ¼ l0 d~ S; Q S; ð~ Bn ~ Bn 1 Þ e ð10Þ where M-dimensional vectors, ~ Enx and ~ Eny , are den n T n ~ fined by Ex ½ex E 1 ; . . . ; ex E M and ~ Eny n n T ~ ½ey E 1 ; . . . ; ey E M . In addition, k is contained in H ðoXÞ and all of its components are the Lagrange multipliers. Furthermore, d~ k and d~ S are test functions, whereas ~ k and ~ S denote trial ones. Let f~ e1 ;~ e2 ; . . . ;~ eM g be an orthonormal system of the M-dimensional vector space and let /j ðxÞ (j ¼ 1; 2; . . . ; N ) be a shape function of the moving least-squares (MLS) approximation [2–4]. Since ~ Sn is not only the M-dimensional vector but a function of x, the product of ~ ep and /j ðxÞ becomes the basis function of the functional space containing ~ S n . The functional space spanned by /j ðxÞ~ ep ðj ¼ 1; 2; . . . ; N ; p ¼ 1; 2; . . . ; MÞ is denoted by V and the space spanned by Nl ðsÞ~ ep ðl ¼ 1; 2; . . . ; K; p ¼ 1; 2; . . . ; MÞ is denoted by U . Here, Nl ðsÞ is a Lagrange interpolant and is a function of arclength s along oX. In the EFG method, both d~ S and ~ S are assumed to be contained in V . Similarly, both d~ k and ~ k are assumed to be contained in U . Under the above assumptions, Eq. (10) can be discretized as W B s eðsÞ u Gðs; kÞ þ Dt
¼0 BT O k 0 0 ð11Þ
725
where the nodal vectors, s and k, correspond to ~ Sn and ~ k, respectively, and eðsÞ originates from the fourth term in the left-hand side of Eq. (10). In addition, u is a nodal vector irrelevant to s. Furthermore, the real symmetric matrix W is calcub þ ð1=eÞbI , lated from the Hermitian operator Q whereas the matrix B is evaluated from the shape functions and the Lagrange interpolants. As is apparent from Eq. (11), the nodal vectors, s and k, are the solution of the nonlinear system Gðs; kÞ ¼ 0. In this way, the initial-boundary-value problem of Eq. (5) is reduced to the problem in which the nonlinear system Gðs; kÞ ¼ 0 is solved at each time step. In order to solve the system, the adaptively underrelaxated Newton method [8] is employed.
4. Numerical results On the basis of the method explained in the previous section, we have developed the numerical code for solving the initial-boundary-value problem of Eq. (5). By means of the code, the spatial distribution of the shielding current density can be easily determined. In the present paper, the frequency of the applied magnetic field is so low that the displacement current can be neglected. Hence, the magnetic flux density can be calculated from the shielding current density through the Biot–Savart law. In this way, we can obtain the time sequence data, Bz0 ; Bz1 ; . . . ; BzJ 1 , of the magnetic field at an arbitrary location. Here, Bzj denotes the z-component of the magnetic flux density generated by the shielding current density at time t ¼ ðj þ J Þ=ðfJ Þ. By means of the fast Fourier transform (FFT) algorithm, the time sequence data are Fourier transformed to get e 0; B e 1; . . . ; B e J 1 . As the measure of the intensity of B the nth harmonics, we adopt the normalized power spectrum defined by , J 1 X 2 e nj e m j2 Pn j B jB ð12Þ m¼0
The HTS is assumed to have a square cross section of side length 2a and its thickness is denoted by D. The geometrical and physical
726
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
Table 1 Values of geometrical and physical parameters Representative length Thickness Temperature Critical current density at B ¼ 0 T Critical electric field Pinning potential Flow resistivity at B ¼ 1 T Parameter in Eq. (3) Frequency
a D T jC0 EC U0 qf0 Bq f
20 mm 2 mm 77 K 1.3 · 107 A/m2 0.1 mV/m 96 meV 7 · 10 9 X m 1T 100 Hz
parameters used in the present paper are tabulated in Table 1. Here, the values of parameters other than the cross-section shape are the same as those in Ref. [8]. The axisymmetric shielding current analysis is performed by using the finite element method in Ref. [8], whereas the axial symmetry is not assumed in the present study. The same values of parameters are adopted only to check the numerical results obtained in the present study. In addition, the time sequence data of the magnetic flux density are calculated at ðx=a; y=a; z=aÞ ¼ ð0; 0; 0:06Þ. The normalized power spectra Pn are evaluated as functions of the layer number M and are depicted in Fig. 1. This figure indicates that the spectra hardly change if M J 6. Hence, the layer number is fixed as M ¼ 6 in the following. Let us first investigate how the amplitude B0 of the applied magnetic field affects the magnetic field generated by the shielding current density. In Fig.
Fig. 1. Dependence of the normalized power spectrum Pn on the layer number M for the case with d ¼ 1:0, B0 =BK ¼ 10 2 and B0 ¼ 1:0 T. ðNÞ n ¼ 0, ðMÞ n ¼ 1, ð Þ n ¼ 2, ð Þ n ¼ 3, and ð Þ n ¼ 4.
Fig. 2. The time evolution of the magnetic flux density Bz at ðx=a; y=a; z=aÞ ¼ ð0; 0; 0:06Þ for the case with d ¼ 1:0 and B0 =BK ¼ 10 2 . (A) B0 ¼ 0:2 T, (B) B0 ¼ 0:4 T, (C) B0 ¼ 0:6 T, (D) B0 ¼ 0:8 T, and (E) B0 ¼ 1:0 T.
2, we show the time evolution of the generated magnetic flux density. This figure indicates that an increase in the amplitude B0 will raise the phase shift and make the wave form deformed from the sinusoidal wave. These tendencies suggest that the harmonics are excited with an increase in B0 . For the purpose of investigating the tendencies in detail, the power spectra Pn of the generated magnetic flux density are calculated as functions of B0 and are depicted in Fig. 3. Since the n ¼ 1 mode means the linear response, P1 is almost constant regardless of B0 . On the other hand, the spectra Pn
Fig. 3. Dependence of the normalized power spectra Pn on the applied magnetic flux density B0 for the case with d ¼ 1:0 and B0 =BK ¼ 10 2 . ðNÞ n ¼ 0, ðMÞ n ¼ 1, ð Þ n ¼ 2, ð Þ n ¼ 3, and ð Þ n ¼ 4.
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
with n ¼ 0; 2; 4 decrease with B0 after a slight increase. In this contrast, the n ¼ 3 mode grows monotonously with B0 until it becomes the second most dominant mode for B0 J 0:5 T. In other words, an increase in the amplitude of the applied magnetic field will induce only the third harmonics. This result agrees qualitatively with the experimental result obtained by Claassen et al. Next, we investigate the effect of the B-dependence of the critical current density on the third harmonics. To this end, P3 is calculated by assuming the flow resistivity as a constant. In Fig. 4, we show the dependence of P3 on B0 and B0 =BK As is apparent from this figure, P3 diminishes with
727
an increase in B0 and the decreasing tendency is not affected by the value of B0 =BK . From this result, we can conclude that the B-dependence of the critical current density does not exert an influence on the amplification of the third harmonics. Finally, we investigate the influence of the Bdependence of the flow resistivity on the third harmonics. The dependence of P3 on B0 and d is numerically determined and is depicted in Fig. 5. We see from this figure that, for d K 0:25, P3 decreases monotonously with B0 . On the other hand, P3 becomes a monotonously increasing function of B0 for the case with d J 0:75. This result implies that the B-dependence of the flow resistivity is closely related to the onset of third harmonics.
5. Conclusions We have developed the numerical code for analyzing the time evolution of the shielding current density in the HTS. By use of the code, we have calculated the magnetic flux density and have performed its spectral analysis. Conclusions obtained in the present study are summarized as follows: Fig. 4. Dependence of the normalized power spectrum P3 on B0 and B0 =BK for the case with d ¼ 0.
(1) An increase in the amplitude of the applied ac magnetic field will cause the onset of the third harmonics of the magnetic field. Such an onset of the third harmonics was also observed in the experiment by Claassen et al. (2) The appearance of the third harmonics is attributable not to the B-dependence of the critical current density but to that of the flow resistivity.
References
Fig. 5. Dependence of the normalized power spectrum P3 on B0 and d for the case with B0 =BK ¼ 10 2 .
[1] J.H. Claassen, M.E. Reeves, R.J. Soulen Jr., Rev. Sci. Instrum. 62 (1991) 996. [2] T. Belytschko, Y.Y. Lu, L. Gu, Int. J. Numer. Methods Eng. 37 (1994) 229. [3] A. Kamitani, A. Saitoh, S. Ikuno, T. Yokono, Physica C 392–396 (2003) 118. [4] A. Kamitani, T. Yokono, K. Hasegawa, IEEE Trans. Appl. Supercond. 13 (2003) 1660.
728
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
[5] M. Murakami, S. Gotoh, H. Fujimoto, K. Yamaguchi, N. Koshizuka, S. Tanaka, Supercond. Sci. Technol. 4 (1991) S43. [6] Y. Yoshida, M. Uesaka, K. Miya, Int. J. Appl. Electromagn. Mater. 5 (1994) 83.
[7] A. Kamitani, S. Ohshima, IEICE Trans. Electron. E82-C (1999) 766. [8] A. Kamitani, K. Hasegawa, T. Yokono, IEEE Trans. Magn. 39 (2003) 1511.
Numerical investigations on influence of B-dependent flow resistivity on third harmonics of generated magnetic field Atsushi Kamitani a
a,*
, Ayumu Saitoh a, Soichiro Ikuno
b
Department of Informatics, Faculty of Engineering, Yamagata University, 4-3-16, Johnan, Yonezawa, Yamagata 992-8510, Japan b School of Computer Science, Tokyo University of Technology, 1404-1, Katakura, Hachioji, Tokyo 192-0982, Japan Received 29 October 2003; accepted 5 January 2004 Available online 6 May 2004
Abstract The numerical code for simulating the time evolution of the shielding current density in the high-temperature superconductor has been developed on the basis of the element-free Galerkin method. The magnetic flux density generated by the shielding current density is calculated by use of the code and its spectral analysis is performed. The results of computations show that an increase in the amplitude of the applied ac magnetic field will cause the appearance of the third harmonics of the magnetic flux density. Furthermore, it is found that the rapid growth of the third harmonics arises not from the B-dependence of the critical current density but from that of the flow resistivity. 2004 Elsevier B.V. All rights reserved. PACS: 02.70.D; 07.05.T; 74.60.J; 74.72.B Keywords: Galerkin method; Computer-modeling and simulation; Critical currents; YBCO
1. Introduction Claassen et al. have developed the new method for measuring the critical current density of the high-temperature superconductor (HTS) without making contact to the sample or modifying it in any way [1]. After placing a multi-turn coil near the HTS surface and driving it with an audio-frequency ac current, the shielding current density flows in the HTS sample. By observing the harmonics of the generated magnetic field, they found that the third harmonics develop rapidly with an increasing coil current. On the basis of the results,
they have developed the novel method for measuring the critical current density. The purpose of the present study is to develop the numerical code for analyzing the time evolution of the shielding current density on the basis of the element-free Galerkin (EFG) method [2–4] and to numerically reproduce Claassen’s experimental results by use of the code. In particular, the influence of the B-dependence of the critical current density and the flow resistivity on the third harmonics is investigated in detail. 2. Governing equations
*
Corresponding author. Tel./fax: +81-238-26-3331. E-mail address: [email protected] (A. Kamitani).
We assume for simplicity that the HTS has the same cross section over the thickness and that it is
0921-4534/$ - see front matter 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.physc.2004.01.097
724
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
placed in the homogeneous magnetic flux density parallel to its thickness direction. By taking z-axis along the thickness direction, the applied magnetic flux density is written as B 0 ¼ B0 sin 2pft ez . Throughout the present paper, the unit vectors along the x-, y-, and z-directions are denoted by ex ; ey , and ez , respectively. In order to reflect the crystallographic anisotropy [5] in the critical current density jC , we further assume the multiplethin-layer approximation [6,7]: the HTS is composed of M pieces of thin layers of thickness 2e and the shielding current density never flows across the interface between every two adjacent layers. In the following, X and oX are the cross section of the HTS and its boundary, respectively. In the pth layer of the HTS, the electric field E p and the shielding current density j p are closely related by the J –E constitutive relation, E p ¼ Eðjj p j; BÞ j p =jj p j
ð1Þ
where B denotes a magnitude of the applied magnetic flux density. As the relation, we adopt the flux-flow and flux-creep model. For the model, the explicit form of Eðj; BÞ is given in Ref. [8]. The B-dependence of the critical current density jC and that of the flow resistivity qf are also assumed as follows: jC ðBÞ ¼ jC0 ½BK =ðB þ BK Þ
ð2Þ
qf ðBÞ ¼ qf0 ðB=Bq Þd
ð3Þ
where qf0 , jC0 , Bq , BK , and d are all constants. Eq. (3) becomes the Bardeen–Stephan model for d ¼ 1, whereas qf is a constant for d ¼ 0. Thus, the Bdependence of the flow resistivity is controlled by the parameter d. On the other hand, the Bdependence of the critical current density is characterized by the parameter B0 =BK . Under the above assumptions, the state of the shielding current density can be expressed by the M-dimensional vector-valued function: ~ Sðx; tÞ ½S1 ðx; tÞ; . . . ; SM ðx; tÞ
T
B and ~ ER , are dewhere M-dimensional vectors, ~ ~ fined by B ½ez B 0 ; . . . ; ez B 0 T and ~ ER ½ez ðr E 1 Þ; . . . ; ez ðr E M ÞT . In addition, bI b is defined by denotes an identity operator and Q Z b~ Qðjx x0 jÞ~ Sðx0 ; tÞ dx0 : ð6Þ Q S X
Here, QðrÞ is the M M matrix whose (i; j) element qij ðrÞ is given in Ref. [7]. The initial and the boundary conditions to Eq. (5) are assumed as follows: ~ S ¼~ 0 at t ¼ 0 and ~ S ¼~ 0 on oX. By solving Eq. (5) together with the conditions, we can investigate the time evolution of the shielding current density.
3. Numerical methods By applying the backward Euler method to Eq. (5), we get b þ 1 bI ~ l0 Q S n þ Dt~ EnR e b þ 1 bI ~ S n 1 ~ ¼ l0 Q Bn ~ Bn 1 ð7Þ e where Dt is the time step size and the superscript n indicates the value at the nth time step t ¼ tn ð nDtÞ. Apparently, the boundary condition to Eq. (7) is written as ~ Sn ¼ ~ 0 on oX. Before the derivation of the weak form equivalent to the boundary-value problem of Eq. (7), let us first define two functional spaces as H ðXÞ f½S1 ; . . . ; SM T jSp 2 L2 ðXÞð1 6 p 6 MÞg T
where Sp ðx; tÞ=e is the z-component of the current vector potential defined by j p ¼ r ðSp ez =eÞ:
Throughout the present study, both x and x0 denote the position vectors of points on the xy plane. The time evolution of ~ Sðx; tÞ is governed by the following integro-differential equation: o b 1b ~ o ~ ~ l0 Q þ I S þ B þ ER ¼ ~ 0 ð5Þ ot e ot
ð4Þ
H ðoXÞ f½S1 ; . . . ; SM jSp 2 L2 ðoXÞð1 6 p 6 MÞg where L2 ðXÞ and L2 ðoXÞ are Hilbert spaces. In addition, the inner products in H ðXÞ and H ðoXÞ are defined by
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
h~ S; ~ S0i
M Z X p¼1
s~ T;~ T 0t
Sp ðxÞSp0 ðxÞ dx
ð8Þ
Tp ðsÞTp0 ðsÞ ds
ð9Þ
X
M I X p¼1
oX
T;~ T 0 2 H ðoXÞ. for 8~ S; ~ S 0 2 H ðXÞ and 8~ After the straightforward calculations, we can derive the weak form that is equivalent to the boundary-value problem of Eq. (7). The resulting weak form can be written in the form
1 b ~n ~ b l0 dS; Q þ I S þ sd~ kt S; ~ kt þ s~ S n ; d~ e
o ~ ~n o ~ ~n dS; Ey þ dS; Ex þ Dt ox oy
D E b þ 1 bI ~ S n 1 d~ ¼ l0 d~ S; Q S; ð~ Bn ~ Bn 1 Þ e ð10Þ where M-dimensional vectors, ~ Enx and ~ Eny , are den n T n ~ fined by Ex ½ex E 1 ; . . . ; ex E M and ~ Eny n n T ~ ½ey E 1 ; . . . ; ey E M . In addition, k is contained in H ðoXÞ and all of its components are the Lagrange multipliers. Furthermore, d~ k and d~ S are test functions, whereas ~ k and ~ S denote trial ones. Let f~ e1 ;~ e2 ; . . . ;~ eM g be an orthonormal system of the M-dimensional vector space and let /j ðxÞ (j ¼ 1; 2; . . . ; N ) be a shape function of the moving least-squares (MLS) approximation [2–4]. Since ~ Sn is not only the M-dimensional vector but a function of x, the product of ~ ep and /j ðxÞ becomes the basis function of the functional space containing ~ S n . The functional space spanned by /j ðxÞ~ ep ðj ¼ 1; 2; . . . ; N ; p ¼ 1; 2; . . . ; MÞ is denoted by V and the space spanned by Nl ðsÞ~ ep ðl ¼ 1; 2; . . . ; K; p ¼ 1; 2; . . . ; MÞ is denoted by U . Here, Nl ðsÞ is a Lagrange interpolant and is a function of arclength s along oX. In the EFG method, both d~ S and ~ S are assumed to be contained in V . Similarly, both d~ k and ~ k are assumed to be contained in U . Under the above assumptions, Eq. (10) can be discretized as W B s eðsÞ u Gðs; kÞ þ Dt
¼0 BT O k 0 0 ð11Þ
725
where the nodal vectors, s and k, correspond to ~ Sn and ~ k, respectively, and eðsÞ originates from the fourth term in the left-hand side of Eq. (10). In addition, u is a nodal vector irrelevant to s. Furthermore, the real symmetric matrix W is calcub þ ð1=eÞbI , lated from the Hermitian operator Q whereas the matrix B is evaluated from the shape functions and the Lagrange interpolants. As is apparent from Eq. (11), the nodal vectors, s and k, are the solution of the nonlinear system Gðs; kÞ ¼ 0. In this way, the initial-boundary-value problem of Eq. (5) is reduced to the problem in which the nonlinear system Gðs; kÞ ¼ 0 is solved at each time step. In order to solve the system, the adaptively underrelaxated Newton method [8] is employed.
4. Numerical results On the basis of the method explained in the previous section, we have developed the numerical code for solving the initial-boundary-value problem of Eq. (5). By means of the code, the spatial distribution of the shielding current density can be easily determined. In the present paper, the frequency of the applied magnetic field is so low that the displacement current can be neglected. Hence, the magnetic flux density can be calculated from the shielding current density through the Biot–Savart law. In this way, we can obtain the time sequence data, Bz0 ; Bz1 ; . . . ; BzJ 1 , of the magnetic field at an arbitrary location. Here, Bzj denotes the z-component of the magnetic flux density generated by the shielding current density at time t ¼ ðj þ J Þ=ðfJ Þ. By means of the fast Fourier transform (FFT) algorithm, the time sequence data are Fourier transformed to get e 0; B e 1; . . . ; B e J 1 . As the measure of the intensity of B the nth harmonics, we adopt the normalized power spectrum defined by , J 1 X 2 e nj e m j2 Pn j B jB ð12Þ m¼0
The HTS is assumed to have a square cross section of side length 2a and its thickness is denoted by D. The geometrical and physical
726
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
Table 1 Values of geometrical and physical parameters Representative length Thickness Temperature Critical current density at B ¼ 0 T Critical electric field Pinning potential Flow resistivity at B ¼ 1 T Parameter in Eq. (3) Frequency
a D T jC0 EC U0 qf0 Bq f
20 mm 2 mm 77 K 1.3 · 107 A/m2 0.1 mV/m 96 meV 7 · 10 9 X m 1T 100 Hz
parameters used in the present paper are tabulated in Table 1. Here, the values of parameters other than the cross-section shape are the same as those in Ref. [8]. The axisymmetric shielding current analysis is performed by using the finite element method in Ref. [8], whereas the axial symmetry is not assumed in the present study. The same values of parameters are adopted only to check the numerical results obtained in the present study. In addition, the time sequence data of the magnetic flux density are calculated at ðx=a; y=a; z=aÞ ¼ ð0; 0; 0:06Þ. The normalized power spectra Pn are evaluated as functions of the layer number M and are depicted in Fig. 1. This figure indicates that the spectra hardly change if M J 6. Hence, the layer number is fixed as M ¼ 6 in the following. Let us first investigate how the amplitude B0 of the applied magnetic field affects the magnetic field generated by the shielding current density. In Fig.
Fig. 1. Dependence of the normalized power spectrum Pn on the layer number M for the case with d ¼ 1:0, B0 =BK ¼ 10 2 and B0 ¼ 1:0 T. ðNÞ n ¼ 0, ðMÞ n ¼ 1, ð Þ n ¼ 2, ð Þ n ¼ 3, and ð Þ n ¼ 4.
Fig. 2. The time evolution of the magnetic flux density Bz at ðx=a; y=a; z=aÞ ¼ ð0; 0; 0:06Þ for the case with d ¼ 1:0 and B0 =BK ¼ 10 2 . (A) B0 ¼ 0:2 T, (B) B0 ¼ 0:4 T, (C) B0 ¼ 0:6 T, (D) B0 ¼ 0:8 T, and (E) B0 ¼ 1:0 T.
2, we show the time evolution of the generated magnetic flux density. This figure indicates that an increase in the amplitude B0 will raise the phase shift and make the wave form deformed from the sinusoidal wave. These tendencies suggest that the harmonics are excited with an increase in B0 . For the purpose of investigating the tendencies in detail, the power spectra Pn of the generated magnetic flux density are calculated as functions of B0 and are depicted in Fig. 3. Since the n ¼ 1 mode means the linear response, P1 is almost constant regardless of B0 . On the other hand, the spectra Pn
Fig. 3. Dependence of the normalized power spectra Pn on the applied magnetic flux density B0 for the case with d ¼ 1:0 and B0 =BK ¼ 10 2 . ðNÞ n ¼ 0, ðMÞ n ¼ 1, ð Þ n ¼ 2, ð Þ n ¼ 3, and ð Þ n ¼ 4.
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
with n ¼ 0; 2; 4 decrease with B0 after a slight increase. In this contrast, the n ¼ 3 mode grows monotonously with B0 until it becomes the second most dominant mode for B0 J 0:5 T. In other words, an increase in the amplitude of the applied magnetic field will induce only the third harmonics. This result agrees qualitatively with the experimental result obtained by Claassen et al. Next, we investigate the effect of the B-dependence of the critical current density on the third harmonics. To this end, P3 is calculated by assuming the flow resistivity as a constant. In Fig. 4, we show the dependence of P3 on B0 and B0 =BK As is apparent from this figure, P3 diminishes with
727
an increase in B0 and the decreasing tendency is not affected by the value of B0 =BK . From this result, we can conclude that the B-dependence of the critical current density does not exert an influence on the amplification of the third harmonics. Finally, we investigate the influence of the Bdependence of the flow resistivity on the third harmonics. The dependence of P3 on B0 and d is numerically determined and is depicted in Fig. 5. We see from this figure that, for d K 0:25, P3 decreases monotonously with B0 . On the other hand, P3 becomes a monotonously increasing function of B0 for the case with d J 0:75. This result implies that the B-dependence of the flow resistivity is closely related to the onset of third harmonics.
5. Conclusions We have developed the numerical code for analyzing the time evolution of the shielding current density in the HTS. By use of the code, we have calculated the magnetic flux density and have performed its spectral analysis. Conclusions obtained in the present study are summarized as follows: Fig. 4. Dependence of the normalized power spectrum P3 on B0 and B0 =BK for the case with d ¼ 0.
(1) An increase in the amplitude of the applied ac magnetic field will cause the onset of the third harmonics of the magnetic field. Such an onset of the third harmonics was also observed in the experiment by Claassen et al. (2) The appearance of the third harmonics is attributable not to the B-dependence of the critical current density but to that of the flow resistivity.
References
Fig. 5. Dependence of the normalized power spectrum P3 on B0 and d for the case with B0 =BK ¼ 10 2 .
[1] J.H. Claassen, M.E. Reeves, R.J. Soulen Jr., Rev. Sci. Instrum. 62 (1991) 996. [2] T. Belytschko, Y.Y. Lu, L. Gu, Int. J. Numer. Methods Eng. 37 (1994) 229. [3] A. Kamitani, A. Saitoh, S. Ikuno, T. Yokono, Physica C 392–396 (2003) 118. [4] A. Kamitani, T. Yokono, K. Hasegawa, IEEE Trans. Appl. Supercond. 13 (2003) 1660.
728
A. Kamitani et al. / Physica C 412–414 (2004) 723–728
[5] M. Murakami, S. Gotoh, H. Fujimoto, K. Yamaguchi, N. Koshizuka, S. Tanaka, Supercond. Sci. Technol. 4 (1991) S43. [6] Y. Yoshida, M. Uesaka, K. Miya, Int. J. Appl. Electromagn. Mater. 5 (1994) 83.
[7] A. Kamitani, S. Ohshima, IEICE Trans. Electron. E82-C (1999) 766. [8] A. Kamitani, K. Hasegawa, T. Yokono, IEEE Trans. Magn. 39 (2003) 1511.