1

Chapter 1
Introduction
This chapter presents a detailed introduction in the subject of coupled field
analysis, which describes particularly the application of coupled field analysis in
electromagnetic devices, with a section dedicated to the history of electromagnetic
forming.
1.1 History and importance of coupled field analysis
A coupled-field analysis is an analysis that takes into account the interaction
(coupling) between two or more disciplines (fields) of engineering, such as
electromagnetic field and mechanical field in electromechanical devices,
electromagnetic field and thermal field in induction-heating applications etc, and can
be used to investigate many physical problems.
When electricity was first introduced, man discovered the effects of strong magnetic
field on conducting bodies: the conducting bodies can be moved and deformed.
These caused by electromagnetic force and both effects (movement and deformation)
have been studied for their positive and negative aspects, ever since the first electrical
machine was designed and brought into operation.
During operation of electrical machines the electromagnetic force produces
sometimes, besides the motion, the undesired deformation of some parts of the
electrical machines. On the other hand, in other applications certain is the desired
effect of electromagnetic force, accompanying the deformation as in electromagnetic
forming.
In the Fig.(1.1), all the fields, that can exist in a transformer, are shown. The coupling
among these fields is also shown in the same figure. It is useful in analyze the
performance and predicting accurate response of the device in different operating
conditions.
2


Fig. 1.1: Coupling between different fields in transformer
Another typical procedure of coupled field analysis is the electromagnetic forming, a
high-velocity forming process aiming at the shaping of conducting objects (the so-
called workpieces) in strong electromagnetic fields. The electromagnetic forming is
the process of giving shape to a metallic object using electromagnetic forces. A high
frequency damped sinusoidal current obtained by L-C resonant circuit is made pass
through a coil kept in front of a metallic object. This damped sinusoidal current,
flowing through the coil, produces a time-varying magnetic field which induces the
eddy currents in the metallic object. This leads to an electromagnetic force, also
known as the Lorentz force, which is given by JB, where J is the eddy current
density and B is the magnetic flux density, acting on the eddy current region. This
force gets transferred to the metal object leading to its deformation in the desired
manner as shown in fig. (1.2)
3


Fig. 1.2: Deformed work piece with electromagnetic forming
1.2 Motivation and scope of work
A coupled field analysis is a study of interaction of multiple fields with at least
a fraction of their domain being common. In coupled magneto-structural analysis, the
electromagnetic fields, produced by current or voltage-fed structure in the
electromagnetic devices, interact with the materials present in that field. This
interaction will result in the production of some types of forces, may be magnetization
or Lorentz forces depending on the materials present. The structural field analysis in
the coupled magneto-structural field analysis investigates the effects of these forces
on the mechanics of the materials. Therefore, the final output of a coupled magneto-
structural system will be the deformation and stresses in the structural model. In this
report, the analysis is carried out over the electromagnetic structure involves non-
ferromagnetic materials and various methodologies are discussed for the computation
of forces for different types of excitations. In subsequent structural analysis these
forces are applied as the input and output will be the deformation and stresses in the
structure. The coupled analysis, discussed here, is carried out by using a numerical
technique, known as Finite Element Method. The main advantage of the finite
4

element method in the coupled field analysis is that the same geometrical
discretization information can be used in the different subsequent field analyses. The
structural analysis, in this work, is restricted to the linear-elastic range. This analysis
may be used to investigate the deformations and vibrations in the winding of electrical
machines. Noise is also one of the main problems to be considered while designing
electrical machines. The present analysis may also be useful for the noise problems in
the electrical machines. This work describes a way to analyze and understand
deformation phenomena of electrical conductors under an electromagnetic field. A
coupled model of magneto-elastic analysis is presented, in which sequential coupling
is used among the different fields involved. The cases, described above, investigate
the negative aspects of the effects of the electromagnetic forces. The same analysis
can be used for the positive aspects in which the deformation in the metals is desired.
It is known as the electromagnetic forming.
Electromagnetic forming process has many advantages over traditional forming
process like improved surface quality, a noncontact type of process, repeatable and
can be controllable by varying electrical parameters. It is useful if we can estimate the
amount and distribution of electromagnetic force acting on the object to be formed.
This work deals with the computation of electromagnetic forces in different operating
conditions as in case of steady flow of direct current, sinusoidal alternating current
and damped-pulse current. Codes in MATLAB are developed which implement force
computation methods through the application of FEM. Further, estimation of the
forces and their distribution on the metallic object (to be formed using
electromagnetic forming process) is done. The results from the code are compared
with the results obtained using ANSYS (a commercial FEM software)
1.3 Outline of the thesis
In chapter 2, a brief literature survey is given on the electromagnetic field
computation, coupled magneto-structural analysis and electromagnetic forming
process by finite element method.
Chapter 3 presents, a general theory about finite element method and electromagnetic
field theory. This chapter describes the electromagnetic field theory, with potential
5

function formulation, results in boundary value problem. The finite element method is
efficient approximation method to solve such type of problems.
Chapter 4 describes the finite element formulation of electromagnetic fields and it
gives linear system of algebraic equations which can be solved by any standard
method for such systems. In this chapter simulation and results are given for linear
static analysis, harmonic and transient case for electromagnetics field analysis.
Chapter 5, this chapter presents a brief theory about elastic and plastic theory and fem
formulation of elastic and elastoplastic range of material. And the various coupling
schemes for magneto-structural problems are explained. The simulation of coupled
magneto-structural by sequential coupling method is presented in the same chapter.
Chapter 6, in this chapter some conclusions are drawn from the analysis of proposed
problem. Some ideas for future research on coupled field analysis of electromagnetic
devices are also given.














6

Chapter 2
Literature review
The first effort to solve the magnetic field problems by finite element method
were made in the late 1960s. In the analysis of electrical machines by finite element
method was first applied to synchronous machines in [1] and to dc machines in [2]
because the operation of these machine types can be approximate modeled by
stationary fields. Even nowadays, most of research work concerning the magnetic
field analysis of the electric machines modeling of the synchronous machines. They
have been analyzed by using step by step methods to solve the time dependence in [3]
and by three dimensional finite element formulations in [4].
The field of induction motor has to be solved with a method takes the time-
independence in to account. Probably the difficulties connected with the solution of
time dependent nonlinear fields have postponed the numerical analysis of induction
motors. The first publication dealing with this problem appeared at the beginning of
the 1980s. In the present work, only two-dimensional formulations have been used.
Ito [5] computed the nonlinear field of an induction motor using an eddy current
formulation and assuming sinusoidal time variation.
Bouillaut & Razek [6] and Brunelli [7] used time-stepping methods to calculate the
time variation of magnetic fields in induction motors.
Shen [8] presents a formulation for the sinusoid ally varying field quantities. In the
1990s, as the computing resources progressed, the electromagnetic forming process
could be simulated better. This progress has been followed by a large number of
papers dealing with numerical simulations of electromagnetic forming. In papers by
Bendjima et al. [9, 10], Azzouz et al. [11], and Meriched and Feliachi [12], numerical
results obtained with a finite-element method using simplified models of the process
were presented and compared with the experimental ones. The numerical results
showed quite a good agreement with the experimental ones. The finite element
formulations for electromagnetic devices are described in [13]. In [14], [15] and [16],
7

various numerical techniques, such as finite difference method, finite element method,
are explained for various electromagnetic cases. The main concern of the presented
work in this report is the coupled-magneto structural analysis, in which the
electromagnetic forces play an important role in the analysis. Some papers on
electromagnetic force computations are investigated. In [17], Maxwell stress method
for the magnetic force computation in the magnetized core is described. In [18]
authors described application of the Maxwell stress method in more accurate way in
the magnetized force computation. T. Kabashima

et al. [19] and

G. Henneberger et
al. [20] investigate an alternative method known as equivalent magnetizing current
source for the magnetizing force computation.

Z. Ren [21] compares various
formulations for the magnetized force computation. In [22], some aspects of magnetic
force computation and mechanical behaviors of the ferromagnetic materials. J. L.
Coulomb and G. Meunier [23] describe the Coulombs virtual work method for the
force and torque computation for 2-D and 3-D cases. In [24], the method for local
magnetic force distribution by using coulombs virtual work method by direct
differentiating the magnetic energy or co-energy with respect to the virtual
displacement of each node.

A. V. Kank and S. V. Kulkarni [25] investigate the
computation of Lorentz forces and magnetic forces by using FEM in the rotting
electrical machines. The magnetic force computation using finite element method
technique in global quantities and local force distribution is investigated in [26].
The mechanical behavior of a metal object, under the applied forces, is described by
elastic, elasto-plastic and visco-elsto-plastic range of the characteristics of the metal
object [27], [28] and [29].
C. Karch and K. Roll [30] investigate the electromagnetic forming process by finite
element method. The numerical model, presented in [30], predicts the electromagnetic
field, temperature, stress, and deformation properties that occur during the forming
process. The numerical results of the tube deformation are compared with available
experimental data.
The finite element formulations for the structural analysis are explained in [31], [32],
and [35]. The electromagnetic forming of aluminum tube in the form of compression
and expansion of the tube is investigated by using a constitutive finite element model
8

in [33] and the simulation of electromagnetic forming is carried out in the commercial
FEM software ANSYS. The simulation results are shown in good agreement with
analytical solutions [34].
The visco-plastic behavior of the metal in metal forming process is analyzed in [36].
The coupled magneto-elastic analysis by finite element method for ferromagnetic
materials is presented in [37]. The magnetostriction effect is also taken in to account
in this analysis.
The mechanical behavior of a solenoidal coil under the electromagnetic forces is
analyzed in [38]. The input data, available in this paper, is also used in the presented
linear sequentially coupled magneto-structural analysis in this report.











9

Chapter 3
FEM: Theory and Application
The finite element method has become a well established method in many
fields of computer aided engineering, such as, electromagnetic field computation,
fluid dynamics and structural analysis.
The finite element method is basically, an efficient approximation method to solve the
partial differential equations or boundary-value problems, which are frequently
occurred in different areas of engineering.
3.1 Finite element method- Basic concept
There are three main steps during the solution of partial differential equation
(PDE) with the finite element method. First the domain on which the PDE should be
solved, descritized into finite elements. Depending on the dimension of the problem,
this can be triangles, squares, rectangles or tetrahedrons, cubes, or hexahedrons. The
solution of PDE is approximated by piecewise continuous polynomials and the PDE
hereby descritized and split into finite number of algebraic equations. Thus, the aim is
to determine the unknown coefficient of these polynomials in such a way, that
distance (which is defined by the norm in a suitable vector space) from the exact
solution becomes a minimum. Therefore, the finite element is essentially a variational
minimization technique.
Since the number of elements is finite, we have reduced the problem of finding a
continuous solution to our PDEs to calculating the finite number of coefficients of the
polynomial. The solution of the Poissons equation, which is required to calculate the
magnetic vector potential, has to be solved for a given current density distribution.
We write the Poissons equation in more general form

In order to apply the finite element method, we have to find variational formulation.
The Galerkins methods leads to the week formulation of the problem: we multiply
Poisson equation by the test function u and integrate over the solution domain
10


( ) ( ) ( ) ( )
2
. . (3.2) u r v r dr f r v r dr
O O
V =
} }



Integration by parts gives


( ) ( ) ( ) ( ) ( ) ( )
(3.3) n u r v r dr u r v r dr f r v r dr
O I O
V V + V =
} } }


Where dr

denotes the surface normal on the boundary I. If appropriate boundary

conditions defines the value of u (Dirichlet boundary condition) or of its derivative
(Neumann boundary conditions) on the boundary, we can simplify (since v
vanishes, where Dirichlet boundary conditions apply

( ) ( ) ( ) ( ) ( )
(3.4)
n
n u r v r dr gv r dr f r v r dr
O I O
V V + =
} } }



The exact solution u
i
shall be approximated by a linear combination of trial functions

( ) ( )
0
(3.5)
n
h i i
i
u r u r
=
=

And we use a finite set of test functions

.
If we insert this expansion into the equation (4) and assume only Dirichlet
boundary conditions
( ) ( ) ( ) ( )
0
(3.6)
n
i i i i
i
u r v r dr f r v r dr
=
V V =

} }


We get a system of algebraic equations.
This can be solved with any standard method for the solution of a system of
algebraic equations, such as Gauss method, the Cholesky decomposition or iterative
scheme like the conjugate gradient method.

11

3.1.1 Boundary conditions
For the solution of partial differential equation like Maxwells equations, we
need boundary conditions to find a unique solution. There are three conditions


3.1.1.1 Dirichlet boundary conditions
The value of the solution is explicitly defined on the boundary (or on a part of
it). The value of solution can be zero or non-zero on the boundary. In
electromagnetics, the magnetic vector potential is usually set to zero along a
boundary, which should not be crossed by magnetic flux.
3.1.1.2 Neumann boundary conditions
The normal derivative of solution is defined on the boundary. If we set the
normal derivative of magnetic vector potential to zero, the boundary can be
interpreted can be interpreted as an interface with a highly permeable metal. Then, the
magnetic flux passes the interface at an angle of 90 degree to the plane of interface. In
order to find a unique solution, a Dirichlet condition must be defined somewhere on
the boundary of the domain.

Fig. 3.1: Problem Domain for boundary value problems
12

3.1.1.3 Mixed boundary condition
A combination of above two boundary conditions is called a mixed boundary
condition or also known as Cauchys boundary condition. In this case the normal
derivative of the solution and the value of the solution itself on the boundary are
connected by a function.
3.2 Application of FEM in electromagnetic field computation
The theory of electromagnetics took a long time to be established, and it can
be understood by the fact that the electromagnetics quantities are abstract or in
other words, cannot be seen or touched (such as mechanical and thermal
quantity). Electromagnetics can be described by the Maxwells equations and
constitutive equations. Actually, the majority of the electromagnetics phenomena
were established by other scientists before Maxwell, such as Ampere (1775-1836),
Gauss (1777-1855), Faraday (1791-1867), Lenz (1804-1865) among others. However,
there was some incompatibility on the formulation and Maxwell (1831-1879), by
introducing an additional term (displacement current) to Amperes law, could
synthesize the electromagnetics in four equations. The genius of this man brought the
electromagnetics to a very simple formalism, kept mainly by only four equations. The
consistency of these equations (along with constitutive ones) is so high that very
distinct phenomena (like microwaves and permanent magnet fields) can be precisely
described by these. While the formalism and the basic concepts of the
electromagnetics relatively simple, realistic problems can be very complicated and
difficult to solve. In fact, when, complicated geometries, non-linearity, many non-
static field sources, etc, appear together (or sometimes even alone) it s virtually
impossible to find analytical solutions for such problems and that is the main reason
why numerical methods have become widely used tools in electrical engineering
nowadays.
In this work, low-frequency phenomenas are of main interest. Electromagnetic can
be divided in different categories as shown in the fig. 3.2


13













3.2.1 Maxwells equations
The general, time dependent, Maxwells equations in differential form (also
called point or local form) are
(3.7)
t
c
V =
c
B
E

(3.8)
t
c
V = +
c
D
H J

(3.9) V = D

0 (3.10) V = B


Equation (3.7) is known as Faradays law of induction. This equation
describes the basic law of induction and, is often associated with the eddy current
applications. Equation (3.8) is the Amperes law in its general form, in which
displacement current t c c D is taken into account. Equation (3.9) is the Gausss law
which, sometimes associated with electrostatic application but which is also required
in other applications. Equation (3.10) is not associated with a particular law and
simply states the nonexistence of the isolated magnetic pole. It is also known as the
magnetic form of the Gausss law [13].
Electromagnetics (Maxwells equations)
Electromagnetics Low
Frequency
Electromagnetics High
Frequency (Waves)
Electrostatic Magnetics
Magnetostatic Magnetodynamics
Fig. 3.2: Classification of electromagnetics
study
14

The integral form of the Maxwells equation is given in following.

(3.11)
c s
dl ds
t
c
=
c
} }
B
E

(3.12)
c s
dl ds
t
c | |
= +
|
c
\ .
} }
D
H J

0 (3.13)
s
ds =
}
B

. (3.14)
s V
ds dv =
} }
D

The differential form of Maxwells equation will be extensively in this work. The
differential form is more convenient in calculations using methods such as finite
element method or finite difference method, while integral form is more convenient in
analytic calculation of fields and in various integral methods of numerical calculations
such as the method of moments and boundary element methods [14].
3.2.2 Constitutive relations
The two additional relations needed to complete the system of equations are the
material constitutive relations:

(3.15) = B H
(3.16) c = D E

These two vector equations are equivalent two six scalar equations

In addition, a constitutive relation involving current densities and the electric field
intensity
= E (3.17) o J
3.2.3 Potential functions
Potential functions are viewed as alternative representation of the
electromagnetic field. These are the simpler, more useful to describe the field
properties rather than to use an abstract field variable like magnetic flux density,
magnetic field intensity etc.
15

These may allow simplification in the computation of fields, better understanding of
results, which we obtained. Potential functions are defined based on properties of
magnetic fields.
3.2.3.1 Electric scalar potential
This potential function is defined based on the irrotational property of the
static electric field. That is
0 (3.18) V = E
Any function that is irrotational can be written as the gradient of a scalar function.

(3.19) V = V E T
his scalar function is called the electric scalar potential function.

2.2.3.2 Magnetic scalar potential
The electric scalar potential is defined based on the irrotational nature of the
electric field. Then it follows that when the magnetic field intensity is irrotational, it
can also be defined in terms of the magnetic scalar potential as

0 (3.20) V = H
(3.21) | = V H
| , is called the magnetic scalar potential.
3.2.3.3 Magnetic vector potential
Scalar potentials are defined based on the irrotational nature of the electric and
magnetic fields. Another type of potential function can also be defined based on the
nature of solenoidal nature of the electric and magnetic fields.
The magnetic flux density B is solenoidal in nature (i.e. 0 V = B ), it can be
derived as the curl of another vector
(3.22) =V B A
Here, A, is called the magnetic vector potential.

16

3.2.3.4 Electric vector potential
Another vector formulation employs the electric vector potential (T) to
describe the current distribution in conducting regions [14]. The electric vector
potential can be regarded as a variation for H as it has same curl relation but,
generally, a different divergence. The magnetic vector potential is defined as

(3.23) V = T J

Since the divergence of H is zero it can be redefined as


(3.24) = VO H T
Where, O, is a magnetic scalar potential.
3.2.4 Gauge conditions
The use of vector potential functions requires the specification of curl and
divergence of the potential functions for unique representation of field quantities. The
curl is normally specified based on the properties of the field (e.g. in the case of
magnetic vector potential the solenoidal nature of B allows the definition of A from
=V B A ). The divergence must then be specified to be consistent with the field
equations. There are normally two gauge conditions. One is the condition:

0 (3.25) V = A
is known as Coulomb gauge condition. The Second is
V
0 (3.26)
t
c
c
V + =
c
A

is known as the Lorentz gauge condition. This relation is consistent with the field
equations because it leads to continuity equation [15].
3.2.5 Field equations in terms of potential functions
The electromagnetic field equations can be written in terms potential
functions. The potential function formulation of field equations is particularly suited
the numerical methods. For this, the potential functions can be simply substituted in
17

the Maxwells equations. Here the magnetic vector potential formulation
(3.27)
t
c
V =
c
B
E
By the definition of magnetic vector potential i.e,
(3.28) =V B A

Substitute this value in the equation (3.27)
( )
(3.29)
t
c V
V =
c
A
E
And therefore, the field intensity can be written as
(3.30)
t
c
=
c
A
E
The electric field intensity is only the time dependent part of the electric field
intensity. In addition the static term given
V (3.31) = V E
Total field intensity is the sum of two electric field intensities
(3.32) V
t
c
= V
c
A
E
In the case of the magnetic vector potential
1
V (3.33)
t t
c

| | c c | |
V V = + V
| |
c c
\ .
\ .
S
A
A J
Where the relation = B Hand c = D E were used, there scalar form and the current
density J was in
S
J (applied or source current density) and
e
J (induced or eddy
current density). Assuming linear media, using vector identity
( ) ( )
2
(3.34) V V = V V V A A A
Equation (3.33) becomes
( )
2
2
2
V
(3.35)
t t
c c
| | c c | |
V V V = V
| |
c c
\ .
\ .
A
A A J
In the case of law frequency (static and quasi-static cases ), last term can be neglected
and if there is no electric scalar potential is present, then the resultant equation will
be
( )
2
(3.36) V V V = A A J

18

Now, using Coulomb gauge, i.e.
0 (3.37) V = A
The resultant equation will be partial differential equation in terms of magnetic vector
potential
2
(3.38) V = A J

this equation can be solved with knowledge of magnetic vector potential on the
boundary of the solution domain [16]. This is now basically a boundary value
problem.















19

Chapter 4
FEM Formulation for Electromagnetic fields and
Force computation
In coupled field analysis for magneto-structural problems, electromagnetic
forming and for other electromagnetic device, the analysis of different fields should
be carried out by an efficient method. The finite element method has been well
established in field computation. In this report, the finite element method is used for
different field analysis. The formulation for electromagnetic field by finite element
method is discussed in following section.
4.1 FEM Formulation for electromagnetic field computation
We start the formulation of magnetic field computation with the governing
partial differential equation in terms of Magnetic vector potential, which is,
2
(4.1) V = A J
Where, J is the total current, which can be decompose in two components
(4.2)
e s
J = J +J
J
e
, is the eddy current density and J
s
, is the source current density. Now above
equation can be written as,

2
( ) (4.3) V =
e s
A J +J

Eddy current density can be rewritten as given in equation (3.4) and (3.5)
(4.4) o =
e
J E
V (4.5)
t
c
= V
c
A
E

The above equation will be the result in equation (4.6)

2
V (4.6)
t

c | |
V = V
|
c
\ .
s
A
A J
20

The term, V V , in magneto-static and quasi-static cases can be assumed to be zero, so
the final equation
2
(4.7)
t

c
V = +
c
s
A
A J
2-D (two dimensional) equation for transient magnetic problems
2 2
2 2
1 1
(4.8)
d d
dx dy t
o

c
+ = +
c
s
A A
J
Let us consider, for an approximated solution A, the residual r is
2 2
2 2
1 1
(4.9)
d d
r
dx dy t
c
= + +
c
s
A A A
J
The weighted residual for element

(4.10)
e
i
N rdxdy
O
=
}}
e
i
R
Where, i=1, 2, 3

2 2
2 2
(4.11)
e e e
e e e
i i i
d d
N dxdy N dxdy N dxdy
dx dy t
o
O O O
| | c
+ +
|
c
\ .
}} }} }}
s
A A A
J
Since udv vdu uv = +
} }

Using 2-D Galerkins Finite element formulation, First integral left in Eqn. (4.11) is

1
. (4.12)
e
e e e
e e
e
e e e i i
i i i
dN dN d d
dxdy N n d N dxdy N dxdy
dx dx dy dy t
o

O O O I
| | c
= + I+
|
c
\ .
}} } }} }}
e
i s
A A A
R D J

Where

1
x y
x y
| | | | c c
= +
| |
c c
\ . \ .
A A
D and considering linear triangular interpolation
function
A=a + bx + cy (4.13)
( ) ( )
3
1
, , A (4.14)
e e
j j
j
x y N x y
=
=

e
A
21


( )
A
, , A 4.15
A
e
i
e e e e
i j k j
e
k
N N N
(
(
( =
(

(

e
A
Where
1
(4.16)
2
e e e e
i i i i
e e e e e
j i j j
e e e e
k k k k
N a b x c y
N a b x c y
N a b x c y
( ( + +
( (
= = + +
( (
( (
+ +

N

A
1
, , A (4.17)
2
A
i
e e e
i j k j
k
d
b b b
dx
(
(
( =
(
A
(

A


A
1
, , A (4.18)
2
A
i
e e e
i j k j
k
d
c c c
dy
(
(
( =
(
A
(

A


1
(4.19)
2
e
i
e
j
e
k
b
d
b
dx
b
(
(
=
(
A
(

e
N


1
(4.20)
2
e
i
e
j
e
k
c
d
c
dy
c
(
(
=
(
A
(

e
N

First integral on the left of equation (3.12) is

1
. (4.21)
e
e
e e
e
e i i
i
dN dN d d
dxdy N n d
dx dx dy dy
O I
| |
= + I
|
\ .
}} }
e
i
A A
R D

Line integral is set to zero, which gives natural boundary condition.
22


1
(4.22)
e
e e
i i
e
dN dN d d
dxdy
dx dx dy dy
O
| |
+
|
\ .
}}
A A

But
(4.23)
e
dxdy
O
= A
}}

Then by using the derivatives of shape function and interpolation functions, which
will be same in the isoparametric elements, which are very frequently used in the
electromagnetics problems,

2 2
2 2
2 2
A
1
A (4.24)
4
A
i i i i j i j i k i k
i j i j j j j k j k j
i k i k j k j k k k
k
b c bb c c bb c c
bb c c b c b b c c
bb c c b b c c b c

( | | + + +
| (
= + + +
| (
A
| (
+ + +
\ .

The second integral of equation becomes,

A
, , A (4.25)
A
e e
e
i
i
e
e e e e e e
j i j k j
e
k
k
N
N dxdy N N N N dxdy
t t
N
o o
O O
(
(
(
c c (
( =
(
(
c c
(
(


}} }}
e
A

Using formula (4.26), equation (4.25) becomes

( ) ( ) ( )
1 2 3

2 (4.26)
( 2)
e
l m n
e e e
l m n
N N N dxdy
l m n
O
= A
+ + +
}}
! ! !
!


A 2 1 1
1 2 1 A (4.27)
12
1 1 2
A
i
e
j
k
t
o
(
| |
c (
|
=
(
|
c
|
(
\ .

The forcing function (third term of equation (4.12) becomes)
23



1
= 1 (4.28)
3
1
e
e
i
N dxdy
O
(
A
(
(
(

}}
s
s
J
J
Element level of equations is,

2 2
2 2
2 2
A
1
A
4
A
i i i i j i j i k i k
i j i j j j j k j k j
i k i k j k j k k k
k
b c bb c c bb c c
bb c c b c b b c c
bb c c b b c c b c

( | | + + +
| (
= + + +
| (
A
| (
+ + +
\ .
-
A 2 1 1
1 2 1 A
12
1 1 2
A
i
e
j
k
t
o
(
| |
c (
|
(
|
c
|
(
\ .

=


1
1 (4.29)
3
1
(
A
(

(
(

s
J

Assembling all elements, global set of equation is given as,

| |{ } | | { } { } (4.30)
t
c
+ =
c
S A T A J

Now, next step is to incorporate the Dirichlets boundary conditions in the above
global equation and since, the Neumanns boundary conditions are implicit in the
formulation so only Dirichlets boundary conditions need special treatment [13].
After incorporation of boundary conditions, the resultant algebraic system of equation
can be solved by any standard method of solving linear algebraic equations.
There are generally three types of analysis in the electromagnetics, which are
described below as,
4.1.1 Magneto-static analysis
In the magneto static analysis, the quantities are not time dependent. This
analysis is generally required for the problems of steady flow of dc electric current,
permanent magnets, an applied external field, moving conductors etc. [14].



24

The governing equation in the magneto static case is given as,

2 2
2 2
1 1
(4.31)
d d
dx dy
+ =
s
A A
J
By using above described finite element formulation, it will result in the algebraic
equation of the form

| |{ } { } (4.32) = S A J
In this case the time-dependent effects (such as eddy currents) are not considered. It
can also be non linear in the presence of saturable, permanent magnet etc. The flux
density can be calculated as,

(4.33) =V B A

4.1.1.1 Algorithm for solving magneto-static problems using FEM
1. Prepare the mathematical model of the problems by considering symmetry of
the problems in terms of geometry and loading (if exist).
2. Discrete the geometry in terms of the linear elements or quadratic elements.
3. By using governing equation of the field as given in equation (3.31), compute
the elemental equations for each element by using weight functions and shape
functions (using Galerkins method or Variational methods) which will be
same in case of isoparametric elements.
4. Assemble all elemental equations in the global system of equation by using
connectivity information which is available from the discretization of the
geometry. It leads to algebraic system of equations.
5. The next step is incorporation of boundary conditions (Dirichlet boundary
condition) in the global system of equations. Neumann boundary conditions
are implicitly used in the formulation hence only Dirichlet boundary
conditions requires special treatment.
6. The resultant algebraic set of equations after incorporation of boundary
conditions is solved using any standard method for solving linear system of
25

equations. The solution of global system of equation gives the values of
unknown variable at each node of the problem domain.
7. The next step is the post processing of the solution. In this phase the secondary
(or derived) quantities can be derived from the solution. These quantities may
be flux density, intensity, magnetic forces etc.

The flow chart for solving magneto-static problems using FEM is given in Fig.
(4.1).














26

4.1.1.2 Flow chart for solving magneto-static problems using FEM



























Fig. 4.1: Flow Chart for Magneto Static Case

Prepare Mathematical Model of the Problem domain
Descritized the problem domain in terms of the finite elements (linear
or quadratic)
Formulate the elemental equations by using Galerkins or variational
formulation for each element of the descritized geometry
Assemble all the elemental equations in a global system of equation
using the connectivity information available in the discrete geometry
data and incorporate boundary conditions in the global system of
equation
Compute the Secondary or derived quantities using the solution of
global system of equation
Solve the resultant Global system of equation by using standard method for
solving linear system of equations
Start
End
27

A linear magneto-static case is analyzed using the algorithm discussed above. This
problem is solved by MATLAB and results obtained with MATLAB are Validated by
a commercial FEM software ANSYS.
4.1.1.3 Linear magneto-static case-study
Case study involves: Determination of the magnetic flux densities, force and
deformation in the conductor
4.1.1.3.1 Geometric data and material properties
Outer radius R = 50 mm,
Inner radius r = 40 mm
Height h = 90 mm and
Input current density J = 510
6
amp/m
2
;
Young Modulus E = 12.2 10
10
N/m
2
;
Poissons ratio = 0.3.




4.1.1.3.2 Assumptions
A long thick solenoid carries a uniform current distribution.
The turns of the solenoid can be modeled as a homogeneous, isotropic
material with modulus of elasticity E and Poissons ratio v.
The symmetries of the problem allow us to consider just a one-half or one-
quarter of the domain of analysis.
In this case study, the coil material is the nonmagnetic, only Lorentz force will
exist.

J
r
h
Fig. 4.2: Solenoidal coil
R
28

4.1.1.3.3 Problem model There are two analysis domains in this model.



4.1.1.3.4 Analysis domain of fields
1. EM analysis domain
Coil and air (both) regions
2. Structural analysis domain
The problem domain is descritized with linear triangular elements and by using
Galerkins formulation, the linear system of equations is obtained. The solution of
these equations gives magnetic vector potential at all the nodes in the problem
domain. The magnetic vector potential seems to be a bundle of information. Other
interested quantities like magnetic flux density, magnetic intensity, magnetic forces
etc.
The plot of magnetic flux lines are shown in fig. (4.4) and the linkage of this flux
with the coil are shown in the fig. (4.5). Local forces distribution in the coil is shown
in fig. (4.6).



Fig. 4.3: Problem model-1
1. Air region
2. Coil region






29

4.1.1.3.5 Magnetic flux lines plot with ANSYS
The flux lines due to current flow in the coil are shown in the fig. (4.4)



4.1.1.3.6 Flux linkage with coil
The flux produced by the current flowing in the coil is linked with the current
carrying part of the coil. This linkage of flux with the coil is shown in Fig. (4.5)


Fig. 4.4: 2-D flux lines plot with ANSYS
Fig. 4.5: Flux linkage with coil
30


4.1.1.3.7 Force distribution in the coil
The force distribution on the plate, due to flux linkage with the coil, is shown in the
Fig. (4.6).



4.1.1.4 Comparison of Results with ANSYS
The values of maximum flux density and maximum force density is compared
with the values obtained by an FEM software ANSYS shown in the table (4.1)
Table 4.1: Comparison of Results with ANSYS
Maximum flux density (Tesla) Maximum force density (N/m
2
)
ANSYS 0.0502 2.4918 10
5

MATLAB 0.0501 2.4866 10
5




Fig. 4.6: Force distribution in coil

31

4.1.2 Harmonic analysis
If one is the interesting in knowing the value of peak force in the first cycle,
the problem can be approximated as a time-harmonic one simplifying the
computational burden. In such a case the governing partial governing differential
equation (PDE) is:

2 2
2 2
1 1
j (4.34)
d d
dx dy
=o

+ = +
s
A A
J A

Where, A has only z-component in case of 2-D problem. The corresponding set of
equations to be solved after FEM discretization is:

| |{ } | |{ } { } (4.35)
r i r i r i
A jA j A jA J jJ e + + + = + S T

Separating

| | | |
| | | |
(4.36)
r r
i i
A J
A J
e
e
| | ( (
=
|
( (
\ .
S T
T S


In harmonic analysis, eddy currents can be calculated with the harmonic magnetic
field, which is in terms of real and imaginary terms. The calculated eddy current will
also be in terms of real and imaginary terms such as,


(4.37) jeo = e J A

e J and

A are the complex quantities.

Magnetic flux density will also be in the complex form, as


(4.38) = V B A


32


4.1.2.1 Algorithm for time-harmonic analysis
Algorithm for time-harmonic analysis is similar to magneto static analysis as
described in [4.1.1]. In comparison to magneto-static analysis, all the quantities in
time-harmonic analysis are in complex form. Due to presence of eddy currents in
time-harmonic analysis, there will be an additional term in the formulation for time-
harmonic analysis as given in equation (4.35).
A time-harmonic problem is analyzed with formulation and algorithm given above.
This problem is solved by MATLAB and results are verified with commercial FEM
software ANSYS.
4.1.2.2 Electromagnetic Analysis in Time-harmonic case
Involves analysis of determination of the magnetic flux densities, force


Table 4.2 Electrical parameters
Input current (J)

116260sin (Amp.) t e
Frequency ()
8761.2 (Hz)

Resistivity of workpiece ()
3.33310
-8
(-m)

Fig .4.7: Electromagnetic forming system
33


Table 4.3 Geometrical parameters of workpiece
Inner radius (r
w
) 6.5mm
Outer radius (R
w
) 7.5mm
Length (h
w
) 50mm
Material Aluminium

Table 4.4 Geometrical parameters of coil
Inner radius (r
c
) 9.25mm
Outer radius (Rc) 34.25mm
Length (hc) 21mm
Number of turns (N) 6
Material Copper

4.1.2.3 Assumptions
A long thick solenoid carries a uniform current distribution.
The turns of the solenoid can be modeled as a homogeneous, isotropic
material.
The symmetries of the problem allow us to consider just a one-half or one-
quarter of the domain of analysis.
In this case study, the coil material is the nonmagnetic, only Lorentz force will
exist.
Analysis is carried out by both formulations i.e. by complex and by time-
stepping for same problem.
4.1.2.4 Problem model
There are three analysis domains work piece, coil and air. All the three regions
are shown in fig. (4.8)
34




The problem domain is descritized with 8-node quadratic rectangular
isoperametric elements and by using Galerkins formulation, the linear system of
equations is obtained. The solution of these equations gives magnetic vector potential
in complex form at all the nodes in the problem domain. The magnetic vector
potential seems to be a bundle of information. Other interested quantities like
magnetic flux density, magnetic intensity, magnetic forces etc.
4.1.2.5 Simulation with complex Magnetic vector potential
formulation (MATLAB)
The simulation of a time-harmonic problem is carried out with complex vector
potential formulation by developing code in MATLAB. The main advantage of
complex vector potential formulation is that it requires computation in single step in
contrast of time-stepping formulation. The main precaution, which must be taken in
the time-harmonic analysis, is the consideration of the skin effect phenomena. To
consider this effect, the number of element in the one skin depth should be sufficient.
The plots of distribution of various field quantities in analysis domain are given in this
section.


Fig. 4.8: Problem model-2
35

4.1.2.5.1 Magnetic vector potential plot in work piece
The magnetic vector potential plot in the work piece is shown in fig. (4.9). It is
clear from figure that it has the maximum value at the midpoint of work piece.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.005
0.01
0.015
Element no.across hiegth of workpeice
m
a
g
n
e
t
i
c

v
e
c
t
o
r

p
o
t
e
n
t
i
a
l
Plot of magnetic vector potential in work peice


4.1.2.5.2 Axial flux density on the work piece
The axial flux density plot is shown in fig. (4.10). It has the maximum value at
the midpoint in the work piece.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
5
10
15
20
25
Element no. in the plate
A
x
i
a
l

f
l
u
x

d
e
n
s
i
t
y

i
n

p
l
a
t
e
Axial flux density in the plate



Fig. 4.9: Magnetic vector potential on the work piece
Fig. 4.10: Axial magnetic flux density on the work piece
36

4.1.2.5.3 Radial flux density on the work piece
The radial flux density is shown in fig. (4.11). It will be symmetric across the
height of work piece.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
-4
-3
-2
-1
0
1
2
3
4
Element no.across the workpeice
r
a
d
i
a
l

f
l
u
x

d
e
n
s
i
t
y

i
n

w
o
r
k

p
e
i
c
e


4.1.2.5.4 Eddy current Density on the work piece
As it is time-harmonic case, eddy current will exist in this case. The plot
of eddy current in the work piece is shown in Fig. (4.12).

0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.5
1
1.5
2
2.5
x 10
10
Element no. across the hieght of plate
E
d
d
y

c
u
r
r
e
n
t

d
e
n
s
i
t
y


Fig. 4.11: Radial flux density on the work piece
Fig. 4.12: Eddy current density in the work piece
37

4.1.2.5.5 Radial Force density on the work piece
The plot of radial force density is shown in fig. (4.13). It has maximum value
at the mid-height of the work piece.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
-0.5
0
0.5
1
1.5
2
2.5
x 10
11
Element no. on the workpeice
R
a
d
i
a
l

f
o
r
c
e

d
e
n
s
i
t
y

o
n

w
o
r
k

p
e
i
c
e
Radial force density plot on workpeice


4.1.2.5.6 Axial force density on the work piece
The plot axial force density across the height of work piece is shown in Fig.
(4.14). It shows symmetric characteristic across the height of work piece.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
-1.5
-1
-0.5
0
0.5
1
1.5
2
x 10
10
A
x
i
a
l

f
o
r
c
e

d
e
n
s
i
t
y

o
n

w
o
r
k

p
e
i
c
e
Element no. on the work peice
Axial force density plot on work peice

Fig. 4.13: Radial Force density on the work piece
Fig. 4.14: axial force density on the work piece
38

4.1.2.6 Simulation with complex Magnetic vector potential
formulation (ANSYS)
The case study of time-harmonic analysis analyzed above by developing own
MATLAB code, here the same case is investigated by using a commercial FEM
software ANSYS. The PLANE 53 elements are used to Descritise the problem
domain in ANSYS. These elements are 8-noded and rectangular in shape. The plots of
distribution of field quantities, obtained with ANSYS, are shown in this section.
4.1.2.6.1 2-D flux lines plot (ANSYS)
2-D lines plot of magnetic flux is shown in fig. (4.15). Flux linkage with the work
piece is also shown in same figure. The flux linkage with work piece is depending on
the input frequency due to skin effect phenomena.






Fig. 4.15: 2-D flux lines plot in electromagnetic domain
39

4.1.2.6.2 Flux density in the work piece
The plot of magnetic flux density in the work piece is shown in the fig. (4.16)



4.1.2.6.3 Nodal force distribution on the work piece
The plot of nodal force in the work piece is shown in fig. (4.17). It shows the
local force distribution on the work piece.



Fig. 4.16: Flux density on the work piece
Fig. 4.17: Nodal force distribution on the
work piece
40

4.1.2.7 Comparison of Results of Time-harmonic analysis with
ANSYS
The results obtained with MATLAB are compared with the ANSYS results, as
shown in table (4.5). The results are shown in good agreement.
Table 4.5: Comparison of Results of Time-harmonic analysis with ANSYS

Max.
Magnetic
vector
potential
(A)

Max. flux
density
(B)
2
(wb/m )

Max. Eddy
current
(J )
e

2
(Amp./m )

Max. force
density
(F)
2
(Newton/m )

Max.
Nodal
force
(F)
(Newton)

MATLAB

0.1496


22.1237


10
2.4063 10

11
2.1511 10

37.89


ANSYS

0.1498


22.759

10
2.4675 10

11
2.3523 10

38.233

4.1.3 Transient magnetic analysis
This type of analysis is required for time-dependent cases, such as short circuit
forces in transformer, pulsed eddy currents, power losses due surges in voltage or
current or pulsed external field.
In this analysis calculating magnetic fields vary with time. A transient analysis can be
either linear or non-linear.
The governing equation (two-dimensional) for transient magnetic problem is given in
equation (4.38), it is

1
(4.39)
t
o

c
A A = +
c
s
A
A J

Using FEM formulation, equation (4.39) will result in the set of algebraic equation as
given as following


| | ( ) { } | | ( ) { } ( ) { }
(4.40) t t t
t
c
+ =
c
S A T A J
41

Various times stepping schemes can be used to generate linear system from equation
(4.26), the linear system of equations to be solved at each time of step is:


| | | | | | | | | | | | | | | |
1 1 1 1
J J (4.41)
t t t t t t
t t
| |
| | | |
+A +A
( ( | | | |
+ = + +
( | | (
A
\ . \ .
S T A T S A

Where, S and T are the metrics which depends upon geometry and material
parameters. The value of the constant | determine whether the algorithm is of the
forward difference type (|=0) backward difference type (|=1) or some intermediate
type (0<|<1). If |=1/2, we have the crank Nicholson method [30]. The linear system
of equation can be solved by using methods such as LU factorization, Cholesky
decomposition, etc.
The eddy current can be calculated with time dependent-magnetic field

( ) ( )
( ) (4.42)
t t t
t t
o o
+ A
c
= =
c A
e
A A
A
J
Magnetic flux density can be calculated at different time instants but at a particular
instant it will be the function of space only.

( , , ) ( , , ) (4.43) x y t x y t =V B A

Time-stepping method has some advantages and disadvantages against complex
formulation are
The advantages are
(1) It can handle easily nonlinearities of the problem
(2) Transient problems can only be solved by this method
The disadvantages are
(1) A lot of computation work is required in this method.
(2) If this method is used for analysis for steady-state sinusoidal excitation, it will
introduce some transient effect in the output. Then it is necessary to take sufficient
time steps to disappear this transient effect.


42

4.1.3.1 Algorithm for transient magnetic analysis
1. The steps 1-5 are same as given in algorithm for magneto-static analysis.
2. After incorporation of boundary condition in step no. 5, the resultant equation (4.39)
is solved by Crank-Nicholson method.
3. Initialize t=0, A(t)=0, |=0.5, t=time-step and n= total no. of time steps
4. Solve equation (4.41) at time instant t+t and get A(t+t). By using equations (4.42)
and (4.43), calculate the values of J
e
and B.
5. Replace A(t) in equation (4.41) by A(t+t), solve for time instant t+2t and get the
values of J
e
and B.
6. Repeat step no 5 up to time instant t+nt.





















43

4.1.3.2 Flow chart for solving transient magnetic problems using
FEM

























Yes
START
Prepare mathematical model of the problem domain.
Descritise the problem domain in terms of the finite elements
(linear or quadratic)

Formulate the elemental equations by using Galerkins or variation
formulation for each element of the descritized geometry.
Assemble all the elemental equations in a global system of equation
using the connectivity information available in the discrete
geometry data and incorporate boundary conditions in the global
system of equation.
Initialize t=0, A(t)=0, |=0.5, t=time-step and n= total no. of
time steps
Solve equation (4.41) at time instant t+t and get A(t+t). By
using equations (4.42) and (4.43), calculate the values of J
e
, B.
Replace A(t) in equation (4.41) by A(t+t), solve for time
instant t+kt, where k=2,3,4,5.. and get the values of J
e
, B.

Is k=< n?
END
Fig. 4.18: Flow-chart for linear transient magneto-structural analysis
44

4.1.3.3 Simulation with Time-stepping (MATLAB)
The case study, investigated above with complex magnetic vector potential, is
analyzed with the time stepping formulation. This type of formulation is well suited
for transient problems such as pulsed eddy current problems, electromagnetic forming
etc.
4.1.3.3.1 Max. Magnetic vector potential in time domain
The plot of magnetic vector potential wave form is shown in fig. (4.19). The
excitation, in this case study, is steady state sinusoidal and some transient effects are
appeared as shown in the same figure. To solve a problem, with sinusoidal excitation,
number of time-steps should be sufficient up to disappear this transient phenomenon.

0 50 100 150 200 250 300 350 400 450
-0.015
-0.01
-0.005
0
0.005
0.01
0.015
0.02
0.025
time steps
M
a
x
,

m
a
g
n
e
t
i
c

v
e
c
t
o
r

p
o
t
e
n
t
i
a
l
Magnetic vector potential plot



Fig. 4.19: Max. Magnetic Vector potential in time
domain
45

4.1.3.3.2 Max. Axial flux density on the work piece
The plot of axial flux density wave form is shown in fig. (4.20)
0 50 100 150 200 250 300 350 400
-25
-20
-15
-10
-5
0
5
10
15
20
25
time steps
f
l
u
x

d
e
n
s
i
t
y

i
n

a
x
i
a
l

d
i
r
e
c
t
i
o
n
Axial flux density


4.1.3.3.3 Radial flux density in the work piece
The plot of radial flux density in the work piece is shown in fig. (4.21)
0 2000 4000 6000 8000 10000 12000
-4
-3
-2
-1
0
1
2
3
M
a
g
n
e
t
i
c

f
l
u
x

i
n

r
a
d
i
a
l

d
i
r
e
c
t
i
o
n

a
t

p
l
a
t
e
Hieght of plate
Radial flux density across plate hight




Fig. 4.20: Max. Axial flux density in time
domain
Fig. 4.21: Radial flux density across the workpiece
46

4.1.3.3.4 Eddy current density on the work piece
The plot of eddy current wave form is shown in the fig. (4.22)

0 50 100 150 200 250 300 350 400
-3
-2
-1
0
1
2
3
x 10
10
Time steps
M
a
x
.
e
d
d
y

c
y
r
r
e
n

d
e
n
s
i
t
y
Eddy current plot


4.1.3.4 Results of Time stepping formulation with MATLAB
The results of time stepping analysis with MATLAB are given in table 4.6. The
results given in table 4.6 can be compared with results for time-harmonic analysis for
the same case.
Table 4.6: Results of Time stepping formulation with MATLAB

Max.
Magnetic
vector
potential (A)

Max. flux
density
(B)

2
(wb/m )

Max. Eddy
current
(J )
e

2
(Amp./m )

Max. force
density
(F)
2
(Newton/m )

Max. Nodal
force
(F)
(Newton)

0.1499

22.7592


10
2.4876 10

11
2.5722 10

38.8933

Fig. 4.22: Max. Eddy current density on the work piece
47

4.2 FEM formulation for force computation
In the analysis of electromagnetic device, electromagnetic force is actually
couple with another field (mechanical field), the accuracy of analysis is strongly
depend upon force computation. There are some methods given in literature for
computation of force field by finite element method. The Lorentz force method is
used for computation of forces in conductive region. For magnetization forces, there
are some methods given in literature. Here in analysis of electromagnetic forming
problem by coupled field theory, the Lorentz forces are essential. Some important
methods, for Lorentz forces and magnetization forces computation, are discussed
following.
The various methods of force computation are-
1-Lorentz force method
2- Maxwell stress tensor method
3-Virtual work method
3-equivalent magnetizing current method
4-equivelent magnetic charges method
Maxwell stress tensor based method and its related aspects such as error of
performance are discussed in [17] and [18]. Theory and implementation of method of
equivalent magnetizing currents for calculation of electromagnetic force acting on a
ferromagnetic material is discussed in [19] and [20]. A comparative study of these
methods of force computation is given [21] and [22].
Practical implementation of virtual work method by determining the derivatives of the
coordinates of the nodes with the respect of the virtual displacement is dealt with [23].
In the virtual work approach, the force is derived at each node by direct differentiation
of the stored energy of the finite elements surrounding the node with respect to virtual
displacements. The magnetic flux distribution should be kept constant while
displacing the node [24].
4.2.1 Lorentz force method
In this method total force on a body is obtained by integrating the forces due
to magnetic field acting on each deferential current element,
48


( ) (4.44)
v
v
dv =
}
F J B

Where, JB

is the force density in the conducting region.
4.2.2 Maxwell stress tensor method
Maxwell stress tensor is widely used for electromagnetic force computations.
Maxwell stress tensor based method and it is related to aspects such as error
performances are discussed in [18] and [19]. A quantity called stress is defined in this
method whose divergence is actually the force density throughout the volume of the
body on which the force is to be determined. Applying the divergence theorem to the
stress tensor, we can consider Maxwell stress as surface force density which when
integrated over a surface enclosing the body gives total force acting on it. The surface
can be chosen so as to satisfy certain performance criterion and to improve the
accuracy of the result. The expression for stress tensor is,

2
1 1
B B B (4.45)
2
i j ij
o

| |
=
|
\ .
ij
T

Where , i j take values (x ,y, z).
ij
o is 1 if i=j, and zero otherwise. The same can be
written in term of force density vector as,

( )
1 1
n n (4.46)
2
ds

| |
=
|
\ .
}}
2
F B B B
Where, n is the normal unit vector to surface under consideration. The expression of
stress tensor in equation (4.36) results in to 33 tensor matrix,
( ) ( )
( ) ( )
( ) ( )
2
2
2
1 1 1 1
2
1 1 1 1
(4.47)
2
1 1 1 1
2
x x x y x z
xx xy xz
yx yy yz y x y y y z
zx zy zz
z x z y z z
B B B B B B B
T T T
T T T B B B B B B B
T T T
B B B B B B B



| |
| |

| |
\ .
|
| |
|
| | |
=
| | |
\ .
|
|
\ .
|
| |

|
|
\ .
\ .

49

In this method body on which the force is calculated, is surrounded by some
imaginary surface. After meshing the region, the edge of the some elements (faces of
the element in case of 3D) will lie on the chosen imaginary surface. Now these tiny
elemental surfaces can be resolved in to three components parallel to three coordinate
surfaces. The x, y and z component of the elemental surface are parallel to y-z, x-z
x-y plans respectively. Now T
xx
gives surface force density on x component (as i = x)
of the surface in x direction (as j=x), T
xy
gives surface force density on x- component
of the surface in y direction and so on. Similarly, T
yx
,T
yy
and T
yz
gives surface force
densities on y- component in x, y and z directions respectively.
T
zx
,T
zy
and T
zz
gives surface force densities on z-component of the surface in x, y and
z direction respectively. One care that needs to be taken while using this method is
that the above procedure gives force densities in positive -x or positive-y or positive-z
direction. But, we need to use surface force densities in the direction away from the
enclosed region. Hence all the x values have to be negated if the direction away from
the enclosed region is along negative coordinate axis. The above procedure equivalent
to evaluating the integral, given in Equation (4.45) over the enclosing the Maxwells
surface.
The Maxwell stress tensor gives both type of force, viz JB forces and magnetization
forces. If the Maxwell surface dose not enclosed any current caring part then it gives
the magnetization force acting on the enclosed body. If the surfaces enclose the
current caring part along with the same permeable material then it gives the addition
of the JB forces acting on the current caring part and the magnetization forces acting
on the permeable material [25].
4.2.3 Virtual work method
The method of virtual work for electromagnetic force calculation based on
generalized principal of virtual displacement. The movable part is assumed to be
displaced and the change of stored magnetic energy divided by the displacement gives
the force acting on the body, as the displacement tends to be infinitesimal. The
displacement is not the actual physical displacement of the body; hence, it is called
virtual displacement. One precaution has to be taken while virtually displacing the
50

body that the flux linkage has to be kept constant throughout the motion. The
implementation can be both at the level of the displacement of nodes. If the nodes are
displaced then the method is called local virtual work method. The expression for the
magnetic energy stored in the field is

0
1
(4.48)
2
B
V
d dV
(
=
(

} }
W H B
Where V is the volume of the field region, B is the flux density and H is the magnetic
flux intensity. The force on the node k, virtually displaced is given by,

(4.49)
q
c
=
c
k
W
F
4.2.3.1 Derivation and description of local virtual work method
In FEM, total domain is divided in to a number of small elements. Elements
used here are triangular in case of two dimensional analyses and tetrahedral in case of
three dimensional analyses. The total energy content of the domain will be the
addition of energy of its all elements. Hence

1
0
1
(4.50)
2
e
e
B
N
e
e
V
d dV
=
(
=
(
(

} }
e e
W H B
Where
e
V the volume of element e is,
e
B is the flux density in the element and
e
H is
the magnetic flux density of that element. N is the total number of elements. If we
want to find the total energy of a particular domain of our problem, the summation
above equation will run over all the elements in that particular domain [26]. If first
order models are used (triangular and tetrahedral are first order models), then the
equation is simplifies to,

2
1
(4.51)
2
N
e e
e
v V
=
=

B
W
Where
e
v is the reluctivity of element e. it can be expressed in term of flux density
squared,
51


2
(4.52)
e e
v v ( =

B
Now, differentiating the expression for the energy with respect to virtual displacement
q, we get,

2 2 2
2
2
1
(4.53)
2 2 2
M
e e e e
Kq e e
e
V V v V
v v
q B q q
=
( c c c c
= + +
(
c c c c

B B B
F
In the above expression, the summation is over M elements which surrounds node K.
Here q is the coordinate of the nodes of this particular element. For triangular
elements, q will successively takes the values
1 2 3 1 2 3
, , , , , x x x y y y which are the x and y
coordinates of local nodes 1, 2 and 3 respectively [24].







52

Chapter 5
Structural field analysis
5.1 Elasticity theory
The theories of elasticity and plasticity describe the mechanics of deformation
of engineering solids. Both theories, as applied to metals and alloys, are based on
experimental relations between stress and strain in a polycrystalline aggregate under
simple loading conditions [27].
Strains, as defined by theory of elasticity are given by (U is displacement along i and j
co-ordinates)

1

2
(5.1)
j
i
ij
j i
x x
c
c
= +
c c
| |
|
\ .
U
U

In the simplest form, the stress relationship between strain and stress is governed by a
linear equation known as Hooks Law.

( )
(5.2)
i i j k
v = + E
Where, E is Youngs Modulus, is the stress, and is the Poissons ratio. This
relation remains valid while the stresses are below the yield stress and deformation is
in the elastic limit. When, the load exceeds the yield stress, the coefficient of strain in
the above equation (E) no longer remains constant, and becomes a function of strain.
The model is then called a quasi-static (strain-rate independent) elastoplastic model.
Some materials exhibit a behavior where the stress is not just a function of strain but
also of strain-rate, i.e. the rate of deformation. Such behavior can be modeled using a
strain-rate dependent model or by using a state equation such as Johnson-Cook law
[28].
Irrespective of the model used to represent the behavior of metal, the laws of
equilibrium must be obeyed at all points of time [29].
53


0 (5.3)
ij
i
i
x

c
+
c
=

f

Where, f
i
are the body forces per unit mass, and is the mass density.
There are three equilibrium equations (equation 5.3), six stress-train relations, one
hydrostatic stress-strain relation and one yield criterion (Trescas or Von mises) in the
elasticity theory. With the help of above 11 equations, all the unknowns can be
simultaneously solved.
Using equations 5.1, 5.2 and 5.3, a second order partial differential equation is
obtained which relates the applied forces to the displacement of the body. Such an
equation is then numerically descritized and solved. In this report finite element
method is used for elastic analysis in the coupled field analysis of magnet-forming
system





Fig. 5.1: Stress-strain curve
54

5.2 Plasticity theory
In the case of plastic deformation, the strains are not proportional to stress as
in the case of elastic deformation. The stress may vary exponentially or by any other
function as shown in Fig. (5.1). The plastic flow rule gives the relation between
incremental stress and incremental strain is

( ) (5.4)
ij p
d d c =
ij
S

Where, d = proportionality constant, S
ij
= Deviatoric stresses and (d
ij
)
p
= Plastic
strain increments.


The proportionality constant d is a scalar quantity and known as plasticity
modulus. By differentiating equation (5.4) with respect to time, we get
(5.5)
ij
d
d
dt dt
c

=
ij
S
Or it can be written as
(5.6)
ij
c =
ij
S


Fig. 5.2: Relation between incremental stress and incremental strain
55

Equations (5.4) and (5.5) are referred to as differential and finite forms of Levy-
Mises equations. The constant

in equation (5.6) is the function of material

constant and strain rate. There are ten equations in the plastic range as given
following:
(1) One continuity equation: 0 =
ij

(2) Three equilibrium equation: 0
ij
i
x
c
c
=

(3) Six Levy-Mises equations: ( )
ij p
d d c =
ij
S
Using above ten equations and the ten equations can be solved simultaneously.
Invariably solving these 10 equations is not easy and several methods are therefore
resorted to. The method of solving to these 10 unknowns is referred as the theory of
plasticity [29].
5.3 Elastoplastic Analysis
The elastoplastic properties of a material can be described by its mass density,
Youngs modulus, Poissons ratio, Lames coefficient, yield strength, ultimate tensile
strength, rupture stress and maximum elongation or maximum strain. Besides, various
non-linear models of stress and strain relations in the plastic deformation use specific
coefficient to account for temperature rise and for the deformation velocity.
The elastoplastic properties of materials are mostly given as experimental data of the
material and depend on the working methods (annealing, heat treatment etc.). The
influence of the deformation velocity on the elastoplastic is taken into account via the
model of elastoplastic behavior.
The temperature rise has a certain influence on the elastoplastic properties of the of
materials, since each material has its own so-called critical temperature, obtained at
the intersection of curves of the yielding and fracture, but the influence of
temperature rise is not taken into account in the analysis reported here.
The pure elastic and pure plastic condition are described above. In some cases the
elastic and plastic deformations are of comparable magnitude, like example of temper
56

rolling of sheets. In such cases elastic deformation cannot be neglected. Prandtle
Reuss equations and Henckys equation describe the eleastoplastic deformation [29].
In the elastoplastic range, it is assumed that the total deformation is the sum of elastic
and plastic deformation, i. e.

( )
N
(5.7)
i
y y
i
u
y
c
c
=
c

(5.8) =
tot el pl
+

Differentiating with respect to time, gives

(5.9) = +
tot el pl



(5.10)
2
ij
G
= +
ij
ij
S
e S

And according to Hencky, in this region total strain is given by

1
( ) (5.11)
2G
| = = +
ij ij ij
e S
Here | is the comparable strain.
5.4 FEM Formulation for linear structural analysis
For modeling and numerical analysis of electromagnetic forming, finite
element method is extensively used worldwide due to its flexibility and ease for use
for various conditions. Finite element analysis is a numerical tool for piecewise
solution of partial or ordinary differential equations. In the simplest form of finite
element method, the entire domain is split in to a finite number of entities, called
elements defined by points in a space or area, called nodes. The differential equations
are reduced to simultaneous linear equations in each element and all the simultaneous
equation across all elements are solved at common nodes, to give approximation
solutions.
Shape functions, which interpolate a given quantity in terms of nodal values of that
quantity, are used to descritized the specific quantity, such as displacement. For a 2-D
case, using Einsteins notation (varies from 1to number of node in an element)
57

( ) (5.12)
x i x
i
u N u =
( )
(5.13)
y i y
i
u N u =
The normal strains and shear strain due to displacement at a given point are,
(5.14)
x
x
u
x
c
c
=
c


( ) (5.15)
i
x x
i
N
u
x
c
c
=
c

(5.16)
y
y
u
y
c
c
=
c

( )
(5.17)
i
y y
i
N
u
y
c
c
=
c

(5.18)
y
x
xy
u
u
x y

c
c
= +
c c

( ) ( )
(5.19)
i i
xy x y
i
i
N N
u u
x y

c c
= +
c c

Therefore the relation between the nodal displacements and strain at any point, are
related by the derivatives of the shape function in respective direction for a triangular
element the above relations in the matrix form are expressed as follows

3 1 2
3 1 2
3 1 1 2 2
0 0 0
0 0 0

x
x
y
y
xy
y
x
u
N N N
x
x x x
u
N N N
y y y y
u N N N N N
u
y x y x y
y x
c
c
c

c
c c c

c
c c c

c
c c c
= = =
` `
c c c c

)
c c c c c c
c
+
c c c c c
c c
)

1
1
2
2
3
3
3
(5.20)

x
y
x
y
x
y
u
u
u
u
N
u
x
u

(

(

(

(

`
(

(

( c

(
c


)


{ } | |{ } (5.21) = B U
The equation which relates the stress to strain or strain rate to any factors on which
the stress is dependent is called constitutive relation. For linear elastic behavior the
58

equation which relates stress to strain is given by the following equation,

{ } | | ( ) { } { } (5.22)
0 0
= D - +
Where
0
c and
0
o are initial strain and stress respectively.
For plane stress condition,

0 0 2
1 0
1 0 (5.23)
1
1
0 0
2
x x
y y
xy xy
E
o c v
o v c c o
v
v
t
| |
( |
( |
= +
( ` `
|

(
|

(
) ) |
\ .

Where E is the Youngs modulus and v is Poissons ratio. In the above equation


| |
2
1 0
1 0 (5.24)
1
1
0 0
2
E
v
v
v
v
| |
|
|
=
|

|
\ .
D


D is a matrix of constant values.
For plastic condition, that is, the behavior of material after yielding, the relation
depend on the strain hence [D] dose not remain constant. The constitutive relations
are mostly based on extensive experimentation of material behavior before and after
yielding. Many constitutive equations have suggested in recent years, some of them
generally applicable, while other for specific uses. One of the most extensively used
constitutive relations is modified Johnson cook law, which suggests that the stress is a
function of strain, strain rate of the temperature [30].
To obtain the integral statement and nodal equivalent forces, the principal of the
virtual work is used in [31]. Consider an infinitesimal virtual nodal displacement.
After the equating the external work done by applied surface tractions and point
forces, and internal work done against stress to causes further strain i,e, strain energy,
we get a statement relating external force (q) to internal body forces (b) and stress (
=[D] [B] {u}) due to deformation.
59


{ } | | { } | |{ } (5.25)
e e
T
e
e
d d
O O
= O O
} }
q B N b
This is valid for any stress-strain relationship. Using the linear stress strain
relationship in eqn. (5.23)

{ } | | { } { } (5.26)
e e e e
= + q K u f
Where the first term of R.H.S represents the forces developed due to displacement,
while the second term stand for the body forces and initial stress and strain if present.



| | | | | || | B D B (5.27)
e
T
e
e
d
O
= O
}
K


{ } | | { } | | | |{ } | | { }
0
(5.28)
e e e
T T
T
e e e
e
d d d c o
O O O
= O O O
} } }
e
f N b B D B
5.5 Quasi-static Elastoplastic Analysis
In quasi-static elastoplastic analysis of the work piece, the strain rate
hardening effects are neglected, and so are the temperature effects. The material is
assumed to be linear elastic till the effective stress or Von Mises stress reaches the
uniaxial yield stress of the material. Thereafter, the material is considered plastic and
linear hardening with the increase in strain. A tangent modulus is calculated at each
iteration defining the relation between Von-mises effective stress and strain
increment. The Von mises stress criterion is used for further yielding.
It is observed through the experimentation, that the volume of workpiece remains
constant while deforming plastically [28],
0 (5.29)
x y z
c c c + + =
Hence all further calculation is done in term of Deviatoric stresses and the normal
component of stresses is neglected, as it does not affect the plastic behavior. Isotropic
60

hardening is assumed, since the bauschinger effect is not significant. The elastic-
plastic stress strain relations, after yielding, are given by prandtl-reuss equations as
given by equation (4.4) in term of the plastic strain increments and Deviatoric stresses
[27]
( ) (5.30)
ij p
d d c =
ij
S
Where, d is the proportionality constant which changes in the course of iterations
and
ij
= Deviatoric stresses and ( )
ij p
dc =plastic strain increment
The algorithm can be summarized as [32]
1. Get the data defining geometry and material properties
2. Evaluate the equivalent nodal forces for applied body forces and surface
traction
3. Load increment loop begins
A. Increments the applied loads according to specified load factors
B. The inner loop begins
a. Calculate element stiffness for elastic and elastoplastic material
behavior
b. Solve the simultaneous equations. calculate if the elements is
-yielding the previous load increment and
-the load is increasing
-yielding in the current load increment
- Not yielding
c. Determine the flow vector and elastoplastic matrix
C. See if the solution coverages
If yes, exit the inner loop
If no, continue the loop
4. Continue the load increment loop until the load reaches final value.
61


5.6 High strain- Rate analysis
While in quasi-static analysis, the strain-rate hardening effects in material
behavior are ignored, in high strain-rate formulation, efforts are made to incorporate
this behavior is obtain a much more realistic simulation of the electroforming
phenomenon as mentioned in [34] introducing mass matrix in transient dynamic
analysis [35].

| |
{ }
| |
{ }
| |{ } { } (5.31) M U + C U + K U = F


And visco-plastic formulation by Perzyna model [36]

(5.32)
e p
= +


1 1- 2
(5.33)
2G E
ij kk ij
v
o = +
Where, G is shear modulus and E is the Youngs modulus.

( ) ( )
(5.34)
p
ij
F
F
v
c |
o
c
=
c

Where, F is the flow function, is the material constant and | (F) is a function which
can be experimentally obtain for each material. Again, the two main loops used in
quasi-static formulation are used here, but incorporating the inertia force and perzyna
model.
5.7 Coupling schemes
The procedure of coupled field analysis depends on which fields are being
coupled, but two distinct methods can be identified in the finite element analysis of
coupled field: sequential coupling (weak coupling) and direct coupling (strong
coupling).
62

5.7.1 Strong or Direct coupling scheme
In a direct coupling, the two field analyses are coupled in such a way that they
result in a single system of equations. The solution of this system gives all the degree
of freedoms. This is useful in the analysis of highly non-linear coupled fields. The
block diagram of direct or strong coupling is given in fig. (5.3).



5.7.2 Sequential or weak coupling scheme
In the sequential coupling, the two fields are solved in sequence and the result
of one analysis is given to subsequent field analysis as an input (fig.5.4). The
sequential coupling is useful in the analysis of linear fields. In the present report, the
sequential coupled model is used in magneto-structural analysis. Firstly, an
electromagnetic analysis is carried out, and subsequently a structural analysis is
performed. A fraction of analysis domain must be common for various fields in
coupled field analysis. Then the same geometrical information can be used in the
various analyses of fields in the coupled field problem. Such an analysis is a step
towards more complex non-linear analyses [15].
Fig. 5.3: Strong coupling scheme
63



In next section, simulation of magneto-structural coupled analysis is carried out
in the following way; firstly the simulation of electromagnetic part of coupled field is
presented for linear static case and then the structural part of coupled field problem is
simulated.
5.8 Coupled linear magneto-structural case
Case study involves: Determination of the magnetic flux densities, force and
deformation in the conductor
5.8.1 Geometric data and material properties
Outer radius R = 50 mm,
Inner radius r = 40 mm
Height h = 90 mm and
Input current density J = 510
6
amp/m
2
;
Young Modulus E = 12.2 10
10
N/m
2
;
Poissons ratio = 0.3.

Field No.1
Analysis
Geometry Information

Field No .2
Analysis
Analysis # 1
Analysis #2
J
R
r
h
Fig. 5.5: Solenoidal coil
Fig. 5.4: Weak coupling scheme
64

5.8.2 Assumptions
A long thick solenoid carries a uniform current distribution.
The turns of the solenoid can be modeled as a homogeneous, isotropic
material with modulus of elasticity E and Poissons ratio v.
The symmetries of the problem allow us to consider just a one-half or one-
quarter of the domain of analysis.
In this case study, the coil material is the nonmagnetic, only Lorentz force will
exist.
5.8.3 Problem model
There are two analysis domains in this model.


5.8.4 Analysis domain of fields
1. EM analysis domain
Coil and air (both) regions
2. Structural analysis domain
Coil region
5.8.5 Coupling scheme:
Sequential coupled model
5.8.6 Electromagnetic analysis
Fig. 5.6: Problem model-1
1. Air region
2. Coil region

coil region




65

5.8.6.1 Magnetic flux lines plot with ANSYS
The flux lines plot due to current flow in the coil is shown in fig. (5.7).


5.8.6.2 Flux linkage with coil
The flux linkage with current-carrying coil is shown in fig.(5.8)


5.8.6.3 Force distribution in the coil
Fig. 5.7: 2-D flux lines plot with ANSYS
Fig. 5.8: Flux linkage with coil
66

The local force distribution in the coil region is shown in fig. (5.9)



5.8.7 Structural analysis for Linear Elastic case
The linear structural analysis is carried out in the elastic region of the stress-strain
curve of the material. The input for structural analysis is the local force distribution
obtained from the electromagnetic part of this coupled problem. The analysis domain
for the structural analysis is only the coil region.
5.8.7.1 Deformation in the coil:
The resultant deformation in the coil region is shown in the fig. (5.10). The
deformation shown in this fig (5.10) is obtained with ANSYS.
Fig. 5.9: Force distribution in coil

67



5.8.7.2 Deformation plot across the height of the coil (MATLAB)
The plot of deformation across the height of plate is shown in fig. (5.11). The plot of
deformation is obtained with MATLAB.



Fig. 5.11: Plot of deformation across the height of coil
Fig .5.10: Deformation in the coil
68


5.8.8 Results of coupled field analysis
The results obtained with the MATLAB are compared with results obtained with
ANSYS.
Table 5.1: Comparison of Results of Coupled magneto-structural analysis with
ANSYS
Maximum flux
density (Tesla)
Maximum force
density (N/m
2
)
Maximum
deformation (m)
ANSYS 0.0502 2.4918 10
5
1.6623 10
-9

MATLAB 0.0501 2.4866 10
5
1.6251 10
-9

Ref. [38] 0.04877 1.946 10
5
1.923 10
-9








69

Chapter 6
Conclusion and Future work
6.1 Summary of work done
There are two cases have been analyzed in this report.
1. In first case, the sequential coupled linear magneto-structural analysis has
been done by developing MATLAB code. The analysis electromagnetic part
of the coupled magneto-structural analysis is investigated for the static
electromagnetic analysis due static excitation in the solenoidal coil. The
results obtained by developing MATLAB code are verified by the results
obtained with a commercial FEM software ANSYS. This analysis is carried
out further for coupled magneto-structural analysis in which the subsequent
structural analysis is coupled with electromagnetic analysis by means of
electromagnetic forces. These electromagnetic forces are applied as input in
the subsequent structural analysis which gives the final deformations, stresses
in the solenoidal coil. The coupled magneto-structural analysis is also carried
out by developing MATLAB code. These results are validated with
commercial FEM package ANSYS.
2. In the second case, various formulations for force computation in the
conductive region have been done. In first, it is analyzed by complex magnetic
vector formulation and in second, same case analyzed by time-stepping
formulation. The time-stepping formulation is particularly suited for transient
and non-linear problems. The time-stepping formulation is used for analyzed
for the same problem as that for in the complex vector potential formulation.
70

Therefore, the results obtained with time-stepping formulation may be
compared with the results obtained with complex vector potential formulation.
All results are validated by commercial FEM package ANSYS.
6.2 Summary of Results
1. In first case, the coupled linear magneto-structural analysis describes the
deformation in the conductor. Results obtained by MATLAB code are in good
agreement with commercial FEM package ANSYS. This analysis may be
useful in the analysis of deformation and vibration problem in electrical
machines. The analysis can be used to explain the deformation in the rotor
conductor in electrical machines. The result obtained with coupled magneto-
structural analysis shows the deformation in the current carrying conductor. In
the case described above can also be used in the calculation of deformation in
the winding conductor under the short-circuit conditions by using the peak
value of the short-circuit current as static excitation.
2. The accuracy of the coupled field analysis is strongly depends on the accuracy
of computation of the electromagnetic forces. The deformation and stress field
patterns in the metal-forming are depend also on the patterns of applied forces
on the metal object. In the second case, the force computation by various
formulations is analyzed. These formulations can be useful for the analysis for
steady-state excitation (sinusoidal input), particularly time-harmonic or
complex magnetic vector potential formulation, and time-stepping formulation
is particularly for transient excitation. Force computation with time-stepping
formulation is the first step for electromagnetic forming case.
71

6.3 Future work
The linear coupled magneto-structural analysis may be useful in the analysis of
deformation and vibration analysis in electrical machines. This analysis can extend to
nonlinear structural analysis, which can be much close to real world values. The
thermal field can also be taken into account to approaching accurate analysis. The
coupled analysis can extend for magnetic materials by computing force with one of
the magnetization force technique. In case of ferromagnetic materials, there is only
modification that can be made in the analysis presented above is in force computation
technique.










72


References
1. M.V.K. Chari, and P. Silvester, Analysis of Turboalternator Magnetic Fields
by Finite Elements, IEEE Transaction on Power Apparatus and Systems,
Vol. 90, no. 2, pp. 454-464, 1971
2. M.V.K. Chari, and P. Silvester, Finite Element Analysis of Magnetically
Saturated D-C Machines, IEEE Transaction on Power Apparatus and
Systems, Vol. 90, no. 5, pp. 2362-2372, 1971
3. S.C. Tandon, A.F. Armor, and M.V.K. Chari, Nonlinear Transient Finite
Element Field Computation for Electrical machines and Devices, IEEE
Transaction on Power Apparatus and Systems, Vol. 102, no. 5, pp. 1089-
1096, 1971
4. M.V.K. Chari, Three-dimensional Vector Potential Analysis for Machine
Field Problem, IEEE Transaction on Magnetics, Vol. 18, no. 2, pp.436-446,
1982.
5. M. Ito, Analytical Model for Magnetic Field Analysis of Induction Motor
Performance, IEEE Transaction on Power Apparatus and Systems, Vol. 100,
no. 11, pp.4582-4590, 1980.
6. F. Bouillault, and A. Rajek, Dynamic Model for Eddy Current Calculation in
Saturated Electrical Machines, IEEE Transaction on Magnetics, Vol. 19, no.
6, pp.2639-2642, 1983.
73

7. B. Brunelli, Transient and Steady-State behavior of Solid Rotor Induction
Machines, IEEE Transaction on Magnetics, Vol. 19, no. 6, pp.2650-2654,
1983.
8. D. Shen, Solution of Magnetic Fields and Electrical Circuits Combined
problems, IEEE Transaction on Magnetics, Vol. 21, no. 6, pp. 2288-2291,
1985.
9. B. Bendjima and M. Feliachi, Finite Element Analysis of Transient
Phenomena in Electromagnetic Forming System, IEE Computation in
Electromagnetics, Vol. pp. 113-114, April 1996.
10. B. Bendjima and M. Feliachi, A Coupling Model for Analyzing Dynamical
Behaviors of an Electromagnetic Forming System, IEEE Transactions on
Electromagnetics, Vol. 33, no. 2, pp. 1638-1641, March 1997.
11. F. Azzouz, B. Bendjima, and M. Feliachi, Applicatin of Macro-Element and
Finite Element Coupling for the Behaviors of Analysis of Magneto Forming
system, IEEE Transactions on Magnetics, Vol. 35, No. 3, pp. 1845-1848,
May 1999.
12. A. Meriched and M. Feliachi, Electromagnetic forming modeling of thin
sheets, IEEE Transactions on Magnetics, Vol. 36, no. 4, pp. 852-853, July
2000.
13. J. P. A. Bastos and N. Sadowski, Electromagnetic modelling by finite element
method, Marcel Dekker, New York, 2003.
74

14. M. N. O. Sadiku, Numerical Techniques in Electromagnetics, CRC Press,
Florida 1992
15. M. V. K. Chari and S. J. Salon, Numerical Methods in Electromagnetism,
Academic Press, San Diego 2000.
16. N. ida, Numerical Modeling for Electromagnetic Non-Destructive Evaluation,
Chapman & Hall, London 1995.
17. AN. Wignall, A. J. Gilbert, and S. J. Yang, Calculation of Forces on
Magnetised Ferrous Cores using the Maxwell Stress Tensor Method, IEEE
Transactions on Magnetics, Vol. 24, no. 1, pp. 459-462, 1988.
18. Z. W. Shi and C. B. Rajanathan, A Method of Reducing Errors in the
Calculation of Electromagnetic Forces by Maxwell Stress Summation
Method, IEEE Transactions on Magnetics, Vol. 35, no. 3, pp. 1358-1369,
May 1999.
19. T. Kabashima, A. Kawahara, and T. Goto, Force Calculation using
Magnetizing Current, IEEE Transactions on Magnetics, Vol. 24, no. 1, pp.
451-454, Jan. 1988.
20. G. Henneberger, P. K. Sattler, and D. Shen, Nature of the Equivalant
Magnetising Current for the Force Calculation, IEEE Transactions on
Magnetics, Vol. 28, no. 2, pp. 1068-1071, Mar. 1992.
21. Z. Ren, Comparison of Different Force Calculation Methods in 3-D Finite
Element Modeling, IEEE Transactions on Magnetics, Vol. 30, no. 5, pp.
3471-3474, 1994
75

22. G. Reyne, J. C. Sabonnadiere, J. L. Coulomb, and P. Brissonneau, A Survey
of the Main Aspects of the Magnetic Forces and Mechanical Behaviour of the
Ferromagnetic Materials under Magnetisation, IEEE Transactions on
Magnetics, Vol. 23, no. 5, pp. 3765-3767, Sept. 1987.
23. J. L. Coulomb and G. Meunier, Finite Element Implementation of Virtual
Work Principle for Magnetic or Electric Force and Torque Computation,
IEEE Transactions on Magnetics, Vol. 20, no. 5, pp. 1894-1896, Sept. 1984.
24. A. Benhama, A. C. Williamson, and A. B. J. Reece, Virtual Work Approach
to the Computation of Magnetic Force Distribution from Finite Element Field
Solution, IEE Proc.-Electr. Power Appl., Vol. 147, no. 6, pp. 437-442, Nov.
2000.
25. A.V. Kank and S. V. Kulkarni, Analysis of Lorentz and Saliency Forces in
Rotating Machines, Proceeding of COMSOL Users Conference, Bangalore,
pp. 18-21, Nov. 2006.
26. J.L. Coulomb, A Methodology for Determination of Global
Electromechanical Quantities from a Finite Element Analysis and Its
Application to the Evaluation of Magnetic Forces, Torques, and Stiffness,
IEEE Transactions on Magnetics, Vol. 20, no. 5, pp. 2514-2519, 1983.
27. W Jonhson and P. B. Mellor, Plasticity for Mechanical Engineers, D. Van
Nostrand Company, London, 1962.
28. R. Hill, The Mathematical Theory of Plasticity, Oxford University Press, New
York, 1950.
76

29. V. Gopinathan, Plasticity Theory and Its Applications in Metal Forming, John
Wiley and Sons, New York, 1982.
30. C. Karch and K. Roll, Transient Simulation of Electromagnetic Forming of
Aluminium tubes, Advanced Materials Research, Vol.6, no.8, pp. 639-648,
2005.
31. O. C. Zienkewicz, R. Taylor, and J. Z. Zhu, The finite element method, its
basis and fundamentals, Butterworth-Heinemann, Elsevier, New Delhi, 2000.
32. D. R. J. Owen and E. Hinton, Finite elements in Plasticity, Pineridge Press
Limited, Swansea, 1950.
33. A. G. Mamalis, D. E. Manolakos, A. G. Kladas, and A. K. Koumoutos,
Physical Principles of Electromagnetic Forming Process: A Constitutive
Finite Element Model, Journal of Material Processing Technology, Vol.
161, pp. 294-299, 2005
34. M. Kleiner, A. Brosius, H. Blum, and M. Steimer, Benchmark Simulation for
Coupled Electromagnetic Mechanical Metal Forming Processes, Journal of
Materials Processing Technology, Vol. 177, no. 1-3, pp. 270-273, 2006.
35. K. J. Bath, Finite Element Procedures, Prentice Hall, New Jersey, 1996.
36. P. Perzyna, Fundamental Problem in Viscoplasticity, Advances in applied
mechanics, Vol. 9, 1966.
37. L. Hirsinger and R. Billardon, Magneto-Elastic Finite Element Analysis
including Magnetic Forces and Magnetostriction Effects, IEEE Transactions
on Magnetics, Vol. 31, no. 3, pp. 1877-1880, May 1995.
77

38. E. Euxibie, J. L. Coulomb, G. Meunier, and J. C. Sabonnadiere, Mechanical
Deformation of a Conductor under the Electromagnetic Stresses, IEEE
Tranactions on Magnetics, Vol. 22, no. 5, pp. 828-830, 1986.