Discussion and Results.

The analytical solution is formulated and solved in terms of the current density vector, which in the presented case has only the z component. The numerical solution to the eddy-current problem is formulated in terms of the magnetic vector potential and the current density is then derived from these calculations by using Maxwell's first equation. The current density is plotted as a function of the radius for different time steps and different radial grid sizes.

The finite element program was run on a cylindrically symmetric grid with as nearly equilateral triangular elements as possible. The mesh is illustrated below in figure (6.1). Along a radius, there are bands of interlaced triangles. These values are averaged along each band for the triangles whose vertex pointed towards the center and for the triangles whose vertex pointed radial outward.

Figure (6.1) Finite element grid with 9 nodes along the radius.

Thus, the values of the current density as a function of the radius are plotted. The finite difference program is formulated in polar coordinates so that the same type of problem could be compared.

It is desirable to estimate how accurate these results are. One could try to estimate the order of convergence of these methods by plotting the error of the current density as a function of the radial grid size and the time step. Unfortunately, it does not make much sense to estimate the order of convergence with no theoretical justification. The plots presented clearly show the convergence of the finite element and finite difference solution to the analytical solution when it is available in the case of non-varying conductivity.

For the case of the finite element method, Zlamal [6.1] gives an error estimate on the magnetic vector potential. This is the closest applicable literature because other authors do not discuss the source term in the diffusion problem. It should be stated that an analysis which disregards this source term is not justifiable and moreover, the results are not valid when the transient effects are large, as in the present problem. The drawback to the above-stated error estimate is that it is only formulated for a time periodic problem and the summation is made over the period. Most of the error estimates given in the literature are global estimates. This is not surprising because this is obtained from the residuals. However, there are many applications in which one is very interested in knowing what the error estimate is at some particular region. Some research has been devoted to this area for the equilibrium and periodic steady-state problems, see Pinchuk and Silvester [6.2]. Bietermann and Babuska [6.3] have investigated error estimates for time-varying problems, however most of the work was oriented towards fluid dynamics problems which lead to slightly different equations. These differences are significant enough not to allow one to use his results without further study and possibly a rederivation under the formulation presented in this thesis.

Electromagnetic accelerators are currently being researched, and the interested reader is referred to Perry [6.4] and Garg and Weiss [6..5] for the engineering considerations. The parameters for the test data were chosen because it is in the range of that being studied practically and could be compared and criticized. A huge current, 0.5 megaamps, with a rise time of 0.5 milliseconds is passed along a narrow copper conductor with a radius of 0.03 meters. The results are displayed at 0.5 milliseconds. and at 2.0 milliseconds. The calculations were all done on a VAX 11/780 under the VMS operating system. The programs were all written in single precision FORTRAN.

In this short period of time, it is possible to completely neglect the thermal diffusion because this effect occurs on a much slower time scale then the one which describes the diffusion of the current density, see Perry [6.4] or Mocanu [6..6] or Miller [6.7].

In the first plots, figures (6.2) the conductivity of copper does not change as a function of the temperature. The results from the finite element and finite difference program are compared against each other and against the analytical results. It is seen that the finite element and finite difference results are very close to each other and approach the calculated analytical solution as the space and time discretization is made finer.

Fig. 6.2 Curves for different space grids for finite difference and finite element methods with temperature independent conductivity diffusion

( s = 0.5 x 108 ), Dt = 0.125 msec.

1) Analytical solution

2) FD Dr = R/8

3) FD Dr = R/16

4) FD Dr = R/32

5) FD Dr = R/64

6) FE Dr = R/8

7) FE Dr = R/16

8) FE Dr = R/32

9) FE Dr = R/64

In the second plots, figure (6.3), the effect of temperature upon the conductivity starting at 273 °K is included in the calculations. The conductivity does not change appreciably in this temperature rang (see appendix 5) and therefore the computed results do not change much from those discussed above.

In the third plot, figures (6.4), the temperature of the conductor is taken to be 80 °K because in this range, the conductivity varies significantly with temperature. Moreover, some of the recent experiments were performed in this range. The change in the current density is significant for this case because the diffusion constant changes as does the temperaure of the element in the conductor. This is the non-linear behavior of the solution. The method used to obtain this effect was simple in that the source current is nothing more than a constant multiple of the conductivity lattice (a lattice

of elements arranged in a grid representing the conductor) added to the "free-running" solution.

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

Figure 6. Comparisons of the diffusion process at various temperatures.

Current density (mega-amp/m2)vs. radial distance (meters)

at time = . msec

1) Analytical solution. (----- ----- -----)

2) Finite element solution using constant material properties s=5.0x10-8 (- - - - - -)

3) Finite element solution using material properties at 273 °K.

(--- --- ---)

4) Finite element solution using material properties at 80 °K.

(----- - ----- - ----- -)

It is worth noting, that when dealing with short rise times,

0.5 x10-6 seconds, the finite difference program gave better results than the finite element program with 16 or 32 nodes on the radius. Interestingly however, decreasing the times steps from .25 x 10-6 seconds to .025 x 10-6 seconds did not noticeably change the calculated values for both methods. It is not clear what causes this, but it suggests that small elements are required for problems with fast time changes. This however disagrees with the Von Neumann Fourier stability arguement that the time steps should be reduced for fixed space grids or given fixed time spacing, the space grid should be made corse in order to keep the quantatity q = n0 Dt / s Dr2 as small as possible. However, since we are using implicit (backward difference) scheme, stability should be independent of the space or time spacing.

Figure 6. The current density is plotted against the radius. The data was obtained using finite differences at time = 0.002 msec, Dt = 0.0005 msec, Dr=R/32.

Figure 6. The current density is plotted against the radius. The data was obtained using finite differences at time = 0.002 msec, Dt = 0.0005 msec, Dr=R/32.

Figure 6. The current density is plotted against the radius. The data was obtained using finite differences at time = 0.002 msec, Dt = 0.0005 msec, Dr=R/32.

Figure 6. The current density is plotted against the radius. The data was obtained using finite elements at time = 0.002 msec., Dt = 0.00005 msec.

(40 time steps), Dr=R/32.

A direct comparison is not done since the finite difference program used the cylindrical symmetry to reduce the dimension of the problem while the finite element program did not.

A one-dimensional finite difference program is presented in the appendix. The first program FD was written to solve the diffusion equation with constant coefficients. The second program FDT does take into account non-linear material properties. A two space dimensional program was not written because a one space dimensional program would work for this problem. A one-dimensional finite element program was not written because a packaged two-dimensional program, TDTEMP, was available. Moreover, one of the earlier motives for this thesis was to verify that this program gives good results. References are made to the literature where it is possible to find a finite element program which could be reconstructed to work on our problem. TDTEMP is the property of the General Electric Company and they have not given permission to reproduce it here. However, the important features of the program are described. TDTEMP might be obtainable through GE CR&D in Schnectady, NY or GE Advanced Electrical Systems Division in Pittsfield, MA. at a later date.

It appears that as applied in practice, the finite element program works well provided that the change in time is limited to be greater than 10-4 seconds. This suggests that there should be a further theoretical and practical study of the program with varying space and true meshes. If the time mesh is small and the relative space mesh is large, then the time change of variable is less than the errors due to the finite element approximation errors. The space errors in finite differences are well behaved, more regularly than the space errors in finite elements.

Conclusion.

The problem of the penetration of the transient electromagnetic field into a conducting metallic cylindrical conductor may be solved by with at least three methods: Laplace transforms, the finite difference method, and the finite element method. The relationships derived in this thesis permit the computation of the magnetic vector potential and the current density in an infinitely long conductor. For simplicity, we consider a one and a two dimensional diffusion equation with given initial and boundary conditions. The solutions are given for a sudden increase in the current flowing through the conductor over a short time interval. Detailed tables and graphs facilitate the use of these methods in practical computation. It is seen that the current is concentrated on the outer surface of this cylinder and creates a strongly nonuniform heating from the outer surface to the interior of the conductor. For temperature dependent material conductivity, the current diffuses faster into the interior of the conductor because the source current in the interior of the conductor is scaled much larger. The comparison of the finite difference and the finite element quantaties gives good results. We see that as both the time and space mesh become finer, the solutions tend to converge to the analytic solution.

The numerical techniques used not only show promise of providing the general mathematical modeling capacity needed for the development of current losses estimation in the design of electromagnetic accelerators, but also might serve to attract the interest of mathematicians to this difficult problem area.

Considerable work remains to be done in extending this formulation to conductors with three space dimensions with non-linear material properties. Moreover, much more work is necessary in the development of error estimation of finite difference and finite element methods applied to a transient diffusion problem with a source term.

References:

[6.1] Zlamal, M., "Finite Element Solution of Quasistationary Nonlinear Magnetic Fields," *R.A.I.R.O.-Numerical Analysis* Vol. 16. No. 2 (1982) pp.161-191.

[6.2] Pinchuk and Silvester, "An Adaptive Mesh Generation Program" *IEEE Transactions on Magnetics*. (1984).

[6.3] Bieterman, M. and Babuska, I., "The Finite Element Method for Parabolic Equations: I. A Posteriori Error Estimations," Numerische Mathematik 40 pp.339-371 (1982).

Bieterman, M. and Babuska, I., "The Finite Element Method for Parabolic Equations: II. A Posteriori Error Estimations," Numerische Mathematik 40 pp.372-399 (1982)

[6.4] Perry, M., *Workshop on Electromagnetic Field Computation*. IEEE Region1, Schnectady, NY, Oct. 1986.

[6.5] Garg, Weiss, and Del Vecchio, R. M. "A Parametric Study of Electomagnetic Launcher Rails Using the WEMAP System", *Workshop on Electromagnetic Field Computation.* IEEE Region1, Schnectady, NY, Oct. 1986.

[6.6] Mocanu, C.I., "Equivalent Schemes of the Circular Conductor at Transient State With Skin Effect," *Rev. Roum. Sci. Techn.- Electrotechnique et Energetique*, Bucarest, vol. 16 (1971), No. 2, pp.235-254.

[6.7] Miller, K. W., "Diffusion of Electric Current into Rods, Tubes, and Flat Surfaces," *AIEE Transactions* (1947), 66, pp. 1496-1504.