COMSOL复现PRL二次谐波产生

COMSOL复现PRL二次谐波产生搞非线性光学模拟的同行应该都懂二次谐波产生(SHG)这玩意儿在COMSOL里有多折腾。最近被PRL某篇论文刺激到了决定自己动手复现他们的微环谐振腔SHG效果。说干就干结果发现这坑比想象中还深...COMSOL复现PRL二次谐波产生首先得搞定材料非线性设置。别傻乎乎直接用内置的χ²非线性选项实测发现对于环形谐振腔这种特殊结构手动定义非线性极化项更靠谱。代码里关键是要把基频波和倍频波的耦合关系写清楚model.param.set(chi2, 2e-12[m/V]); //典型LN晶体的二阶非线性系数 model.physics(ewfd).feature(nlin1).set(Pnl, {0, chi2*ewfd.Ey^2, 0}); //y方向极化分量这里有个魔鬼细节——坐标系方向必须和实际微环结构匹配。之前翻车三次才发现论文用的环形波导是切向极化Y轴对应环的切线方向这个设置直接影响了非线性极化的空间分布模式。频率耦合设置才是真地狱。得在波动方程模块手动添加两个相互作用的频域研究study1 model.study.create(study1); study1.feature.create(freq, Frequency); study1.feature(freq).set(plist, {1.55e14}); //基频 study1.feature(freq).set(notlistsol, on); study2 model.study.create(study2); study2.feature.create(freq, Frequency); study2.feature(freq).set(plist, {3.1e14}); //倍频 study2.feature(freq).set(notlistsol, on);重点来了必须开启频域扰动选项让两个频率产生耦合否则就是在做无用功。这一步的GUI藏得比三星堆文物还深建议直接调参数对象model.sol(sol1).feature(freq).set(perturb, on); model.sol(sol1).feature(freq).set(kspType, fgmres); //迭代求解必选网格划分要特别注意环形波导的弯曲处这里电场强度可能突变。建议用边界层网格曲率修正不然计算出的转换效率能差两个数量级。有个取巧办法是用极坐标系剖分mesh1.feature(size).set(customize, on); mesh1.feature(size).set(hgrad, 1.5); //曲率自适应 mesh1.feature(dis1).set(type, cylindrical); //极坐标剖分最后验证能量守恒时发现倍频波功率异常震荡排查发现是边界条件相位不匹配。解决方法是在端口边界添加相位补偿项model.physics(ewfd).bnd(port1).set(PhaseOffset, exp(-i*k0*x)); //行波相位修正折腾完这些终于看到那个优美的二次谐波共振峰了。不过说真的COMSOL这非线性耦合的调试过程简直比让猫走直线还难。下次再碰SHG模拟可能得先备两瓶速效救心丸...