/* 
 * WMM-Application:
 * single raised strip: dispersion curves,
 * modes and propagation constants for varying width
 */

/*
 * WMM
 * Wave-matching method for mode analysis of dielectric waveguides
 * Manfred Hammer 
 * (2002)
 */

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

/* waveguide parameters */
#define Wgpnb 1.45   // substrate refractive index
#define Wgpnr 1.596  // core refractive index 
#define Wgpna 1.0    // cladding refractive index 
#define Wgpt  1.0    // core thickness
#define Wgpl  1.55   // vacuum wavelength 

Polarization Mpol=QTE; // QTE or QTM, semivectorial simulation    

#define Neffmin 1.45  // refractive index interval that is investigated for
#define Neffmax 1.53  // prospective effective mode indices;
                      // suitable choice: Neffmin = Wgpnb
                      // Neffmax = effective index of the fundamental 
                      // properly polarized mode of a slab of thickness Wgpt 

/* range of core widths [Minw, Maxw], sampled at Numw equidistant points */
#define Wmin 0.5 
#define Wmax 3.5 
#define Numw 30  

/* waveguide definition */
Waveguide wgdef(double w)
{
	Waveguide g(1, 1);

	g.hx(0) = 0.0;
	g.hx(1) = Wgpt;

	g.hy(0) = -w/2.0;
	g.hy(1) =  w/2.0;

	g.n(0,0) = Wgpnb;
	g.n(0,1) = Wgpnb;
	g.n(0,2) = Wgpnb;

	g.n(1,0) = Wgpna;
	g.n(1,1) = Wgpnr;
	g.n(1,2) = Wgpna;

	g.n(2,0) = Wgpna;
	g.n(2,1) = Wgpna;
	g.n(2,2) = Wgpna;

	g.lambda = Wgpl;

	return g;
}

/* analysis parameters */
WMM_Parameters pardef()
{
	WMM_Parameters p;

	p.vform = HXHY;
	p.vnorm = NRMMH;
	p.ccomp = CCALL;

	p.ini_d_alpha   = 0.01;
	p.ini_N_alpha   = 20;
	p.ini_alpha_max = 2.0;

	p.ini_steps = 50;
	p.ref_num   = 5;
	p.ref_exp   = 4.0;
	p.ref_sdf   = 0.5;

	p.fin_d_alpha   = 0.01;
	p.fin_N_alpha   = 30;
	p.fin_alpha_max = 2.5;

	p.btol   = 1.0e-7;
	p.mshift = 1.0e-8;

	return p;
}

/* analysis: mode properties versus core width */
int main()
{
	Fcomp fc;	
	WMM_Parameters par = pardef();
	WMM_ModeArray ma;
	WMM_ModeArray modes;
	Waveguide wg;
	int nfms, nfma;
	// files for propagation constants
	char pcnam[20] = "__pc__";
	// ... and effective indices
	char efnam[20] = "__neff__";
	double width;

	fc = principalcomp(Mpol);
	pcnam[0] = poltochar1(Mpol);
	pcnam[1] = poltochar2(Mpol);
	efnam[0] = poltochar1(Mpol);
	efnam[1] = poltochar2(Mpol);

	for(int w=0; w<=Numw; w+=1)
	{
// define the waveguide 
		width = Wmin+((double) w)*(Wmax-Wmin)/((double) Numw);
		wg = wgdef(width);

		modes.clear();

// compute symmetrical modes
		nfms = WMM_modeanalysis(wg, Mpol, SYM, Neffmin, Neffmax, 
					par, dig10(w), dig1(w), ma);
// append propagation constants to appropriate file
		for(int m=0; m<=nfms-1; ++m)
		{
			pcnam[4] = 's';
			pcnam[5] = '0'+m;
			apptoxyf(pcnam, width, ma(m).beta);
			efnam[6] = 's';
			efnam[7] = '0'+m;
			apptoxyf(efnam, width, ma(m).neff);
			fprintf(stderr, "\n -> sm%d(%g)  beta: %g  neff: %g\n", 
			                m, width, ma(m).beta, ma(m).neff);
		}
		ma.clear();

// compute antisymmetrical modes
		nfma = WMM_modeanalysis(wg, Mpol, ASY, Neffmin, Neffmax, 
					par, dig10(w), dig1(w), ma);
// append propagation constants to appropriate file
		for(int m=0; m<=nfma-1; ++m)
		{
			pcnam[4] = 'a';
			pcnam[5] = '0'+m;
			apptoxyf(pcnam, width, ma(m).beta);
			efnam[6] = 'a';
			efnam[7] = '0'+m;
			apptoxyf(efnam, width, ma(m).neff);
			fprintf(stderr, "\n -> am%d(%g)  beta: %g  neff: %g\n", 
			                m, width, ma(m).beta, ma(m).neff);
		}
		ma.clear();
	}
	return 0;
}