Math4FEM is a finite element method (FEM) program for two dimensional static fields that runs on Mathematica 4 (Math4, Wolfram Research, Inc.) on MS Windows.

[Japanese]

This program was written by Wataru Mizutani based on the codes in 
`Mathematica ni yoru denjikai simulation nyuumon' by Dr. Kiyoshi Yoda, 
(Introduction to electromagnetic field simulation by Mathematica)
Morikita Shuppan Co., Ltd., Tokyo, 1997, Chapter2. [ISBN4-627-71520-XC3354]
For further information of the book, refer to the following URL:
http://www.morikita.co.jp/math.html (in Japanese only). 

This code can be freely used for non-commercial use; 
for commercial use, you need to contact Dr. Yoda through W. Mizutani. 
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

Math4FEM has no good GUI so far. 
Users of Math4FEM need to know basic skills of Mathematica. 
You have to input the coordinates and boundary conditions by hand. 
Any suggestions or tools for improving mesh generation are appreciated. 

Usage:

0. Copy math4fem.nb and fem-init.mx into the same folder. 
Run Mathematica 4, select the menu (file->open) and 
open math4fem.nb[Notebook file]. 
(This would be safer for Math4 to recongnize the NotebookDirectory 
than to run the program by double-clicking the nb file.) 

Initialize the program by excuting the cell starting with 
(* THIS CELL MUST BE EXCUTED ONCE IN A SESSION FOR INTIALIZATION *)
If the initialization does not work well, try to specify the name of 
the directory containing the initializing file, 
eg. if you copied the initializing file in c:\tmp, excute 
<< "c:/temp/fem-init.mx";


1. Input X-Y coordinate in pts and excute the cell (Shift-Return). 

eg. pts = Join[{{x,y}},{{x2,y2}},{x3,y3}}, ... ,{{xn,yn}}];

Don't put "," at the end of the list.

2. Excute (* Make Triangular Elements *) cell. 
ele is the triangular elements generated by pts. 

3. Set the boundary conditions in the cell starting with 
(* Set Boundary Conditions *)
For the electrode surfaces, set the voltages. 
v[n] is the voltages at the point pts[n]. 

eg. bc = {v[1] -> 0.0, v[2] -> 0.0, ..., v[n] -> 1.0};

Don't put "," at the end of the list.

4. Set relative dielectric constants for the elements, and 
excute the cell for calculation. 
It takes minutes with my pentium3, 880 MHz CPU. 

5. Enlarge the region of your interest. 

6. Read out the field at {{x, y},  {Ex, Ey}} 
Field {Ex, Ey} is calculated for the triangular element ele. 
{x, y} is the center of the ele. 

References:

[1] K. Yoda, `Mathematica ni yoru denjikai simulation nyuumon'
(Introduction to electromagnetic field simulation by Mathematica)
Morikita Shuppan Co., Ltd., Tokyo, 1997, Chap. 2 [in Japanese].

[2] T. Onuki and K. Nemoto, `Personal computer ni yoru yuugenyousohou nyuumon'
(Introduction to finite element method by personal computer)
Ohm sha, Tokyo, 1983 [in Japanese].

Download

Math4FEM は、有限要素法による二次元電界分布の計算プログラムです。 中核部分は 依田先生の「Mathematicaによる電磁界シミュレーション入門」(森北出版)の第二章がもとになっています。

水谷が拡張した主な点は、

1. 三角要素を追加モジュール(DelaunayVertices)を使わずに作製します (CDからインストールしたままのMathematicaで動きます)。

2. 節点データを入力し、境界条件を決めれば、境界条件の決まっていない節点を自動的に判別して計算します。

3. 三角要素の節点の値を元にして等電位線を描きます(メッシュの切り直し不要、計算領域全体の形状に制限はありません)。

有限要素法の計算では、節点と要素を使って材料の構造や材料定数を指定します。 二次元の電界分布を計算するので、節点は二次元座標で与えます。 各節点には節点番号を割り振ります。 次に近隣の三節点を結んで三角形の要素を作製します。 三角形の作り方、節点と要素の番号は(* Make Triangular Elements *)で始まるセルを実行すると表示されます。 例えば ele[[1]] を見ると {1, 2, 13} と3個の節点番号がリストされています。 不都合があれば、ele を直接変更して別の三角要素を作ることもできます。
要素の中心座標は centg[要素番号]、節点の電位は V[節点番号] で表しています。 各要素における電界ベクトルは varray[[要素番号]] として求められます。

テストに使っているデータでは、二枚の平行平板電極に±0.5Vをかけたときの電界分布を 計算しています(y=30 の電極に -0.5V、y=70 の電極に +0.5V)。
使っている境界条件では x=0 で左右対称になっています。 一番上の要素(y＞90)にはε = 2.2 (有機物)、 一番下の要素(y＜10)にはε = 3.82 (シリコン酸化膜) のように物質を設定していますが電極から離れているのであまり効果はみえていません。

【図を描くための記号の説明】
g1は、節点番号です。
g2は、三角要素で描かれたメッシュです。
g3は、要素番号ですが、単独では表示されないので、他のグラフィックオブジェクトと同時に表示してください。その場合、一番最後に指定します。

Show[g2,g3];

g6 = potential[19];(* 19 は描く線の本数+1 *)

【計算された電界強度について】
電界強度の絶対値は、長さの単位に対してリニアなので、座標をm(メートル)と考えれば、V/m となり、nmと思えば V/nm です。 電極間の距離を 40m と考えると、電圧差 1V なので、電極の間では電界強度は 0.025V/m になるはずです。計算値では 0.0249-0.0253 程度の値になっています。

不具合、御要望があれば、メールでお知らせください。