/* 
 * METRIC-Application:
 * A three layer waveguide, mode spectrum discetization
 */

#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include"metric.h"

#define Pol     TE     // light polarization
#define Wgpns   1.0    // substrate refractive index
#define Wgpnf   2.0    // core refractive index
#define Wgpnc   1.0    // cover refractive index
#define Wgpt    1.2    // slab thickness 
#define Wavel   1.55   // vacuum wavelength 

#define NumM    150    // number of terms in the spectral discretizations
#define CW      5.0    // computational window: 
                       // distance from the waveguide axis (core center)

/* waveguide definition */
Waveguide wgdef()
{
	// a waveguide with one inner layer
	Waveguide g(1);

	// positions of the dielectric interfaces
	g.hx(0) = -Wgpt/2.0;
	g.hx(1) = Wgpt/2.0;

	// refractive index values 
	g.n(0) = Wgpns;
	g.n(1) = Wgpnf;
	g.n(2) = Wgpnc;

	// the vacuum wavelength
	g.lambda = Wavel;

	return g;
}

/* calculate the discretized eigenmodes, write profile data;
   output files: squared pairs (mode order, squared propagation constant) */
int main()
{
	ModeArray ma;
	Fcomp cp = principalcomp(Pol);
	Dvector bq(NumM);
	Dvector od(NumM);
	char name[10] = "t__bq";
		
	// define the waveguide
	Waveguide wg = wgdef();
	// the computational window
	Interval cw(-CW, CW);

	// find modes defined by Dirichlet boundary conditions
	modespectrum(wg, Pol, cw, LIMD, NumM, ma); 
	// ma.orthonormality();
	for(int j=0; j<=ma.num-1; ++j)
	{
		bq(j) = ma(j).betaq;
		od(j) = (double) j; 
	}
	name[1] = polchr(Pol); 
	name[2] = 'D'; 
	toxyf(name, od, bq);

	// plots and data files for a few low order profiles 
	for(int i=0; i<=9; ++i)
	{
		char t = 'L'; if(i>=1) t = 'V';
		ma(i).plot(cp, ORG, cw, 500, 'D', dig1(i), t);
		ma(i).writeprofile(cp, ORG, cw, 500, 'D', dig1(i));
	}

	// find modes, defined by Neumann boundary conditions
	modespectrum(wg, Pol, cw, LIMN, NumM, ma); 
	// ma.orthonormality();
	for(int j=0; j<=ma.num-1; ++j)
	{
		bq(j) = ma(j).betaq;
		od(j) = (double) j; 
	}
	name[1] = polchr(Pol); 
	name[2] = 'N'; 
	toxyf(name, od, bq);

	// plots and data files for a few low order profiles 
	for(int i=0; i<=9; ++i)
	{
		char t = 'L'; if(i>=1) t = 'V';
		ma(i).plot(cp, ORG, cw, 500, 'N', dig1(i), t);
		ma(i).writeprofile(cp, ORG, cw, 500, 'N', dig1(i));
	}

	return 0;
}