/*
Metaline.c / OpenFDTD / ch.19
*/

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

int main(void)
{
	const double fstart = 2.0e9;
	const double fstop  = 5.0e9;
	const int fdiv = 120;
	const double w = 2.0e-3;
	const double B = 1.6e-3;
	const double er = 2.6;
	const double p = 10e-3;
	const double p0 = 4e-3;
	const double dg = 1e-3;
	//const double rvia = 0.5e-3;
	//const double hvia = 0.6e-3;
	const double Cz = 0.53e-12;
	const double Ly = 3.83e-9;
	const double RB = 80;
	const int M = 20;
	const double GPx = 240e-3;
	const double GPy = 50e-3;
	const double d = 1e-3;   // 基本セルサイズ
	double x;
	char title[BUFSIZ], fout[BUFSIZ];

	sprintf(title, "Metaline Antenna M=%d", M);
	sprintf(fout, "Metaline_M%d.ofd", M);

	// initialize

	ofd_init();

	// title

	ofd_title(title);

	// mesh

	const double x0 = -GPx / 2;
	const double x1 = +GPx / 2;
	ofd_xsection(2, x0, x1);
	ofd_xdivision(1, NINT(x1 - x0, d));

	const double y0 = -GPy / 2;
	const double y1 = -w / 2;
	const double y2 = -y1;
	const double y3 = -y0;
	ofd_ysection(4, y0, y1, y2, y3);
	ofd_ydivision(3, NINT(y1 - y0, d), NINT(y2 - y1, d), NINT(y3 - y2, d));

	const double z0 = -2 * d;
	const double z1 = 0;
	const double z2 = z1 + B;
	const double z3 = 20e-3;
	ofd_zsection(4, z0, z1, z2, z3);
	ofd_zdivision(3, NINT(z1 - z0, d), NINT(z2 - z1, d), NINT(z3 - z2, d));

	// material

	ofd_material(er, 0.0, 1.0, 0.0, "基板");

	// geometry

	// substrate
	ofd_geometry(2, 1, x0, x1, y0, y3, z1, z2);

	// ground
	ofd_geometry(1, 1, x0, x1, y0, y3, z0, z1);

	// feed
	x = -(M / 2.0) * p;
	ofd_geometry(1, 1, x, x, 0, 0, z1, z2);
	ofd_feed('Z', x, 0, z1 + 0.1 * d, 1, 0, RB);
	ofd_rfeed(10);

	// metaline
	for (int i = 0; i < M; i++) {
		// patch
		x = (-M / 2.0 + i) * p;
		ofd_geometry(1, 1, x, x + p0 / 2, -w / 2, +w / 2, z2, z2);
		// C
		x += p0 / 2;
		ofd_geometry(1, 1, x, x + dg, 0, 0, z2, z2);
		ofd_load('X', x + 0.1 * d, 0, z2, 'C', 2 * Cz);
		// patch
		x += dg;
		ofd_geometry(1, 1, x, x + p0,     -w / 2, +w / 2, z2, z2);
		// via + L
		ofd_geometry(1, 1, x + p0 / 2, x + p0 / 2, 0, 0, z1, z2);
		ofd_load('Z', x + p0 / 2, 0, z1 + 0.1 * d, 'L', Ly);
		// C
		x += p0;
		ofd_geometry(1, 1, x, x + dg, 0, 0, z2, z2);
		ofd_load('X', x + 0.1 * d, 0, z2, 'C', 2 * Cz);
		// patch
		x += dg;
		ofd_geometry(1, 1, x, x + p0 / 2, -w / 2, +w / 2, z2, z2);
	}

	// load
	x = +(M / 2.0) * p;
	ofd_geometry(1, 1, x, x, 0, 0, z1, z2);
	ofd_load('Z', x, 0, 0.1 * d, 'R', RB);

	// ABC

	ofd_pml(8, 3, 1e-8);

	// frequency

	ofd_frequency1(fstart, fstop, fdiv);

	ofd_frequency2(2.4e9, 4.2e9, 6);

	// solver

	ofd_solver(200000, 200, 1e-3);

	// iteration

	ofd_plotiter(1);

	// frequency

	//ofd_plotsmith();
	ofd_plotzin(1, 0, 0, 0);
	ofd_plotyin(1, 0, 0, 0);
	ofd_plotref(1, 0, 0, 0);

	// far-1d

	ofd_plotfar1d('Y', 180, 0);
	ofd_far1dcomponent(1, 0, 0);
	ofd_far1dscale(-20, +10, 6);

	// window

	ofd_window3d(400, 400, 12, 60, 30);

	// output

	ofd_outdata(fout);

	return 0;
}