phc13.cpp

基于时域有限差分的遗传算法主程序

#include <stdio.h>
#include <math.h>
#include <string.h>
#include <vector>
#include <iostream>
#include <fstream>
using namespace std;

/* 全局变量 */
double p=0.3975e-6;
vector<vector<double> > E;
vector<vector<double> > ER;
double pgrid = 50;
double dx=p/pgrid;
double dz=p/pgrid;
double Epsilon=1;
double EpsilonR=2.8*2.8;
double Lstraight=5*p;


/* 函数定义 */
int round(double,double);
double shape(double,double,double);
void shape2(int,int,double,int);
void shape3(int,int,double,int);

double shape(double x,double y,double width)
{
	double z;
    if (x<10)             z=0.5*width*((9-x)+1)*p+0.5*p;
	if ((x>9)&&(x<20))    z=0.5*p;
    if ((x>19)&&(x<30))   z=0.5*width*((x-20)+1)*p+0.5*p;
	return z;
}

void shape2(int x, int y, double width, int txstraight)
{
	int shapecosin=round(10*p,dx);
	if (x<=shapecosin-1)
	{
		if (y<=round(0.5*(width)*p*((shapecosin-x)/shapecosin),dz))
			E[y][x+txstraight]=EpsilonR;
	}
}

void shape3(int x, int y, double width, int txstraight)
{
	int shapecosin=round(10*p,dx);
	if (x<=shapecosin-1)
	{
		if (y<=round(0.5*(width)*p*((shapecosin-1-x)/shapecosin+1),dz))
			E[y][x+txstraight+2*shapecosin]=EpsilonR;
	}
}

int round(double x,double y)
{
	int z;
    z=int (x/y+0.5);
	return z;
}
int txstraight=round(Lstraight,dx);


void main()
{
	int i,j;
    double r=0.3*p;
    double L=40*p;
    double H=7*p;
	double dist;
    double icenter;
    double jcenter;
	double offset;
	double width=0.2;
	double Lcouple=30*p;
	int txcouple=round(Lcouple,dx);
    int tx=round(L,dx);
	int tz=round(H,dz);
	E = vector<vector<double> > (tz,vector<double>(tx,Epsilon));
    ER = vector<vector<double> > (2*tz,vector<double>(tx,Epsilon));

	/*char file1[] = "E:\\GA\\SGA\\data\\file1.txt";
	ifstream filedata(file1,ios::in);
	filedata >> width;
	filedata.close();*/

	double Hstraight=(5*width)*p;//
	int tzstraight=round(Hstraight,dz);

	for (i=0; i<=txstraight-1; i++)
	{
		for(j=0; j<=tz-1; j++)
		{
			if (j<=tzstraight-1)
				E[j][i]=EpsilonR;
		}
	}
	
    for (i=0; i<=txcouple-1; i++)
	{
		for (j=0; j<=tz-1; j++)
		{
			icenter = floor((i+0.5)*dx/p);
		    if(icenter < 10) 
			{
				jcenter=floor((j+0.5)*dz/p-offset/p);
                offset=shape(icenter,jcenter,width);
			    if((jcenter > -1) && (jcenter < 5))
					dist=sqrt(pow(((i+0.5)*dx-(icenter+0.5)*p),2.0)+pow(((j+0.5)*dz-((jcenter+0.5)*p+offset)),2.0));
                else
                    dist=10000;
				if (dist<=r)
					E[j][i+txstraight]=EpsilonR;
			}
			shape2(i,j,width,txstraight);
			
            if((icenter > 9) && (icenter < 20))
			{
                jcenter=floor((j+0.5)*dz/p-offset/p);
				offset=shape(icenter,jcenter,width);
                if((jcenter > -1) && (jcenter < 5))
                    dist=sqrt(pow(((i+0.5)*dx-(icenter+0.5)*p),2.0)+pow(((j+0.5)*dz-((jcenter+0.5)*p+offset)),2.0)); 
                else
                    dist=10000;
                if (dist<=r)
					E[j][i+txstraight]=EpsilonR;
			}
                    
            if((icenter > 19) && (icenter < 30))
			{
                jcenter=floor((j+0.5)*dz/p-offset/p);
				offset=shape(icenter,jcenter,width);
                if((jcenter > -1) && (jcenter < 5))
                    dist=sqrt(pow(((i+0.5)*dx-(icenter+0.5)*p),2.0)+pow(((j+0.5)*dz-((jcenter+0.5)*p+offset)),2.0));
                else
                    dist=10000;
                if (dist<=r)
                    E[j][i+txstraight]=EpsilonR;
			}
			shape3(i,j,width,txstraight);
		}		
	}
	
	for (i=0; i<=txstraight-1; i++)
	{
		for(j=0; j<=tz-1; j++)
		{
			if (j<=tzstraight-1)
				E[j][i+txstraight+txcouple]=EpsilonR;
		}
	}
	
	for (i=0; i<=tx-1; i++)
	{
		for (j=0; j<=tz-1; j++)
		{
			ER[j+tz][i]=E[j][i];
	        ER[tz-j-1][i]=E[j][i];
		}
	}
	
	char file[100];
    sprintf(file,"E:\\PhC\\PhC data\\file13.txt");
	ofstream eps(file,ios::out);
	for(j=0; j<=2*tz-1; j++)
	{
		for(i=0; i<=tx-1; i++)
			eps<< ER[j][i] << '\t';
		eps<< endl;
	}
	eps.close();
}