色散:金屬波導(2DTM)

設定參數

1. 主要設定 (Main-Setup)

二維 2D TM (ExEzHy)

設定總大小

2. 波源與觀測 (Wave & Ob.)

設定為softsource ，使用載入檔案做為激發源，設計結構中輸入以下函數

LoadSource(ib/2,1:jb,1:kb)=1;

LoadSource(ib/2+3,1:jb,1:kb)=1;

3. 頻譜 & 場 (Spe.&Field)

4. 邊界條件 (B.C.)

因為是要計算單一圓柱體散射吸收頻譜，因此四周圍皆為使用吸收層PML

5. 結構設計 (Geo. Design)

勾選 ☑ 設計結構 (Design Structures) 勾選 ☑ E,H,V,I (FFT) : 量測電場、磁場、電壓、電流

【指令輸入方塊】輸入以下程式碼，然後按【執行】，確定好結構無誤的話後再按【輸出結構】

d=200e-9;

%=================矩形 Brick===========================

xstart=ib/2*gdx+d/2;

xend=ib*gdx;

ystart=1*dy;

yend=1*dy;

zstart=1*gdz;

zend=1*gdz;

nindex=1^2;

sigma=0;

choice='E_Model1'; %E_Iso,PEC,M_Iso,PMC,E_Model1,M_Model1,EM_Model1

gridtype=-1; %

Iso_Brick(choice,gridtype,nindex,sigma,xstart,xend,ystart,yend,zstart,zend)

%=================矩形 Brick===========================

xstart=1*gdx;

xend=ib/2*gdx-d/2;

ystart=1*dy;

yend=1*dy;

zstart=1*gdz;

zend=1*gdz;

nindex=1^2;

sigma=0;

choice='E_Model1'; %E_Iso,PEC,M_Iso,PMC,E_Model1,M_Model1,EM_Model1

gridtype=-1; %

Iso_Brick(choice,gridtype,nindex,sigma,xstart,xend,ystart,yend,zstart,zend)

%=================矩形 Brick===========================

%------------探測磁場 1 Probe Magnetic 1------------

% direction 1(Hx),2(Hy),3(Hz),4(|H|)

direction=2;

x1=(ib/2+1)*gdx;

x2=(ib/2+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_magnetic(1,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%------------探測磁場 2 Probe Magnetic 1------------

% direction 1(Hx),2(Hy),3(Hz),4(|H|)

direction=2;

x1=(ib/3+1)*gdx;

x2=(ib/3+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_magnetic(2,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%------------探測磁場 3 Probe Magnetic 1------------

% direction 1(Hx),2(Hy),3(Hz),4(|H|)

direction=2;

x1=(ib/4+1)*gdx;

x2=(ib/4+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_magnetic(3,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%------------探測電場 1 Probe Electric 1------------

% direction 1(Ex),2(Ey),3(Ez),4(|E|)

direction=3;

x1=(ib/2+1)*gdx;

x2=(ib/2+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_electric(1,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%------------探測電場 2 Probe Electric 1------------

% direction 1(Ex),2(Ey),3(Ez),4(|E|)

direction=3;

x1=(ib/3+1)*gdx;

x2=(ib/3+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_electric(2,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%------------探測電場 3 Probe Electric 1------------

% direction 1(Ex),2(Ey),3(Ez),4(|E|)

direction=3;

x1=(ib/4+1)*gdx;

x2=(ib/4+1)*gdx;

y1=(jcenter)*gdy;

y2=(jcenter)*gdy;

z1=(kcenter)*gdz;

z2=(kcenter)*gdz;

probe_electric(1,:)=[x1,x2,y1,y2,z1,z2,direction];

%----------------------------------

%看結構是否填入

figure(1);plot(geox)

LoadSource(ib/2,1:jb,1:kb)=1;

LoadSource(ib/2+3,1:jb,1:kb)=1;

6. 電,磁性材料模型(E,M Models)

前面範例有金屬材料交界面的反射率範例(交界面interface)，在這裡我們試著解說如何自行輸入電性金屬材料參數、及如何利用CompactFDTD求得色散曲線。

例如若有一組材料是eV單位，要在風行裡使用必須經過單位轉換，參考

S.-H. Chang, S. Gray, and G. Schatz, “Surface plasmon generation and light transmission by isolated nanoholes and arrays of nanoholes in thin metal films,” Opt. Express 13, 3150-3165 (2005).

http://www.opticsinfobase.org/oe/abstract.cfm?uri=oe-13-8-3150

裡有提到模擬金材料的Drude參數 500nm-1000nm 是

epsilon_infinite = 11.4577;

omegad =9.40274 eV;

gamad =0.0831375 ev;

轉成風行可用單位

ehbar = 1.51926751447914e+015; % e/hbar where hbar=h/(2*pi) and e=1.6e-19

omegad= 9.40274/27.2214/(2.41888*10^-17)=9.40274*1.51926751447914e+015= 1.428527742909359e+016 rad/s

gamad =0.0831375/27.2214/(2.41888*10^-17)=0.0831375*1.51926751447914e+015=1.263081029850095e+014 rad/s

7. 模擬計算 (Calculation)

因為使用探測場Probe，所以不需要勾選穩態場輸出 (Output Steady State - Fields) -(Poynting DFT)。

運算結束後

分析結果(Result Analysis)

N_Probed_E_Td(觀測位置編號,計算步數,分割數) %量測電場(時域) Probed Electric (Time-domain)

N_Probed_H_Td(觀測位置編號,計算步數,分割數) %量測磁場(時域) Probed Magnetic (Time-domain)

N_Probed_E_Fd(觀測位置編號,計算步數,分割數) %量測電場(頻域) Probed Electric (Frequency-domain)

N_Probed_H_Fd(觀測位置編號,計算步數,分割數) %量測磁場(頻域) Probed Magnetic (Frequency-domain)

★這裡您只需要調整兩個參數即可，第一個是觀測場的輸入，我們前面設定四個觀測場，其中三個磁場位置，一個電場位置

★要做模態判定，一次只能輸入電場或磁場其中一種。

★Probes=N_Probed_H_Td; 這個方式，是將三個觀測磁場都輸入，Mode_Analysis會將三個觀測的磁場做平均化來分析尋找模態

★Probes=N_Probed_H_Td(3,:,:); 這個方式，是將第三個觀測磁場都輸入，Mode_Analysis只會用第三個觀測的磁場來做分析尋找模態

★FieldMax 設定臨界值，高於此值的將視為找到的模態，小於此值的將過濾掉。若找到太多雜訊模態，調整前面的係數0.01變大來過濾

★max_band_number 在一輸入波數k(wave number)位置上，最多找幾個模態

★lambda0 中心入射波長

【指令輸入方塊】輸入以下程式碼，然後按【執行】

Probes=N_Probed_E_Td;

%figure(1);plot(Probes(1,:,10))

FieldMax=0.05*max(Probes(:))

max_band_number=6;

Mode_Analysis(Probes,FieldMax,max_band_number)

%==理論值對比 Analytical Solution======

clear RunNumber

RunNumber=100;

% lambda0=500e-9;

% cc=3e8;

% freq0=3e8/lambda0;

% omegas=linspace(0.01,1.0,RunNumber)*2*pi*freq0;

w1=2*pi*cc/lambda_end;

w2=2*pi*cc/lambda_start*1.2; %x1.2是將理論值觀察範圍延伸至更高頻的位置

omega_input=linspace(w1,w2,RunNumber);

epsr_insulator= 1;

% 自行輸入Drude模型

% epsInf =11.4577;

% WDp = 9.40274/27.2214/(2.41888*10^-17);

% gama = 0.0831375/27.2214/(2.41888*10^-17);

% epsr_metal= epsInf-WDp^2./(omegas.^2+omegas*gama*j);

%使用風行內建函數

epsr_cal(omega_input) ;

epsr_metal= epsr_model_1;

figure(1);hold on;

plot(1./omega_input,real(epsr_model_1),'r','linewidth',3);

plot(1./omega_input,imag(epsr_model_1),'b','linewidth',3);

xlabel('\lambda (nm)');

ylabel('\epsilon r');

hold off;

% 交界面色散關係式 Metal-insulator interface dispersion relation

kxx=(omega_input/cc).*((epsr_insulator*epsr_metal)./(epsr_metal+epsr_insulator) ).^0.5;

f_para=kxx; %尋找 k 的參考值 reference value

d=200e-9;

%風行計算金屬波導色散的函數(MIM-Waveguide Dispersion)，輸入方式如下

MIM_Waveguide(RunNumber,omega_input,d,epsr_metal,epsr_insulator,f_para);

figure(2);

hold on;

plot(omega_input/cc,omega_input,'k','LineWidth',2);

plot(real(kxx),omega_input,'g','LineWidth',2);

plot(real(MIM_Tanh),omega_input,'b' ,'LineWidth',2);

plot(real(MIM_Coth),omega_input,'m' ,'LineWidth',2);

plot(betaz,Eigen_omega(:,:),'r*');

legend('light-line','Analytical-interface','Analytical-MIM(Tanh)','Analytical-MIM(Coth)','FDTD')

title('金屬波導色散關係曲線')

ylabel('ω');

xlabel('k');

axis([0 2.e7 0 5e15])

hold off;

關於金屬波導(Metal-insulator-Metal)色散關係，請參考

H. Raether, Surface Plasmons (Springer, New York, 1998).

或者直接google查詢

MIM dispersion relation

就會有許多相關資料