fdtd.c

来自「用C++或C平台实现的一维时域有限差分算法的实现」· C语言 代码 · 共 72 行

/*FD1D_1.1.c  1D FDTD simulation in free space*/
#include <math.h>
#include <stdlib.h>
#include <stdio.h>

#define KE 200                /*KE is the number of cells to be used */

main ()
{
  float ex[KE],hy[KE];
  int n,k,kc,ke,NSTEPS;
  float T;
  float t0,spread,pulse;
  FILE *fp , *fopen();

  /*initialize*/
  for (k=1; k<KE;k++)
  {
    ex[k]=0;
    hy[k]=0;
    }
  kc=KE/2;	/*center of the problem space*/
  t0=40.0;	/*center of the incident pulse*/
  spread=12;	/*width of the incident pulse*/
  T=0;
  NSTEPS=1;

  while (NSTEPS>0) {
  printf("NETEPS --> ");	/*NETEPS is the center of times the main loop has excuted*/
  scanf ("%d",&NSTEPS);
  printf("%d\n",NSTEPS);
  n=0;

  for (n=1;n<=NSTEPS;n++)
  {
    
    T=T+1;    /*T keeps track of the totle number of times the main loop is excuted*/
 /*main FDTD loop*/
 /*calculate the EX filed*/
 for (k=1;k<KE;k++)
 {
   ex[k]=ex[k]+0.5*(hy[k-1]-hy[k]);}
   /*put a gaussian pulse in the middle*/
   pulse=exp(-0.5*(pow((t0-T)/spread,2.0)));
   ex[kc]=pulse;
   printf("%6.2f %6.3f\n",t0-T,ex[kc]);
   /*calculate the Hy field */
  for (k=0;k<KE-1;k++)
 {
   hy[k]=hy[k]+0.5*(ex[k]-ex[k+1]);
 }
 /*end of the main FDTD loop*/
 /*at the end of the calculation,print out the ex and hy field*/
 for (k=1;k<KE;k++)
 {
  printf("%3d %6.3f %6.3f\n",k,ex[k],hy[k]); }
  
  /*write the E filed out to a file "EX"*/
  fp=fopen("Ex","w");
  for (k=1;k<=KE;k++)
  {
    fprintf (fp,"%6.3f\n",ex[k]);}
    fclose(fp);
    fp=fopen("Hy","w");
  for (k=1;k<=KE;k++)
  {
    fprintf (fp,"%6.3f\n",hy[k]);}
    fclose(fp);
    printf("T=%5.0f\n",T);
  }
}
}