826阅读
1回复

[求助]Taflove的FDTD书附带程序的问题请教

UID ：74059

注册：2011-03-19
登录：2013-01-14
发帖：64
等级：仿真一级

0楼发表于: 2011-08-15 22:34:28

请教高手指点一下，我有几个地方不太明白
1，一维中不是只有EX和HY吗，程序中怎么会有EZ呢？
2，movie初始化中，    rect=get(gcf,'Position');
                                     rect(1:2)=[0 0];
                                     M=moviein(nmax/2,gcf,rect);
这几句说的是什么意思？看不大懂
3，ez(2:ie)=ca*ez(2:ie)+cb*(hy(2:ie)-hy(1:ie-1))    hy(1:ie)=hy(1:ie)+ez(2:ib)-ez(1:ie) 这两个迭代公式怎么和葛德彪书上的不太一样啊，还有括号中的是时间还是空间？如果是空间的话，HY不应该是k+1/2吗？
4，if mod(n,2)==0; 这句为什么呢？只要时间的2的整数倍？

源程序：
clear
%***********************************************************************
% Fundamental constants(基本常量)
%***********************************************************************
cc=2.99792458e8; %speed of light in free space
muz=4.0*pi*1.0e-7; %permeability of free space(真空中磁导率)
epsz=1.0/(cc*cc*muz); %permittivity of free space（真空中介电常数）
freq=1.0e+9; %frequency of source excitation（激励源的周期，频率）
lambda=cc/freq; %wavelength of source excitation（激励源的波长）
omega=2.0*pi*freq;%（角频率 =欧米噶）
%***********************************************************************
% Grid parameters (网格参量)
%***********************************************************************
ie=200; %number of grid cells in x-direction（X轴向的格点数）
ib=ie+1;
dx=lambda/20.0; %space increment of 1-D lattice（空间步长条件：dx=<lambda/12）
dt=dx/cc; %time step(时间步长条件：dt=<sqir(2)*dx/2*cc)
omegadt=omega*dt;
nmax=round(12.0e-9/dt); %total number of time steps(12.0e-9计算可得,时间步数)
%***********************************************************************
% Material parameters（材料背景参量）
%*********************************************************************** .,:
eps=1.0;%（介质电导率）
sig=5.0e-3;%（电导率） .
%***********************************************************************
% Updating coefficients for space region with nonpermeable media（非渗透性材料）
%***********************************************************************
scfact=dt/muz/dx;
ca=(1.0-(dt*sig)/(2.0*epsz*eps))/(1.0+(dt*sig)/(2.0*epsz*eps));
cb=scfact*(dt/epsz/eps/dx)/(1.0+(dt*sig)/(2.0*epsz*eps));
%***********************************************************************
% Field arrays （场列阵）
%***********************************************************************
ez(1:ib)=0.0;
hy(1:ie)=0.0;
%***********************************************************************
% Movie initialization （初值化）
%***********************************************************************
x=linspace(dx,ie*dx,ie); '[M2Q"X
subplot(2,1,1),plot(x,ez(1:ie)/scfact,'r'),axis([0 3 -1 1]); /
ylabel('EZ'); ;Z0cD*Jb
subplot(2,1,2),plot(x,hy,'b'),axis([0 3 -3.0e-3 3.0e-3]);
xlabel('x (meters)');ylabel('HY');
rect=get(gcf,'Position'); J};
rect(1:2)=[0 0]; k{&E}:
M=moviein(nmax/2,gcf,rect);
%***********************************************************************
% BEGIN TIME-STEPPING LOOP（时间步循环）
%***********************************************************************
for n=1:nmax
%***********************************************************************
% Update electric fields
%***********************************************************************
ez(1)=scfact*sin(omegadt*n);
rbc=ez(ie); =C}
ez(2:ie)=ca*ez(2:ie)+cb*(hy(2:ie)-hy(1:ie-1));
ez(ib)=rbc;
%***********************************************************************
% Update magnetic fields
%***********************************************************************
hy(1:ie)=hy(1:ie)+ez(2:ib)-ez(1:ie);
%***********************************************************************
% Visualize fields
%***********************************************************************
if mod(n,2)==0;
rtime=num2str(round(n*dt/1.0e-9));
subplot(2,1,1),plot(x,ez(1:ie)/scfact,'r'),axis([0 3 -1 1]);
title(['time = ',rtime,' ns']);
ylabel( ..

未注册仅能浏览部分内容，查看全部内容及附件请先登录或注册

离线103595669

UID ：74059

注册：2011-03-19
登录：2013-01-14
发帖：64
等级：仿真一级

1楼发表于: 2011-08-16 00:12:49

有没有好心人指教一下呀。。？？？？？搞算法太累了。。。。

发帖回复