1052阅读
3回复

[求助]新手请教问题，谢谢高手指教呀

离线thiker_2009

UID ：58085

注册：2010-04-26
登录：2010-06-02
发帖：2
等级：旁观者

0楼发表于: 2010-05-18 17:06:53

下面是一个一维的FDTD算法，其中的scfact=dt/muz/dx是什么意思和用途呀
谢谢高手指点。
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);
subplot(2,1,1),plot(x,ez(1:ie)/scfact,'r'),axis([0 3 -1 1]);
ylabel('EZ');
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');
rect(1:2)=[0 0];
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);
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('EZ');
subplot(2,1 ..

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