⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 patch1.cpp

📁 一个用C++计算微带天线的完整的3D例子
💻 CPP
📖 第 1 页 / 共 2 页
字号:
#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <time.h>
#include <complex>
#include <iostream>
#include <fstream>
#include <string.h>
typedef std::complex<double> Complex;
using namespace std;

#define This_file "Patch.txt"

#define Z0 Complex(50,0)
#define nsteps 10000
#define f0 7.45e9

#define patch_L 12.45e-3
#define patch_W 16e-3
#define patch_h 0.794e-3
#define sub_eps 2.2

#define microstrip_st 2.09e-3
#define microstrip_W 2.46e-3

#define i_patch_L 62
#define j_patch_W 80
#define k_top 4
#define surround 15

#define n_pml 6
#define R0 0.0001

#define pi 3.14159
#define eps0 (8.85e-12)
#define mu0 (4*pi*1.0e-7)
#define cc (2.9986e+008)

#define IE (n_pml*2+i_patch_L+surround*2+1)	/* the total FDTD grID length (including pml) */
#define JE (n_pml*2+j_patch_W+surround+50+1)
#define KE (k_top+n_pml+3+1)

#define 	ddx (patch_L/i_patch_L)
#define 	ddy (patch_W/j_patch_W)
#define 	ddz (patch_h/k_top)
#define 	ra_x (ddz/ddx)
#define 	ra_y (ddz/ddy)
#define 	dt (ddz/(2*cc))
#define 	k0 (2*pi*f0/cc)
#define 	W0 (2*pi*f0*dt)

#define 	i_center (IE/2)
#define 	j_center (JE/2)
#define 	i_max (IE-n_pml-1)
#define 	j_max (JE-n_pml-1)
#define 	k_max (KE-n_pml-1)

#define 	i_patch_st (n_pml+surround)
#define 	i_patch_end (i_patch_st+i_patch_L)
#define 	j_patch_end (j_max-surround)
#define 	j_patch_st (j_patch_end-j_patch_W)
#define 	j_ref (n_pml+15)
#define 	k_ref (k_top-1)
#define 	i_mic_start (i_patch_st+(int)(microstrip_st/ddx+.5))
#define 	i_mic_end (i_mic_start+(int)(microstrip_W/ddx+.5))
#define 	i_ref ((i_mic_start+i_mic_end)/2)

float Dx[IE][JE][KE],Dy[IE][JE][KE],Dz[IE][JE][KE];
float Ex[IE][JE][KE],Ey[IE][JE][KE],Ez[IE][JE][KE];
float Hx[IE][JE][KE],Hy[IE][JE][KE],Hz[IE][JE][KE];
float gax[IE][JE][KE],gay[IE][JE][KE],gaz[IE][JE][KE];

float Dxy[IE][JE][KE],Dxz[IE][JE][KE];
float Dyz[IE][JE][KE],Dyx[IE][JE][KE];
float Dzx[IE][JE][KE],Dzy[IE][JE][KE];
float Hxy[IE][JE][KE],Hxz[IE][JE][KE];
float Hyz[IE][JE][KE],Hyx[IE][JE][KE];
float Hzx[IE][JE][KE],Hzy[IE][JE][KE];

float CAx[n_pml],CBx[n_pml],CPx[n_pml],CQx[n_pml];
float CAy[n_pml],CBy[n_pml],CPy[n_pml],CQy[n_pml];
float CAz[n_pml],CBz[n_pml],CPz[n_pml],CQz[n_pml];

const Complex _J=Complex(0,1);
Complex Voltage,Current,Zin,S11;
Complex fx[181][2],fy[181][2],fmx[181][2],fmy[181][2];
Complex E_theta[181][2],E_phi[181][2],E_total[181][2];
Complex Jx[IE][JE],Jy[IE][JE],Jmx[IE][JE],Jmy[IE][JE];


int n,i,j,k;
int ix,jy,kz;

char Ez_file[80]={"Ez_f0_"};
char Pattern_file[80]={"Pattern_f0_"};
char time_file[80]={"Time_f0_"};
char S11_file[80]={"S11_f0_"};

FILE *fp;
fstream outfile;
void init_parameter();
void init_conductor();
void show_parameter();
void cal_pattern();

void main()
{
    clock_t time_start, time_finish;
    double duration;
	time_start=clock();
	float aaa;
	init_parameter();
	init_conductor();
	fp=fopen(time_file,"w");
	/* Start of tJ_Jm main FDTD loop */
	for (n=1;n<=nsteps;n++){
		/* FDTD (not PML) E field */
		for (i=0;i<(IE-1);i++){
			for (j=(n_pml+1);j<j_max;j++){
				for (k=1;k<k_max;k++){
					Dx[i][j][k]=Dx[i][j][k]+0.5*(ra_y*(Hz[i][j][k]-Hz[i][j-1][k])-Hy[i][j][k]+Hy[i][j][k-1]);
		}}}
		for (i=(n_pml+1);i<i_max;i++){
			for (j=0;j<(JE-1);j++){
				for (k=1;k<k_max;k++){
					Dy[i][j][k]=Dy[i][j][k]+.5*(Hx[i][j][k]-Hx[i][j][k-1]-ra_x*(Hz[i][j][k]-Hz[i-1][j][k]));
		}}}
		for (i=(n_pml+1);i<i_max;i++){
			for (j=(n_pml+1);j<j_max;j++){
				for (k=0;k<(KE-1);k++){
					Dz[i][j][k]=Dz[i][j][k]+.5*(ra_x*(Hy[i][j][k]-Hy[i-1][j][k])-ra_y*(Hx[i][j][k]-Hx[i][j-1][k]));
		}}}
		/* PML field */
		/* k=k_max:KE */
		for (k=k_max;k<(KE-1);k++){
			kz=(k-k_max);
			for (i=0;i<(IE-1);i++){
				for (j=1;j<(JE-1);j++){
					Dxz[i][j][k]=CAz[kz]*Dxz[i][j][k]-CBz[kz]*(Hy[i][j][k]-Hy[i][j][k-1]);
				}
				for (j=(n_pml+1);j<j_max;j++){
					Dxy[i][j][k]=Dxy[i][j][k]+.5*ra_y*(Hz[i][j][k]-Hz[i][j-1][k]);
					Dx[i][j][k]=Dxy[i][j][k]+Dxz[i][j][k];
			}}
			for (j=0;j<(JE-1);j++){
				for (i=1;i<(IE-1);i++){
					Dyz[i][j][k]=CAz[kz]*Dyz[i][j][k]+CBz[kz]*(Hx[i][j][k]-Hx[i][j][k-1]);
				}
				for (i=(n_pml+1);i<i_max;i++){
					Dyx[i][j][k]=Dyx[i][j][k]-.5*ra_x*(Hz[i][j][k]-Hz[i-1][j][k]);
					Dy[i][j][k]=Dyz[i][j][k]+Dyx[i][j][k];
		}}}
		/* I=1:i_min */
		for (i=1;i<=n_pml;i++){
			ix=n_pml-i;
			for (j=0;j<(JE-1);j++){
				for (k=1;k<k_max;k++){
					Dyz[i][j][k]=Dyz[i][j][k]+.5*(Hx[i][j][k]-Hx[i][j][k-1]);
				}
				for (k=1;k<(KE-1);k++){
					Dyx[i][j][k]=CAx[ix]*Dyx[i][j][k]-CBx[ix]*(Hz[i][j][k]-Hz[i-1][j][k]);
					Dy[i][j][k]=Dyz[i][j][k]+Dyx[i][j][k];
			}}
			for (k=0;k<(KE-1);k++){
				for (j=1;j<(JE-1);j++){
					Dzx[i][j][k]=CAx[ix]*Dzx[i][j][k]+CBx[ix]*(Hy[i][j][k]-Hy[i-1][j][k]);
				}
				for (j=(n_pml+1);j<j_max;j++){
					Dzy[i][j][k]=Dzy[i][j][k]-.5*ra_y*(Hx[i][j][k]-Hx[i][j-1][k]);
					Dz[i][j][k]=Dzx[i][j][k]+Dzy[i][j][k];
		}}}
		/* I=i_max:(IE-1) */
		for (i=i_max;i<(IE-1);i++){
			ix=i-i_max;
			for (j=0;j<(JE-1);j++){
				for (k=1;k<k_max;k++){
					Dyz[i][j][k]=Dyz[i][j][k]+.5*(Hx[i][j][k]-Hx[i][j][k-1]);
				}
				for (k=1;k<(KE-1);k++){
					Dyx[i][j][k]=CAx[ix]*Dyx[i][j][k]-CBx[ix]*(Hz[i][j][k]-Hz[i-1][j][k]);
					Dy[i][j][k]=Dyz[i][j][k]+Dyx[i][j][k];
			}}
			for (k=0;k<(KE-1);k++){
				for (j=1;j<(JE-1);j++){
					Dzx[i][j][k]=CAx[ix]*Dzx[i][j][k]+CBx[ix]*(Hy[i][j][k]-Hy[i-1][j][k]);
				}
				for (j=(n_pml+1);j<j_max;j++){
					Dzy[i][j][k]=Dzy[i][j][k]-.5*ra_y*(Hx[i][j][k]-Hx[i][j-1][k]);
					Dz[i][j][k]=Dzx[i][j][k]+Dzy[i][j][k];
		}}}
		/* j=1:j_min */
		for (j=1;j<=n_pml;j++){
			jy=n_pml-j;
			for (i=0;i<(IE-1);i++){
				for (k=1;k<k_max;k++){
					Dxz[i][j][k]=Dxz[i][j][k]-.5*(Hy[i][j][k]-Hy[i][j][k-1]);
				}
				for (k=1;k<(KE-1);k++){
					Dxy[i][j][k]=CAy[jy]*Dxy[i][j][k]+CBy[jy]*(Hz[i][j][k]-Hz[i][j-1][k]);
					Dx[i][j][k]=Dxy[i][j][k]+Dxz[i][j][k];
			}}
			for (k=0;k<(KE-1);k++){
				for (i=(n_pml+1);i<i_max;i++){
					Dzx[i][j][k]=Dzx[i][j][k]+.5*ra_x*(Hy[i][j][k]-Hy[i-1][j][k]);
				}
				for (i=1;i<(IE-1);i++){
					Dzy[i][j][k]=CAy[jy]*Dzy[i][j][k]-CBy[jy]*(Hx[i][j][k]-Hx[i][j-1][k]);
					Dz[i][j][k]=Dzx[i][j][k]+Dzy[i][j][k];
		}}}
		/* j=j_max:(JE-1) */
		for (j=j_max;j<(JE-1);j++){
			jy=j-j_max;
			for (i=0;i<(IE-1);i++){
				for (k=1;k<k_max;k++){
					Dxz[i][j][k]=Dxz[i][j][k]-.5*(Hy[i][j][k]-Hy[i][j][k-1]);
				}
				for (k=1;k<(KE-1);k++){
					Dxy[i][j][k]=CAy[jy]*Dxy[i][j][k]+CBy[jy]*(Hz[i][j][k]-Hz[i][j-1][k]);
					Dx[i][j][k]=Dxy[i][j][k]+Dxz[i][j][k];
			}}
			for (k=0;k<(KE-1);k++){
				for (i=(n_pml+1);i<i_max;i++){
					Dzx[i][j][k]=Dzx[i][j][k]+.5*ra_x*(Hy[i][j][k]-Hy[i-1][j][k]);
				}
				for (i=1;i<(IE-1);i++){
					Dzy[i][j][k]=CAy[jy]*Dzy[i][j][k]-CBy[jy]*(Hx[i][j][k]-Hx[i][j-1][k]);
					Dz[i][j][k]=Dzx[i][j][k]+Dzy[i][j][k];
		}}}

		/* Source */
		aaa=sin(W0*n);
		if (W0*n<(5*pi)){
			aaa=aaa*.5*(1-cos(W0*n/5.));
		}
		for (i=i_mic_start;i<=i_mic_end;i++){
			for (k=0;k<k_top;k++){
				Dz[i][n_pml+2][k]=Dz[i][n_pml+2][k]+aaa;
		}}
		/* Calculate the E from D field */
		/* Ex and Ey are zero at k=0 */
		for (i=0;i<IE;i++){
			for (j=0;j<JE;j++){
				for (k=0;k<=k_top;k++){
					Ex[i][j][k]=gax[i][j][k]*Dx[i][j][k];
					Ey[i][j][k]=gay[i][j][k]*Dy[i][j][k];
					Ez[i][j][k]=gaz[i][j][k]*Dz[i][j][k];
				}
				for (k=k_top+1;k<KE;k++){
					Ex[i][j][k]=Dx[i][j][k];
					Ey[i][j][k]=Dy[i][j][k];
					Ez[i][j][k]=Dz[i][j][k];
		}}}
		/* Calculate the Hx field */
		for (i=1;i<(IE-1);i++){
			for (j=n_pml;j<j_max;j++){
				for (k=0;k<k_max;k++){
					Hx[i][j][k]=Hx[i][j][k]+.5*(Ey[i][j][k+1]-Ey[i][j][k]-ra_y*(Ez[i][j+1][k]-Ez[i][j][k]));
		}}}
		for (i=n_pml;i<i_max;i++){
			for (j=1;j<(JE-1);j++){
				for (k=0;k<k_max;k++){
					Hy[i][j][k]=Hy[i][j][k]+.5*(ra_x*(Ez[i+1][j][k]-Ez[i][j][k])-Ex[i][j][k+1]+Ex[i][j][k]);
		}}}
		for (i=n_pml;i<i_max;i++){
			for (j=n_pml;j<j_max;j++){
				for (k=1;k<(KE-1);k++){
					Hz[i][j][k]=Hz[i][j][k]+.5*(ra_y*(Ex[i][j+1][k]-Ex[i][j][k])-ra_x*(Ey[i+1][j][k]-Ey[i][j][k]));
		}}}
		/* PML field */
		/* k=k_max:KE */
		for (k=k_max;k<(KE-1);k++){
			kz=k-k_max;
			for (i=1;i<(IE-1);i++){
				for (j=0;j<(JE-1);j++){
					Hxz[i][j][k]=CPz[kz]*Hxz[i][j][k]+CQz[kz]*(Ey[i][j][k+1]-Ey[i][j][k]);
				}
				for (j=n_pml;j<j_max;j++){
					Hxy[i][j][k]=Hxy[i][j][k]-.5*ra_y*(Ez[i][j+1][k]-Ez[i][j][k]);
					Hx[i][j][k]=Hxy[i][j][k]+Hxz[i][j][k];
			}}
			for (j=1;j<(JE-1);j++){
				for (i=0;i<(IE-1);i++){
					Hyz[i][j][k]=CPz[kz]*Hyz[i][j][k]-CQz[kz]*(Ex[i][j][k+1]-Ex[i][j][k]);
				}
				for (i=n_pml;i<i_max;i++){
					Hyx[i][j][k]=Hyx[i][j][k]+.5*ra_x*(Ez[i+1][j][k]-Ez[i][j][k]);
					Hy[i][j][k]=Hyz[i][j][k]+Hyx[i][j][k];
		}}}
		/* i=0:i_min-1 */
		for (i=0;i<n_pml;i++){
			ix=n_pml-i-1;
			for (j=1;j<(JE-1);j++){
				for (k=0;k<k_max;k++){
					Hyz[i][j][k]=Hyz[i][j][k]-.5*(Ex[i][j][k+1]-Ex[i][j][k]);
				}
				for (k=0;k<(KE-1);k++){
					Hyx[i][j][k]=CPx[ix]*Hyx[i][j][k]+CQx[ix]*(Ez[i+1][j][k]-Ez[i][j][k]);
					Hy[i][j][k]=Hyz[i][j][k]+Hyx[i][j][k];
			}}
			for (k=1;k<(KE-1);k++){
				for (j=0;j<(JE-1);j++){
					Hzx[i][j][k]=CPx[ix]*Hzx[i][j][k]-CQx[ix]*(Ey[i+1][j][k]-Ey[i][j][k]);
				}
				for (j=n_pml;j<j_max;j++){
					Hzy[i][j][k]=Hzy[i][j][k]+.5*ra_y*(Ex[i][j+1][k]-Ex[i][j][k]);
					Hz[i][j][k]=Hzx[i][j][k]+Hzy[i][j][k];
		}}}
		/* i=i_max:(IE-1) */
		for (i=i_max;i<(IE-1);i++){
			ix=i-i_max;
			for (j=1;j<(JE-1);j++){
				for (k=0;k<k_max;k++){
					Hyz[i][j][k]=Hyz[i][j][k]-.5*(Ex[i][j][k+1]-Ex[i][j][k]);
				}
				for (k=0;k<(KE-1);k++){
					Hyx[i][j][k]=CPx[ix]*Hyx[i][j][k]+CQx[ix]*(Ez[i+1][j][k]-Ez[i][j][k]);
					Hy[i][j][k]=Hyz[i][j][k]+Hyx[i][j][k];
			}}
			for (k=1;k<(KE-1);k++){
				for (j=0;j<(JE-1);j++){
					Hzx[i][j][k]=CPx[ix]*Hzx[i][j][k]-CQx[ix]*(Ey[i+1][j][k]-Ey[i][j][k]);
				}
				for (j=n_pml;j<j_max;j++){
					Hzy[i][j][k]=Hzy[i][j][k]+.5*ra_y*(Ex[i][j+1][k]-Ex[i][j][k]);
					Hz[i][j][k]=Hzx[i][j][k]+Hzy[i][j][k];
		}}}
		/* j=0:j_min-1 */
		for (j=0;j<n_pml;j++){
			jy=n_pml-j-1;
			for (i=1;i<(IE-1);i++){
				for (k=0;k<k_max;k++){
					Hxz[i][j][k]=Hxz[i][j][k]+.5*(Ey[i][j][k+1]-Ey[i][j][k]);
				}

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -