2dimensionFDTDprogram.rar_FDTD.c_MATLAB 波

// 2 dimension FDTD program (TM mode) presented by T. Nakanishi. (Ver. 1.0) // All boundaries are perfect conductors. // For Visualization of the field, I use PGPLOT which is available as port // (for FreeBSD) or binary package (for Debian Linux). // (1) In FreeBSD, compile with // gcc -o fdtd fdtd.c -I/usr/local/include -L/usr/local/lib -L/usr/X11R6/lib // -lcpgplot -lpgplot -lg2c -lX11 -lm // (2) In debian Linux, complile with // gcc -o fdtd fdtd.c -lcpgplot -lpgplot -L/usr/X11R6/lib -lX11 -lpng // -lz -lm #include "cpgplot.h" #include <math.h> /* define dimension of cells */ int xdim=100,ydim=100; int main() { int i,j,k; float t,e=2,m=2,mat[100][100]; float ez[100][100],hx[100][100],hy[100][100]; cpgopen("/xwindow"); /* Draw on Xwindow */ cpgpap(5.0,1.0); /* Size (5') and ratio (1:1) of Window */ cpgenv(1,100,1,100,1,0); /* Plot erea (100x100) */ for(t=0;t<400;t++){ /* calculate Ez */ for(i=1;i<xdim-1;i++){ for(j=1;j<ydim-1;j++){ ez[i][j] = ez[i][j] + 1/e*((hy[i][j]-hy[i-1][j])-(hx[i][j]-hx[i][j-1])); } } /* Source of Ez */ for(i=0;i<xdim-1;i++){ ez[i][0] = sin(t/4)*sin(i*3.141592/(xdim-1)); } /* calculate Hx and Hy */ for(i=1;i<xdim-1;i++){ for(j=0;j<ydim-1;j++){ hx[i][j] = hx[i][j] - 1/m*(ez[i][j+1]-ez[i][j]); } } for(i=0;i<xdim-1;i++){ for(j=1;j<ydim-1;j++){ hy[i][j] = hy[i][j] + 1/m*(ez[i+1][j]-ez[i][j]); } } /* call plot part */ plot(ez); } cpgend(); /* end of pgplot */ } /* call PGPLOT for illustrating the field */ plot(float mat[]) { int i,j,i1=1,i2=100,j1=1,j2=100; float a1=1,a2=-1,tr[6]={0,0,1.0,0,1.0,0}; float RL[9]={0,0.5, 1}; float RR[9]={1.0,0,0.0}; float RG[9]={0.0,0,0.0}; float RB[9]={0.0,0,1.0}; cpgwedg("BI", 2.0, 5.0, -1.0, 1.0, ""); cpgctab(RL, RR, RG, RB, 9, 1.0, 0.5); cpgimag(mat,xdim, ydim, i1, i2, j1, j2, a1, a2, tr); return; }