Sample Scripts from GB Books
GB003: 2次元における電導

ここで紹介するスクリプトはGunnar Backstrom氏の承諾のもと、書籍 “Simple Fields of Physics by Finite Element Analysis” に記されている多数のFlexPDE適用事例の中からその一部を紹介するものです。

PDF版 (1353KB)

GB002のケースでは金属オブジェクト全体が等電位の状態に置かれていました。しかしオブジェクトの異なる部分の間に電位がかけられた場合には電流が流れることになります。普通の金属の場合には電流密度（electrical current density）ベクトル J と電場ベクトル E、ポテンシャル U の間には次の関係が成立します。

(1)

ここに σ は電導率（electric conductivity）を表すパラメータです。一方、連続の式から

(2)

が導かれます。ここに ρ は電荷密度を表します。定常状態ではこの第2項が0となるため、(1)と(2)を結合することにより

(3)

というラプラス型の偏微分方程式が導かれます。ただし σ は空間的に変化する可能性があるため括弧の外に出すことはできません。しかしこのような場合でもFlexPDEでの扱いは容易です。なお、FlexPDEに対し方程式を定義する際、従属変数 U を明示的に含んだ(3)の形の数式を入力する必要はありません。div(J) = 0 と記述するだけで(1)により U に関する2階の偏微分方程式が自動的に誘導されます。この U に関する解から E や J が算出されることになります。

1. 矩形プレート中の電導

最初に最も簡単な例として右のような銅製の矩形プレートについて考えることにします。この上辺に 1.0V の電位をかけ、下辺は 0V に維持、一方、左右の境界は絶縁状態に置かれているものとします。今、下辺の位置を y = 0 におき、矩形のy方向の長さを Ly とすると U の厳密解は

という単純な数式で表されることになるので、有限要素法で算出される U の値との対比も行ってみます。
なお、ここでは2次元の問題として扱いますが、対象のオブジェクトに z 軸方向の厚みがあっても構いません。U(x, y, z) が z に依存しないと仮定できるなら、結局2次元の問題に帰着されるからです。

1.1 Problem descriptor [ econduction01a.pde ]

まずタイトルを設定します。
TITLE
    'Conduction in a Rectangular Plate'    { econduction01a.pde }

次に演算精度に関するセレクタをセットします。デフォルトは 0.002 なのですが、ここでは問題が単純なので精度を 10^(-5) とします。
SELECT
    Errlim = 1e-5

従属変数を定義します。
VARIABLES
    U                     { Electric potential }

偏微分方程式の定義に先立ち、パラメータ類を定義します。これらは境界の定義、及び境界条件の設定に際して使用されます。ここでは銅の電導率として 5.99e7 という値を使用します。
DEFINITIONS
    Lx = 0.5 Ly = 1.0
    cond = 5.99e7         { Conductivity of Cu }
    Ex = -dx(U) Ey = -dy(U) E = -grad(U) Em = magnitude(E)
    Jx = cond*Ex Jy = cond*Ey J = cond*E Jm = magnitude(J)
                          { Electrical current density }
    U_ex = y/Ly*1.0       { Exact solution for 1.0V }

方程式を定義します。従属変数 U の数式表現になっていませんが、DEFINITIONSセクションの情報に基づき、U に関する2階のPDEが自動的に導出されます。
EQUATIONS
    div(J)=0              { 2nd order PDE in U }

次に境界の形状と境界条件を設定します。上辺、下辺に対してはValue文を、左右の絶縁境界に対してはNatural文を使用して境界条件を定義します。
BOUNDARIES
    Region 1
      Start (-Lx, 0)
        Value(U) = 0    Line to (Lx, 0)
        Natural(U) = 0 Line to (Lx, Ly)    { Insulated, Ex = Jx = 0 }
        Value(U) = 1.0 Line to (-Lx, Ly)
        Natural(U) = 0 Line to Close       { Insulated }

最後に出力すべき情報を指定します。
PLOTS
    Grid(x, y)
    Contour(U)    Surface(U)
    Vector(E)
    Contour(Jx)    Contour(Jy)    Contour(Jm)
    Contour(U - U_ex)

END