📄 fdtd3d.txt
字号:
for ( i=1; i < IE; i++ ) {
for ( j=ja; j <= jb; j++ ) {
for ( k=1; k < KE; k++ ) {
curl_h = ( hx[i][j][k] - hx[i][j][k-1]
- hz[i][j][k] + hz[i-1][j][k]) ;
dy[i][j][k] = gi3[i]*gk3[k]*dy[i][j][k]
+ gi2[i]*gk2[k]*.5* curl_h ;
} } }
for ( i=1; i < IE; i++ ) {
for ( j=jb+1; j < JE; j++ ) {
jyh = j - jb - 1;
for ( k=1; k < KE; k++ ) {
curl_h = ( hx[i][j][k] - hx[i][j][k-1]
- hz[i][j][k] + hz[i-1][j][k]) ;
idyh[i][jyh][k] = idyh[i][jyh][k] + gj0[j]*curl_h;
dy[i][j][k] = gi3[i]*gk3[k]*dy[i][j][k]
+ gi2[i]*gk2[k]*.5*( curl_h + gj1[j]*idyh[i][jyh][k]);
} } }
/* Calculate the Dz field */
for ( i=1; i < IE; i++ ) {
for ( j=1; j < JE; j++ ) {
for ( k=0; k < ka; k++ ) {
curl_h = ( hy[i][j][k] - hy[i-1][j][k]
- hx[i][j][k] + hx[i][j-1][k]) ;
idzl[i][j][k] = idzl[i][j][k] + gk0[k]*curl_h;
dz[i][j][k] = gi3[i]*gj3[j]*dz[i][j][k]
+ .5*gi2[i]*gj2[j]*( curl_h + gk1[k]*idzl[i][j][k] );
} } }
for ( i=1; i < IE; i++ ) {
for ( j=1; j < JE; j++ ) {
for ( k=ka; k <= kb; k++ ) {
curl_h = ( hy[i][j][k] - hy[i-1][j][k]
- hx[i][j][k] + hx[i][j-1][k]) ;
dz[i][j][k] = gi3[i]*gj3[j]*dz[i][j][k]
+ gi2[i]*gj2[j]*.5* curl_h ;
} } }
for ( i=1; i < IE; i++ ) {
for ( j=1; j < JE; j++ ) {
for ( k=kb+1; k < KE; k++ ) {
kzh = k - kb - 1;
curl_h = ( hy[i][j][k] - hy[i-1][j][k]
- hx[i][j][k] + hx[i][j-1][k]) ;
idzh[i][j][kzh] = idzh[i][j][kzh] + gk0[k]*curl_h;
dz[i][j][k] = gi3[i]*gj3[j]*dz[i][j][k]
+ .5*gi2[i]*gj2[j]*( curl_h + gk1[k]*idzh[i][j][kzh] );
} } }
/* Source */
pulse = exp(-.5*(pow((t0-T)/spread,2.0) ));
dz[ic][jc][kc] = pulse;
printf("%4.0f %6.2f \n ",T,pulse);
/* Calculate the E from D field */
/* Remember: part of the PML is E=0 at the edges */
for ( i=1; i < IE-1; i++ ) {
for ( j=1; j < JE-1; j++ ) {
for ( k=1; k < KE-1; 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] ;
}
}
}
/* Calculate the Hx field */
for ( i=0; i < ia; i++ ) {
for ( j=0; j < JE-1; j++ ) {
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ey[i][j][k+1] - ey[i][j][k]
- ez[i][j+1][k] + ez[i][j][k]) ;
ihxl[i][j][k] = ihxl[i][j][k] + fi0[i]*curl_e;
hx[i][j][k] = fj3[j]*fk3[k]*hx[i][j][k]
+ fj2[j]*fk2[k]*.5*( curl_e + fi1[i]*ihxl[i][j][k] );
} } }
for ( i=ia; i <= ib; i++ ) {
for ( j=0; j < JE-1; j++ ) {
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ey[i][j][k+1] - ey[i][j][k]
- ez[i][j+1][k] + ez[i][j][k]) ;
hx[i][j][k] = fj3[j]*fk3[k]*hx[i][j][k]
+ fj2[j]*fk2[k]*.5* curl_e ;
} } }
for ( i=ib+1; i < IE; i++ ) {
ixh = i - ib-1;
for ( j=0; j < JE-1; j++ ) {
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ey[i][j][k+1] - ey[i][j][k]
- ez[i][j+1][k] + ez[i][j][k]) ;
ihxh[ixh][j][k] = ihxh[ixh][j][k] + fi0[i]*curl_e;
hx[i][j][k] = fj3[j]*fk3[k]*hx[i][j][k]
+ fj2[j]*fk2[k]*.5*( curl_e + fi1[i]*ihxh[ixh][j][k] );
} } }
/* Calculate the Hy field */
for ( i=0; i < IE-1; i++ ) {
for ( j=0; j < ja; j++ ) {
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ez[i+1][j][k] - ez[i][j][k]
- ex[i][j][k+1] + ex[i][j][k]) ;
ihyl[i][j][k] = ihyl[i][j][k] + fj0[j]*curl_e ;
hy[i][j][k] = fi3[i]*fk3[k]*hy[i][j][k]
+ fi2[i]*fk3[k]*.5*( curl_e + fj1[j]*ihyl[i][j][k] );
} } }
for ( i=0; i < IE-1; i++ ) {
for ( j=ja; j <= jb; j++ ) {
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ez[i+1][j][k] - ez[i][j][k]
- ex[i][j][k+1] + ex[i][j][k]) ;
hy[i][j][k] = fi3[i]*fk3[k]*hy[i][j][k]
+ fi2[i]*fk2[k]*.5*curl_e ;
} } }
for ( i=0; i < IE-1; i++ ) {
for ( j=jb+1; j < JE; j++ ) {
jyh = j - jb-1;
for ( k=0; k < KE-1; k++ ) {
curl_e = ( ez[i+1][j][k] - ez[i][j][k]
- ex[i][j][k+1] + ex[i][j][k]) ;
ihyh[i][jyh][k] = ihyh[i][jyh][k] + fj0[j]*curl_e ;
hy[i][j][k] = fi3[i]*fk3[k]*hy[i][j][k]
+ fi2[i]*fk3[k]*.5*( curl_e + fj1[j]*ihyh[i][jyh][k] );
} } }
/* Calculate the Hz field */
for ( i=0; i < IE-1; i++ ) {
for ( j=0; j < JE-1; j++ ) {
for ( k=0; k < ka; k++ ) {
curl_e = ( ex[i][j+1][k] - ex[i][j][k]
- ey[i+1][j][k] + ey[i][j][k] );
ihzl[i][j][k] = ihzl[i][j][k] + fk0[k]*curl_e;
hz[i][j][k] = fi3[i]*fj3[j]*hz[i][j][k]
+ .5*fi2[i]*fj2[j]*( curl_e + fk1[k]*ihzl[i][j][k] );
} } }
for ( i=0; i < IE-1; i++ ) {
for ( j=0; j < JE-1; j++ ) {
for ( k=ka; k <= kb; k++ ) {
curl_e = ( ex[i][j+1][k] - ex[i][j][k]
- ey[i+1][j][k] + ey[i][j][k] );
hz[i][j][k] = fi3[i]*fj3[j]*hz[i][j][k]
+ fi2[i]*fj2[j]*.5* curl_e ;
} } }
for ( i=0; i < IE-1; i++ ) {
for ( j=0; j < JE-1; j++ ) {
for ( k=kb+1; k < KE; k++ ) {
kzh = k - kb - 1;
curl_e = ( ex[i][j+1][k] - ex[i][j][k]
- ey[i+1][j][k] + ey[i][j][k] );
ihzh[i][j][kzh] = ihzh[i][j][kzh] + fk0[k]*curl_e;
hz[i][j][k] = fi3[i]*fj3[j]*hz[i][j][k]
+ .5*fi2[i]*fj2[j]*( curl_e + fk1[k]*ihzh[i][j][kzh] );
} } }
}
/* ---- End of the main FDTD loop ---- */
/* printf( "idxl \n");
for ( i=1; i < ia; i++ ) {
printf( "%2d ",i);
for ( j=0; j < JE; j++ ) {
printf( "%6.3f",idxl[i][j][kc]);
}
printf( " \n");
}
printf( "idxh \n");
for ( i=ib+1; i < IE; i++ ) {
ixh = i - ib - 1;
printf( "%2d ",i);
for ( j=0; j < JE; j++ ) {
printf( "%6.3f",idxh[ixh][j][kc]);
}
printf( " \n");
}
printf( "ihxl \n");
for ( i=1; i < ia; i++ ) {
printf( "%2d ",i);
for ( j=0; j < JE; j++ ) {
printf( "%6.3f",ihxl[i][j][kc]);
}
printf( " \n");
}
printf( "ihxh \n");
for ( i=ib+1; i < IE; i++ ) {
ixh = i - ib - 1;
printf( "%2d ",i);
for ( j=0; j < JE; j++ ) {
printf( "%6.3f",ihxh[ixh][j][kc]);
}
printf( " \n");
}
printf( "idyl \n");
for ( j=1; j < ja; j++ ) {
printf( "%2d ",j);
for ( i=0; i < IE; i++ ) {
printf( "%6.3f",idyl[i][j][kc]);
}
printf( " \n");
}
printf( "idyh \n");
for ( j=jb+1; j < JE; j++ ) {
jyh = j - jb - 1;
printf( "%2d ",j);
for ( i=0; i < IE; i++ ) {
printf( "%6.3f",idyh[i][jyh][kc]);
}
printf( " \n");
}
printf( "ihyl \n");
for ( j=1; j < ja; j++ ) {
printf( "%2d ",j);
for ( i=0; i < IE; i++ ) {
printf( "%6.3f",ihyl[i][j][kc]);
}
printf( " \n");
}
printf( "ihyh \n");
for ( j=jb+1; j < JE; j++ ) {
jyh = j - jb - 1;
printf( "%2d ",j);
for ( i=0; i < IE; i++ ) {
printf( "%6.3f",ihyh[i][jyh][kc]);
}
printf( " \n");
} */
/* Write the E field out to a file "Ez" */
fp = fopen( "Ez","w");
for ( j=0; j < JE; j++ ) {
for ( i=0; i < IE; i++ ) {
fprintf( fp,"%9.6f ",ez[i][j][kc]);
}
fprintf( fp," \n");
}
fclose(fp);
/* Write the E field out to a file "Ezk" */
fp = fopen( "Ezk","w");
for ( k=0; k < KE; k++ ) {
for ( i=0; i < IE; i++ ) {
Eav = ( ez[i][jc][k-1] + ez[i][jc][k])/2.;
fprintf( fp,"%7.4f ",Eav);
}
fprintf( fp," \n");
}
fclose(fp);
/* Write the E field out to a file "Exk" */
fp = fopen( "Exk","w");
for ( k=0; k < KE; k++ ) {
for ( i=0; i < IE; i++ ) {
Eav = ( ex[i-1][jc][k] + ex[i][jc][k])/2.;
fprintf( fp,"%7.4f ",Eav);
}
fprintf( fp," \n");
}
fclose(fp);
/* Write the Hx field out to a file "Hx" */
fp = fopen( "Hx","w");
for ( k=0; k < KE; k++ ) {
for ( i=0; i < IE; i++ ) {
fprintf( fp,"%7.4f ",hx[i][jc][k]);
}
fprintf( fp," \n");
}
fclose(fp);
printf( "JC Plane \n");
printf( "Ez \n");
for ( k=0; k < KE; k++ ) {
printf( "%2d ",k);
for ( i=0; i < IE; i++ ) {
printf( "%5.2f",ez[i][jc][k]);
}
printf( " \n");
}
printf( "KC Plane \n");
printf( "Ez \n");
for ( j=0; j < JE; j++ ) {
printf( "%2d ",j);
for ( i=0; i < IE; i++ ) {
printf( "%5.2f",ez[i][j][kc]);
}
printf( " \n");
}
printf("T = %4.0f \n",T);
}
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -