BELA
Finite Element
Electrostatics Solver
Version 1.0
User’s Manual
March 29, 2003
David Meeker
http://femm.berlios.de/bela.htm
©2003
Many thanks to those
people who offered feedback on the beta version of Bela:
· Mark Smith
BELA is a suite of
programs for solving linear electrostatic problems. The programs currently
address problems on twodimensional planar and axisymmetric domains. BELA is
divided into three parts:
Two additional
programs are also called to perform specialized tasks. These are:
The Lua scripting
language (see http://www.lua.org) is also
integrated into the pre and postprocessors. Lua allow “batch” runs to be
performed without user interaction. In addition, all edit boxes in the user
interface are parsed by Lua, allowing equations or mathematical expressions to
be entered into any edit box in lieu of a numerical value. In any edit box in
BELA, a selected piece of text can be evaluated by Lua via a selection on the
right mouse button click menu.
The purpose of this
document is to give a brief explanation of the kind of problems solved by BELA
and to provide a fairly detailed documentation of the programs’ use.
The goal of this
section is to give the user a brief description of the problems that BELA
solves. This information is not really crucial if you are not particularly
interested in the approach that BELA takes to formulating the problems. You can
skip most of Overview, but take a look at Section 2.2. This section
contains some important pointers about assigning enough boundary conditions to
get a solvable problem.
Some familiarity with
magnetism and Maxwell’s equations is assumed, since a review of this material
is beyond the scope of this manual. However, the author has found several
references that have proved useful in understanding the derivation and solution
of Maxwell’s equations in various situations. A very good introductorylevel
text is Plonus’s Applied electromagnetics
[1]. A good intermediatelevel review of Maxwell’s equations, as well as a
useful analogy of electromagnetism to similar problems in other disciplines is
contained in Hoole’s Computeraided analysis and design of electromagnetic
devices [2]. For an advanced treatment, the reader has no recourse but to refer
to Jackson’s Classical electrodynamics [3].
For the electrostatics
problems addressed by BELA, only a subset of Maxwell’s equations are
required. These problems describe the
behavior of electric field intensity, E,
and electric flux density (alternatively electric displacement), D.
There are two conditions that these quantities must obey. The first condition is the integral form of
Gauss’ Law, which says that the flux out of any closed volume is equal to the
charge contained within the volume:
_{} (1)
where r represents charge density. The second is the integral form of Ampere’s
loop law:
_{} (2)
Displacement and field
intensity are also related to one another via the constitutive relationship:
_{} (3)
where e is the electrical permittivity. Although some electrostatics problems might
have a nonlinear constitutive relationship between D and E, BELA only
considers linear problems.
To simplify the
computation of fields which satisfy these conditions, BELA employs the electric
scalar potential, V, defined by its
relation to E as:
_{} (4)
Because of the vector
identity _{} for any scalar y, Ampere’s loop law is automatically
satisfied. Substituting into Gauss’ Law
and applying the constitutive relationship yields the secondorder partial
differential equation:
_{} (5)
which applies over
regions of homogeneous e. This is the equation that BELA solves—the
program solves for voltage V over a
userdefined domain with userdefined sources and boundary conditions.
Some discussion of
boundary conditions is necessary so that the user will be sure to define an adequate
number of boundary conditions to guarantee a unique solution. Boundary
conditions for BELA come in three flavors:
_{} (6)
This boundary condition is most often used by BELA as a computationally
cheap way of emulating unbounded solution domain.
If no boundary
conditions are explicitly defined, each boundary defaults to a _{} Neumann boundary
condition. However, a nonderivative boundary condition must be defined
somewhere so that the problem has a unique solution.
Although the
differential equations that describe V
appear relatively compact, it is very difficult to get closedform solutions
for all but the simplest geometries. That’s where finite element analysis comes
in. The idea of finite elements is to break the problem down into a large
number regions, each with a simple geometry (e.g. triangles). For example,
Figure 1 shows a map of the Massachusetts broken down into triangles.
Figure 1.
Triangulation of Massachusetts
Over these simple
regions, the “true” solution for V is
approximated by a very simple function. If enough small regions are used, the
approximate V closely matches the
exact V.
The advantage of
breaking the domain down into a number of small elements is that the problem
becomes transformed from a small but difficult to solve problem into a big but
relatively easy to solve problem. Specifically, triangulation of the problem
results in a linear algebra problem with perhaps tens of thousands of unknowns.
However, techniques exist that allow the computer to solve for all the unknowns
in only seconds.
BELA uses triangular
elements. Over each element, the solution is approximated by a linear
interpolation of the values of V at
the three vertices of the triangle. The linear algebra problem is formed by
choosing V on the basis of minimizing
the total energy of the problem.
The preprocessor is
used for drawing the problems geometry, defining materials, and defining
boundary conditions. Drawing a valid
geometry usually consists of four (though not necessarily sequential) tasks:
·
Drawing the
endpoints of the lines and arc segments that make up a drawing.
·
Connecting the
endpoints with either line segments or arc segments
·
Adding “Block
Label” markers into each section of the model to define material properties and
mesh sizing for each section.
·
Specifying
boundary conditions on the outer edges of the geometry.
This section will
describe exactly how one goes about performing these tasks and creating a
problem that can be solved.
The key to using the
preprocessor is that the preprocessor is always in one of five modes: the Point
mode, the Segment mode, Arc Segment mode, the Block mode, or the Group mode.
The first four of these modes correspond to the four types of entities that
define the problems geometry: nodes that define all corners in the solution
geometry, line segments and arc segments that connect the nodes to form
boundaries and interfaces, and block labels that denote what material
properties and mesh size are associated with each solution region. When the
preprocessor is in a one of the first four drawing modes, editing operations
take place only upon the selected type of entity. The fifth mode, the group
mode, is meant to glue different objects together into parts so that entire
parts can be manipulated more easily.
One can switch between
drawing modes by clicking the appropriate button on the Drawing Mode potion of
the toolbar. This section of the toolbar is pictured in Figure 2.
Figure 2.
Drawing Mode toolbar buttons.
The buttons correspond
to Point, Line Segment, Arc Segment, Block Label, and Group modes respectively.
The default drawing mode when the program begins is the Point mode.
Although most of the tasks that need to be
performed are available via the toolbar, some important functions are invoked
only through the use of “hot” keys. A summary of these keys and their
associated functions is contained in Table 1.
Point
Mode Keys 

Key 
Function 
Space 
Edit the
properties of selected point(s) 
Tab 
Display
dialog for the numerical entry of coordinates for a new point 
Escape 
Unselect
all points 
Delete 
Delete
selected points 
Line/Arc
Segment Mode Keys 

Key 
Function 
Space 
Edit the
properties of selected segment(s) 
Escape 
Unselect
all segments and line starting points 
Delete 
Delete
selected segment(s) 
Block Label Mode Keys 

Key 
Function 
Space 
Edit the properties of selected block
labels(s) 
Tab 
Display dialog for the numerical entry of
coordinates for a new label 
Escape 
Unselect all block labels 
Delete 
Delete selected block label(s) 
F3 
Scale mesh size in all block labels by ½ 
F4 
Scale mesh size in all block by 2 
Group
Mode Keys 

Key 
Function 
Space 
Edit the
properties of selected objects 
Escape 
Unselect
all 
Delete 
Delete
selected objects(s) 
View Manipulation Keys 

Key 
Function 
Left Arrow 
Pan left 
Right Arrow 
Pan right 
Up Arrow 
Pan up 
Down Arrow 
Pan down 
Page Up 
Zoom in 
Page Down 
Zoom out 
Home” 
Zoom “natural 
Table 1:
BELADRAW hot keys
Likewise, specific functions are associated
with mouse button input. The user employs the mouse to create new object,
select obects that have already been created, and inquire about object
properties. Table 2 is a summary of the mouse button click actions.
Point Mode 

Action 
Function 
Left Button Click 
Create a new point
at the current mouse pointer location 
Right Button Click 
Select the nearest
point 
Right Button
DblClick 
Display coordinates
of the nearest point 
Line/Arc Segment
Mode 

Action 
Function 
Left Button Click 
Select a start/end
point for a new segment 
Right Button Click 
Select the nearest
line/arc segment 
Right Button
DblClick 
Display length of
the nearest arc/line segment 
Block Label Mode 

Action 
Function 
Left Button Click 
Create a new block
label at the current mouse pointer location 
Right Button Click 
Select the nearest
block label 
Right Button
DblClick 
Display coordinates
of the nearest block label 
Group Mode 

Action 
Function 
Right Button Click 
Select the group
associated with the nearest object 
Table 2. BELADRAW
Mouse button actions
Generally, the user
needs to size or move the view of the problem geometry displayed on the screen.
Most of the view manipulation commands are available via buttons on the
preprocessor toolbar. The functionality of can generally also be accessed via
the ‘View Manipulation Keys’ listed in Table 1. The View Manipulation
toolbar buttons are pictured in Figure 3.
Figure 3. View
Manipulation toolbar buttons.
The meaning of the
View Manipulation toolbar buttons are:
The arrows on the toolbar correspond to moving the view in the direction
of the arrow approximately _{} of the current screen
width. The “blank page” button scales the screen to the smallest possible view
that displays the entire problem geometry. The “+” and “” buttons zoom the
current view in and out, respectively. The “page with magnifying glass” button
allows the view to be zoomed in on a userspecified part of the screen. To use
this tool, first push the toolbar button. Then, move the mouse pointer to one
of the desired corners of the “new” view. Press and hold the left mouse button.
Drag the mouse pointer to the opposite diagonal corner of the desired “new”
view. Last, release the left mouse button. The view will zoom in to a window
that best fits the user’s desired window.
Some infrequently used
view commands are also available, but only as options off of the Zoom selection of the main menu. This menu contains
all of the manipulations available from the toolbar buttons, plus the options Keyboard, Status Bar,
and Toolbar.
The Keyboard selection allows the user to zoom in to a
window in which the window’s corners are explicitly specified by the user via
keyboard entry of the corners’ coordinates. When this selection is chosen, a
dialog pops up prompting for the locations of the window corners. Enter the
desired window coordinates and hit “OK”. The view will then zoom to the
smallest possible window that bounds the desired window corners. Typically,
this view manipulation is only done as a new drawing is begun, to initially
size the view window to convenient boundaries.
The Status Bar selection can be used to hide or show the
oneline status bar at the bottom of the Beladraw window. Generally, it is
desirable for the toolbar to be displayed, since the current location of the
mouse pointer is displayed on the status line.
The Toolbar selection can be used to hide or show the
toolbar buttons. The toolbar is not fundamentally necessary to running
Beladraw, because any selection on the toolbar is also available via selections
off of the main menu. If more space on the screen is desired, this option can
be chosen to hide the toolbar. Selecting it a second time will show the toolbar
again. It may be useful to note that the toolbar can be undocked from the main
screen and made to “float” at a userdefined location on screen. This is done
by pushing the left mouse button down on an area of the toolbar that is not
actually a button, and then dragging the toolbar to its desired location. The
toolbar can be docked again by moving it back to its original position.
To aid in drawing your
geometry, a useful tool is the Grid. When the grid is on, a grid of light blue
pixels will be displayed on the screen. The spacing between grid points can be
specified by the user, and the mouse pointer can be made to “snap” to the
closest grid point.
The easiest way to
manipulate the grid is through the used of the Grid Manipulation toolbar
buttons. These buttons are pictured in Figure 4.
Figure 4. Grid
Manipulation toolbar buttons.
; The leftmost button in Figure 0 shows
and hides the grid. The default is that the button is pushed in, showing the
current grid. The second button, with an icon of an arrow pointing to a grid
point, is the “snap to grid” button. When this button is pushed in, the
location of the mouse pointer is rounded to the nearest grid point location. By
default, the “snap to grid” button is not pressed. The rightmost button brings
up the Grid Properties dialog. This dialog is shown in Figure 5.
Figure 5. Grid
Properties dialog.
The Grid Properties
dialog has an edit box for the user to enter the desired grid sizing. When the
box appears, the number in this edit box is the current grid size. The edit box
also contains a drop list that allows the user to select between Cartesian and
Polar coordinates. If Cartesian is selected, points are specified by their
(x,y) coordinates for a planar problem, or by their (r,z) coordinates for an
axisymmetric problem. If Polar is selected, points are specified by an angle
and a radial distance from the origin. The default is Cartesian coordinates.
Several useful tasks
can be performed via the Edit menu
off of the main menu.
Perhaps the most
frequently used is the Undo
command. Choosing this selection undoes the last addition or deletion that the
user has made to the model’s geometry.
For selecting many
objects quickly, the Select Group command is
useful. This command allows the user to select objects of the current type
located in an arbitrary rectangular box. When this command is selected, move
the mouse pointer to one corner of the region that is to be selected. Press and
hold the left mouse button. Then, drag the mouse pointer to the opposite
diagonal corner of the region. A red box will appear, outlining the region to
be selected. When the desired region has been specified, release the left mouse
button. All objects of the current type completely contained within the box
will become selected.
Any objects that are
currently selected can be moved, copied, or pasted. To move or copy selected
objects, simply choose the corresponding selection off of the main menu’s Edit menu. A dialog will appear prompting for an
amount of displacement or rotation.
The definition of
problem type is specified by choosing the Problem
selection off of the main menu. Selecting this option brings up the Problem
Definition dialog, shown in Figure 6.
Figure 6.
Problem Definition dialog.
The first selection is
the Problem
Type drop list. This drop box
allows the user to choose from a 2D planar problem (the Planar selection), or an axisymmetric problem (the Axisymmetric selection).
Next is the Length Units drop list. This box identifies what unit is
associated with the dimensions prescribed in the model’s geometry. Currently,
the program supports inches, millimeters, centimeters, meters, mils, and _{}meters.
The first edit box is
the Depth specification. If a Planar problem is
selected, this edit box becomes enabled. This value is the length of the
geometry in the “into the page” direction. This value is used for scaling
integral results in the post processor (e.g. force, inductance, etc.) to the
appropriate length. The units of the Depth selection are the same as the
selected length units.
The second edit box is
the Solver
Precision edit box. The number in
this edit box specifies the stopping criteria for the linear solver. The linear
algebra problem could be represented by:
Mx=b (7)
where _{} is a square matrix, _{} is a vector, and _{} is a vector of
unknowns to be determined. The solver precision value determines the maximum
allowable value for _{}. The default value is _{}.
Lastly, there is an
optional Comment edit box. The user can enter in a few lines of
text that give a brief description of the problem that is being solved. This is
useful if the user is running several small variations on a given geometry. The
comment can then be used to identify the relevant features for a particular
geometry.
To make a solvable
problem definition, the user must identify boundary conditions, block materials
properties, and so on. The different types of properties defined for a given
problem are defined via the Properties
selection off of the main menu.
When the Properties selection is chosen, a drop menu appears that
has selections for Materials, Boundary, Point, and Conductors. When any one of
these selections is chosen, the dialog pictured in Figure 7 appears.
Figure 7.
Property Definition dialog box
This dialog is the
manager for a particular type of properties. All currently defined properties
are displayed in the Property Name drop list at
the top of the dialog. At the beginning of a new model definition, the box will
be blank, since no properties have yet been defined. Pushing the Add Property button allows the user to define a new
property type. The Delete Property button
removes the definition of the property currently in view in the Property Name box. The Modify Property button allows the user to view and edit the property currently selected
in the Property
Name box. Specifics for defining
Point, Segment, and Block properties are addressed in the following
subsections.
If a new point
property is added or an existing point property modified, the Nodal Property dialog box appears. This dialog box is
pictured in Figure 8.
Figure 8. Nodal
Property dialog.
The first selection is
the Name edit box. The default name is “New Point
Property,” but this name should be changed to something that describes the
property that you are defining.
Next are edit boxes
for defining the voltage at a given point, or prescribing a point charge
density at a given point. The type of point property is chosen via the radio
buttons, and the value is entered in the enabled edit box.
The Boundary Property dialog box is used to specify the properties
of line segments or arc segments that are to be boundaries of the solution
domain. When a new boundary property is added or an existing property modified,
the Boundary
Property dialog pictured in
Figure 8 appears.
Figure 9.
Boundary Property dialog.
The first selection in
the dialog is the Name of the property. The
default name is “New Boundary,” but you should change this name to something
more descriptive of the boundary that is being defined.
The next selection is
the BC
Type drop list. This specifies the
boundary condition type. Currently, BELA supports the following types of
boundaries:
·
Fixed Voltage With this
type of boundary condition, potential V
is set to a prescribed along a given boundary
·
Mixed This denotes a boundary
condition of the form:
_{} (8)
The parameters for this class of boundary condition are specified in the
Mixed BC
parameters box in the dialog. By
the choice of coefficients, this boundary condition can either be a Robin or a Neumann
boundary condition. By carefully
selecting the _{} coefficient and
specifying _{}, this boundary condition can be applied to the outer
boundary of your geometry to approximate an unbounded solution region. For more
information on open boundary problems, refer to the Appendix.
·
Surface Charge Density This selection is used to apply distributions
of line charge to segments or arc segments in the problem. Unlike all of the other boundary conditions,
this BC type is often used on internal boundaries between materials or on
isolated segments. Typically, other BCs
are only used on exterior boundaries.
·
Periodic This type of boundary condition is applied
to either two segments or two arcs to force the potential to be identical along
each boundary. This sort of boundary is useful in exploiting the symmetry
inherent in some problems to reduce the size of the domain which must be
modeled. The domain merely needs to be periodic, as opposed to obeying more
restrictive _{} or _{} line of symmetry
conditions. Another useful application of periodic boundary conditions is for
the modeling of “open boundary” problems, as discussed in Appendix 3.
Often, a periodic boundary is made up of several different line or arc
segments. A different periodic condition must be defined for each section of
the boundary, since each periodic BC can only be applied to a line or arc and a
corresponding line or arc on the remote periodic boundary.
·
Antiperiodic The
antiperiodic boundary condition is applied in a similar way as the periodic
boundary condition, but its effect is to force two boundaries to be the
negative of one another. This type of boundary is also typically used to reduce
the domain which must be modeled, e.g. so that an electric machine might be
modeled for the purposes of a finite element analysis with just one pole.
The Block Property dialog box is used to specify the properties
to be associated with block labels. The properties specified in this dialog
have to do with the material of which the block is composed. When a new
material property is added or an existing property modified, the Block Property dialog pictured in Figure 10 appears.
Figure 10. Block
Property dialog.
As with Point and
Boundary properties, the first step is to choose a descriptive name for the
material that is being described. Enter it in the Name edit box in lieu of “New Material.”
Next, permittivity for
the material needs to be specified. BELA allows you to specify different
relative permittivities in the vertical and horizontal directions (e_{x} for the x or horizontal direction, and e_{y} for the y or vertical .
A volume charge
density (r) can also be prescribed by filling in the
appropriate box in the material properties dialog.
Since one kind of
material might be needed in several different models, BELA has a builtin
library of Block Property definitions. The user can access and maintain this
library through the Properties  Materials Library selection off of the main menu. When this option is selected, the Materials Library dialog pictured in Figure 11
appears.
Figure 11.
Materials Library dialog.
The options on this
dialog allow the user to exchange Block Property definitions between the
current model and the materials library.
The materials library
should be located in the same directory as the BELA executable files, under the
filename statlib.dat. If you move the materials library, Beladraw will
not be able to find it.
The purpose of the
conductor properties is mainly to allow the user to apply constraints on the
total amount of charge carried on a conductor.
Alternatively, conductors with a fixed voltage can be defined, and the
program will compute the total charge carried on the conductor during the
solution process.
For fixed voltages,
one could alternatively apply a Fixed Voltage boundary condition. However,
applying a fixed voltage as a conductor allows the user to group together
several physically disjoint surfaces into one conductor upon which the total
net charge is automatically computed.
The dialog for
entering conductor properties is pictured in Figure 12.
Figure 12.
Conductor Property dialog
To actually mesh the
model, analyze the model, and view the results, the Beladraw editor must spawn
external programs. These tasks are most easily performed by the toolbar buttons
pictured in Figure 13.
Figure 13.
Toolbar buttons for spawning external tasks.
The first of these
buttons (with the “yellow mesh” icon) runs the mesh generator. The solver
actually automatically calls the mesh generator to make sure that the mesh is
up to date, so you never have to call the mesher from within Beladraw. However,
it is almost always important to get a look at the mesh and see that it “looks
right.” When the mesh generation button is pressed, the mesher is called. While
the mesher is running, an entry labeled “triangle” will appear on the Windows taskbar.
Triangle is actually a console application that runs in a minimized window.
After the geometry is triangulated, the finite element mesh is loaded into
memory and displayed underneath the defined nodes, segments, and block labels
as a set of yellow lines.
If you have a very
large model, just keeping all of the mesh information in core can take up a
significant amount of memory. If you are about to analyze a very large problem,
it might be a good idea to choose the Mesh  Purge Mesh option off of the main menu. When this option is selected, the mesh is
removed from memory, and the memory that it occupied is freed for other uses.
The second button,
with the “handcrank” icon, executes the solver, Belasolv.exe. Before Belasolve is actually run, the Triangle
is called to make sure the mesh is up to date. Then, Belasolve is called. When
Belasolve runs, it opens up a console window to display status information to
the user. However, Belasolve requires no user interaction while it is running.
When Belasolve is finished analyzing your problem, the console window will
disappear. The time that Belasolve requires is highly dependent on the problem
being solved. Solution times are typically on the order of 1 to 10 seconds,
depending upon the size and complexity of the problem and the speed of the
machine analyzing the problem.
The “big magnifying
glass” icon is used to run the postprocessor once the analysis is finished. A
detailed description of the postprocessor addressed in Section 4.
For interfacing with
CAD programs and other finite element packages, Beladraw supports the import
and export of the AutoCAD dxf file format. Specifically, the dxf interpreter in
Beladraw was written to the dxf revision 13 standards. Only 2D dxf files can be
imported in a meaningful way.
To import a dxf file,
select Import
DXF off of the File menu. A dialog will appear after the file is
seleted asking for a tolerance. This tolerance is the maximum distance between
two points at which the program considers two points to be the same. The
default value is usually sufficient. For some files, however, the tolerance
needs to be increased to import the file correctly. Beladraw does not
understand all the possible tags that can be included in a dxf file; instead, it
simply strips out the commands involved with drawing lines, circles, and arcs.
All other information is simply ignored. Generally, dxf import is a useful
feature. It allows the user to draw an initial geometry using their favorite
CAD package. Once the geometry is laid out, the geometry can be imported into
Beladraw and detailed for materials properties and boundary conditions.
Do not despair if
Beladraw takes a while to import dxf files (especially large dxf files). The
reason that Beladraw can take a long time to import dxf files is that a lot of
consistency checking must be performed to turn the dxf file into a valid BELA
geometry. For example, large dxf files might take up to a minute or two to
import.
The current Beladraw
geometry can be exported in dxf format by selecting the Export DXF option off of the File menu. The dxf files generated from Beladraw
can then be imported into CAD programs to aid in the mechanical detailing of a
finalized design, or imported into other finite element or boundary element
programs.
The executable Belaview.exe is the postprocessor used to view solutions
generated by the Belasolve solver. This
program can either be run on its own, from the “Start” menu (to view previously
solved problems), or spawned from within Beladraw
to view a newly generated solution. Data files for Belaview have the .ans prefix.
Similar to the
preprocessor, the postprocessor always operates in one of three modes,
depending upon the task to be performed. These modes are: Point Values Mode In this mode, the user can
click on various points in the solution region. Local field values are then
listed in the Belaview Output window.
Contour Mode This mode allows the user to define arbitrary contours in the solution
region. Once a contour is defined, plots of field quantities can be produced
along the contour, and various line integrals can be evaluated along the
contour. Block Mode The Block Mode lets the user define a subdomain in the
solution region. Once the block has been defined, a variety of area and volume
integrals can be taken over the defined subdomain. Integrals include stored
energy (inductance), various kinds of losses, total current in the block, and
so on.
The current
postprocessor mode is controlled via the Analysis Mode toolbar buttons, shown
in Figure 14.
Figure 14.
Analysis Mode toolbar buttons
The buttons denote,
respectively, Point Values mode, Contour Mode, and Block Mode. The depressed
button denotes the current mode. The default mode when Belaview starts is Point
Values mode.
The aspects of the
current view and of the current grid are regulated via the use of toolbar
buttons. The view is manipulated by the following toolbar buttons:
Figure 15. View
Manipulation toolbar buttons.
and the grid settings
are manipulated by these grid manipulation toolbar buttons:
Figure 16. Grid
Manipulation toolbar buttons.
The grid and view
manipulation work in exactly the same fashion as these same features in the preprocessor.
Refer to Section 3.4 for a detailed description of grid manipulation, and
to Section 3.3 for view manipulation.
Unlike the
preprocessor, Belaview is not very dependent on keyboard commands.
In the Point Values
mode, there is only one relevant key press. In this mode, the Tab key allows
the user to enter the coordinates of a specific point. After the coordinates of
a point are entered, the field values at that point are displayed in the Belaview Output window.
In the Contours, there
are three relevant keys. Pressing the Escape key wipes out any current contour
or block definition. Pressing the Delete removes the last point added to the
current contour or block edge. Last, pressing the Tab keys allows the user to
numerically enter the coordinates of a point to be included in the current
contour.
In the Block mode, the
Escape and Delete keys have the same definition as in Contours mode. In the
Block mode, the Tab key does not do anything, since all points on the contour
must also be Points defining the model’s geometry.
In contrast, the
operation of the postprocessor is very dependent upon input from the mouse.
In the Point Values
mode, a Left Button Click is used to display field values at the current mouse
location. If Snap to Grid is enabled, values are displayed at the closes grid
point instead.
In the Contours mode,
mouse clicks are used to define the contour. A Left Button Click adds the
closest Point in the model’s geometry. Via a Right Button Click, the current
mouse pointer location is added to the contour. A contour appears as a red line
on the screen.
Blocks are defined in
Block mode in a fashion very similar to the way in which contours are defined.
A block is defined by drawing a contour around the region of interest. The
contour appears as a green line on the Belaview screen. When the ends of the
contour meet, the block is defined. All elements enclosed by the contour (all
elements that form the block) then turn green in the Belaview window.
A Left Button Click
attempts to add the nearest Point in the input geometry to the Block’s contour.
However, a block can only be defined along line and arc segments from the input
geometry. Each node on the boundary of the block must be selected in order to
form the Block boundary. In Block mode, the Right mouse button has no function.
One of the most useful
ways to get a subjective feel for a solution is by plotting the eqipotentials
of voltage. By default, a set of 19 equipotential lines are plotted when a
solution is initially loaded into Belaview. The number and type of
equipotential lines to be plotted can be altered using the Contours Plot icon
in the Graph Mode section of the toolbar (see Figure 17). The Contour Plot
icon is the icon with the black contours.
Figure 17. Graph
Mode toolbar buttons.
When this button is
pressed, a dialog pops up, allowing the choice of the number of contours.
In the contour plot
dialog, a check box is also present titled “Show stress tensor mask”. If this
box is checked, the contour lines associated with the last Weighted Stress
Tensor integration are also displayed, by default as orange flux lines.
Density plots are also
a useful way to get a quick feel for the flux density in various parts of the
model. By default, a flux density plot is not displayed when Belaview first
starts. However, the plot can be displayed by pressing the middle button in the
Graph Mode section of the toolbar (see Figure 17). A dialog the pops up that
allows the user to turn density plotting on.
The user can select
between density plots of Voltage (V) or the magnitude of Electric Field Intensity (E) or Electric
Flux Density (D). The field at each point is classified into one of 12 contours
distributed evenly between either the minimum and maximum flux densities or
userspecified bounds.
When Belaview is in
Contours Mode, various field values of interest can be plotted along the
defined contour. A plot of a field value defined contour is performed by
pressing the “graphed function” icon in the Plot and Integration group of
toolbar buttons, shown in Figure 18.
Figure 18. Line
Plot and Integration toolbar buttons
When this button is
pressed, the XY Plot dialog (see
Figure 19) appears with a drop list containing the types of line plots
available. Choose the desired type of plot and press “OK.”
Figure 19. XY
Plot dialog.
After “OK” is pressed,
the program computes the desired values along the defined contour. These values
are then plotted using the femmplot
program, which is called automatically to display the plot.
By default, the Write data to text
file box is not checked. If the
user selects this option, the file selection dialog will appear and prompt for
a filename to which to write the data. The data is written in twocolumn text
format. If Write data to text file is
selected, a femmplot window will not appear.
Currently, the type of
line plots supported are: Potential
along the contour; Magnitude of the flux density along the contour; Component
of flux density normal to the contour; Component of flux density tangential to
the contour; Magnitude of the field intensity along the contour; Component of
field intensity normal to the contour; Component of field intensity tangential
to the contour;
In all of these plots,
the direction of the normal is understood to be as shown in Figure20. The
tangential direction is understood to be the direction in which the contour was
defined.
Figure 20. When
in doubt plots and integrals taken on this side of a contour.
In certain cases, the
quantity to be plotted can be ambiguous. This can occur, for example, if a plot
of the tangential field intensity is requested on a contour running along an
interface between two materials of differing permittivity. In this case, there
is a discontinuity in the tangential field intensity, and the value of this
quantity is different on each side of the interface. Belaview resolves the
conflict by always evaluating the plots at a differentially small distance to
the “normal” side of the line. Therefore, by defining the same contour but
reversing the order in which the points are specified, plots of the quantity of
interest on each side of a boundary can be obtained.
Once a contour has been
specified in Contours mode, Line Integrals can be performed along the specified
contour. These integrals are performed by evaluating a large number of points
at evenly spaced along the contour and integrating using a simple
trapezoidaltype integration scheme.
To perform an
integration, press the “integral” icon on the toolbar (as shown in
Figure 0). A small dialog will appear with a drop list. Choose the desired
integral from the drop list and press OK. The
amount of time required to perform the integral will be virtually instantaneous
for some types of integrals; however, some types may require several seconds to
evaluate. When the evaluation of the integral is completed, the answer appears
on the screen in a popup box.
The line integrals
currently supported are:
·
E.t This integral returns
the voltage drop along the defined contour
·
D.n. This integral returns
the total electrix flux passing through a volume defined by extruding or
sweeping the defined contour. If this
integral is performed over a closed contour, the resulting quantity is equal to
the charge contained inside the contour.
·
Contour Length/Area. This
integral returns the length of the defined contour in meters, as well as the
area of the extruded or swept volume associated with the defined contour.
·
Force from stress tensor.
This integral totals the force produced on the contour derived from Maxwell’s
stress tensor. Deriving meaningful force results requires some care in the
choice of integration path; refer to Section 4.10 for a detailed
discussion of force and torque calculation.
·
Torque from stress tensor.
This selection integrates the torque about the point (0,0) inferred from
Maxwell’s stress tensor. Again, some guidelines must be followed to get
accurate torque results (see Section 4.10).
To select the regions
over which a block integral is to be performed, leftclick with the mouse in
the desired region. The selected region
will appear highlighted in green. For
some block integrals (i.e. weighted
stress tensor force and torque), one desires to select a conductor composed of
lines and points, rather than a block.
In this case, the desired conductor can be selected by clicking on it
with the right mouse button. The
selected conductor will appear in red.
To perform an
integration, press the “integral” icon on the toolbar (as shown in
Figure 18), and a dialog will appear with a drop list. Choose the desired
integral from the drop list and press OK. The
integral is then performed by analytically integrating the specified kernel
over each element in the defined region, and summing the results for all
elements. Volume integrals may take
several seconds to evaluate, especially on dense meshes. Be patient. When the
evaluation of the integral is completed, the answer appears on the screen in a
popup box.
The block integrals
currently supported are:
If the user is interested in the contours along
which the integral was performed, the "stress tensor mask" box can be
checked in the contour plot dialog. A set of orange (by default) lines will be
displayed that.
Often, the calculation
of electrostatic and torques is a goal of finite element analysis. This section
discusses some of the different methods of deducing forces and torques using
BELA.
If one is attempting
to compute the force on a collection of currents in a region containing only
materials with a unit relative permeability, the volume integral of Lorentz
torque is always the method to employ. Lorentz force results tend to be very
accurate. However, again, they are only applicable for the forces on conductors
of with unit permeability (e.g. coils in a voice coil actuator).
This volume integral
greatly simplifies the computation of forces and torques. Merely select the
blocks or conductors upon which force or torque are to be computed and evaluate
the integral. No particular “art” is required in getting good force or torque
results (as opposed to the Stress tensor line integral), although results tend
to be more accurate with finer meshing around the region upon which the force
or torque is to be computed.
One limitation of the
Weighted Stress Tensor integral is that the regions upon which the force is
being computed must be entirely surrounded by air and/or abutting a boundary.
In cases in which the desired region abuts a nonair region, force results may
be deduced from differentiation of stored electric field energy.
Generally, you should
not use the Stress Tensor line integral to compute forces and torques if it can
be avoided (i.e. use the volume integral version instead). The indiscriminate
use Maxwell’s Stress Tensor can result in bad predictions forces and torques.
Maxwell’s stress
tensor prescribes a force per unit area produced by the electric field on a
surface. The differential force produced is:
(9) 
where n denotes
the direction normal to the surface at the point of interest. The net force on
an object is obtained by creating a surface totally enclosing the object of
interest and integrating the stress over that surface.
For best results,
never integrate the stress tensor along an interface between materials. Always
define the integration contour as a closed path around the object of interest
with the contour displaced several elements (at least two elements) away from
any interfaces or boundaries.
Always use as fine a
mesh as possible in problems where force results are desired. A good way to proceed in finding a mesh
that is “dense enough” is to solve the problem on progressively finer meshes,
evaluating the force on each mesh. By comparing the results from different mesh
densities, you can get and idea of the level of accuracy (by looking at what
digits in the answer that change between various mesh densities). You then pick
the smallest mesh density that gives convergence to the desired digit of
accuracy.
For torque computations,
all the same rules apply as for force computations (i.e. define integration
contours away from boundaries and interfaces, and use a dense mesh).
Ultimately, you may
want to export graphics from BELA for inclusion in reports and so on. It is
possible to get what you are seeing on the screen onto disk in several
different graphics formats.
Probably the easiest
way to get graphics out of Belaview is to use the Copy as Bitmap or Copy as Metafile selections off of the main menu’s Edit
list. These command takes whatever is currently in the Belaview window and
copies it to the clipboard as a Device Independent Bitmap (.bmp format) and Extended Metafile (.emf format), respectively. The clipboard data can
then be pasted directly into most applications (e.g. Word, MS Paint, etc).
L_{A}T_{E}X
afficionados typically find PostScript to be the most useful type of graphics
output. BELA does not support postscript output directly, but it is still
relatively easy to create postcript figures with BELA. To obtain a postscript
version of the current view, you first must set up a postscript printer driver
that outputs to File:. This is done via the
following steps: Choose
+Settings/Printers+ off of the Windows Start menu. A window containing the list
of currently defined printers will appear. Double click on the Add Printer icon in this list. The Add Printer Wizard will appear on the screen. Choose Local Printer; hit Next; A
list of printers will appear. Choose a postscript printer off of this list. The
Apple Laserwriter II NT is a good choice. Select FILE: as the port which will be used with this printer. Accept the defaults
for all remaining questions. Now, when
you want a postscript picture of the currently displayed screen, just choose
+File/Print+ off of Belaview’s main menu. As the printer, choose the postscript
printer that you have previously defined. When you print to this printer, you
will be prompted for a file name, and graphic will be written as a postscript
figure to the specified file name.
If conductor
properties are used to specify the excitation, a useful byproduct is ready
access to the voltage and charge on the conductor. To view the conductor
results, select ViewConductor Props off of
the Belaview main menu. A dialog, as pictured in Figure 21 will appear. There
is a drop list on the dialog, from which the user selects the conductor for
which results are desired. When a conductor is selected, the voltage and charge
associated with that conductor are displayed.
Figure 21.
Conductor results dialog.
There are some
additional entries on the Belaview View menu
that might be useful to you from time to time. These are:
·
Smoothing By default, a smoothing algorithm is applied
to the flux density solution. Because firstorder triangles are used as trial
functions for potential, the resulting flux density and field intensity
distributions are piecewise constant in each element. The smoothing algorithm
uses a nearest neighbor interpolation to obtain linear D and E distributions
over each element. The smoothed solution generally looks better on the screen,
and somewhat increases the accuracy of D
and E near the vertices of each
element. However, if you want to toggle smoothing, this can be done by
selecting the Smoothing option.
·
Show Points Especially when making graphics for reports,
presentations, etc, it may be desirable to hide the small boxes on the screen
that denote input node points. The Show Points
option allows the user to toggle whether or not the input points are shown.
·
ToolBar Use this toggle to hide and show the
floating toolbar.
·
Point Props Use this toggle
to hide and show the floating dialog box used to display point property
information.
The Lua extension
language has been used to add scripting/batch processing facilities. The
preprocessor and postprocessor can either run Lua scripts through a selection
on the Files menu, or Lua commands can be entered in directly to the Lua
Console Window in either program.
Lua is a complete,
opensource scripting language. Source code for Lua, in addition to detailed
documentation about programming in Lua, can be obtained from the Lua homepage
at http://www.lua.org.
Because the scripting files are text, they can be edited with any text editor
(e.g. notepad).
In addition to the
standard Lua command set, a number of BELAspecific functions have been added
for manipulating files in both the pre and postprocessor. These commands are
described in the following sections.
Lua scripts are
invoked by selecting the Open Lua Script selection off of the File menu of either the pre or postprocessor. A
file selection dialog then appears, and the selected Lua script file is
executed.
A number of different
commands are available in the preprocessor. Two naming conventions can be used:
one which separates words in the command names by underscores, and one that
eliminates the underscores. A list of alternate, equivalent preprocessor
scripting function names is shown in Table 3.
open_bela_file openbelafile zoom_out zoomout add_material addmaterial 
add_point_prop addpointprop 
Table 3.
Alternate preprocessor scripting function names
run_post("c:\\myluascript.lua","luavar=filename=myfilename")
All variables are passed as strings so lua
internal conversion routines must be employed to obtain numbers from the
strings. The postprocessor window can also be minimized by including the
parameter "windowhide", e.g.:
run_post("c:\\myluascrip.lua","windowhide")
runpost("\"c:\\program
files\\femm30\\bin\\testpost.lua\"")
This command will affect all subsequent uses of the other editing
commands, if they are used WITHOUT the editaction
parameter.
propnum Symbol Description
0 BlockName Name of
the material
1 ex x
(or r) direction relative permittivity
2 ey y
(or z) direction relative permittivity
3 qs Volume
charge
propnum Symbol Description
0 BdryName Name
of boundary property
1 Vs Fixed
Voltage
2 qs Prescribed
charge density
3 c0 Mixed
BC parameter
4 c1 Mixed BC
parameter
5 BdryFormat Type
of boundary condition:
0 = Prescribed V
1 = Mixed
2 = Surface charge density
3 = Periodic
4 = Antiperiodic
propnum Symbol Description
0 PointName Name of
the point property
1 Vp Prescribed nodal voltage
2 qp Point charge density in C/m
propnum Symbol Description
0 ConductorName Name
of the conductor property
1 Vc Conductor voltage
2 qc Total conductor
charge
3 ConductorType 0 =
Prescribed charge,
1 = Prescribed voltage
There are a number of
Lua scripting commands designed to operate in the postprocessor. As with the
preprocessor commands, these commands can be used with either the underscore
naming or with the nounderscore naming convention. The equivalent function
names in the two conventions are shown in Table 4.
get_point_values getpointvalues 
hide_density_plot hidedensityplot 
Table 4.
Alternate postprocessor scripting function names
type Integral
0 _{}
1 _{}
2 Contour
length
3 Force
from stress tensor
4 Torque
from stress tensor
This integral returns either 1 or 2
values, depending on the integral type, e.g.:
Fx, Fy =
lineintegral(3)
type Integral
0 Stored Energy
1 Block Volume
2 Average D over
the block
3 Average E over
the block
4 Weighted Stress Tensor
Force
5 Weighted Stress Tensor
Torque
Returns one or two floating point values as results, e.g.:
Fx, Fy = blockintegral(4)
Symbol Definition
V Voltage
Dx x or r direction component of displacement
Dy y or z direction component
of displacement
Ex x or r direction component
of electric field intensity
Ey y or z direction component
of electric field intensity
ex x or r direction component
of permittivity
ey y or z direction component
of permittivity
nrg electric field energy density
Example: To catch all values at (0.01,0) use
V,Dx,Dy,Ex,Ey,ex,ey,nrg= getpointvalues(0.01,0)
PlotType Definition
0 V
(Voltage)
1 D (Magnitude of
flux density)
2 D . n (Normal flux density)
3 D . t (Tangential flux density)
4 E (Magnitude of
field intensity)
5 E . n (Normal field intensity)
6 E . t (Tangential field intensity)
Valid file formats are
FileFormat Definition
0 Multicolumn text with legend
1 Multicolumn text with no legend
2 Mathematicastyle formatting
For example, if one wanted to plot V
to the screen with 200 points
evaluated to make the graph, the command would be:
makeplot(0,200)
If this plot were to be written to disk as a metafile, the command would be:
makeplot(0,200,"c:temp.emf")
To write data instead of a plot to disk, the command would be of the form:
makeplot(0,200,"c:temp.txt",0)
–
legend Set to 0 to hide the
plot legend or 1 to show the plot legend.
–
gscale Set to 0 for a
colour density plot or 1 for a grey scale density plot.
–
upper_D Sets the upper
display limit for the density plot.
–
lower_D Sets the lower
display limit for the density plot.
–
type Sets the type of
density plot. A value of 0 plots
voltage; 1 plots the magnitude of D, and 2 plots the magnitude of E
–
numcontours Number of
equipotential lines to be plotted.
–
upper_V Upper limit for
contours.
–
lower_V Lower limit for
contours.
If numcontours
is 1 all parameters are ignored and default values are used, e.g. show_contour_plot(1)
save_bitmap("\"c:\\temp\\femm30\\bin\\screenshot.bmp\"")
For those of you
interested in what's going on behind the scenes in the Belasolve solver, this
section is meant as a brief description of the methods and techniques used by
BELA. References are cited as applicable.
All elements were
derived using variational formulations (based on minimizing energy, as opposed
to Galerkin, least squares residual, and so on). Explanations of the
variational approach for 2D planar problems with firstorder triangle elements
are widely available in the literature ([6] in particular was referred to
during the creation of BELA).
For all problems,
variations of the iterative Conjugate Gradient solver are used. This technique
is appropriate for the sort of problem that BELA solves, because the matrices
are symmetric and very sparse. A rowbased storage scheme is used in which only
the nonzero elements of the diagonal and upper triangular part of the matrix
are solved.
The preconditioned
conjugate gradient (PCG) code is based on the discussion in [7]. Minor
modifications are made to this algorithm to avoid computing certain quantities
more than once per iteration. Although Silvester promotes the use of the
Incomplete Cholesky preconditioner, it is not used in BELA, because it nearly
doubles the storage requirementsfor each element of the matrix stored, a
corresponding element of the preconditioner must also be stored. Instead, the
Symmetric Successive OverRelaxation (SSOR) preconditioner, as described in
[8], is used. The advantage of this preconditioner is that it is built on the
fly in a simple way using only the matrix elements that are already in storage.
In general, the speed of PCG using SSOR is said to be comparable to the speed
of PCG with Incomplete Cholesky.
Before a problem is
formulated, a node renumbering algorithm is employed. Although the conjugate
gradient schemes work well without renumbering, the renumbering seems to
roughly halve the solution time. There is an overall advantage to using the
renumbering because the time required to perform the renumbering is small
compared to the time required to run the Conjugate Gradient solver. Although
there are many possible approaches to renumbering, BELA uses the CuthillMcKee
method as described in [2]. Although there are newer schemes that yield a
tighter profile, CuthillMcKee does a relatively good job and requires very
little time to execute. The renumbering seems to speeds up the CG algorithm by
reducing the error between the SSOR approximation of the matrix inverse and the exact matrix inverse. An
interesting paper on the effect of the ordering of the unknowns on convergence
in conjugate gradient methods is [9].
Since firstorder triangles are used by BELA, the resulting solution for E obtained by differentiating V is constant over each element. If the raw D and E are used by the postprocessor, density plots of D and 2D plots of field quantities along userdefined contours look “crusty”. Also, the values of D and E are not so accurate at points in an element away from the element's centroid. ^{}
The use of smoothing
to recover the accuracy lost by differentiating V is known as “superconvergence”. There are quite a few researchers
actively pursuing this area. Of the greatest interest to BELA are socalled
“patch recovery” techniques. The basic idea is the the solutions for D are most accurate at the centroid of
the triangular element (known as its Gauss Point). One desires a continuous
profile of D that can be interpolated
from nodal values, in the same way that voltage V can be represented. However, the “raw” solution of D is multivalued at any node point,
those values being the different constant values of D in each element surrounding the node point. The general approach
to estimating the “true” value of D
at any node point is to fit a leastsquares plane through the values of D at the Gauss points of all elements
that surround a node of interest, and to take the value of the plane at the
node point's location as its smoothed value of D [10].
However, this approach
to patch recovery has a lot of shortcomings. For the rather irregular meshes
that can arise in finite elements, the leastsquares fit problem can be
illcondition, or even singular, at some nodes in the finite element mesh.
Furthermore, the superconvergence solution can actually be less accurate than
the piecewise constant solution in the neighborhood of boundaries and
interfaces.
One can note that the
patch recovery method is merely a weighted average of the flux densities in all
of the elements surrounding a given node. Instead of a leastsquares fit, BELA
simply weights the values of flux density in each adjacent element's Gauss
point with a value inversely proportional to the distance from the Gauss point
to the node point of interest. Away from boundaries, the results seem to be
nearly as good as a leastsquares fit. At boundaries and interfaces, the
smoothed solution is no worse than the unsmoothed solution.
[1] M. Plonus, Applied electromagnetics.
McGrawHill, 1978.
[2] S. R. Hoole, Computeraided analysis and
design of electromagnetic devices, Elsevier, 1989.
[3] J. D. Jackson, Classical electrodynamics, 2^{nd}
ed, Wiley, 1975.
[4] S. McFee, J. P. Webb, and D. A. Lowther, “A
tunable volume integration formulation for force calculation in finiteelement
based computational magnetostatics,” IEEE Transactions on Magnetics,
24(1):439442, January 1988.
[5] F. Henrotte, G. Deliege, and K. Hameyer,
“The eggshell method for the computation of electromagnetic forces on rigid
bodies in 2D and 3D,” CEFC 2002, Perugia, Italy, April 1618, 2002.
[6] P. E. Allaire, Basics of the finite element
method, 1985.
[7] P. P. Silvester, Finite elements for
electrical engineers, Cambridge University Press, 1990.
[8] C. A. Fletcher, Computational techniques for
fluid dynamics, SpringerVerlag, 1988.
[9] E. F. D'Azevedo, P. A. Forsyth, and W. Tang,
“Ordering methods for preconditioned conjugate gradient methods applied to
unstructured grid problems,” SIAM J. Matrix Anal. Appl., 12(4), July 1992.
[10] O. C. Zienkiewicz and J. Z. Zhu, “ The
superconvergent patch recovery and a posteriori estimates, part 1: the recovery
technique,” International Journal for Numerical Methods in Engineering, 33:13311364,
1992.
[11] Q. Chen and A. Konrad, “A review of finite
element open boundary techniques for static and quasistatic electromagnetic
field problems,” IEEE Transactions on Magnetics, 33(1):663676, January 1997.
[12] E. M. Freeman and D. A. Lowther, “A novel
mapping technique for open boundary finite element solutions to Poissons
equation,” IEEE Transactions on Magnetics, 24(6):29342936, November 1988.
[13] D. A. Lowther, E. M. Freeman, and B.
Forghani, “A sparse matrix open boundary method for finite element analysis,”
IEEE Transactions on Magnetics, 25(4)28102812, July 1989.
Typically, finite
element methods are best suited to problems with welldefined, closed solution
regions. However, a large number of problems that one might like to address
have no natural outer boundary. A prime example is a parallel plate capacitor
in air. The boundary condition that one would like to apply is V=0 at
infinity. However, finite element methods, by nature, imply a finite domain.
Fortunately, there are methods that can be applied to get solutions that
closely approximate the “open boundary” solution using finite element methods.
The simplest, but least accurate, way to
proceed is to pick an arbitrary boundary “far enough” away from the area of
interest and declare V=0 on this
boundary. According to [11], a rule of thumb is that the distance from the
center of the problem to the outer boundary should be at least five times the
distance from the center to the outside of the objects of interest. Truncation
is the method employed by most finite element programs, because it requires no
additional effort to implement.
The down side to
truncation is that get an accurate solution in the region of interest, a volume
of air much larger than the region of interest must also be modeled. Usually,
this large region exterior to the area of interest can be modeled with a
relatively coarse mesh to keep solution times to a minimum. However, some extra
time and space is still required to solve for a region in which one has little
interest.
A thorough review of open boundary techniques
is contained in [11]. Perhaps the simple way to approximate an “open” boundary
(other than truncation) described in [11] is to use asymptotic boundary
conditions. The result is that by carefully specifying the parameters for the
“mixed” boundary condition, and then applying this boundary condition to a
circular outer boundary, the unbounded solution can be closely approximated.
Consider a 2D planar problem in polar coordinates. The domain is a circular shell of radius r_{o} in an unbounded region. As
r®¥, potential V goes to zero. On the surface of the circle, the vector is a prescribed function of q. This problem has an analytical solution, which is:

(10) 
where the v_{m} and a_{m} parameters are chosen so that the solution matches the prescribed potential on the surface of the circle.
One could think of this solution as describing the solution exterior to a finite element problem with a circular outer boundary. The solution is described inside the circle via a finite element solution. The trick is to knit together the analytical solution outside the circle to the finite element solution inside the circle.
From inspecting (10), one can see that the highernumbered harmonic, the faster the magnitude of the harmonic decays with respect to increasing r. After only a short distance, the highernumbered harmonics decay to the extent that almost all of the openspace solution is described by only the leading harmonic. If n is the number of the leading harmonic, the openfield solution for large, but not infinite, r is closely described by:

(11) 
Differentiating with respect to r yields:

(12) 
If (12) is solved for a_{n} and substituted into (11), the result is:

(13) 
Now, (13) is a very useful result. This is the same form as the “mixed” boundary condition supported by FEMM. If the outer edge of the solution domain is circular, and the outer finite element boundary is somewhat removed from the area of primary interest, the open domain solution can be closely approximated by applying (13) the circular boundary.
To apply the Asymptotic Boundary Condition, define a new, mixedtype boundary condition. Then, pick the parameters so that:

where r_{o} is the outer radius of the region in meters (regardless of the working length units), and e_{o} = 8.85418781762e012. Note that e_{o} is defined in the Lua implementation in both the pre and postprocessors as the global variable eo, which can be used in any script or edit box in the program.
Although the above derivation was specifically for 2D problems, it turns out that when a similar derivation is performed for the axisymmetric case, the definition of the mixed boundary condition coefficients are exactly the same as (14). One subtle difference, however, is that n = 1 in the axisymmetric case corresponds to the case in which there is a net charge with (i.e. the geometry looks like a point charge when viewed from a distance), whereas the n = 1 case corresponds to a dipole charge distribution in 2D planar problems. If charge is conserved in the geometry of interest in the axisymmetric case, one needs to use n = 2 in Eq. (14). It should be noted that this is a departure from the magnetostatic case with the vector potential formulation (e.g. as implemented in FEMM) in which n = 1 corresponds to a dipole arrangement in the axisymmetric case.
Some care must be used in applying this boundary condition. Most of the time, it is sufficient to take n = 1 (i.e the objects in the solution region look like a dipole when viewed from a large distance). However, there are other in which the leading harmonic is something other than n = 1. You need to use your insight into your specific problem to pick the appropriate n for the leading harmonic. You also must put the objects of interest roughly in the center of the circular finite element domain to minimize the magnitude higherorder field components at the outer boundary.
Although the application of this boundary condition requires some thought on the part of the user, the results can be quite good. Applying the absorbing boundary condition imposes no additional computing cost on the problem. The ABC is computationally no more timeconsuming to apply than enforcing V = 0 at the outer boundary. Solution times for the PCG solver are equivalent in either case.
A particularly good approach to “open boundary” problems is the Kelvin Transformation, a technique first discussed in the context of computational electromagnetics in [12] and [13]. The strengths of this technique are:
· the effects of the exterior region are, in theory, exactly modeled by this approach;
· a sparse matrix representation of the problem is retained (unlike FEMBEM methods, which give the same “exact solution” but densely couples together the boundary nodes).
· requires no “special” features in the finite element solver to implement the technique, other than the ability to apply periodic boundary conditions.
The purposes of this note are to explain what the Kelvin transformation is
derived and to show how it is implemented in the context of the FEMM finite
element program.
In the “far field” region, the material is typically homogeneous (e.g. air and free of sources. In this case, the differential equation that describes potential V is the Laplace equation:

(15) 
If we write (15) in polar notation, A is described by:

(16) 
Assume that the “near field” region of the problem can be contained in a circle of radius r_{o} centered at the origin. The farfield region is then everything outside the circle.
One approach to unbounded problems is to attempt to map the unbounded region onto a bounded region, wherein problems can more easily be solved. Specifically, we desire a way to transform the unbounded region outside the circle into a bounded region. One simple way to make such a mapping is to define another variable, R, that is related to r by:

(17) 
By inspecting (17), it can be seen that this relationship maps the exterior region onto a circle of radius r_{o}.
The next step is to transform (15), the differential equation that the field must satisfy, into the mapped space. That is, (15) must be written in terms of R and q rather than r and q. We can evaluate derivatives in terms of R instead of r by employing the chain rule:

(18) 
Now, we can note that at r = R = r_{o},

(19) 
and we can substitute (18) into (16) to yield, after some algebraic manipulation:

(20) 
Eq. (20), the transformed equation for the outer region, has exactly the same form as inner region, only in terms of R rather than r. The implication is that for the 2D planar problem, the exterior can be modeled simply by creating a problem domain consisting of two circular regions: on circular region containing the items of interest, and an additional circular region to represent the “far field.” Then, periodic boundary conditions must be applied to corresponding edges of the two circles to enforce the continuity of V at the edges of the two regions. For a finite element formulation consisting of firstorder triangles, (19) is enforced automatically at the boundaries of the two regions as long as the values of V at the boundaries are slaved together. The second circular region exactly models the infinite space solution, but does it on a bounded domainone could always back out the field for any point in space by applying the inverse of (17).