Maxewllの方程式の定型化が終わったので、電界、磁界を計算するプログラムをコーディングした。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 |
private int nx, ny, nz; //セル数 private float[,,] Ex, Ey, Ez, Hx, Hy, Hz; //電界,磁界 private float[] Cex, Cexy, Cexz, Cey, Ceyx, Ceyz, Cez, Cezx, Cezy, Chxy, Chxz, Chyx, Chyz, Chzx, Chzy; //定数 private int[,,] idex, idey, idez, idhx, idhy, idhz; //媒質インデックス //電界の計算 private void calc_erectric_field() { for(int k = 2; k < nz; k++) { for(int j = 2; j < ny; j++) { for(int i = 1; i < nx; i++) { int id = idex[i, j, k]; Ex[i, j, k] = Cex[id] * Ex[i, j, k] + Cexy[id] * (Hz[i, j, k] - Hz[i, j - 1, k]) - Cexz[id] * (Hy[i, j, k] - Hy[i, j, k - 1]); } } } for (int k = 2; k < nz; k++) { for (int j = 1; j < ny; j++) { for (int i = 2; i < nx; i++) { int id = idey[i, j, k]; Ey[i, j, k] = Cey[id] * Ey[i, j, k] + Ceyz[id] * (Hx[i, j, k] - Hx[i, j, k - 1]) - Ceyx[id] * (Hz[i, j, k] - Hz[i - 1, j, k]); } } } for (int k = 1; k < nz; k++) { for (int j = 2; j < ny; j++) { for (int i = 2; i < nx; i++) { int id = idez[i, j, k]; Ez[i, j, k] = Cez[id] * Ex[i, j, k] + Cezx[id] * (Hy[i, j, k] - Hy[i -1, j, k]) - Cezy[id] * (Hy[i, j, k] - Hy[i, j - 1, k]); } } } } //磁界の計算 private void calc_magnetic_field() { for (int k = 1; k < nz; k++) { for (int j = 1; j < ny; j++) { for (int i = 2; i < nx; i++) { int id = idhx[i, j, k]; Hx[i, j, k] = Hx[i, j, k] - Chxy[id] * (Ez[i, j + 1, k] - Ez[i, j, k]) - Chxz[id] * (Ey[i, j, k + 1] - Ey[i, j, k]); } } } for (int k = 1; k < nz; k++) { for (int j = 2; j < ny; j++) { for (int i = 1; i < nx; i++) { int id = idhy[i, j, k]; Hy[i, j, k] = Hy[i, j, k] + Chyz[id] * (Ex[i, j, k + 1] - Ex[i, j, k]) - Chyx[id] * (Ez[i + 1, j, k] - Ez[i, j, k]); } } } for (int k = 2; k < nz; k++) { for (int j = 1; j < ny; j++) { for (int i = 1; i < nx; i++) { int id = idhz[i, j, k]; Hz[i, j, k] = Hz[i, j, k] + Chzx[id] * (Ey[i + 1, j, k] - Ey[i, j, k]) - Chzy[id] * (Ex[i, j + 1, k] - Ex[i, j, k]); } } } } |
Cex等の係数は、全てのセルで異なった値になることはほぼ無いと考えられるので、媒質の種類を格納する配列を設け、その値を使う。
また、セル数の初期値と終値が解析領域の内側になっているのは、解析領域の境界は、吸収境界条件が適用され、別の処理を行う必要があるため。