Gmsh
今回、軸対称の形状が多いので、押し出し(Extrude)を用いた。 ここでは、コイル部分の扱いについて紹介する。
コイルの形状は電線を個別に扱うのでは無く、トロイダル電流がある領域で流れているとする。 そのコイルはチューブ上になっており、その横断面は矩形である。 そこで、まず矩形平面を定義し、それを回転押し出しをすることで、コイル領域を定義していく。 以下にその部分のGEOファイルを示す。
// // inductor coil //+ by = -ht-gap; p = newp; Point(p) = {rt1, by, 0, lc2}; Point(p+1) = {rt2, by, 0, lc2}; Point(p+2) = {rt2, by+ht, 0, lc2}; Point(p+3) = {rt1, by+ht, 0, lc2}; l = newl; Line(l) = {p, p+1}; Line(l+1) = {p+1,p+2}; Line(l+2) = {p+2,p+3}; Line(l+3) = {p+3,p}; ll = newll; Line Loop(ll) = {l,l+1,l+2,l+3}; s = news; Plane Surface(s) = {ll}; coil0[] = Extrude {{0, 1, 0}, {0, 0, 0}, Pi/2} { Surface{s}; }; coil1[] = Extrude {{0, 1, 0}, {0, 0, 0}, Pi/2} { Surface{coil0[0]}; }; coil2[] = Extrude {{0, 1, 0}, {0, 0, 0}, Pi/2} { Surface{coil1[0]}; }; coil3[] = Extrude {{0, 1, 0}, {0, 0, 0}, Pi/2} { Surface{coil2[0]}; }; vol_coil = {coil0[1],coil1[1],coil2[1],coil3[1]}; coil_crs = s; //+ Surface Loop(4) = {coil0[{2,3,4,5}],coil1[{2,3,4,5}], coil2[{2,3,4,5}],coil3[{2,3,4,5}] };
いつものように、Pointから定義し、Line, Surface, Volumeと定義していく。 中頃の、Plane Surface s が横断面の矩形である。 これを用いてY軸の周りに1回転させれば形状が作成できる。 ところが、回転押し出しは、$[0,\pi/2]$までしか許されない。 そのため、$\pi/2$ずつ4回押し出している。 ここでExtrude関数は返り値としてリジョンの配列を返してくれている。 その配列の意味は、0が押し出された元の平面の変換面がはいり、つぎには押し出しされた体積、そのあとには押し出された時に生じる側面が必要に応じた数だけ含まれている。 ここでは、押し出された面を再帰的に用い、コイル全体のVolume(vol_coil)を定義し、またAir領域からこの領域を除外するために必要なSurface Loop(4)を定義するするために用いている。
Gmshでは、簡単に体の定義がGUIを用いて出来るが、生の数字が用いられるため、少し変更したりすると番号が全部狂ってしまうので、なるべくこのような方法をとると思いもしないミスを防ぐのに一役買ってくれる。
GetDP
導かれた弱形式方程式(Weak Formulation)を表記するのは今まで通りである。 電磁方程式、熱方程式を別々のFormulationで表記する。 この時、変数名は共通にしなければならない。 以下、GetDPスクリプトを示す。
//--------------------------------------------------------------------------------------------- Formulation { { Name MagStaDyn_av_js0_3D ; Type FemEquation ; Quantity { { Name a ; Type Local ; NameOfSpace HSpace ; } { Name v ; Type Local ; NameOfSpace USpace ; } } Equation { Galerkin { [ nu[] * Dof{d a} , {d a} ] ; In Domain ; Jacobian Vol ; Integration II ; } Galerkin { DtDof[ sigma[] * Dof{a} , {a} ] ; In DomainC ; Jacobian Vol ; Integration II ; } Galerkin { [ sigma[] * Dof{d v} , {a} ] ; In DomainC ; Jacobian Vol ; Integration II ; } Galerkin { [ -js0[], {a} ] ; In DomainS ; Jacobian Vol ; Integration II ; } Galerkin{ DtDof[ sigma[] * Dof{a}, {d v} ] ; In DomainC ; Jacobian Vol ; Integration II ; } Galerkin{ [ sigma[] * Dof{d v} , {d v} ] ; In DomainC ; Jacobian Vol ; Integration II ; } } } { Name Thermal; Type FemEquation; Quantity { { Name t; Type Local; NameOfSpace TSpace; } { Name a; Type Local; NameOfSpace HSpace; } { Name v; Type Local; NameOfSpace USpace; } } Equation { Galerkin { [ K[] * Dof{d t}, {d t} ]; In DomainC; Integration II; Jacobian Vol; } Galerkin { DtDof [ rho[]*Cp[] * Dof{t}, {t} ]; In DomainC; Integration II; Jacobian Vol; } Galerkin { [ -0.5*sigma[]*<a>[SquNorm[Dt[{a}]+{d v}] ], {t} ]; In DomainC; Integration II; Jacobian Vol; } Galerkin { [ Ht[]*Dof{t}, {t} ]; In Skin_ECore; Jacobian Sur ; Integration II ; } Galerkin { [-Ht[]*Tamb[], {t} ]; In Skin_ECore; Jacobian Sur ; Integration II ; } Galerkin { [ sigma_sb*Ep[]*({t})^4, {t} ]; In Skin_ECore; Jacobian Sur ; Integration II ; } Galerkin { [ -sigma_sb*Ep[]*(Tamb[])^4, {t} ]; In Skin_ECore; Jacobian Sur ; Integration II ; } } } } Resolution { { Name Analysis ; System { { Name Thermal ; NameOfFormulation Thermal; } { Name Sys ; NameOfFormulation MagStaDyn_av_js0_3D ; Type ComplexValue ; Frequency Freq ; } } Operation { CreateDir[Dir] ; InitSolution[Sys] ; Generate[Sys] ; Solve[Sys] ; SaveSolution[Sys] ; InitSolution[Thermal]; TimeLoopTheta[ time0, time1, dt[], theta ] { IterativeLoop[ NL_NbrMax, NL_Eps, NL_Relax ] { GenerateJac[Thermal]; SolveJac[Thermal]; } SaveSolution[Thermal]; } PostOperation[MagDyn] ; PostOperation[Thermal] ; } } } //-----------------------------------------------------------------------------------------------
Resolution内では、二つのFormulationがそれぞれ、Thermal, Sysとして定義されている。この時、各々のFormulation内の変数名が共有されるらしい。
FormulationのSysで、渦電流方程式の解を求め、そこで複素解:a, vが得られる。 Therma Formulation内でそれらの値を,
Galerkin { [ -0.5*sigma[]*<a>[SquNorm[Dt[{a}]+{d v}] ], {t} ];
として扱っている。ここで、0.5は振動するエネルギーの平均値のための係数である。 また、Thermal Formulationでは、基本的に「実数」なので、複素数であるaとvを扱うときに、複素数である事を明示し、ある変数(複数も可)を複素数として計算する領域を指示する関数、<
変数名>
[式]を用いている。 今回のように複素変数が複数あっても、一つを指示し、領域を[]で指示すれば内部の変数は全て複素数と認識してくれる。 この関数を使わないと、基本的に実部しか用いない計算結果が得られることになる。 また、時間微分演算子Dtは自動的に複素変数に対応してくれる。