Lumerical官方案例、FDTD时域有限差分法仿真学习（九）——布拉格光栅（Bragg gratings）

最新推荐文章于 2023-07-24 23:19:39 发布

子非鱼icon

最新推荐文章于 2023-07-24 23:19:39 发布

阅读量1.1w

点赞数 17

分类专栏： FDTD学习笔记文章标签： FDTD 时域有限差分法 Lumerical 布拉格光栅硅光器件仿真

本文链接：https://blog.csdn.net/weixin_45719141/article/details/123843250

版权

FDTD学习笔记专栏收录该内容

34 篇文章 663 订阅

订阅专栏

一、使用 FDTD 的布拉格光栅初始设计

在本例中，将使用 3D FDTD 模拟来了解布拉格光栅的性能如何受几何参数（如波纹深度和未对准）的影响。

在这里插入图片描述
1.1 背景

波导布拉格光栅是一维光子带隙结构的示例，其中对直波导的周期性扰动形成特定波长的介电镜。这些器件通常用作实现波长选择功能的滤光片。

1.1 仿真设置

在本例中，我们将使用光栅的单个晶胞的 3D FDTD 模拟来找到无限周期器件的中心波长和带宽。在 Bragg_FDTD_unit_cell.fsp 中，模拟区域恰好包含光栅的一个单位单元。 Bloch 边界用于 x 最小/最大边界，这允许我们为无限周期设备设置 kx 值。使用模式源作为激发，使用能带结构分析组（ bandstructure analysis group）计算光谱。对于此模拟，我们对带边缘 kx = π/a （文献里面有这个）处的光谱感兴趣，这将为我们提供光栅带隙的大小和位置。

1.2 仿真结果

1.2.1 带宽和中心波长

一旦模拟完成运行，谱带结构分析组将返回 kx = π/a 处的谱，如下所示。从光谱中的两个峰值我们可以找到带隙的大小和位置，它们对应于布拉格光栅的带宽和中心波长。请注意，峰的宽度随着模拟时间的增加而减小；因此，模拟时间应该足够长以产生清晰的单独峰值。附件 Bragg_FDTD_unit_cell.fsp 中的原始仿真时间为 1250fs。

在这里插入图片描述

不同模拟时间的频谱结果。

1.2.2 耦合系数作为光栅深度的函数

即使这种 FDTD 方法适用于无限周期光栅，它对于设计有限长度光栅仍然非常有用。通过这种方法，人们可以快速了解各种设计参数如何改变光栅的性能。例如，我们可以运行波纹深度的参数扫描（参见 Bragg_FDTD_unit_cell.fsp 中的优化和扫描选项卡）并计算光栅耦合系数

$\kappa=\frac{\pi n_g \Delta \lambda}{\lambda_0^2}$

其中 Δλ 是带宽，λo 是中心波长，ng 是 λo 处的群折射率。

在这里插入图片描述上述结果需要将模拟时间增加到 3750fs，因为对于较弱的光栅（光栅深度较小），光谱中的峰值更难区分。

1.2.3 耦合系数作为侧壁波纹中未对准的函数

在参考文献 [1] 中，应用相同的 FDTD 方法来显示如何通过改变两个侧壁上的波纹之间的错位来调整光栅强度。
在这里插入图片描述

1.3 仿真结构设置

波导材料采用洛伦兹模型，无色散损耗（名字是这个）。
在这里插入图片描述

光栅结构参数设置

在这里插入图片描述
参考文献中的结构设置：

在这里插入图片描述

光源范围：
1.52-1.56um

可以选择注入基本模式（fundamental mode），即具有最高有效指数的模式。或基本TE或基本TM模式，他们对应于具有最高有效折射率的形式，并具有所需要的类TE（TE-like）或类TM（TM-like）极化。
在FDTD求解中，类TE模态对应于沿x轴注入的源电场极化主要沿y方向的模态，而沿y轴或z轴注入的源电场极化主要沿x方向的模态。

在这里插入图片描述

仿真区域：

仿真时间设置为1250fs
在这里插入图片描述
Bloch边界条件类似于周期边界，但他们考虑了在给定的光源注入角度下，器件各周期之间的相位差。Bloch边界使每个单元之间的相位差不受波长的影响。

在这里插入图片描述
能带结构分析组：

Setup里的脚本

deleteall;
##############################################
# Bandstructure
# This script creates n randomly positioned monitors
# 此脚本创建 n 个随机定位的监视器
# Input properties
# x,y,z span:     The span of the random cloud of monitors
# 监视器随机云的跨度
# n monitors:     The number of monitors
# 监视器随机云的跨度
# Tags: bandstructure time monitors band gap photonic crystal resonant frequencies
# 能带结构时间监控器带隙光子晶体谐振频率
# Copyright 2012 Lumerical Solutions Inc
##############################################

# simplify variable names by removing spaces
x_span = %x span%;
y_span = %y span%;
z_span = %z span%;
n_monitors = %n monitors%;

# Set the random number generator seed.  Change this values
# to get different random orientation of objects
# 设置随机数生成器种子。 更改此值以获得对象的不同随机方向
randreset(123);

for(i=1:n_monitors) {  
  addtime;
  set("name","m"+num2str(i));
  set("x",(rand-0.5)*x_span);# rand是一个在0-1之间的随机数
  set("y",(rand-0.5)*y_span);
  set("z",(rand-0.5)*z_span);
}

1.3.1 randreset
在这里插入图片描述
1.3.2 addtime

1.3.3 rand

在这里插入图片描述
Analysis里的脚本

在这里插入图片描述

##########################################################################
# Bandstructure
# This script is to be used in bandstructure calculations. 
#   此脚本用于能带结构计算
# Input properties
# f1, f2: min, max frequency of interest 感兴趣频率的最大最小值
# apod width: apodization width as a fraction of total simulation time 
# 变迹宽度占总模拟时间的部分
# apod center: center location of the apodization filter as a fraction of total simulation time
# 变迹滤波器的中心位置占总模拟时间的一小部分
# Output properties
# spectrum: result from the summation of ffts of time signals from each 
#     monitor - this is a spectrum with peaks at resonant frequencies
# 频谱：来自每个监视器的时间信号 ffts 总和的结果 - 这是在谐振频率处具有峰值的频谱
# Tags: bandstructure time monitors band gap photonic crystal resonant frequencies
#
# Copyright 2012 Lumerical Solutions Inc
###########################################################################

# simplify variable names by removing spaces
apod_width = %apod width%;
apod_center = %apod center%;
n_monitors = %n monitors%;

if (f1<=0) { f1 = f2/10000; }  # don't allow f1=0, as this causes a minor problem:  f=0 corresponds to lambda=infinity, and the infinity causes problems in the data visualizer.
# # 不允许 f1=0，因为这会导致一个小问题：f=0 对应于 lambda=infinity，而无穷大会导致数据可视化器出现问题。



t = getdata("m1","t");
fs = 0;

# loop over each time monitor # 循环每个时间监视器
for(j=1:n_monitors) {
  for (component_i=1:3) { # loop over field components 遍历场分量
    if (component_i==1) { component="Ex"; }
    if (component_i==2) { component="Ey"; }
    if (component_i==3) { component="Ez"; }
	
    signal = pinch(getdata("m"+num2str(j),component)); # collect the field data  收集场数据

    #apodize the signal with a gaussian, to ignore start and end effects
    # 用高斯变迹信号，忽略开始和结束的影响
    signal = signal*exp( - 0.5*(t-max(t)*apod_center)^2/(apod_width*max(t))^2);

    # fft the time signal and add it to the other fft's
    # 将时间信号fft时间信号并将其添加到其他 fft并将其添加到其他 fft
    f = linspace(f1,f2,5000);
    fs = fs + abs(czt(signal,t,2*pi*f))^2; # czt调频z变换
  }
}

# create data set
spectrum = matrixdataset("spectrum"); # 创建了名为spectrum的空数据集
spectrum.addparameter("f",f,"lambda",c/f);   # in most cases, lambda is before f in datasets.  However, for this particular application, frequency is far more common, so it is first.
# # 在大多数情况下，lambda 在数据集中位于 f 之前。 然而，对于这个特定的应用，频率要普遍得多，所以它是第一位的。
spectrum.addattribute("fs",fs);

注：这个程序有关高斯变迹信号的部分，其实没咋看懂，需要去补理论知识。

1.3.4 czt
在这里插入图片描述
参数扫描：

扫描的结果bandwidth和kappa是从哪儿来的？

在model的分析组里。右键model，然后点edit，会发现analysis选项卡下会有这两个结果。

在这里插入图片描述
下面分析其中的脚本：

Setup的Variables：
在这里插入图片描述
Setup的Script：

# set period and bloch vector of simulation
# 设置仿真周期和bloch向量
name="FDTD";
setnamed(name,"x span",ax); # FDTD的x span确实为0.32
setnamed(name,"bloch units","SI");        # we generally use bandstructure units, but for the simulation region, it's 
setnamed(name,"kx",kx*2*pi/ax);             # most convenient to use SI units becuase kx and ky can have different 
# # 我们通常使用能带结构单位，但是对于模拟区域，使用 SI 单位最方便，因为 kx 和 ky 可以有不同的
# 下面有解释

# set period and apodization settings for banstructure analysis object
# 为banstructure分析对象设置周期和变迹设置
name="bandstructure";
setnamed(name,"x span",ax);
setnamed(name,"apod center",apod_center);# 这里相当于是把model里的参数值传到分析组去
setnamed(name,"apod width",apod_width);
setnamed(name,"f1",f1);
setnamed(name,"f2",f2);

select('w1');
set('y span',0.5e-6-w);  # w就是那个波纹宽度

select('w2');
set('y span',0.5e-6+w);

setnamed(name,“bloch units”,“SI”); 这一句对应这里的设置
在这里插入图片描述
setnamed(name,“kx”,kx2pi/ax);

Variables里面定义了ax=0.32，其实就是光栅周期；kx=0.5，kx2pi/ax = π/ax 对应前文

当布拉格光栅周期约为半波导波长的整数倍时，变回抑制透射，形成禁带

在论文里的解释：
在这里插入图片描述
两个共振峰对应于 kx =π∕Λ 点处的能带边缘，两个峰之间的波长范围对应于没有传播解的带隙。

类比光纤布拉格光栅：
在这里插入图片描述
参考

Analysis的Variables：

在这里插入图片描述
Analysis的Script：

spectrum = getresult('bandstructure','spectrum');
fs = spectrum.fs;
lambda = spectrum.lambda;

bandwidth = matrix(length(w));# 扫描起来w就不只有一个值了
center_wavelength = matrix(length(w));

temp = fs(1:length(lambda));
ind = findpeaks(temp,2); # 找到两个峰值
bandwidth = abs(lambda(ind(2))-lambda(ind(1)));
center_wavelength = lambda(ind(1)) + bandwidth/2;

ng = 4.3874;
kappa = pi*ng*bandwidth/(center_wavelength^2); # 依据了kappa计算公式

1.4 仿真结果

设置了x范围为（1.5-1.6um之间），不然就是分析组脚本里的1.2-2.15um之间
在这里插入图片描述
扫描结果

论文中的结果

在这里插入图片描述
二、布拉格光栅全器件仿真与 EME

在此示例中，我们将使用 MODE 中的本征模扩展 (EME) 求解器来计算具有任意周期数的波导布拉格光栅的完整透射光谱。请注意，在许多情况下，光栅的中心波长和带宽将提供有关光栅性能的足够信息，并不总是需要完整的透射光谱。在这种情况下，建议使用布拉格光栅初始设计和 FDTD 中所示的更简单的 FDTD 方法。

2.1 背景

波导布拉格光栅是一维光子带隙结构的示例，其中对直波导的周期性扰动形成特定波长的介电镜。这些器件通常用作实现波长选择功能的滤光片。

在这里插入图片描述

2.2 仿真设置

使用 EME 模拟周期性结构时，只需要定义几何的一个单位单元。在 Bragg_EME.lms 中，EME 求解器覆盖光栅的单个单元格，如下所示。

在这里插入图片描述
求解器的每一端都设置了一个端口，用于计算基本 TE 模式的透射和反射。在 EME 设置选项卡下，我们为 EME 求解器定义了 2 个单元组，一组用于波导宽度较大的区域，另一组用于波导宽度较小的区域。由于每个细胞组内的折射率和几何形状是一致的，我们只需为每个细胞组使用 1 个细胞。对于初始模拟，我们将在 EME 计算中为每个单元组使用十种模式。请注意，此处使用对称性来减少此计算所需的模式数量。

为了设置光栅的周期，我们将在 EME 设置选项卡右侧的“周期组定义”（periodic group definition）表下定义一个周期组。起始单元组和结束单元组分别设置为 1 和 2，周期数设置为 500。这意味着单元单元（由 2 个单元组组成）将传播 500 次，设备的最终长度将是 160 微米。
在这里插入图片描述

由于 EME 是一种频域方法，原则上我们需要为每个感兴趣的波长运行一次模拟。这可能非常耗时，因为需要为每个波长重新计算所有模式，并且通常需要大量模拟来准确描述完整的透射光谱。但是，我们可以使用微扰方法覆盖波长来快速扫描波长，如果需要，可以通过全波长扫描来验证。

2.3 仿真结果

波导布拉格光栅的 500 个周期的透射光谱如下所示，使用 EME 分析窗口中的波长扫描功能计算得出，具有 400 个点。或者，也可以通过使用参数扫描直接扫描波长来获得相同的透射/反射光谱，这对于最终验证来说是一个好主意，但这种方法比使用我们的微扰方法直接扫描波长更耗时，并且只应在波长数量较少时使用。

在这里插入图片描述
EME 方法非常适合模拟有限长度波导布拉格光栅的透射光谱，因为由于所需的计算时间和内存量，整个设备对于基于 FDTD 的方法可能具有挑战性。

注：微扰法（Perturbative method）

波长扫描也可以使用波导布拉格光栅和谐振器中描述的方法有效地完成。这种方法在本质上与这里使用的内置波长扫描工具相似；但是，内置工具更通用，不需要额外的脚本即可运行。

2.4 结构设置

波导结构的绘制：

（比较简单粗暴，用polygon，输入了一个个点坐标，从左下角逆时针画到了左上角，连起来。）

在这里插入图片描述

在这里插入图片描述

材料还是第一章案例的洛伦兹模型的Si.

EME求解器设置：
在这里插入图片描述

同学习笔记（八）一样，此处由于关于y轴对称，为了节省计算量，也在此处使用了anti-symetric条件。
在这里插入图片描述
观察端口设置，已选择为基本 TE 模式

在这里插入图片描述
2.5 仿真自运行结果

（区别于上面的“2.3仿真结果”，这个主要是自己运行一下看看结果）

EME 分析窗口中的波长扫描功能位置：

在这里插入图片描述
扫描结果包含其他S11,S22等，不展示可将其remove。

通过sweep也可以进行波长扫描，如下图所示，

扫描结果如下：

这里所用点数为100，而EME自带的波长扫描那里点数设置为400，所以结果上有出入。

参考文献：
[1] X. Wang, et al., “Precise control of the coupling coefficient through destructive interference in silicon waveguide Bragg gratings”, Opt. Lett. 39, 5519-5522 (2014).
[2]X. Wang, “Silicon photonic waveguide Bragg gratings,” Ph.D. Thesis, University of British Columbia (2013).