静電気学
目次
はじめに
このモノグラフでは偏微分方程式を使って,時間非依存の電場,電位,電荷密度分布をモデル化して解析する.静電気シミュレーションにおける一番の前提は,電場と磁場が結合しないということである.
静電方程式は,静電荷によって生成される電場をモデル化するのに使われる.マクスウェルの方程式の簡約化として,これは断熱材,または誘電材料を記述するためにだけ使われる.
偏微分方程式(PDE)を使って電磁気装置をモデル化することだけが電磁気装置のモデル化の方法ではない.他の方法として常微分方程式(ODE)を設定するというものもある.Wolfram System Modelerはこのアプローチに従っている.簡単に言うと,System Modelerのアプローチは大規模な電磁気装置系の相互作用に適しており,偏微分方程式のアプローチは特定の装置のきめ細かい分析に適している.この2つのアプローチを組み合わせて使った方がよい場合もある.
最初のセクションでは,利用可能な静電気解析および機能を紹介するためにコンデンサを使う.その後で,根底にあるアイディアや概念をより理論的に説明する.静電解析の知識があると,理論的背景はかなり理解しやすくなる.その次に利用可能な境界条件について説明する.
この電気装置解析の目標は,指定された制約条件下における電位分布 (単位は)を求めることである.次のステップは電場の強度 等の二次的な場および静電容量 等の値を求めることである.これらの物理量の解析と解釈は考慮する装置の工学的設計をよりよくするのに役立つ.
そのようなモデリング過程は,結果としてNDSolveで解くことのできる偏微分方程式系になる.
電磁モデリングを拡張した例はモデルコレクションで見ることができる.
電磁装置解析は通常段階を追って行われる.まず物体を解析するために幾何学的モデルを作成する必要がある.幾何学的モデルは通常コンピュータ支援設計(CAD)の過程で作成される.CADモデルはインポートすることも製品内で作成することもできる.形状をインポートするために,DXF,STL,STEP等の一般的なファイル形式がサポートされている.これらの形状はImportを使ってインポートできる.別の方法として,例えばOpenCascadeLinkを使って製品内で幾何学的モデルモデルを作成するというものもある.
幾何学的モデルが利用可能になると,どのような種類の解析を行うかを考える必要がある.次のステップは境界条件と制約条件を設定することである.使用される材料がPDEモデルを指定する.PDEモデルが完全に指定されると,引き続き有限要素解析が調査中の装置の所望の数量を計算する.この数量は後処理される.Either by visualizing them or some derived quantities are computed. このノートブックはCADモデル生成以外のすべての必要なステップを示す.
電磁装置は通常三次元物体である.モデルでは結果として簡約された2Dモデルになる特殊な場合が存在する.事実,2Dモデルは電気工学の重要な問題の90%を解決するということが分かっている[Cardoso, 2018].
このチュートリアルを通して使われる記号および対応する単位は用語集セクションにまとめてある.
基本的な例
電磁気学のこのようなサブケースにおける有限要素法の使用法を説明するためには,簡単な例を提示し,設定の概要,さまざまな解析タイプ,可能な後処理のステップを説明するのがよい.
最初の例は簡単な電磁気学PDEモデルを設定するワークフローを導入する方法である.これは静電気解析をカバーする.
2枚の伝導性のある板とその間に誘電体を挟んだ3Dの円形平行コンデンサを使う以下の設定を考える.
3Dコンデンサの図と,高さ ,半径 の誘電体を青で示した横断面.
このモデルを簡単にしておくために,青い誘電体領域だけをモデル化する.
形状
コンデンサ装置をモデル化する場合,2種類の方法がある.一つは装置と周りの体積を考える方法である.このアプローチでは板の端に現れ,装置から周囲の体積に向かって外側に拡張する場をシミュレーションすることができる.この場は漏れ磁場と言われる.2つ目のアプローチは,ここで使うものであるが,アプローチは誘電領域を考える.どちらのアプローチを使ったらよいか,またモデルに対する効果はどのようなものかについての詳細は,静電容量セクションをご参照いただきたい.
覚えておかなければならないのは,幾何学的モデルで使われるスケールである.例えば境界メッシュの長さがメートルという単位ならば,材料パラメータは同じ単位で指定しなければならない.
3D幾何学的モデルの生成およびインポートに関する詳細は「OpenCascadeLinkを使う」というチュートリアルに記載されている.
材料パラメータ
次のステップは材料パラメータの指定である.一般に電磁モデルのパラメータはすべて必要なパラメータ値を含むAssociation pars で集められる.
このモノグラムではデフォルトの電気定数値 を使う.方程式をスケールするために定数が1等の異なる値を持たなければならない場合は"VacuumPermittivity"パラメータも指定することができる.
このモデルでは,コンデンサの誘電体は電気定数を持つ.この値は一般的な誘電材料を表す.
材料を直接Entityとして指定することは現在不可能である.
値がQuantityオブジェクトとして与えられる材料パラメータを使うことができる.必要な材料特性を指定する厳密な特性名はElectrostaticPDEComponentのリファレンスページで見付けることができる.
いくつかの材料からなる形状を使うこともできる.そのような例は「静電気学の複数材料」のセクションに示されている.
単位
形状の単位が材料の単位と異なる場合は材料単位をスケールすることができる.
内部的にすべての材料データの単位は"SIBase"単位に変換される.その結果長さのデフォルト単位は"Meters"である.形状の単位もメートルならば,何も変更する必要はない.形状の単位がメートルでない場合は,PDEか材料特性をその単位にスケールするか,形状を"Meters"にスケールするかする必要がある.PDEおよび材料パラメータの単位をスケールするために,"ScaleUnits"を与えることができる.明示的に記されていない限り,このノートブックの例ではデフォルトの"SIBase"単位を使う.
境界条件
問題を完全に指定するためには,境界条件を物理的境界で定義する必要がある.この場合,装置は平板間の電圧差を適用することによって帯電される.要するにこれは末端でDirichletConditionまたはNeumannValueをそれぞれ適用することで行われる.さまざまな境界条件を使うことができ,その詳細は「静電気学の境界条件」セクションに記載されている.この概要の目的では電位条件で十分である.
このセクションの目的は境界条件が適用される位置を構築することである.
境界条件が適用される位置を求める方法は,先の幾何学セクションで定義した輪郭形状とともに可視化することである.
誘電体の辺ではゼロNeumannValueになる.他の何も指定されていない場合,ゼロNeumannValue境界条件はデフォルトの境界条件なのでその境界条件は省略することができる.
より複雑な形状の場合は,境界条件の述語を指定する別の方法が適していることがある.これについては「球面コンデンサチュートリアル」で見ることができる.
メッシュ生成
有限要素解析を実行するためには,幾何学的領域をメッシュに離散化する必要がある.
メッシュ生成過程についての詳細は「要素メッシュ生成」チュートリアルで見ることができる.あるいはメッシュをインポートすることもできる.FEMAddOnsを使うと一般的なメッシュファイル形式のいくつかをインポートすることができる.
静電解析
静電解析とはモデルに課せられた静電荷分布によって起こる誘電体内部の電位および電場を調べるものである.従属変数は電位 である.
ここで は電気定数, は比誘電率である.この方程式はElectrostaticPDEComponent関数で提供される.
比誘電率 を設定することによって,関数はこの線形構成関係に基づいた方程式を生成する.偏極ベクトル 等の他のパラメータを指定することができる.
結果はInterpolatingFunctionオブジェクトであり,これは電位分布 を与える.
解法およびオプションに関する詳細は「有限要素のためのNDSolveオプション」チュートリアルに記載されている.
後処理
また,関心があるのは電場の強さ と電束密度 であり,これらから容量値 を得ることができる.これらの場は電位から回復される.容量値の計算方法は静電容量セクションに示してある.
3次元モデルの場合,関数Gradはそれぞれ3つの独立変数 ,, を持つ3つのInterpolatingFunctionのリストを返す.
電場のさまざまな成分にはPartを使ってアクセスすることができる.
電場の強さ を計算すると,コンデンサの静電容量を求めることができる:
ここで(単位は)はコンデンサの平板の一つで保持される電荷の大きさ,は平板間の電位差, は平板の一つにおける電束密度の法線成分である.
平板間の電位差を計算するためには,まず上の平板の電圧を計算する必要がある.上の電圧は円板領域上で を積分しそれを同じ領域で割ることで求められる.
この場合,上の平板の内向きの法線成分は負の方向 の 成分である. 成分は で与えられる.
を計算するためには,モデルで使われた電気定数値を抽出する必要がある.計算全体で同じ値が確実に使われるようにするために,デフォルトのモデルパラメータ計算を使う.
取得された数値的静電容量は,解析的な静電容量とほぼ一致する.
上の場合,円筒領域の上部として記号的領域表現Diskを使った.現実世界では,問題となっている境界の記号表現が常に利用可能だとは限らない.境界表現が可能であっても,計算中に使われた数値表現とは異なる可能性がある.これはメッシュがDiskのような記号表現の近似に過ぎないからである.このため,一般に問題となる領域の数値近似を直接使った方がよい.これにはFEMNBoundaryIntegrateを使うことができる.
FEMNBoundaryIntegrateを使う場合,所望の境界表面の要素マーカーを指定する必要がある.この場合2は上の平板の境界を表す.メッシュのマーカーについての詳細は「要素メッシュの生成」チュートリアルの「マーカー」セクションでご覧いただきたい.
数値メッシュの誤差の方が記号的表面領域よりも大きいが,これは許容範囲である.表面の記号的表現が利用できないときでもFEMNBoundaryIntegrateは通常利用可能なワークフローなので,これは便利である.
別の方法として,静電エネルギー を使って静電容量を計算することができる:
静電容量を計算するために使われたメソッドについての詳細は「静電容量」セクションに示してある.
静電方程式
静電気学は,静電荷によって生成される電場を扱う電磁気学のサブフィールドである.したがって,生成された電場の原動力は分布電荷あるいは指定された電位差(単位はボルト)である.マクスウェルの方程式から派生した静電方程式は,断熱材料または誘電材料内の静電荷によって生成された定常電場のモデリングに使われる.ここでは誘電体内の方程式の最も一般的な形式を示す:
この方程式の従属変数は電位 であり,位置 により変化する.この方程式は拡散係数として電気定数 を持つ拡散項,および変数 が偏極ベクトルである導関数項からなる.項 (単位は)は領域内の体積電荷密度を表す.これについては「ソースタイプ」セクションで説明する.
領域内の金属成分はこの方停止いないでは考慮できない.導電性材料では電位差がないため,この方程式は誘電材料または断熱材料だけのためのものである.例えば電極は境界条件としてモデル化しなければならないのである.また,領域内のすべての金属は取り除く必要がある.
静電方程式は1D軸対称モデル,1Dモデル,2Dモデル,2D軸対称モデル,3Dモデルについて解くことができる.
モデルの設定
2Dモデルでは,電位が 方向と 方向だけで変化し 方向では一定である特殊な対称性が存在する.この対称性によって, 方向の厚さを表す厚さ変数 が加えられるという点が異なる同様の方程式を解く:
平面でシミュレーションを行うと対称性があるように,1Dモデルにも対称性がある.1Dの静電モデルでは, 方向と 方向を表す断面積 が方程式に与えられる:
1Dおよび2Dの軸対称モデルは,回転軸を持つ3D立体で見られる別のタイプの対称性を考慮に入れる.これらの場合の詳細は,セクション「軸対称モデル」に記載してある.
このモデルの定義には非アクティブなPDE演算子が使われる.「偏微分方程式の数値解法」には非アクティブな演算子の使い方を説明するいくつかのセクションがある.
静電方程式
静電という名前から分かるように,このサブケースは時間で変化することのない電場,つまり静電気と言われる場を扱う.簡約化のために,真空をシミュレーション領域と仮定する.
静電の場合は2つのマクスウェルの方程式,ガウスの法則で記述される:
ガウスの法則(方程式1)は,空間の任意の点における電場の強さ の発散はその点における体積電荷密度 に比例するというものである.つまり電荷密度はその領域の の挙動に影響を及ぼすのである.電場 の発散は,空間内の任意の点で電場がどのくらい広がっているか,または収束しているかを示す.
その一方,ファラデーの法則(方程式2)は, は回転がないため保存場であるというものである.この方程式を使うとスカラー電位 を導くことができる.
これは任意のスカラー関数 の勾配の回転はゼロと述べている.したがって,両方の方程式を比較し,方程式3に電場 を代入すると, は関数 の勾配として表せるということが推定できる.
ここで方程式4を方程式5に代入すると,自由空間の静電方程式が得られる:
誘電材料と構成方程式
上の方程式は真空中で有効である.真空中の静電方程式は誘電材料を考えるために拡張することができる.誘電材料は自由電子を持たないことで伝導体と区別される.誘電材料の電荷は力で拘束される.外部電場が誘電材料に適用されると,原子の中の電子雲が 場の反対方向に移動され,原子核が平衡位置から 場が指すのと同じ方向に移動されるため,双極子モーメントが誘導される.
この電荷の移動により拘束体積電荷密度 が生じ,その結果永続的あるいは誘導された電気双極子モーメントの密度を表すベクトル場 が生成される.この偏極ベクトルがあるため,誘電材料内部の電場は自由空間の電場と異なるのである[Sadiku, 2011].
左:非極性原子.右:電場 が適用されたときの非極性原子の極性形成.
偏極ベクトル と拘束体積電荷密度 は以下の方程式によって関連付けられる:
これはガウスの法則に似ている.負の記号は,外部 場によって誘導される誘電体の整列した双極子の電荷と逆の符号であるためである.
誘電体のこのような拘束電荷が,その誘電体に適用された外部電場とどのように相互作用するのかが知りたい.自由電荷と誘電体自体を構成する拘束電荷を見分けると役に立つ. で表される自由電荷は,コンデンサの平板のように伝導体に導入するか伝導体から取り去るかすることができる.一方,拘束電荷は誘電体内で移動できない.
したがって,自由電荷と拘束密度電荷の両方を持つシステムでは,ガウスの法則を変更して全電荷密度を含むようにする必要がある:
これは電束密度という新しい数量 を導入することでさらに簡約することができる:
この構成関係は基本的に,外部電場 を誘電材料に適用すると電束密度が自由空間内よりも大きくなるということを述べている.自由空間では である.
線形誘電材料
線形材料では, は電場を使って と記述することができる.ここで無単位の は材料の電気感受率である.したがって線形材料では構成関係は以下になる:
また,は絶対誘電率である.自由空間では である.一方, が と線形に変化しない場合,材料は非線形と言われる.
比誘電率は空間座標 に依存し,考慮される領域内で変化(不均質材料等)したり,方向付きで変化(異方性材料)したりする.
複数材料
シミュレーション領域でさまざまな材料をモデル化することはよくある.例えばコンデンサを周囲空間にまで拡張する電場のシミュレーションを行う場合等がある.この場は漏れ磁場と言われる.コンデンサの漏れ磁場のシミュレーションを行うためには,誘電領域とその周囲の空気を含む必要がある.
ここで,電極が 軸に沿って両側に拡張している( 平面の外)と仮定して,並列コンデンサの横断面をモデル化する.このモデルは後のセクションでエネルギー,静電容量,力を計算するときに役立つ.
また,この例では要素マーカーを使って各部分領域の比誘電率の値および境界条件を指定する.要素マーカーとメッシュ内でのその使用法は「要素メッシュの生成」チュートリアルの「要素マーカー」セクションで説明してある.
次にToElementMeshの"RegionMarker"オプションで領域マーカーを指定する.そのためには誘電体および空気領域の座標を整数マーカーが必要である.オプションとして,部分領域を調整するために追加の最大セル寸法を指定することができる.複数材料領域をもっと正確で効率的に扱う方法は,要素マーカーを使うことである.これでメッシュがその材料に対する特定の部分領域を持つため正確な解が得られるのである.
金属が解析領域内にある場合,それはメッシュから除去しなければならない.
効率的なPDE係数の設定方法はこちらを参照のこと.
このコンデンサでは上の電極の電位が,下の電極は接地されている.これらの境界条件の他,箱の境界にゼロ表面電荷密度が適用される.これは電場が箱内部にとどまるよう制約する.
ゼロ表面電荷密度はノイマンのゼロ値と同等である.これらは自然なデフォルトの境界条件であるため省略することができる.
このモデルでは,境界条件指定でマーカーを使うことができる.これらはDirichletCondition境界なので,点要素マーカーが使われる.
電場 と を正しく計算するために,異なる材料の接触面では成分のうちの一つが非連続になる可能性があることを考慮に入れる必要がある.接触面の境界条件についての詳細は「静電学の境界条件」を参照されたい.
この非連続性を解決するために,関数DiscontinuousInterpolatingFunctionまたはEvaluateOnElementMeshを使う.これらの関数は部分領域の一つからの値に優先順位を付ける.デフォルトでは"MeshElementMarkerUnion"で与えられる順番が使われる.優先順位を指定したい場合はオプション"MarkerPriority"を与える.
誘電体と空気の接触面では,この場合のデフォルトは空気よりもコンデンサの誘電体領域の優先順位の方が高い.この例では の非連続性は接触面で見られる.
電場の強さ は電位 の勾配なので,まず電位関数でDiscontinuousInterpolatingFunctionが使われ,新しい非連続関数から導関数が計算される.
DiscontinuousInterpolatingFunctionは2つの異なる材料の接触面における値を計算し,それからプロットに使ったり他の計算の引数として使ったりするのに便利である.
出力はDiscontinuousInterpolatingFunctionとなる.これはまず領域の範囲を示してマーカーの優先順位を示す.
このベクトルプロットは境界におけるゼロノイマン値が場の線を領域内に留めるよう制約しているのを可視化する.換言すると,場は外部領域の接線となるよう強制されるのである.この場の制約は場の中央部に直接影響し,結果的に静電容量値にも影響を及ぼす.これは「静電容量」セクションで計算し説明する.
このモデルは比誘電率 に対して2つの異なる値を持つので,EvaluateOnElementMeshを使って異なる材料領域における電束密度 を正しく計算する必要がある.この電束はDiscontinuousInterpolatingFunctionでも表せる.
EvaluateOnElementMeshは非連続の補間関数のリストを返す.これらの補間関数は独立変数 および を持たない.Throughコマンドを使うとこれらを加えることができる.
電束非連続性を見るために,電束密度の大きさのDensityPlotを作成することができる.プロットされた変数の最大値と最小値を抽出する必要がある場合がある.まず,場の"ValuesOnGrid"を抽出することを許可する規則を作成してからMinMaxを適用する.
異方性誘電材料
誘電体は等方性でも異方性でもあり得る.等方性材料はすべての方向で同じ電気特性を持つ材料である.等方性の では, と は という関係にあるため,常に同じ方向を持つ.異方性材料の場合,比誘電率はテンソル になり,電気特性は方向とともに変化する.例えば結晶材料は異方性である[Sadiku, 2011].
異方性材料では または はテンソルになる.3Dの場合,比誘電率は9成分を持つ3×3のテンソルである.
ここで は比誘電率テンソルである.およびはそれぞれ主誘電率係数,非対角誘電率係数と言われる.
比誘電率がテンソルなので, と の関係は に変化する.このため と はもはや並列ではなくなる.
次の例では電気光学変調器(EOM)をモデル化する.このEOMはニオブ酸リチウム()のコア,シリカクラッド(),コンデンサの役割を果たす2本の金電極からできており,コアに電場を誘導する.EOMは多数の光学装置で広く使われている.
ニオブ酸リチウムのコアは光が通り抜ける光導波路である.この材料は直交異方性材料である.これは主方向 および について対称であり,非対角誘電率係数はゼロである.
このEOMの設計では台形の導波路がの平板上にあり,2本の金電極がコアの上にあり,これらがシリカクラッドに埋め込まれている.この設計は[Loucas, 2023]に基づいている.
電気光学変調器の断面図.電極は黄色,シリカクラッドは灰色,LiNbO3の平板と導波路は紫で描かれている.
このモデルでは,電極によって誘導される電場 のシミュレーションを行う. 軸を中心とする軸対称性とEOMが十分長いという仮定のため,電場には への依存性はない.このためモデルを簡約し,EOMの横断面だけをモデル化することができる.
2つの異なる材料がモデル化されるので,領域マーカーを使ってその部分領域を定義する.
電極間の電位差はである.右の正の電極の境界が,左の負の電極がである.
このモデルでは,境界条件の指定にマーカーが使われる.これらはDirichletCondition境界なので,点要素マーカーが使われる.
非線形誘電材料
非線形材料の場合,電場の強さ と電束密度 の関係は次のように変化する:
ここで は電場が存在しないとき分極であり,残留分極ベクトルと言われる.
この方程式は Pr によって引き起こされる追加の分極だけでなく線形関係 も考慮に入れている.簡単に言うと,残留分極は外部電場がない場合でも,全体的な電束密度に寄与しているということである.
方程式10で表された構成関係は,強誘電体等のヒステリシス効果を示さない非線形材料に限られる.
軸対称モデル
3D立体が回転軸を持つ場合,軸対称シミュレーションを行うことができる.軸対称モデルは円筒座標の代りに独立変数 を持つ2Dの切頭円筒座標系を使う.この系は 軸について回転対称なので,円筒座標変数 は消失する.
軸対称モデルには,時間とメモリについての計算コストが完全な3Dモデルを解く場合よりも格段に少ないという利点がある.
ElectrostaticPDEComponent関数は静電方程式の軸対称形式を生成することができる.そのためには,パラメータ"RegionSymmetry"を"Axisymmetric"に設定する.
この場合,電位は 方向で一定である.このことは,電場が 平面に対して接線方向であることを示している.2D軸対称モデルについては,「同心球コンデンサモデル 」で見ることができる.
1D軸対称形状では,電位は 方向だけではなく 方向でも一定である.線形材料に対する1D軸対称静電方程式は以下で与えられる:
円筒コンデンサを1D軸対称領域でモデル化する.円筒コンデンサでは,電位は 方向でのみ変化する.このコンデンサは半径 の同軸円筒シェルで囲まれた,半径 ,長さ の立体円筒状導体からなる.立体導体の電位は でありシェルは接地している.
立体とシェルの間の空間は比誘電率1である.コンデンサの長さは 方向なので,この装置の長さは厚さ を指定することに等しい.
モデルを検証するために静電容量を計算し,それを理論的値と比較する.
ソースタイプ
静電気源は電場に影響を及ぼす電荷密度である.静電気源は電荷分布に基づいて分類される.ここでは体積源,線源,点源に分ける.特別な場合,2Dモデルの点源は3Dモデルの線電荷に等しい.
静電方程式の右辺のソース項 は内部電束密度場 を生成するソースとして使われ,のときはソースとして,のときはシンクとして振舞う.
メッシュがソース項 の幾何学的境界に一致していることが重要である.このための最適な方法は明示的にその境界のためのメッシュを生成することである.解法として使われる有限要素法は,個々の要素が材料境界を横断しないメッシュの恩恵を受けている.別の方法としてMeshRefinementFunctionを使うというものがある.この場合,要素は材料境界を横断するがメッシュが十分細かければわずかである.両方の方法を組み合せることもできる.
体積電荷密度 または体積源は,領域内の任意形状のソースまたはシンクをモデル化するために使うことができる.体積電荷密度の値は常にの単位で指定される.この名前は体積電荷の3D版から来ているが,他の次元でも使われている.
である1Dの場合を考えてみよう.体積電荷密度 をの単位で指定すると,その値に単位の横断面積 が掛けられ,最終的な単位はとなる.
である2D領域でも同様に,単位の体積電荷密度 に単位の厚さ が掛けられ,最終的な単位はとなる.
次の2Dの例では,2つの反対の値 を持つ円形体積源 を使って電極をモデル化する.
この例では,ソース領域を指定するために関数RegionMemberを使う.別のオプションとして,要素マーカーを使ってソースを指定する方法があり,こちらの方がより正確な解となる.
効率的なPDE係数の設定についてはこちらを参照のこと.
このプロットで,双極子の等電位面とそれぞれの電場 を見ることができる.電場はどこにおいても等電位線に垂直である.
点電荷 または点源は のときの内部ソースまたは のときのシンクをモデル化するために使うことができる.点電荷は空間拡張がないと考えられる.点電荷は3Dモデルまたは2D軸対称モデルで使うことができる.
2D軸対称モデルでは,点電荷は対称軸上に指定されたときだけ存在することができる.軸対称モデルの回転の様子が考慮されると,対称軸の1点が線電荷になる.
点電荷は平面外の線電荷を意味するので,2Dモデルでは使われない.これについての詳細は「線電荷」セクションに記述してある.
2D軸対称領域の対称軸(灰色)上に指定された点電荷(赤).軸対称領域を回転しても点電荷が線電荷に変換されないことが分かる.
点電荷 の単位は常にであり,体積電荷密度の単位はである.通常以下になる:
ここで はディラックのデルタ関数である. はソースの位置Xのそれぞれの方向で単位を持つ.これによって使うことのできる電荷密度 の式が導ける:
しかし,ディラックのデルタ関数は離散化された空間領域で解かれないため,数値シミュレーションにおいて問題を引き起こす.これはディラックのデルタ関数がソース位置 において特異だからである.2つ目の問題は係数の評価が常にメッシュ要素内で起こり,決して辺では起こらないという点である.したがって,ディラックのデルタ関数の近似が必要になる.ディラックのデルタ関数を近似する過程は正則化と呼ばれる.
さまざまな正則化デルタ関数 が利用可能である[Peskin, 2002][Bilbao & Hamilton, 2017].このチュートリアルでは,次を選ぶ:
ここで は正則化デルタ関数 のサポートを制御する正規化パラメータである.通常 はメッシュ間隔 に相当する大きさである. は差分 を表す.
この概念は例を使うとよく分かる.その前にヘルパー関数を作成し,正則化ディラックのデルタ関数を計算する.
正則化デルタ関数 を使うために,正規化パラメータがメッシュ間隔の値の9倍になるよう選ぶ().
正則化パラメータ は,積分点がソース点の領域にあるように,十分大きくなければならない.次の可視化は,メッシュおよび積分点内の点源の位置における半径 の球を示している.
上の図はに位置する半径 の赤い球で表される点源を示している.青い点は各メッシュ要素の積分点を表す.積分点はPDE係数が評価される座標である.点源は,数値積分ができるよう,十分な積分点が内部に含まれるような十分な拡張 が必要である.
軸対称モデルの概念は前と同じであり,電荷密度体積の式は以下で与えられる:
ここでディラックのデルタ関数は対称軸のソース位置 において単位を与える.ただし は 平面上にあり,0から2 の任意の数である.軸対称の場合を考えるので, を指定する必要はない.
線電荷または層ソース は薄すぎてモデル形状の厚さが保てないソースまたはシンクをモデル化する.線電荷は3D領域,2D領域,2D軸対称で使うことができる.
3D領域では,線電荷は形状に任意に浮いている線を含む,形状の辺に指定することができる.
2D軸対称領域では,2つのオプションが使える.対称軸で指定される線電荷と対称軸を中心とするシミュレーション領域の回転を介した長さになる点として指定される線電荷である.
2Dでも,点は平面外の方向におけるモデルの厚さ のため,点は線電荷を表す:
3Dモデルでは線電荷は形状の任意の辺で指定することができる.
ここで はソース位置 におけるディラックのデルタ関数であり,単位を与える. は線電荷の値である.
平面外の線電荷をモデル化するためには点をソースとして指定することが必要である.
点にはすべての方向における空間拡張がないので,、モデリング領域の各次元(例:)にディラックのデルタ関数 を適用する必要がある:
ここで はを単位とするソース位置 のディラックのデルタ関数, は線電荷の値, はを単位とする厚さである.
次の2Dの例では,二次元双極子が位置の 軸上の電荷密度を持つ2つの線で構成されている.線電荷の値はである.
対称なので形状の四分の一だけを考える.左境界では地電位が適用される.地電位 は反対称境界条件を指定することに等しい.この境界では電場 の接成分はゼロである.上と右の境界でも地電位が適用される.下の境界では対称境界条件が適用される.
正則化デルタ関数 を使うために,正則化パラメータがメッシュ間隔の四分の一 になるようにする.
この問題では電位 に対する解析解が存在し,以下で与えられる:
2Dの軸対称モデルでは,線電荷は対称軸 上に指定することができる.体積電荷密度 は以下で与えられる:
ここで はソース位置 における単位のディラックのデルタ関数, は線電荷の値である.ただし は 平面上にあり,0から2 の任意の数である.軸対称の場合を考えるので, を指定する必要はない.
増分量
シミュレーションが完了したら,モデルがより深く理解できる増分量を導くことができる.この増分量には静電容量,力,トルク,エネルギー関係のパラメータが含まれる.これらの数量は,静電気系の挙動をより深く説明し,広範な応用範囲での設計及び解析処理を助ける.
静電エネルギー
静電エネルギー は荷電粒子の電場に保存されるエネルギーであり,以下の方程式で表すことができる:
を持つ線形誘電材料の場合は,静電エネルギーは以下で与えられる:
上の並列コンデンサの例では,場EfieldC,DfieldCを計算した.このセクションでは並列コンデンサのエネルギーを計算する.
静電容量
静電容量 は電荷を保存するための装置の特性であり,装置に見られる電位差に対する電荷の大きさの比率で表される.コンデンサでは,静電容量は以下で与えられる:
ここではコンデンサの平板の一つに保存される電荷の大きさであり,は平板間の電位差である.
「複数材料」セクションでモデル化した並列平板コンデンサの静電容量を計算するためには,まず平板の一つに分布する電荷を知る必要がある.電位差は既知の変数である.
電荷 は表面電荷密度 として分布される. から を得るために面積分を行う:
積分形式のガウスの法則から,表面電荷密度は変位場の法線成分に等しい()ことを導くことができる.したがって以下の積分を解く:
場 は法線方向(この場合 方向)で金属板から出たり入ったりする.場 が金属板を出るところを完全に捉えるためには,からまでの直線上で を積分する.
上の並列コンデンサの例では,場EfieldC,DfieldCを計算した.このセクションでは並列コンデンサのエネルギーを計算する.
例えばコンデンサへ流れてくる電流を指定した場合,電極の電圧 はFEMNBoundaryIntegrateまたはNIntegrateを使って求めることができる.電圧は電極の境界上で積分され,その面積で除算される.
上の場合,記号領域表現を使うこともできた.現実世界のアプリケーションでは,記号表現は常に利用可能だとは限らないし,利用可能だとしても数値表現とは異なる可能性がある.これはメッシュが形状の記号表現への近似に過ぎないからである.このため,問題となる領域の数値近似を直接使った方がよい場合がある.これはFEMNBoundaryIntegrateを使って直接行うことができる.
FEMNBoundaryIntegrateを使うとき,所望の境界面の要素マーカーを指定する必要がある.この場合3は上の平板境界を意味する.メッシュのマーカーについての詳細は「要素メッシュの生成」チュートリアルの「マーカー」セクションをご覧いただきたい.
最初に気が付くのは理論値が数値よりも低いということである.この差は理論方程式が装置を囲む漏れ磁場を考慮していないがモデルではされていたから生じたのである.さらに境界に課された制約条件が静電容量の値に影響を及ぼしている..
装置を囲む空間が考慮されるとき,最も正確な静電容量の値を得るためには無限空間をモデル化するべきである.これを解くためには2つの方法がある.静電容量の値が特定の値に収束するまで領域を拡張する方法と,無限領域を使う方法である.
エネルギー方法
静電容量の値計算する別の方法として,以下で与えられるコンデンサの静電エネルギー方程式を使うというものがある:
ここで は静電容量,は電位差である. について解くと,以下が得られる:
エネルギー法は,どの境界で積分するかを指定する必要がないので,静電容量の値を計算するのにより適した方法である.
静電気力
このセクションでは,物体間の静電気力の計算方法を示す.静電気力の計算はさまざまな工学および科学分野においていろいろな理由で必要なものである.例えば,多くの装置が制御や作動について静電気力に依存している.センサーも静電気力の変化に依存しているため,設計が正確さを確実にするために精密なモデリングが必要となる.
ここでは2つの方法を示す.最初の方法は, から導出した方程式を使って金属面間の力だけを計算するものである.この方法は2つ目の方法よりも少し簡単である.2つ目の方法はより一般的なもので,誘電体や金属の間の力を計算するのに使うことができる.この2つ目の方法にはマクスウェルの応力テンソル を使う.
クーロン力
電荷と電荷の間にクーロン力があるように,連続電荷分布の間にも力がある.例えば表面電荷密度では,単位面積あたりの力は以下で表せる:
このセクションでは並列コンデンサの平板間の静電気力の 成分を計算する.この計算では法線の電束密度場は 成分の になる.
あるいは,問題となっている線の記号表現が利用できない場合FEMNBoundaryIntegrateを使うことができる.
モデルの形状がより複雑な場合は,BoundaryUnitNormalを使って任意の表面の法線成分を計算することができる.
力の計算では,BoundaryUnitNormalとFEMNBoundaryIntegrateを組み合せて使わない方がよい.内部境界における法線の計算と向きが十分定義されていないからである.
ここで は平板の面積である.並列コンデンサでは,方程式は以下に還元される:
ここで は で与えられる. は平板間の距離である.コンデンサの比誘電率は値 を持つ.
マクスウェルの応力テンソル
静電気力を計算する別の方法としてマクスウェルの応力テンソルがある.これは誘電体にも静磁気力にも適用できるより一般的な方法である.
マクスウェルの応力テンソル を使うと,力は一般に以下のように表すことができる:
ここで はマクスウェルの応力テンソル, は単位法線ベクトル, は微分表面, は領域表面境界, は領域体積, は微分体積, はポインティングベクトルである.
ここで は E の大きさ, はクロネッカーのデルタである.クロネッカーのデルタは指標が同じ場合は1,それ以外の場合はゼロである.
静電気の場合,時間依存がないため方程式は以下のように簡約される:
静電気学では,成分ごとのマクスウェルの応力テンソルは以下のように簡約される:
最初の指標 は力の向きを示し,2つ目の指標 は表面の向きを示す.
3Dの場合は空間変数があるので,マクスウェルの応力テンソルには9要素あり,以下のようになる:
マクスウェルの応力テンソル は表面に作用する単位面積当たりの力である.与えられた表面に対する実際の応力 を計算するためには,以下のような所望の表面法線のドット積を取る:
方程式の導出をすべて見たい場合は,[Griffiths, 2017]を参照されたい.同等の定式は[Ida & Bastos, 2013]である:
方程式11は応力 の別の計算方法であり,平板間の静電気力を計算するために使った.力が 方向なので, だけが必要である.
誘電体への力
次の例では誘電体棒および金属棒に働く力を計算する.両方のメソッドを使ってそれらを比較する.これらの棒は真空中にあり横の平行電界に晒されている.ニュートンの第三法則の結果,力は等しく反対向きである.
このセクションでは力は境界で計算されるので,特別なメッシュを作成して力の計算での正確さを期す.領域の内部境界での調整が必要なので,境界メッシュを手動で生成する.さらにNumericalRegionを使うと,最小要素の高品質二次メッシュが生成される.力の式を積分するときNumericalRegionを使うことで他の利点があり,より正確にするためにNIntegrateを使うことができる.
境界メッシュを手動で生成するためには,領域の各部分を別々に生成してから新しい境界メッシュでそれらを組み合せる.
この件についての詳細は「要素メッシュの生成」チュートリアルのセクション「数値領域」で見ることができる.
マーカー1,2,3はそれぞれ金属棒,誘電体棒,ボックスの境界に割り当てられた.点マーカーは同じID番号に従う.
ElectricFluxDensityValueを使って平行電界を設定する.詳細は「電束境界条件」セクションを参照のこと.
誘電体と別の材料の間の接触面で力が計算されるとき,計算される場に注意する必要がある.接触面では電場の強さ と電束密度場 は不連続性を示す可能性がある.
この不連続性を解決するために,DDiscontinuousInterpolatingFunctionまたはEvaluateOnElementMeshを使う.これらの関数は部分領域の一つからの値に優先順位を付ける.デフォルトでは"MeshElementMarkerUnion"で与えられる順序が使われる.優先順位を指定するためにはオプション"MarkerPriority"を与える.
この場合,誘電体領域を優先するデフォルトの優先順位と真空領域を優先する"MarkerPriority"の両方を使うことができる.
まず金属棒にかかる力の 成分を計算する.金属表面を扱っているので,方程式13が使える:
金属棒と真空との間の接触面では不連続性はない.金属棒は領域内にあると考慮されないためどの場が使われても関係ないからである.
あるいはマクスウェルの応力テンソルを使って を計算することもできる:
誘電体に適用される力を計算したいので,マクスウェルの応力テンソル法を使う必要がある.
誘電体と真空の接触面における誘電体への力の 成分を正確に計算するためには,誘電体の値より真空の値を優先させる.
誘電体上でこれらの力を生成する場は,誘電体を取り囲んでいる場である.したがって,力を正確に計算するためには,誘電体を囲む表面上の場の値,つまり領域の真空側の値を知っておく必要がある.
この値は反対の記号であり,予想通り先に計算した値に非常に近い.
誘電体の力を計算するための別の式[Küpfmüller, 1968]が以下である:
ここで下付き文字 は誘電体を, は真空を意味する. は場の接線成分, は法線成分である.
マクスウェルの応力テンソルを使う利点は,これが複数の誘電材料にも,1つの誘電体にも,非誘電体にも使えるという点である.
法線力の 成分だけを計算するので,式には完全な単位法線 ではなく,法線と 軸の間の角度 の方向余弦を掛けるだけである.
静電気学における境界条件
静電気分野の境界条件は2つのカテゴリに分けられる.一つは伝導体の電位を指定するディリクレ条件である.その例として,電極が電位に固定されている伝導体がある.この場合,通常一方が地電位 に固定され,もう一方がゼロ以外の異なる特定の値に固定される.
もう一つは,変位場 の法線成分の値を指定する境界がある.これらの境界はノイマン値境界条件で実現される.
この2つのタイプの境界条件に加えて,周期的境界条件として知られる3つ目のタイプがある.この境界条件は,形状に例えば多極回転電機等の反復部分を持つ場合や領域の拡張を減少させたいときに便利である.このタイプの境界条件は境界の一部で電位 を指定して別の部分と同じようにする.
この3つの一般的な境界条件のタイプでは,以下の静電気学特有の境界条件が導入される:
2つの異なる誘電材料を分ける内部境界では,以下の条件を満足しなければならない:
方程式14と方程式15は, が法線成分を, が接線成分を表すという新しい用語を導入した方が説明しやすい.
これは表面電荷密度が境界に故意に置かれると,電束密度の法線成分は不連続になるということを意味する.置かれなければ電束は境界全体で連続である.
それでも電場は連続ではない. に を代入すると以下が得られる:
これは の法線成分が境界において不連続であることを示している.
これは電場Eの接線成分が内部境界に沿って連続であることを意味する.しかし に を代入すると,以下の方程式が導ける:
電位境界条件
電位境界条件の目的は境界のある部分に特定の電位を設定することである.
電位境界条件は通常電極のような金属接触をモデル化する境界の一部で適用される.もとんどの例において,一つの電極は値 で固定され,もう一方は で接地している.
反対称境界条件は,関数 が軸に対して反対称性を示すと考えられる境界をモデル化する.つまり,片方で調べている関数の値が軸の反対側でその負の値に等しいということである.
静電気学では,電場はこの反対称挙動を見せることがある.例えば双極子源からの電場の場合である.距離 によって分離され 軸に沿っている正と負の電荷があるとする.以下の図から,反対称軸を表す赤い破線について電場が反対称挙動を示していることが分かる.
反対称境界条件を理解するよい方法として,対称線に沿った鏡を考えてみるとよい.鏡を見ることが対称境界条件を表す.反対称の場合,鏡像は転置する.
上の図は双極子からの電場 を示す.ここで青と赤の点はそれぞれ正と負の電荷である.反対称軸の右側と左側には反対称性が存在する.である.
電場の接成分は反対称軸ではゼロとなり,電場は法成分しか持たないということを知ることは重要である.したがって,反対称境界条件は以下で与えられる:
そしてこれは,定電位の勾配がゼロであるため,標準定電位 が境界 で設定されると実現される:
対称境界条件が軸について鏡面対称である境界をモデル化しモデルのサイズを減少させるために使われるように,反対称境界条件も以下の図で示されているようにシミュレーション領域のサイズを減少させるために使うことができる.
電気双極子のシミュレーション領域を減少させ,正電荷だけをモデル化する.地電位は境界 で設定され,反対称挙動 を強制する.
電束密度境界条件
電束密度境界条件の目的は,境界における電束密度場の法線成分を定義することである.換言すると,指定された境界に電束条件を与えるために使われる.
電束密度境界条件は2つの異なる方法で指定することができる.一つは表面電荷密度 を指定することで,もう一つは境界電束密度 を指定することである.
ここで の定義は で与えられ, を代入すると が得られる.方程式18に負の符号を加えて場 が方程式に代入できるようにする.
このノイマン値はこのノイマン値は静電方程式に見られるように, の発散に直接関連している.
正の値は内向きの電束,負の値は外向きの電束を表す.言い換えると,電束は外向きの法線 と逆向きに指定される.
表面電荷密度は,2つの非導電材料間の外側の境界または内側の境界で指定することができる.2つ目の場合,境界条件は以下をモデル化する:
境界電束密度 と境界単位法線 を使うと,境界条件は以下で与えられる:
電束密度がゼロ値のベクトルとして指定される場合,または表面電荷密度がゼロ値として指定される場合は,境界には表面電荷密度はない.これはノイマンゼロ境界条件に等しい.指定された境界で他の境界条件が指定されていない場合,ゼロ電束が暗示的なデフォルトの境界条件となる.
内側境界では,ゼロ電束値は境界を通過できる電束はなく,電位は境界全体で不連続であることを意味する.
以下の例では,上の平板に表面電荷密度,下の平板に地電位を指定することによってコンデンサが帯電する.
次の例では円形断面のガラス棒が外部場に晒されている.電束密度の法線成分が,ガラス棒を取り囲む箱の右および左の境界に適用される.
地電位が後で適用されるため,の点をメッシュに加えた.また異なる比誘電率を割り当てるためにメッシュに領域マーカーを指定した.
対称境界条件
対称境界条件は,計算領域および想定される電位がシミュレーション領域に軸を中心にして鏡面対称の場合に使われる.
境界のある部分で境界条件が設定されていない場合,暗示的なノイマンゼロ境界条件が使われる.対称境界条件はゼロ電荷密度境界条件に等しい.
周期的境界条件
周期的境界条件の目的は,電位を一つの境界から別の境界にマップし,領域の周期性をモデル化することである.
電位 を周期的境界 からターゲット境界にマップする関数 があるとき,周期的境界条件は次のように書ける:
例として,完全真空の領域で円形電極の配列をモデル化する.周期的境界条件を使うことで1対の電極だけをモデル化してシミュレーションを行うことが可能である.配列はの等間隔の電極である.上の電極は電位 を持ち,下の電極は地電位を持つ.
真空内の円形電極の配列.赤で描かれた上の電極の境界は電位 を持ち,青で描かれた下の境界は電位0を持つ.
縮小された領域では,1対の電極だけがモデル化される.箱の左境界では周期的境界条件が適用され,境界電位を箱の右側のターゲット境界にマップする.
電極は領域から取り除かれ,電位境界条件で置き換えることでモデル化される.
周期的条件はディリクレ条件が適用された境界には適用できない.
静電容量を計算するのと同じ手法をここでも使う.ただし最終的な結果にセクションの総数を掛ける必要がある.エネルギー法を使うので,以下の方程式を使う:
静電容量は,解析を完全な領域で行った場合よりも小さい.これは電極配列の左側と右側を考慮しなかったためである.
用語集
参考文献
1. Bilbao, S., & Hamilton, B. (2017, October). Directional source modeling in wave-based room acoustics simulation. In 2017 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics (WASPAA) (pp. 121-125). IEEE.
2. Cardoso, J. R. (2018). Electromagnetics through the finite element method. Boca Raton: CRC Press.
3. Griffiths, D. J. (2017). Introduction to electrodynamics (4th ed.). Cambridge University Press.
4. Ida, N., & Bastos, J. P. (2013). Electromagnetics and calculation of fields. Springer Science & Business Media.
5. Küpfmüller, K. (1968). Einführung in die theoretische Elektrotechnik. Springer.
6. Loucas, A. (2023). Optimization of lithium niobate electro-optical modulators (Master's thesis). Stevens Institute of Technology, New Jersey, NJ.
7. Sadiku, M. N. O. (2011). Elements of electromagnetics (5th ed.). Oxford University Press.
8. Peskin, C. S. (2002). The immersed boundary method. Acta numerica, 11, 479-517.