fdmenu.c

matlab编写的三维FDTD程序 pml边界条件 实现电磁波传播数值模拟和仿真

/* program fdmenu.c  - last modified 7/26/94 */
/*
  This program creates a binary file called "fdmenu.par" which the program fdtd3d
  reads when it is executed.  The main purpose of this program is to modify flexible
  parameters in ftdt3d.  The ascii "fdmenasc.par" contains default parameters
  that are read into fdmenu when fdmenu is executed.
*/
#include <stdio.h>
#include <math.h>
#include <string.h>
#include "findif.h"
FILE *ffdi1,*fasc;
void fdmenu();
void fdmain();
void fd3dmenu();
void outputit();
void antmenu();
void automenu();
void inasc();
void outasc();
void newparms();
main() {
fdmenu();
antcalcs();
if(fd1.obj==1) sphcalc();
if(fd1.obj==2) cylcalc();
if(fd1.patact==1) profpat();
if((fd1.obj==7)||(fd1.obj==8)) pipcalc3();
}
void fdmenu() {
char ch;
short int runflag=0,dum1,ascflag=1;/*when ascflag=1, ascii field fdmenasc.par is opened and read*/
float compstab;
 if(ascflag==0) {
  if ((ffdi1=fopen("fdmenu.par","r+")) == NULL) {
   printf("GENERATING NEW PARAMETER DATASET\n");
   if((ffdi1=fopen("fdmenu.par","wb"))== NULL) {
        printf("CAN'T GENERATE SPHMOD.PAR FILE - ABORT\n");
	exit();
   }
        fd1.filout[0]='a'; fd1.filout[1]='1'; fd1.centfreq=0.0; fd1.amp=0.0;
	fd1.perm=0.0; fd1.cond=0.0; fd1.dist=0.0; fd1.size=0.0; 
	fd1.dx=fd1.dy=fd1.dz=0.02; fd1.dt=0.03125; fd1.modsiz=500.; fd1.pol='y';
	fd1.timnum=1200; fd1.xmax=2.565; fd1.ymax=2.565; fd1.zmax=2.565; fd1.bounswi=3;
	fd1.tlength=3.25; fd1.tcval=4.5; fd1.moddim=0; fd1.indat=0; fd1.infile[0]='n';
        fd1.infile[1]='o'; fd1.infile[2]='n'; fd1.infile[3]='e'; fd1.obj=0;
	fd1.rad=0.041; fd1.xpos=1.285; fd1.ypos=1.285; fd1.zpos=0.545; fd1.xcen=1.285;
 	fd1.ycen=1.285; fd1.xlen=0.481; fd1.xang=60.0; fd1.xrcen=1.685; fd1.yrcen=1.285;
 	fd1.xrlen=0.481; fd1.xrang=60.0; fd1.anthei=0.059; fd1.zinter=0.865;
	fd1.zcond=0.005; fd1.zperm=5.5; fd1.antpat=0; 
	fd1.fullEx=fd1.fullEy=fd1.fullEz=fd1.fullHx=fd1.fullHy=fd1.fullHz=fd1.HzEz=fd1.ExEy=0;
	fd1.HxHy=0; fd1.xlin=fd1.ylin=0.0; fd1.timpro=8; fd1.timx=1.285; fd1.timy=1.285;
	fd1.dump=0; fd1.dumpfile[0]='n'; fd1.dumpfile[1]='a'; fd1.dumptime=0;
 	fd1.filauto[0]='n'; fd1.filauto[1]='o'; fd1.zinc=0.0; fd1.xinc=0.02; fd1.yinc=0.2;
 	fd1.numsteps=1; fd1.anttype=2; fd1.sheet=0; fd1.objlen=2.0; fd1.shtydim=0.80;
	fd1.shtxdim=0.80; fd1.shtthk=0.005; fd1.shtcond=0.0; fd1.shtperm=2.0;
	fd1.shthei=0.041; fd1.patact=0; fd1.patrad=0.5; fd1.patphi=0.0; fd1.cylphi=90.;
	fd1.cyldip=0.; fd1.perconflg=1; fd1.targperm=1.; fd1.targcond=0.0; fd1.zur=1.0;
	fd1.tmproflg=0; fd1.absfreq=125; fd1.pdflag=0; fd1.attcoax=0; fd1.coaxurt=1.13;
	fd1.coaxst=0.0; fd1.coaxert=1.; fd1.coaxurr=1.13; fd1.coaxsr=0.0; fd1.coaxerr=1.;
	fd1.prdepth=0.50; fd1.pvala=-0.055; fd1.pvalb=0.245; fd1.pvalc=0.485;
	fd1.delaya=0.0; fd1.delayb=0.0; fd1.delayc=0.0;
        
 
 }
 else fread(&fd1,sizeof(struct fddat),1,ffdi1);
 fclose(ffdi1);
}
else if(ascflag==1) inasc();
 for( ;;) {
	printf("1 = OUTPUT DATASET NAME               %s\n",fd1.filout);      
	printf("8 = DX:  (m)                          %f\n",fd1.dx);      
	printf("9 = DY:  (m)                          %f\n",fd1.dy);      
	printf("a = DZ:  (m)                          %f\n",fd1.dz);      
	printf("b = DT:  (ns)                         %f\n",fd1.dt);      
	printf("c = PEAK VOLTAGE OF PULSE:            %f\n",fd1.modsiz);      
	printf("e = CHECK STABILITY \n");      
	printf("f = TIME STEP REPITITIONS:            %d\n",fd1.timnum);
        printf("g = BURIED OBJECT PARAMETERS	\n");      
        printf("h = X-DIMENSION OF GRID (m)           %f\n",fd1.xmax);      
        printf("i = Y-DIMENSION OF GRID (m)           %f\n",fd1.ymax);      
        printf("j = Z-DIMENSION OF GRID (m)           %f\n",fd1.zmax);      
        printf("k = BOUNDARY SWITCH (1,2)             %d\n",fd1.bounswi);      
        printf("m = TIMELENGTH OF INPUT PULSE (ns)    %f\n",fd1.tlength);      
        printf("n = TIME DELAY BEFORE PEAK OF PULSE(ns)  %f\n",fd1.tcval);      
        printf("o = DATA OUTPUT OPTIONS     \n");      
        printf("s = ANTENNA AND GROUND PROP. OPTIONS \n");      
        printf("t = AUTO MODE OPTIONS\n");      
        printf("u = INPUT DATA FROM PREVIOUS RUN:     %d\n",fd1.indat);      
        printf("w = ACCESS TO COAXIAL CABLE PARAMETERS                 \n");      
	printf("r = EXIT\n");
	scanf(" %c",&ch);
        switch(ch) {
           case '1':
	      printf("ENTER NEW OUTPUT DATASET NAME: ");
	      scanf("%s",fd1.filout);
	   break;
	   case '8':
	      printf("ENTER DX (m): ");
	      scanf("%f",&fd1.dx);
           break;
	   case '9':
	      printf("ENTER DY (m): ");
	      scanf("%f",&fd1.dy);
           break;
	   case 'a':
	      printf("ENTER DZ (m): ");
	      scanf("%f",&fd1.dz); 
           break;
	   case 'b':
	      printf("ENTER DT (ns): ");
	      scanf("%f",&fd1.dt);
           break;
	   case 'c':
	      printf("ENTER PEAK VOLTAGE OF INPUT PULSE: ");
	      scanf("%f",&fd1.modsiz);
           break;
	   case 'e':
compstab=pow(1./(fd1.dx*fd1.dx)+1./(fd1.dy*fd1.dy)+1./(fd1.dz*fd1.dz),-0.5)-fd1.dt*0.3;
          printf("COMPUTATIONAL STABILITY CHECK = %f\n",compstab);
           break;
	   case 'f':
              printf("ENTER NUMBER OF TIME STEPS: ");
	      scanf("%d",&fd1.timnum); 
           break;
	   case 'g':
             fd3dmenu(); 
           break;
	   case 'h':
              printf("ENTER X-DIMENSION OF GRID (m): ");
              scanf("%f",&fd1.xmax);
           break;
	   case 'i':
              printf("ENTER Y-DIMENSION OF GRID (m): ");
              scanf("%f",&fd1.ymax);
           break;
	   case 'j':
              printf("ENTER Z-DIMENSION OF GRID (m): ");
              scanf("%f",&fd1.zmax);
           break;
	   case 'k':
              printf("ENTER BOUNDARY SWITCH OPTION (1= F.O., 2= F.O.+): ");
              scanf("%d",&fd1.bounswi);
           break;
	   case 'm':
              printf("ENTER TIME LENGTH OF PULSE (ns): ");
              scanf("%f",&fd1.tlength);
           break;
	   case 'n':
              printf("ENTER TIME DELAY BEFORE PEAK OF PULSE (ns) ");
              scanf("%f",&fd1.tcval);
           break;
	   case 'o':
            outputit(); 
           break;
	   case 's':
            antmenu(); 
           break;
	   case 't':
           automenu();
           break;
	   case 'u':
           printf("ENTER '1' TO INPUT FILE FROM PREVIOUS RUN, ELSE '0': ");
           scanf("%d",&fd1.indat);
           break;
	   case 'w':
              newparms();
           break;
	   case 'r':		   
	      runflag=1;
	   break;
       }  
       if (runflag==1) {
          if((ffdi1=fopen("fdmenu.par","wb"))==NULL){
	   printf("CAN'T OPEN fdmenu.par FOR WRITING - ABORTING\n");
	   exit();
	  }
          fwrite(&fd1,sizeof(struct fddat),1,ffdi1);	
          fclose(ffdi1);
          outasc();
          break;
       }
 }
}
void fd3dmenu() {
char ch;
short int xflag=0,dum1;
float compstab;
 for( ;;) {
	printf("1 = OBJECT (0-5)                      %d\n",fd1.obj);      
	printf("2 = RADIUS (m)                        %f\n",fd1.rad);      
	printf("3 = DEPTH (m)                         %f\n",fd1.zpos);
	printf("4 = X POSITION OF CENTER IN GRID (m)  %f\n",fd1.xpos);      
	printf("5 = Y POSITION OF CENTER IN GRID (m)  %f\n",fd1.ypos);      
	printf("6 = LENGTH (for pipe and cylinder)  (m)  %f\n",fd1.objlen); 
        if((fd1.obj==2)||(fd1.obj>=7)) {
           printf("7 = TREND OF CYLINDER (OR PIPE) AXIS FROM X-AXIS (deg.) %f\n",fd1.cylphi);
           printf("8 = DIP OF CYLINDER FROM HORIZONTAL (deg.) %f\n",fd1.cyldip);
           printf("CYLINDER DIP TRENDS IN DIRECTION CHOSEN IN OPTION 7\n");
        }
        printf("9 = PERFECT CONDUCTOR FLAG (0,1)       %d\n",fd1.perconflg);  
        if(fd1.perconflg==0) {
           printf("a = RELATIVE PERMITTIVITY OF TARGET  %f\n",fd1.targperm);  
           printf("b = CONDUCTIVITY OF TARGET (S/m)     %f\n",fd1.targcond);  
        }
	printf("z = BACK TO MAIN MENU\n");
	scanf(" %c",&ch);
        switch(ch) {
           case '1':
	      printf("ENTER OBJECT DESIRED: \n");
              printf("   0 = NONE\n");
              printf("   1 = SPHERE\n");
              printf("   2 = CYLINDER \n");
              printf("   3 = CUBE \n");
              printf("   4 = HORIZONTAL SQUARE METAL SHEET (DIMEN=LENGTH*LENGTH) \n");
              printf("   7 = HORIZONTAL METAL PEOPE(S) ORIENTED AT ANGLE WITH RESPECT TO TRANSMIT ATENNA\n");
              printf("   8 = HORIZONTAL DIELECTRIC PEOPE(S) ORIENTED AT ANGLE WITH RESPECT TO TRANSMIT ATENNA\n");
              printf("  10 = HALFSPACE POSITIONED AT TARGET DEPTH AND POSSESSING TARGET ELEC. PROP.\n");
	      scanf("%d",&fd1.obj);
	   break;
	   case '2':
	      printf("ENTER RADIUS OF TARGET (m)");
	      scanf("%f",&fd1.rad);
           break;
	   case '3':
	      printf("ENTER DEPTH TO CENTER OF OBJECT:   ");
	      scanf("%f",&fd1.zpos);
	   break;
	   case '4':
	      printf("ENTER X POSITION OF CENTER OF OBJECT: ");
	      scanf("%f",&fd1.xpos);
           break;
	   case '5':
	      printf("ENTER Y POSITION OF CENTER OF OBJECT: ");
	      scanf("%f",&fd1.ypos);
           break;
	   case '6':
	      printf("ENTER LENGTH OF TARGET (IF PIPE OR CYLINDER): ");
	      scanf("%f",&fd1.objlen);
           break;
	   case '7':
	      printf("ENTER 'PHI' OF CYLINDER AXIS (DEGREES FROM X-AXIS): ");
	      scanf("%f",&fd1.cylphi);
           break;
	   case '8':
	      printf("ENTER DIP OF CYLINDER AXIS FROM HORIZONTAL (DEG.) ");
	      scanf("%f",&fd1.cyldip);
           break;
	   case '9':
	      printf("ENTER PERFECT CONDUCTOR FLAG ('0' =TARG. NOT PER. COND.,ELSE ENTER '1' ");
	      scanf("%d",&fd1.perconflg);
           break;
	   case 'a':
	      printf("ENTER RELATIVE PERMITTIVITY OF TARGET: ");
	      scanf("%f",&fd1.targperm);
           break;
	   case 'b':
	      printf("ENTER CONDUCTIVITY OF TARGET: ");
	      scanf("%f",&fd1.targcond);
           break;
	   case 'z':
           xflag=1;		   
	   break;
       }
       if(xflag==1) break;  
 }
}
void antmenu() {
char ch;
short int dum1,xflag=0;
 for( ;;) {
	printf("1 = X-POSITION OF TRANSMIT ANTENNA:   %f\n",fd1.xcen);      
	printf("2 = Y-POSITION OF TRANSMIT ANTENNA:   %f\n",fd1.ycen);      
	printf("3 = FLARE ANGLE OF TRANSMIT ANTENNA:  %f\n",fd1.xang);      
	printf("4 = X-POSITION OF RECEIVE ANTENNA:    %f\n",fd1.xrcen);      
	printf("5 = Y-POSITION OF RECEIVE ANTENNA:    %f\n",fd1.yrcen);      
	printf("6 = FLARE ANGLE OF RECEIVE ANTENNA:   %f\n",fd1.xrang);      
	printf("7 = LENGTH OF TRANSMIT ANTENNA        %f\n",fd1.xlen);      
	printf("8 = LENGTH OF RECEIVE  ANTENNA        %f\n",fd1.xrlen);      
	printf("9 = ANTENNA HEIGHT ABOVE INTERFACE:   %f\n",fd1.anthei);      
	printf("a = DEPTH OF INTERFACE FROM GRID TOP  %f\n",fd1.zinter);      
	printf("b = GROUND CONDUCTIVITY:              %f\n",fd1.zcond);      
	printf("c = GROUND PERMITTIVITY:              %f\n",fd1.zperm);      
	printf("d = GROUND MAGNETIC PERMEABILITY:     %f\n",fd1.zur);