New features for automatic program operation

One of the outstanding features of our software is the use of input scripts, a classic feature of early programs that we have updated. The term script refers to a text file of input data for a program run. Figure 1 shows the role of input scripts in Field Precision software. Years ago, the primary input method for technical programs followed the sequence in the middle of the illustration. A user would employ a text editor (or a punch card machine, if we go far enough back) to create the set of instructions passed to the program. Early scripts were rigid in the order and syntax of instructions and documentation was often sparse. As a result, the script method acquired a bad reputation.

Input options for Field Precision technical programs

Figure 1. Input options for Field Precision technical programs.

Users welcomed the appearance of graphical-user interfaces. Here, the programs remembered the instruction syntax and provided setup guidance. Most modern technical programs follow the sequence show by the red arrow in Fig. 1. The user sets controls or fills in fields, and the data is passed directly to the program. In contrast, Field Precision programs follow a two-step process where the interactive routines produce a text script which is then read by the program[1]. There are several reasons why we have chosen this route:

  • Compact scripts are an effective method for archiving run setups and exchanging them with colleagues.
  • Input information is transparent to the user.
  • Scripts provide expanded input options.

We will concentrate on the third advantage in this article.

Field Precision programs offer two alternatives to the graphical-user-interface. Direct editing of the script is useful for checking values and making small changes (e.g., moving objects, modifying material properties,…) without the need to negotiate the GUI. In the past, our programs have had a limited capability for control through Windows batch files and external programs (e.g., Python scripts,…), indicated by the black arrow in Fig. 1. All of our programs had the capability to run from the command line or to be called by an external program. The single pass parameter was the name of the input script.

The change described in this article, represented by the green arrow in Fig. 1, applies to the all 2D and 3D Field Precision programs. Batch files and external programs can now control how the technical program interprets variable quantities in the input script. The easiest way to describe the feature is to follow an application example. Figure 2-top shows the geometry of a test calculation, the magnetic field of a cylindrical coil in an air volume. Finite-element solutions are performed in a finite-volume and require a specified condition on the boundary. The most boundary for magnetic-field calculations is a Dirichlet condition (constant vector potential). In the case, the boundary acts like a perfectly-conducting wall (Fig. 2-bottom). The purpose of the study is to determine how large the boundaries must be to approach the infinite-space result.

Test calculation, magnetic field of a cylndrical coil in space

Figure 2. Test calculation, magnetic field of a cylndrical coil in space. Top: Mesh. Bottom: Calculated field lines.

 

The procedure is to compute several solutions with different boundaries, comparing the values of the field at the center of the coil, Bz(0,0). The mesh input script BATCHCONTROL.MIN looks  like this:

*  -------------------------------------------------------
GLOBAL
 ZMESH 
  %1  -3.5  0.2
  -3.5  3.5  0.1
  3.5  %2   0.2
 END
 RMESH
  0.0  4.5  0.1
  4.5  %3   0.2
 END
END
*  -------------------------------------------------------
REGION  FILL AIR
 L     %1    0.00000    %2    0.00000
 L     %2    0.00000    %2    %3
 L     %2    %3         %1    %3
 L     %1    %3         %1    0.00000
END
*  -------------------------------------------------------
REGION  FILL COIL
 L     -3.00000    2.00000    3.00000    2.00000
 L      3.00000    2.00000    3.00000    4.00000
 L      3.00000    4.00000   -3.00000    4.00000
 L     -3.00000    4.00000   -3.00000    2.00000
END
*  -------------------------------------------------------
REGION BOUNDARY
 L     %1    0.00000    %2    0.00000
 L     %2    0.00000    %2    %3
 L     %2    %3         %1    %3
 L     %1    %3         %1    0.00000
END
*  -------------------------------------------------------
ENDFILE

In the GLOBAL section, notice that there is region of small elements around the coil and then larger elements to the outer boundaries. Further, note the symbolic representation of the boundary limits, a convention familiar to users of Windows batch files. The symbol %1 represents ZMin, %2 represents ZMax and %3 represents RMax.

The calling batch file has the following content:

START /B /WAIT C:\fieldp_pro\tricomp\mesh.exe C:\Simulations\BatchControl  -5.00   5.00   7.50
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.PIN
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.SCR
START /B /WAIT C:\fieldp_pro\tricomp\mesh.exe C:\Simulations\BatchControl  -6.00   6.00   9.00
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.PIN
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.SCR
START /B /WAIT C:\fieldp_pro\tricomp\mesh.exe C:\Simulations\BatchControl  -7.00   7.00  10.50
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.PIN
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.SCR
START /B /WAIT C:\fieldp_pro\tricomp\mesh.exe C:\Simulations\BatchControl  -8.00   8.00  12.00
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.PIN
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.SCR
...

The file seems verbose, but it is mainly copy-and-paste from a template prepared with the Create task button of the TriComp program launcher. The interesting lines are those that call Mesh. The first pass parameter is the input script listed above. The three additional string parameters give numerical values for the variables %1, %2 and &3. There are multiple solutions with expanding boundaries with the constraint RMax = 1.5 ZMax. The two commands that follow run PerMag with the modified mesh and then execute the analysis script BATCHCONTROL.SCR. This file has the following content:

INPUT BatchControl.POU
OUTPUT BatchControl.DAT Append
POINT 0.0 0.0
ENDFILE

The Append specification in the second command ensures that all the data will be added in sequence to a single output file.

The full set of runs is generated in about two seconds by executing the batch file. The data are available as text entries in BATCHCONTROL.DAT. Figure 3 shows a plot of the results. The dashed line is the theoretical result for an infinite space.

Central field of a cylindrical coil as a function of the boundary dimension

Figure 3. Central field of a cylindrical coil as a function of the boundary dimension. The dashed line is the theoretical result for infinite space.

You can imagine the infinite possibilities for automation. For example, suppose we wanted to do a more intensive analysis of each solution and record the results in a separate file. In this case, the calling batch might look like this:

...
START /B /WAIT C:\fieldp_pro\tricomp\mesh.exe C:\Simulations\BatchControl  -6.00   6.00   9.00
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.PIN
START /B /WAIT C:\fieldp_pro\tricomp\permag.exe C:\Simulations\BatchControl.SCR BatchControl02.DAT
...

and the OUTPUT command in the analysis script would become

OUTPUT %1

Footnotes

[1] Field Precision scripts are quite different from those of older programs. They are readable, and the technical programs use a free form input parser that does not require rigid numerical formats or a fixed order of commands.

[2] Contact us : techinfo@fieldp.com.

[3] Field Precision home page: www.fieldp.com.

Test magnetostatic solution: simple coil with boundaries

In this article, we’ll advance to 2D magnetostatic solutions using the programs Mesh and PerMag. The previous articles on electrostatics have given background on the basic concepts of FEM calculations and program operation. Therefore, in this lecture and following ones I’ll concentrate on the special features of magnetic field calculations. In preparation, I suggest you try the prepared walk-through example discussed in the PerMag manual.

To review, the calculation discussed in the previous articles determined the electrostatic potential φ (a scalar quantity) at the node points of a mesh that we created. Components of the electric field vector components at a location could then be determined by collecting local values of phi and taking numerical derivatives:

E = – ∇φ

The components were Ex and Ey for planar solutions and Ez and Er for cylindrical solutions. Things get more involved for magnetic field calculations[1]. The calculated node quantity is the vector potential A. The magnetic flux density B is given by

B = ∇ × A.

Fortunately, in 2D calculations there is only a single component of A. The node quantity is Az for planar solutions (with flux density components Bx and By) and rAθ for cylindrical solutions (with components Bz and Br). The vector potential has useful properties for making plots:

  • In planar solutions, contours of Az separated by a uniform interval ΔAz lie along lines of B separated by equal intervals of magnetic flux per length.
  • In cylindrical solutions, contours of rAθ separated by a uniform interval Δ(rAθ) lie along lines of B separated by equal intervals of magnetic flux.
Dialog to start a Mesh script in the text mode

Figure 1. Dialog to start a Mesh script in the text mode.

Let’s get to work and generate some fields. We’ll start with a cylindrical coil in free space. Because the geometry is simple, we’ll write the boundary specifications directly with an editor. In Mesh, click the New mesh (text) tool to bring up the dialog of Fig. 1. The values define a solution volume that covers the range -10.0 cm ≤ z ≤ 10.0 cm, 0.0 ≤ r ≤ 15.0 cm[2]. When you click OK, the program opens the internal text editor with the default content shown in Fig. 2. The Global section at the top sets the foundation mesh (the element geometry before adjustment of boundary nodes). It covers the range we specified with a default element size of 0.1 cm. Available advanced commands (like TriType) are listed as comment lines. We won’t need them for this calculation, so delete the comments.

Starting a Mesh script in text mode

Figure 2. Starting a Mesh script in text mode, display of the internal editor.

A default region named SolVolume (Region 1) that covers the entire solution volume has been added —the default is appropriate for this calculation[3]. Change the region name to Air. Mesh has also started a default second region. We’ll use it to define the rectangular cross section of a cylindrical coil. The comment lines show the commands that could appear within a region section. Erase the comments and rename the region Coil. The coil will have the dimensions -3.0 cm ≤ z ≤ 3.0 cm and 4.0 cm < r < 6.0 cm. Copy and paste the line vectors of Region 1 and use Search/Replace to modify the dimensions. We’ll also add a third region by copying and pasting Region 1 verbatim. Rename it Boundary and remove the word Fill in the Region command. We’ll discuss the purpose of this extra region later. For now, the modified script should look like this:

GLOBAL
 XMesh
  -10.00000   10.00000    0.10000
 End
 YMesh
  0.00000   15.00000    0.10000
 End
END
*  -------------------------------------------------------
REGION  FILL Air
 L   -10.00000    0.00000   10.00000    0.00000
 L    10.00000    0.00000   10.00000   15.00000
 L    10.00000   15.00000  -10.00000   15.00000
 L   -10.00000   15.00000  -10.00000    0.00000
END
*  -------------------------------------------------------
REGION  FILL Coil
 L    -3.00000    4.00000    3.00000    4.00000
 L     3.00000    4.00000    3.00000    6.00000
 L     3.00000    6.00000   -3.00000    6.00000
 L    -3.00000    6.00000   -3.00000    4.00000
END
*  -------------------------------------------------------
REGION  Boundary
 L   -10.00000    0.00000   10.00000    0.00000
 L    10.00000    0.00000   10.00000   15.00000
 L    10.00000   15.00000  -10.00000   15.00000
 L   -10.00000   15.00000  -10.00000    0.00000
END
*  -------------------------------------------------------
ENDFILE

Save your work and exit the editor.

Click the Load script (MIN) tool to read the boundary specifications in BARECOIL.MIN and then click Process. When complete, click Save mesh (MOU) to create the file BARECOIL.MOU. This file will be used by PerMag. You can go to the Plot/Repair menu to view a cross-section of the resulting mesh, consisting of 60,000 elements.

Run PerMag and click the “1” tool to bring up the dialog of Fig. 3. As in the previous calculation, this dialog creates a script of run and material parameters to control the solution program. We can ignore several of the fields because the values apply yo advanced runs with non-linear materials (e.g., saturable iron). For this run, we need only set the Geometry to Cylindrical and fill in the magnetic properties of the regions. Air has relative permeability μr = 1.0. The coil also has μr = 1.0. In addition, we need to set a total current that flows in the θ direction for a cylindrical solution. If the coil has N = 1000 turns and the drive current is I = 2.0 A, then the total A-turns over the cross section is 2000.0. Equivalently, the region carries a uniform azimuthal current density of (2000)/(0.06*0.02) = 1.667 MA/m2[4]. Finally, we want to see what happens if we make no special provisions for the Boundary region, so don’t set any property. Click OK and save the data in the file BARECOIL.PIN. You can check the file content with the internal PerMag editor.

PerMag dialog to define run parameters and material properties

Figure 3. PerMag dialog to define run parameters and material properties.

Click the “2” button and choose BARECOIL.PIN to find the finite-element solution. Then click “3” and load BARECOIL.POU to analyze the solution. To begin, let’s check the lines of B. Click the Plot type tool and set the option Contour lines. Then click the Plot quantity tool and choose rAthet to produce the plot of Fig. 4a. There are two interesting features:

  • Even though the field inside the coil is fairly uniform, the space between the contours is larger near the axis. This is because the lines are separated by intervals of equal magnetic flux and the area normal to z gets gets smaller near the axis. In other words, the spacing between lines is proportional to |B|/r in cylindrical coordinates.
  • The lines of B are normal to the boundaries, the standard default Neumann condition for a finite-element magnetic solution (i.e., the derivative of vector potential normal to the boundary is zero). With regard to the axial boundaries, the condition would occur if there were an infinite array of coils with alternating positive and negative currently. Certainly, this is not what we intended. Furthermore, the boundary on the outer radius is non-physical in the sense that we could find no set of real coils that would generate the field pattern.
Lines of B for the solution

Figure 4. Lines of B for the solution. a) Default Neumann boundary condition. b) Dirichlet condition, fixed vector potential.

With regard to the second point, we must do some thinking about the boundary conditions. The alternative is to set the vector potential equal to a constant along the boundary. The Dirichlet condition is rAθ = 0.0. In this case, the component of B normal to boundary is zero, so that lines of B are parallel. To do this, we change the specification of Region 3 to:

* Region 3: BOUNDARY
VecPot(3) = 0.0000E+00

The command states that the node quantity has a fixed value.

Figure 4b shows the modified solution. It is better in one respect. The field pattern could be achieved in an actual physical system — a coil inside a superconducting can. At this point, you may protest that you wanted a coil in infinite space, not in a superconducting can. We must recognize that there is a limitation to finite-element calculations. They must be performed in a finite volume. There will always be boundaries, and the conditions on the boundaries must be specified. Without additional information, the choices are Neumann or Dirichlet[5].

But, is it a limitation? Often, issues in a numerical solution mirror concerns in the physical system. When is the last time you operated a magnet in infinite space? In a real system, there is usually an array of nearby objects in unpredictable locations. These objects would strongly alter the far fields and could affect field components in the working volume of the magnet. The inmplication is that the magnet of Fig. 4, with its widely-distributed pattern of fringing fields, is a poor design. An additional disadvantage of the external fields is that they represent wasted field energy. We shall see that the extra energy means that more power must be supplied to the drive  coil. In a following article, we will consider the use of ferromagnetic materials to improve the design. In the next article, we will carry out a numerical study to quantify how much the boundaries in a finite-element solution affect the results.

Footnotes

[1] For a review of magnetic-field physics, see S. Humphries, Finite-element Methods for Electromagnetics (CRC Press, 1997), available at http://www.fieldp.com/femethods.html

[2] Mesh displays an error message if you try to open a cylindrical solution with r < 0.0. As an exercise, can you think of a solution where you would use r > 0.0?s

[3] Although Region 1 often covers the full solution, it is not a necessary condition (as we saw in the previous electrostatic solution).

[4] Coil regions in PerMag always have uniform average current density. To represent a coil with variations of winding density, divide it into multiple regions.

[5] Looking ahead to later articles, Field Precision programs can approximate infinite-space magnetic calculations via the Boundary procedure (computed Dirichlet boundary)

[6] Contact us : techinfo@fieldp.com.

[7] Field Precision home page: www.fieldp.com.

Aether simulation: heating dielectric samples in a waveguide

Recently I had an inquiry from a user interested in microwave heating of dielectric samples in a air-filled waveguide. The application presented a good opportunity to apply the recently-developed capability for coupling Aether solutions to HeatWave, so I wrote and tested a set of template input files.

In the solution, the sample was suspended at the center of a rectangular waveguide that carried a traveling wave in the TE10 mode. The wave was polarized with Ey and the guide had dimensions a = 10.0 cm along x and b = 8.0 cm along y. The sample was a cylinder of radius 2.0 cm and height 3.0 oriented as shown in Fig. 1. Heating occurred via capacitive currents driven in the resistive object. The mesh had cubic elements with sides 0.1 cm for good resolution of the sample volume.

Mesh for the calculation

Figure 1. Mesh for the calculation showing the waveguide walls, the source and absorbing layers and the heated sample.

The first task was to generate a good approximation to a TE10 mode solution. For guidance, I followed Chap. 9 of the Aether Tutorial Manual, TE10 mode in a rectangular waveguide. The cutoff frequency was

f0 = c/2a = 1.5 GHz.

At a working frequency of f = 2.0 GHz, the vacuum wavelength was λ0 = 15.0 cm and the ratio of the group velocity to the speed of light was:

vg/c = √[1 - (λ0/2a)^2] = 0.6614.

I modeled a waveguide section of length 20.2 cm. The mesh was built by defining a metal box with dimensions Wx = 10.2 cm, Wy = 8.2 cm and Wz = 20.2 cm. The interior was created by over-writing with an air box of dimensions Wx = 10.0 cm, Wy = 8.0 cm and Wz = 20.2 cm. Ideal absorbing layers of thickness Wz = 0.1 cm filled the entrance and exit. The entrance layer carried a current density

Jy = (1.0 × 10^4) cos (0.3142 x).

When waves are normally incident, the conductivity for a matched absorbing layer is given by[1]:

σ = 1/η Wz (S/m),

where η is the characteristic impedance of the adjacent medium (377.3Ω for air). As discussed in the Aether tutorial, the value must be adjusted by vg/c for modes in a waveguide. The calculated value is σ = (0.6614)(2.654) = 1.756 S/m. Figure 2 shows a plot of |E| in the plane y = 0.0 cm with the properties of the sample set to those of air. For an ideal traveling wave, the magnitude of the electric field should be constant in z.

Plot of |E| for the TE10 mode with no sample

Figure 2. Plot of |E| for the TE10 mode with no sample.

If the heating is to be effective and the finite-element calculation accurate, the electric field should penetrate the object rather than concentrate in a thin surface layer (i.e., field exclusion). I picked the following material properties: εr = 1.8, μr = 1.0 and σ = 0.5 S/m. Figure 3 shows the altered electric field when the sample is present. A substantial fraction of the electromagnetic power is reflected, resulting in an upstream standing-wave component. About half the power is transmitted as a downstream traveling wave and half is absorbed by the sample (19.2 W). Figure 4 shows the spatial variation of time-averaged power density in the sample in the plane x = 0.0 cm. The peak value is p = 1.644 × 10^6 W/m^3.

Plot of |E| for the TE10 mode with the sample present

Figure 3. Plot of |E| for the TE10 mode with the sample present.

Power-density p in the sample, plot in the plane x = 0.0 cm

Figure 4. Power-density in the sample, plot in the plane x = 0.0 cm.

The next step was to set up a HeatWave solution using the output from Aether as a thermal source. The goal of the calculation was to investigate heating of the sample and the surrounding air. Thermal effects in the virtual source and absorber layers are not of interest. Therefore, these regions along with the metal wall were set to the fixed-temperature condition T = 0.0 °C. The sample had thermal properties k = 0.02 W/m^2-°C, Cp = 2000.0 J/kg-°C and ρ = 1000.0 kg/m^3. The air region has properties k = 0.024 W/m^2-°C, Cp = 1005.0 J/kg-°C and ρ = 1.293 kg/m^3. Both the sample and air were at initial temperature T = 0.0 °C. The dynamic thermal calculation ran for Δt = 10.0 s. In the absence of thermal condition, the peak temperature in the sample is predicted to be

Tmax = p Δt/ρ Cp = 8.22 °C.

Figure 5 shows the end result, the temperature in the sample and surrounding air at t = 10.0 s in the plane x = 0.0 cm. The temperature variation tracked the power-density profile of Fig. 4with a peak value of T = 8.11 °C.

Temperature in the sample and surrounding air

Figure 5. Temperature in the sample and surrounding air at 10 s, plot in the plane x = 0.0 cm.

The complete set of input files is available to Aether/HeatWave users. To make a request, please contact us.

Footnotes

[1] See S. Humphries, Finite-element Methods for Electromagnetics (CRC Press, Boca Raton, 1997) Sect. 13.1.

[2] Contact us : techinfo@fieldp.com.

[3] Field Precision home page: www.fieldp.com.

Test electrostatic solution: meshing and accuracy

In this article (the sixth in the online course Electric and magnetic field calculations with finite-element methods), we’ll make several calculations to understand the effect of element size on solution accuracy. To expand our technique toolbox, we’ll automate steps in the analysis procedure.

We’ll continue with the solution discussed in the previous two articles. The best way to start a project is with a clear statement of the goal. In this case, we want to check the accuracy of numerical estimates of the electric field Er(r) at locations in the dielectric as a function of element size. Specifically, we’ll check field interpolations near the inner electrode (r = 5.0 cm), near the outer electrode (r = 15.0 cm) and at the center of the gap (r = 10.0 cm) for the following choices of element size: 0.1 cm, 0.2 cm, 0.3 cm, 0.4 cm and 0.5 cm.

As mentioned in previous articles, good supporting software (text editors, file managers, image editors,…) is essential for effective technical work. In this case, we’ll use a spreadsheet created with OpenOffice Calc[1]. Figure 1 shows the upper section of the sheet. Parameters at the top include the inner and outer radii, the applied voltage and the theoretical value for the capacitance per length. The section beneath lists the radial positions for the calculations. Note that the inner radius is slightly larger than Ri and the outer radius smaller than Ro. Remember that the electrode boundaries are a set of straight lines that approximate sections of circles. We have to include a little slack to make sure that the test points lie within the dielectric region for all choices of element size. We may as well let the spreadsheet do the mathematical work. The cells for the x-y coordinates of the measurement points on a 45° line contain the formula r/√2. Similarly, the cells for the theoretical values contain the equations for Er(r) and φ(r) discussed in the previous articles.

Top section of a Calc spreadsheet: solution setup

Figure 1. Top section of a Calc spreadsheet: solution setup.

The next section (listed below) demonstrates another use of a spreadsheet. We can use the program to prepare input scripts automatically. In this case, we’ll generate an analysis script — a set of instructions to EStat for standard calculations to carry out for each of the solutions:

* Coaxial cylinders analysis script
INPUT COAXIALCYLINDERS.EOU
OUTPUT COAXIALCYLINDERS.DAT
POINT 3.53624 3.53624
POINT 7.07107 7.07107
POINT 10.60589 10.60589
VOLUMEINT 1
ENDFILE

The spreadsheet has supplied the coordinates. To use the section information, simply copy and paste it into a text file COAXIALCYLINDERS.SCR (SCRipt). The commands listed perform the following actions. Open the EStat output file COAXIALCYLINDERS.EOU and write the results of calculations to the text file COAXIALCYLINDERS.DAT. Make three point interpolations at the desired locations. Also, calculate volume-integral quantities over Region 1 to find the capacitance per length.

We can quickly change the element size by directly editing the Mesh input script. Figure 2 shows the file displayed in the Context editor with syntax highlighting installed for our programs. The red arrows show the quantities to be modified. Given a modified Mesh input script, there are three steps to regenerate the solution and to update values in the data file COAXIALCYLINDERS.DAT:

  • Run Mesh, load COAXIALCYLINDERS.MIN, process the mesh and save the result.
  • Run EStat, load COAXIALCYLINDERS.EIN and run the solution.
  • In EStat, run the analysis script COAXIALCYLINDERS.SCR.

With only five solutions, the procedure wouldn’t take long. On the other hand, it takes only a few seconds to set up a batch file[2] that handles everything automatically.

File COAXIALCYLINDERS displayed in Context

Figure 2. File COAXIALCYLINDERS displayed in Context.

In the Tricomp Program Launcher, press the Create task button to bring up the dialog of Fig. 3. Supply the file prefix COAXIALCYLINDERS and activate the Beep at completion box. We need to define the three steps of the solution. Clicking a cell in the first column raises of popup menu (inset) that shows your installed Field Precision programs as well as useful batch commands. Create the three tasks shown to run Mesh and then to run EStat in the solution and analysis modes. When you press OK, the program creates the file COAXIALCYLINDERS.BAT in the working directory. Here’s the content:

REM TriComp batch file, Field Precision
START /B /WAIT C:\fieldp_basic\tricomp\mesh.exe COAXIALCYLINDERS
START /B /WAIT C:\fieldp_basic\tricomp\estat.exe COAXIALCYLINDERS.EIN
START /B /WAIT C:\fieldp_basic\tricomp\estat.exe COAXIALCYLINDERS.SCR
START /B /WAIT C:\fieldp_basic\tricomp\NOTIFY.EXE
IF EXIST COAXIALCYLINDERS.ACTIVE ERASE COAXIALCYLINDERS.ACTIVE
Create task dialog in tc.exe

Figure 3. Create task dialog in tc.exe.

Note that most Field Precision programs can be run from the command line with the input file as a parameter. The /WAIT option ensures that a program does not start before the previous program completes it action. Under batch file control, a full solution consists of the following operations:

  • Edit and save COAXIALCYLINDERS.MIN with the desired element size.
  • Press Run task in the program launcher and choose COAXIALCYLINDERS.BAT.
  • Open the resulting file COAXIALCYLINDERS.DAT and extract values of interest
Bottom section of a Calc spreadsheet: solution analysis

Figure 4. Bottom section of a Calc spreadsheet: solution analysis.

I performed the solutions and copied/pasted values to create the sections of the spreadsheet shown in Fig. 4. I set up formulas in the spreadsheet to determine the relative errors in electric-field values. For a visual reference, Fig. 5 shows the mesh for element sizes of 0.2 and 0.5 cm. We can now consider some implications of the results:

  • For a size 0.1 cm the mesh had 35336 elements. The mesh with size 0.5 cm had only 1412 elements. Therefore, the solution with the fine mesh requires about 25 times as much computational work.
  • The error in the electric field calculation at the center of the gap (r = 10.0 cm) was infinitesimal for element size 0.1 cm (0.0009%) and acceptable for many applications at 0.5 cm (0.0289%).
  • The errors in electrical field near the electrodes were higher, as we would expect when representing a curved surface with a set of line segments. The relative error was higher on the inner electrode because there were fewer segments. For the 0.5 cm elements, the error of 2.7% on the inner electrode could be an issue for some applications.

 

Mesh geometries for small and large elements

Figure 5. Mesh geometries for small and large elements.

In general, I have found that inexperienced code users often employ far more elements than is necessary. For the determination of electric field values in the dielectric volume (removed from electrodes), the accuracy with the fine mesh of Fig. 5 may not justify the extra computational work. If you have concerns about mesh resolution, the standard procedure is to do two calculations with different mesh sizes and to check whether values in critical regions differ by more than the accuracy requirement. Note that there are two techniques to increase accuracy without greatly increasing the number of elements:

  • Variable mesh resolution.
  • Boundary method to generate a microscopic solution within a macroscopic solution.

These topics will be discussed in future articles.

We’ve covered a lot of ground in the last three articles by following a relatively simple electrostatic solution.

  • Organizing simulation data and using a file manager.
  • Defining region boundaries with the Mesh Drawing Editor.
  • Employing symmetry boundaries to reduce the size of a calculation.
  • Using a spreadsheet to plan and to analyze calculations.
  • Understanding the set of input and output files for an electrostatic solution.
  • Creating program input scripts in interactive dialogs and modifying them with a text editor.
  • Generating and inspecting standard mesh definition files with Mesh.
  • Defining run-control and material parameters for EStat solutions.
  • Running Mesh and EStat interactively or from the Windows command prompt.
  • Automatically controlling EStat with analysis scripts.
  • Creating batch files for automatic run control with the task generator of the TriComp Program Launcher.
  • Gaining insights into the relationship between element size and interpolation accuracy.

In the next article, we’ll turn our attention to 2D magnet-field calculations.

Footnotes

[1] Over the long existence of Open Office (starting from Star Office), I have been amazed that people will pay over $100 per computer for Microsoft Office when they can get a free product that will serve all their needs with considerably less aggravation. I suspect that the computer experience of the most people is governed more by anxiety than logic. The commercial program offers the same sense of security that you’d get from having the Mob protect your candy store.  At any rate, here is a link to get OpenOffice 4.1.1: http://www.openoffice.org/

2) Batch files are one of the most useful features of Windows, but many users do not realize they exist. These articles review how to apply batch files: Discovering your inner DOS, Part 1 and Part 2.

[3] Contact us : techinfo@fieldp.com.

[4] Field Precision home page: www.fieldp.com.

Test electrostatic solution: calculating and analyzing fields

In this fifth section of the online course Electric and magnetic field calculations with finite-element methods, we’ll use the mesh file discussed in the previous section to build a field solution. To get started, run EStat from the TriComp program launcher. Click the “1” tool and choose the Mesh file COAXIALCYLINDERS.MOU. The program determines the defined regions in the mesh and displays the dialog of Fig. 1. The values shown are appropriate to the solution parameters discussed in the previous article. Fill in the fields and click OK, accepting the output file name COAXIALCYLINDERS.EIN (EStat INput).

Dialog to create the EStat input script

Figure 1. Dialog to create the EStat input script.

To understand the action of the dialog, choose the File/Edit script (EIN) menu command and load the file in the editor. Here is the content:

* File: CoaxialCylinders.EIN
Mesh = CoaxialCylinders
Geometry = Rect
DUnit =   1.0000E+02
ResTarget =   5.0000E-08
MaxCycle =    5000
* Region 1: DIELECTRIC
Epsi(1) =   1.0000E+00
* Region 2: INNERELECTRODE
Potential(2) =   1.0000E+02
* Region 3: OUTERBOUNDARY
Potential(3) =   0.0000E+00
EndFile

You’ll recognize many of the entries. The ResTarget and MaxCycle parameters that control the solution accuracy are described in the EStat manual. The default values are fine for this application.

Exit the editor and click on the “2” tool. After you choose COAXIALCYLINDERS.EIN, EStat generates and solves the finite-element equations in less than a second. Larger meshes will take longer, but generally the run time of any practical solution is less than a minute. The program creates the files COAXIALCYLINDERS.ELS (a diagnostic listing file that you can inspect with an editor) and COAXIALCYLINDERS.EOU (a record of the mesh coordinates and potential values at the nodes).

To see the results, press the “3” tool and choose COAXIALCYLINDERS.EOU. You can generate interesting plots in the Analysis menu (Fig. 2). There are also useful tools like the movable field probe shown in the figure. For this discussion, let’s concentrate on hard numbers. First, let’s check the absolute accuracy of the solution with a single point calculation. At radius r = 10.0 cm, the formulas in the previous article give the values:

Er(r) = 910.23923 V/m
φ(r) = 36.90703 V

Numerical interpolations in EStat involve the collection of potential values from surrounding nodes. On symmetry boundaries, there are only half the available nodes, so the accuracy is not optimal. We should use an internal point for a good comparison. On a 45° line, the position r = 10.0 cm corresponds to x = y = 7.071068 cm. Click the Point calculation tool. We can specify the location by moving the mouse inside the solution volume and clicking the left button. This method is not accurate enough for this comparison. Instead, press the F1 key after clicking the Point calculation tool to enter coordinates manually. Set the x-y values as 7.071068 and click OK. Several calculated quantities are displayed at the bottom (Fig. 2). The calculated values of potential (36.90767 V) and electric field (910.23111 V/m) agree with the theoretical value to within a few thousands of a percent.

Filled-contour plot of potential with results of a point calculation displayed

Figure 2. Filled-contour plot of potential with results of a point calculation displayed.

Rather than copy numbers from the screen, why not let EStat write them for us? Click on the Open data record command and accept the default name COAXIALCYLINDERS.DAT. Now, every time you do an interactive calculation, the results are recorded in the text file. For example, choose the menu command Analysis/Volume integrals. You can inspect the resulting file with an external editor. To use the internal EStat editor, click the Close data record command first (two instances of a file cannot be open in the same program). Here is the field-energy calculation (volume integrals of the field-energy-density):

Quantity: Energy
Global integral:   1.708936E-07
RegNo   Integral
==========================
1    1.708936E-07
2    9.881184E-37

The field inside the inner electrode (Region 2) is numerically zero. We can use the field energy in the dielectric (Region 1) to find the capacitance per length. Because the calculation covers only the first quadrant, we need to multiply by a factor of 4.0 to find the field energy per length of coaxial cylinders at 100 V, Ue = 6.835744E-7 J/m. The equation C = 2*Ue/100^2 gives c = 1.367149E-10 F/m, within 0.8% of the theoretical value cited in the previous article (c = 1.3567E-10 F/m).

To understand the effects of element size, we want to compare accuracy at several radii for different choices of element dimensions. The procedure would be tedious if we had to run the Point calculation tool and type coordinates for each datum. In the next article, we’ll automate the analysis procedure.

Footnotes

[1] Contact us : techinfo@fieldp.com.

[2] Field Precision home page: www.fieldp.com.