fdtd2-tm.txt

fdtd的一些一维到三维的程序,通常为pml边界,用txt文本编辑,无需解压密码

// 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;
}