电磁场数值模拟论文
研究生《自然辩证法》
课程论文
题目：道德与利益的辨证分析
研究生姓名：韩友华
学号：1623090026
电话号码: 184********
专业：16级光学工程
单位：云南师范大学物理与电子信息学院授课时间：2016-2017第二学期
授课教师：李际
学生签名：
成绩：教师签名：
评分标准
1．研究内容有价值和意义
2．研究有一定的结果或结论
3．书写规范整齐
简易给定中心热源温度扩散模型
韩友华1，姚斌1
（1 云南师范大学物理与电子信息学院，昆明650500）
摘要：本文通过对中心稳定面热源向四周扩散的分析，在matlab环境下进行模拟，采用时域有限差分法（FDTD）的显示差分策略，对其进行仿真显示。结果表面：模场分布变化较大，随着激励过后时间的增长，中心温度逐渐降低，热量向四周扩散，四周温度逐渐上升；结果还表明，越靠近中心热源，其温度变化越迅速，越远离热源，其温度变化越缓慢。模拟结果与理论讨论结果一致，验证了模拟和理论的正确性。
关键字：温度扩散，时域有限差分，温度变化；
1 引言
本文运用这种方法对给定的中心热源的温度扩散进行研究，使用给定的热源进行激励，以计算机作为数值计算的工具进行理论求解，所得的模场分布图与理论结果一致，证明了时域有限差分法（FDTD）法用于激励模场研究的可行性和有效性。
2 理论
扩散现象（diffusion）是指物质分子从高浓度区域向低浓度区域转移直到均匀分布的现象，速率与物质的浓度梯度成正比。扩散是由于分子热运动而产生的质量迁移现象，主要是由于密度差引起的。
温度不同，水的内能不同,内能是分子运动剧烈程度的标志,那么分子运动剧烈程度就不同。温度越高分子运动越剧烈,扩散能力就越强。
时域有限差分法（FDTD）是用中心差商代替场量对时间和空间的一阶偏微商, 通过在时域的递推模拟波的传播过程, 从而得出场分布；它最早由 K.S.Yee 于1966 年提出，自 80 年代末，时域有限差分法成为电磁场数值计算的重要方法之一。
时域有限差分法（FDTD）的显示差分策略是指已知外界条件（包括激励源和边界条件）下，可以模拟出下一状态的模场分布。
（a）显示差分策略
对于温度扩散模型，满足时域有限差分法（FDTD）的显示差分策略公式：
（1）
化简之后得：
（2）其中：
3 模拟结果与分析
对所分析的温度扩散模型空间进行了等网格划分，二维参数边界为x=0和50，y=0和50，z为随时间变化的温度大小。激励源为：
if (j-25)*(j-25)+(k-25)*(k-25)<=25
u(:,:,1)=4+abs((j-25)*(j-25)+(k-25)*(k-25)-25);
else
u(:,:,1)=100/((j-25)*(j-25)+(k-25)*(k-25));
end
end
所表示的热源，对此激励源下的模场分布进行模拟。模拟结果如下图：
（b）激励源（循环次数40）
（c）循环次数233
（d）循环次数600
如（b）（c）（d）图所示，在前200次循环中，温度下降14.5℃，下降幅度速率明显；后400次循环中，温度下降4.6℃，下降幅度较小，速率降缓，向四周面积扩散的速率也变慢。随着时间的增长，中心温度越来越低，四周温度逐渐越来越高。不仅如此，热量向四周扩散，扩散范围越来越大。
中心热源温度平面向四周扩散趋势如图所示：
（e）循环次数20
（f）循环次数211
（g）循环次数600
如上图（e）（f）（g）所示，前200次循环中扩散面积为100π，扩散速率较快；200到600次循环时，扩散面积位125π，速率相对之前明显降缓。中心温度向四周逐渐扩散，向四周面积扩散速率较快可以推测，如果时间足够长，整个平面区域温度分布会几乎均匀。
本文用标量FDTD 法, 分析了给定中心热源模场分布情况，得到了模场各个时间点的分布图，与解析法和其它方法相比，具有直观、便捷的特点，而且与理论分析所得结果一致。
参考文献
[1] 黄华. 时域有限差分法分析混晌室中场的均匀性(西安电子科技大学电子工程学院，陕
西西安7l0071)
[2]张祺. 加脊型波导系统的二维频域有限差分法分析[J]. 广西科学，2010,17(3):239-241.
[3] 黄重庆，张国平，刘靖. 条形波导光场分布的时域有限差分法模拟[J]. 光电子激光，
2002 ,13(1):20-22.
[4] 张岩，吕善伟. ADI-FDTD + GRT 在波导电路分析中的应用.（北京航空航天大学电子
信息工程学院，北京100083）
[5] 曹世昌. 电磁场的数值计算和微波的计算机辅助设计[M]. 电子工业出版社，1989:1-10.
程序目录：
%thid program is used for the presentation of energy diffusion
clc
clear all
b=1.0e-3;
l=50;
x0=0; xf=l;
y0=0; yf=l;
tmax=600;
syms x y t;
Phi_xy=0;
mu1_yt=0;
mu2_yt=0;
v1_xt=0;
v2_xt=0;
Mx=50; hx=l/Mx; j=0:1:Mx; x=x0+j*hx;
My=50; hy=l/My; k=0:1:My; y=y0+k*hy;
Mt=10; tau=tmax/Mt; n=0:1:Mt;t=n*tau;
for n=1:tmax
for k=0:My
u(k+1,1,n)=0;
end
for k=0:My
u(k+1,Mx+1,n)=0;
end
for j=0:Mx
u(1,j+1,n)=0;
end
for j=0:Mx
u(My+1,j+1,n)=0;
end
h=hx;
c=b*tau/(2*h^2);
for k=2:My
for j=2:Mx
if (j-25)*(j-25)+(k-25)*(k-25)<=25
u(:,:,1)=4+abs((j-25)*(j-25)+(k-25)*(k-25)-25);
else
u(:,:,1)=100/((j-25)*(j-25)+(k-25)*(k-25));
end
u(k,j,n+1)=u(k,j,n)+c*(u(k,j+1,n)-2*u(k,j,n)+u(k,j-1,n))+c*(u(k+1,j,n )-2*u(k,j,n)+u(k-1,j,n));
end
end
u(:,:,n+1)
% surf(x,y,u(:,:,n+1),'linewidth',4)
contour( u(:,:,n+1))
title(n)
shading interp
axis([0 50 0 50 ])
xlabel('x','fontsize',10)
ylabel('y','fontsize',10)
zlabel('z','fontsize',10)
pause(0.0000000000000000000000000001) end