?Gerrit Mur, retired, but not
quite
Summary
I
look back at the
discussion about edge elements for electromagnetics. Not much progress
seems to
have been made since I published ‘the Fallacy’ and the points against
edge
elements I made in that paper remain to be refuted. In this paper I
shortly
describe the situation and propose a new finite-element method for
computing
electromagnetic fields. Unfortunately I do not have the means any more
to put
them to the test.
1.
What were the problems?
a.
The problems with edge elements
In
my paper ‘The fallacy of
edge elements’ I argued against the use of edge elements for computing
electromagnetic fields. The main points against edge elements are
listed below,
for a full discussion I refer to my paper.
1. Contrary to what is stated by the advocates of edge elements these elements do allow spurious solutions, they even do so by definition.
2. The use of edge elements is very inefficient because of generating a unnecessarily large number of unknowns.
3. The
representation of the field by edge elements has a poor condition as
compared with Cartesian bases.
The
good thing about
edge-element methods is that, despite the fact that edge elements do
allow
spurious solutions, methods using them do not seem to suffer very much
from it
in practice and it is inefficiency that is their
main
disadvantage.
a.1
Why do edge elements not usually produce wrong results?
Edge
elements do allow
spurious,
i.e. erroneous, solutions, there is no doubt about that, but methods
using
them do
not usually produce errors, how can we understand that? The
answer
to this question is found by studying the properties of 1: the
resulting
system of
equations together with 2: the solution vector and 3: the method of
solution.
First
we note that the system of
equations is going to be solved iteratively. (The reason for this is,
of
course, that only iterative methods are practical for practical
problems, but
remember that the matrix is highly singular because of which a direct
solution
would be impossible anyhow.) In an iterative solution we have a
starting
solution vector and this vector is iteratively updated until a
sufficiently
accurate result is obtained. About this process there are two
observations we
can make here.
1. For the initial solution all
elements of the solution vector
are set to zero. Note that this initial ‘solution’ is 100% incorrect
but
because all
coefficients in the expansion of the field are zero the normal
continuity
conditions across the interfaces between elements are satisfied exactly
and no
error is made in the normal fluxes. It are the normal fluxes where edge
elements freely allow errors.
2. The iterative solution proceeds by
manipulating the solution
vector using the system matrix and possibly, for speeding up the
process, an approximate
inverse of this matrix. Neither the system matrix, nor its approximate
inverse,
contains any information about the normal continuity and consequently
the
errors in the normal fluxes remain unaltered i.e. zero.
b.
The problem with nodal elements
The
well known nodal
elements, I prefer to refer to them as Cartesian elements, are very
efficient
and provide an optimum condition of the representation of the field.
The
problem with Cartesian elements is, however, that they have three
unknowns (the
three components of either the electric or magnetic field
vector) at each
node whereas four conditions that apply to the simplicial star of this
node (3
for the curl and 1 for the divergence) are imposed upon them. Hence we
have
more independent equations than unknowns and problems must be expected,
and they are
in fact found.
2.
What are the options?
Thinking about our
problem we can start
from two sides:
1. Use
edge elements and try to
remove the problems with normal fluxes by adding additional equations.
In this
way we may be able to exclude the possibility of spurious solutions but
we add
to the complexity of the large system matrix, i.e. we increase the
already
existing
inefficiency of the method. This must be a dead end.
2.
Use
Cartesian elements and add,
at each node, one degree of freedom making the number of conditions
applying at
each node (4) equal to the number of unknowns (4). The question then
becomes: ’How to add a meaningful additional degree of freedom
at each
node?’. It is that question that we want to answer below.
When
deciding how to add an
additional degree of freedom at a node we first note that the curl
equations
already provide sufficient information for computing the 3-component
vector. So
it is the modelling of the divergence that we need additional solution space
for.
A hint to
solve our problem is
given by edge elements that do not have this type of problem, with edge
elements we have a much to
wide solution space since the normal fluxes across
interfaces are left free
to jump.
3. My proposal
Consider,
for
simplicity of the argument, a node in a regular rectangular
3-dimensional mesh
that consists of rectangular bricks, each of these bricks being
subdivided in
tetrahedra. Let us, again for simplicity of the argument,
assume that the
origin of our coordinate system coincides with this node and that we
aim at
computing the electric field strength E using
first order, i.e. linear, elements and that the medium is locally
homogeneous. We now know that we have 4 equations for the 3 unknowns, Ex,
Ey and Ez,
at the centre of this simplicial
star.
For
solving our shortage of
unknowns we now propose to widen our solution space by splitting
one of them, e.g. Ex
but this choice is
arbitrary, in two unknowns viz Ex+,
applying to the
part of the simplicial star above the plane x=0,
and Ex-
applying to the part below this plane. In this way we allow
for a
discontinuity in Ex across
the plane x=0, and
consequently a possible divergence across this ‘interface’. Our fourth
equation, modelling the divergence, will now ensure that the divergence
in
the
simplicial star of this node is set to 0 (or to any
other value
required, an option not available with edge elements). In case we omit
this
fourth
equation we obtain a situation comparable with not prescribing normal
flux conditions at edge-element interfaces and the difference between
the
final
results Ex+
and Ex- will
most probably be 0.
Applying this splitting to
each node of the
mesh we have four independent equations and four unknowns at each node
throughout the domain of computation and reliable results can be
obtained
with an optimum efficiency. With linear elements we have 4, or perhaps
6,
unknowns per node, with equivalent linear edge elements we would have
tens of
unknowns related to each node, the exact number depending on the mesh.
4. Some additional notes
1. Note
that the choice for
splitting up Ex is
arbitrary. Each of the other components
could have been chosen with the same result.
2. For
symmetry reasons one could
opt for splitting up all three components of the electric field
strength. This
gives rise to two more unknowns but we can solve this symmetrically by
adding
the two equations Ex+
- Ex-
= Ey+
- Ey-
=Ez+
- Ez-.
This choice may be useful for the modelling
of the field near corners and points where interfaces intersect.
3. It is easily understood that our proposal can straightforwardly be generalized to the application at general meshes.
4. At
arbitrarily oriented interfaces between different
media the Cartesian basis of the representation of the field may
locally be
rotated so as to easily accommodate the local continuity conditions.
Another
option is to use edge elements along (intersecting) interfaces, but
only there, like
we did in
the FEMAX packages. In this way we combine the best things of both
worlds, simplicity and efficiency in homogeneous regions and flexbility
along interfaces and their intersections. Note that the
number of
elements related to interfaces (2D) is, for practical problems, much
smaller than the number of elements in the homogeneous subdomains (3D)
because of which the efficiency of our method is not seriously reduced
by the introduction of edge elements along interfaces.
5. Final conclusion
We have proposed a new and
very efficient way
of representing the field in a finite element method for the direct
computation
of electric or magnetic field strength distributions i.e. without the
use of
methods like vector potentials and/or methods based on the exclusive
use of edge elements.
On a personal level:
When my health forced me to retire I was not anymore able to deal with the challenges I saw in the finite element modelling of electromagnetic fields. With the above thoughts I have closed the book.
