/*
sinkansen.c
VŠ²ό
OpenFDTDƒf[ƒ^μ¬ƒ‰ƒCƒuƒ‰ƒŠ—p
*/

#include <stdlib.h>
#include <stdio.h>
#include "ofd_datalib.h"

// ƒpƒ‰ƒ[ƒ^[
const double X0 = 0;          // Œγ•Η(-X)
const double X1 = 0.825;      // ‘ζ1Θ
const double X2 = 21.410;     // ‘O•Η(+X)
const int Nx = 20;            // Θ‚̐”
const double Y0 = -1.560;     // ‰E•Η(-Y)
const double Y6 = +1.560;     // Ά•Η(+Y)
const double Dy = 0.510;      // Θ‚Μ•
const int Ny = 5;             // 1—ρ‚̐Ȃ̐”(=2+3)
const double Z0 = 0;          // °(-Z)
const double Z1 = 1.710;      // ‰Χ•¨’uζ’[Z
const double Z2 = 2.000;      // “Vˆδ•Η‘€
const double Z3 = 2.170;      // “Vˆδ’†S
const double Wx = 0.600;      // ‘‹‚Μ‰‘•
const double Wz1 = 0.790;     // ‘‹‚Μ‰Ί•Σ‚‚³
const double Wz2 = 1.280;     // ‘‹‚̏γ•Σ‚‚³
const double Yb = 1.020;      // ‰Χ•¨’uζ’[Y
const double Ds = 0.600;      // ΐΘ‚̐is•ϋŒό’·‚³
const double D1 = 0.090;      // ΐΘ”w’†‚ΜŒϊ‚³
const double H1 = 0.200;      // ΐΘ‚Μ‰Ί–Κ
const double H2 = 0.400;      // ΐΘ‚̏γ–Κ
const double H3 = 1.200;      // ”w’†‚̏γ–Κ
const double Dg = 0.020;      // ‘‹ƒKƒ‰ƒX‚ΜŒϊ‚³
const double Xa = 0.060;      // ƒAƒ“ƒeƒi’†SXΐ•W
const double Ya = 0;          // ƒAƒ“ƒeƒi’†SYΐ•W
const double Za = 2.000;      // ƒAƒ“ƒeƒi’†SZΐ•W
const double Ze = 0.700;      // ŠΟ‘ͺ–ΚZ‚‚³

static int nint(double l, double d) {return (int)((l / d) + 0.5);}

// Šτ‰½ƒf[ƒ^
static void geometry(void)
{
	int m;
	double g[10];
	double yc[5];
	const double dx = ((X2 - X0) - 2 * (X1 - X0)) / (Nx - 1);
	const double y1 = Y0 + Dy;
	const double y2 = y1 + Dy;
	const double y3 = y2 + Dy;
	const double y5 = Y6 - Dy;
	const double y4 = y5 - Dy;

	// ‘OŒγ•Η
	ofd_geometry(1, 1, X0, X0, Y0, Y6, Z0, Z3);
	ofd_geometry(1, 1, X2, X2, Y0, Y6, Z0, Z3);

	// Ά‰E•Η
	ofd_geometry(1, 1, X0, X2, Y0, Y0, Z0, Z3);
	ofd_geometry(1, 1, X0, X2, Y6, Y6, Z0, Z3);

	// °
	ofd_geometry(1, 1, X0, X2, Y0, Y6, Z0, Z0);

	// “Vˆδ1
	ofd_geometry(1, 1, X0, X2, Y0, Y6, Z3, Z3);

	// “Vˆδ2(ŽOŠp’Œ)
	double hy = 0.500;
	for (int i = 0; i < 2; i++) {
		g[0] = X0;
		g[1] = X2;
		g[2] = (i == 0) ? Y0 : Y6;
		g[3] = (i == 0) ? -hy : +hy;
		g[4] = g[3];
		g[5] = g[2];
		g[6] = Z3;
		g[7] = g[6];
		g[8] = g[7];
		g[9] = Z2;
		ofd_geometry_array(1, 21, g);
	}

	// ΐΘA‘‹ƒKƒ‰ƒX
	yc[0] = (Y0 + y1) / 2;
	yc[1] = (y1 + y2) / 2;
	yc[2] = (y2 + y3) / 2;
	yc[3] = (y4 + y5) / 2;
	yc[4] = (y5 + Y6) / 2;
	for (int i = 0; i < Nx; i++) {
		const double xc = X1 + (i * dx);

		// ΐΘ
		m = 2;
		for (int j = 0; j < Ny; j++) {
			// ΐΘ
			ofd_geometry(m, 1, xc - 0.5 * Ds,      xc + 0.5 * Ds, yc[j] - Dy / 2, yc[j] + Dy / 2, H1, H2);
			// ”w’†
			ofd_geometry(m, 1, xc + 0.5 * Ds - D1, xc + 0.5 * Ds, yc[j] - Dy / 2, yc[j] + Dy / 2, H2, H3);
		}

		// ‘‹ƒKƒ‰ƒX
		m = 3;
		ofd_geometry(m, 1, xc - 0.5 * Wx, xc + 0.5 * Wx, Y0 - Dg / 2, Y0 + Dg / 2, Wz1, Wz2);
		ofd_geometry(m, 1, xc - 0.5 * Wx, xc + 0.5 * Wx, Y6 - Dg / 2, Y6 + Dg / 2, Wz1, Wz2);
	}

	// ‰Χ•¨’I
	ofd_geometry(1, 1, X0, X2, Y0, -Yb, Z1, Z1);
	ofd_geometry(1, 1, X0, X2, Y6, +Yb, Z1, Z1);
}

int main(void)
{

	const double freq = 2.5e9;
	const double d = 10e-3;
	const double lant = 60e-3;
	const char fout[] = "sinkansen_2.5GHz.ofd";

/*
	const double freq = 5e9;
	const double d = 5e-3;
	const double lant = 30e-3;
	const char fout[] = "sinkansen_5GHz.ofd";
*/
	// initialize

	ofd_init();

	// title

	ofd_title("VŠ²ό");

	// mesh

	ofd_xsection(2, X0, X2);
	ofd_xdivision(1, nint(X2 - X0, d));

	const int sy = 5;  // ‘‹‚ΜŠO‚Ι5‘w
	ofd_ysection(2, Y0 - sy * d, Y6 + sy * d);
	ofd_ydivision(1, nint(Y6 - Y0, d) + 2 * sy);

	ofd_zsection(2, Z0, Z3);
	ofd_zdivision(1, nint(Z3 - Z0, d));

	// material

	ofd_material(3.0, 0.2, 1, 0, "ΐΘ");       // =2
	ofd_material(7.0, 0.5, 1, 0, "‘‹ƒKƒ‰ƒX");   // =3

	// geometry

	geometry();

	// antenna

	ofd_geometry(1, 1, Xa, Xa, Ya, Ya, Za - lant / 2, Za + lant / 2);
	ofd_geometry_name("ƒAƒ“ƒeƒi");

	// feed

	ofd_feed('Z', Xa, Ya, Za, 100, 0, 50);

	// frequency

	ofd_frequency1(freq, freq, 0);
	ofd_frequency2(freq, freq, 0);

	// solver

	ofd_solver(300000, 200, 1e-2);

	// iteration

	ofd_plotiter(1);

	// near-2d

	ofd_plotnear2d("E", 'Z', Ze);
	ofd_near2ddim(1, 0);
	ofd_near2ddb(1);
	//ofd_near2dscale(-30, 30);
	ofd_near2dcontour(1);
	ofd_near2dobj(1);

	// output

	ofd_outdata(fout);

	return 0;
}