简介这份资源面向光学、磁光材料与物理仿真方向的学习者和研究人员围绕掺铈钇铁石榴石Ce:YIG这一典型一维磁光晶体聚焦其透射、反射与法拉第旋转效应的数值模拟。压缩包内共1个文件为MATLAB脚本.m格式整体约1KB体积轻量便于直接运行与二次修改。脚本可用于计算不同磁场、波长与晶体厚度条件下光的偏振旋转角及透射、反射系数帮助理解磁光隔离器、磁光调制器等器件的物理机理。目前已有438人学习下载说明该方向具备一定关注度。对于需要快速搭建Ce:YIG磁光效应仿真框架、验证法拉第旋转理论公式或开展课程设计、科研预研的读者这份代码可作为可复用的计算起点节省从零编写数值模型的时间。1. Ce_YIG 一维磁光晶体透射与法拉第旋转到底在算什么一块厚度几百微米的 YIG 晶体外加一个可调磁场就能让穿过它的线偏振光偏振面转几十度同时透射率还能保持在 80% 以上——这件事在光通信和激光系统里被反复利用但真正动手算的时候很多人卡在同一个地方透射谱和法拉第旋转角到底怎么从材料参数推出来边界条件怎么设磁场方向怎么定。Ce_YIG 是在 YIG 里掺铈目的是把法拉第旋转角在 1550 nm 附近拉高一个量级代价是吸收边红移、损耗上升。一维磁光晶体的透射和法拉第旋转计算本质上是解一层或多层各向异性介质中的传播问题核心变量是介电张量的非对角元。适合做磁光隔离器、环形器、磁场传感的从业者也适合想从材料参数直接算器件性能、不想只靠仿真软件黑匣子出结果的人。下面从物理量定义一路写到可复现的计算流程和参数边界。2. 从介电张量到透射矩阵一维磁光传播的物理骨架2.1 磁光效应的微观来源与介电张量形式YIG 的磁光效应来自磁化后电子自旋轨道耦合导致的介电张量非对角元。在饱和磁化、磁化方向沿 z 轴时介电张量写成ε ε0 * [ ε1 -i*ε2 0 i*ε2 ε1 0 0 0 ε3 ]其中 ε1 是普通介电常数ε2 是磁光耦合项正比于磁化强度 Mε3 在立方晶系中通常等于 ε1。法拉第旋转角 θ_F 在弱吸收近似下正比于 ε2而 ε2 又正比于 M所以外加磁场通过改变 M 来调旋转角。Ce 掺杂的作用是增强自旋轨道耦合把 ε2 在 1550 nm 附近抬高代价是 ε1 的虚部吸收也增大。这里有一个容易翻车的地方很多教材直接给 θ_F V·B·LV 是费尔德常数但这个公式只在远离吸收带、且磁化饱和时成立。Ce_YIG 在 1550 nm 附近吸收不可忽略必须回到介电张量解 Maxwell 方程否则算出来的旋转角会偏大 20% 以上。2.2 一维传播的 4×4 传输矩阵怎么建一维意味着只考虑沿 z 方向传播、层状结构。对每一层把电场写成四个分量Ex、Ey 及其对应的磁场分量。代入 Maxwell 方程后得到本征值问题解出四个本征模式两个前向、两个后向每个模式有各自的传播常数和偏振态。层内传播用对角矩阵表示层间界面用边界条件匹配。具体步骤对每层材料由介电张量构造 4×4 矩阵 Δ解其特征值得到传播常数 k1~k4 和特征向量。构造层内传播矩阵 P diag(exp(-ik1d), exp(-ik2d), exp(-ik3d), exp(-ik4d))d 是层厚。构造界面矩阵 D把电场和磁场切向分量从一层映射到下一层。总传输矩阵 M D1^-1 * D2 * P2 * D2^-1 * D3 * ... 按层序乘起来。由 M 和入射/出射半空间的边界条件解出反射系数 r 和透射系数 t。透射率 T |t|^2 * (n_out/n_in)法拉第旋转角 θ_F 0.5 * atan(2Re(t_xyt_yy* t_xxt_yx) / (|t_xx|^2 |t_yy|^2 - |t_xy|^2 - |t_yx|^2))椭圆率由虚部给出。注意这里 t 是 2×2 琼斯矩阵不是标量。只取 t_xx 算透射率、忽略交叉项是新手最常见的错误会导致旋转角算出来恒为零。2.3 材料参数从哪来Ce_YIG 的 ε1 和 ε2 取值Ce_YIG 没有像 Si 或 SiO2 那样通用的数据库。常见做法是ε1 的实部由折射率 n ≈ 2.2~2.41550 nm反推虚部由吸收系数 α 换算ε2 由法拉第旋转角实验值反推。如果手头没有自己的椭偏或磁光克尔谱数据可以用文献里 Ce_YIG 的典型值起步ε1 ≈ 5.3 i*0.02ε2 ≈ 0.01~0.03具体随 Ce 浓度和退火条件变化很大。我一般会先固定 ε1 实部扫 ε2 从 0.005 到 0.05看透射率和旋转角怎么变再拿实验值去卡。这样比一上来就拟合所有参数快得多也能看出哪个参数是主导。3. 用 Python 跑通 Ce_YIG 透射谱与法拉第旋转角的最小实现3.1 环境准备与依赖只需要 numpy 和 matplotlib不需要 COMSOL 或 Lumerical。Python 3.9 以上即可。pip install numpy matplotlib不依赖任何磁光专用库因为 4×4 传输矩阵本身只有几十行自己写反而可控。用现成软件包的问题是参数含义不透明出了异常值不知道是物理还是设置问题。3.2 构造介电张量与 4×4 传输矩阵import numpy as np def eps_tensor(eps1, eps2): 构造磁化沿z轴的介电张量 eps np.array([ [eps1, -1j*eps2, 0], [1j*eps2, eps1, 0], [0, 0, eps1] ], dtypecomplex) return eps def layer_matrices(eps, k0, d): 返回单层的传播矩阵P和界面矩阵D # 构造Delta矩阵解本征值 # 这里用简化形式对正入射、磁化沿z本征模式为左右圆偏振 n_r np.sqrt(eps[0,0] eps[0,1]*1j) # 右旋 n_l np.sqrt(eps[0,0] - eps[0,1]*1j) # 左旋 k_r k0 * n_r k_l k0 * n_l P np.diag([np.exp(-1j*k_r*d), np.exp(-1j*k_l*d), np.exp(1j*k_r*d), np.exp(1j*k_l*d)]) # 界面矩阵D把圆偏振基映射到线偏振基 D np.array([ [1, 1, 0, 0], [1j, -1j, 0, 0], [0, 0, 1, 1], [0, 0, 1j, -1j] ], dtypecomplex) / np.sqrt(2) return P, D逻辑说明对正入射、磁化沿传播方向的情况本征模式是左右圆偏振折射率分别为 sqrt(ε1 ± ε2)。这个简化在偏离正入射超过 10 度后误差迅速增大那时必须回到完整 4×4 本征值求解。参数 d 是层厚单位与 k0 一致k0 2π/λ。3.3 计算透射率和法拉第旋转角的主循环def calculate(lam, eps1, eps2, d): k0 2*np.pi / lam eps eps_tensor(eps1, eps2) P, D layer_matrices(eps, k0, d) # 单层M D^-1 * P * D M np.linalg.inv(D) P D # 半空间边界入射n01出射n01 # 提取琼斯矩阵t简化取M的左上2x2块 t M[:2, :2] T np.abs(t[0,0])**2 np.abs(t[1,0])**2 # 法拉第旋转角 theta 0.5 * np.angle(t[0,0] 1j*t[1,0]) - 0.5 * np.angle(t[0,0] - 1j*t[1,0]) return T, np.degrees(theta) # 扫波长 lams np.linspace(1500, 1600, 200) T_list, th_list [], [] for lam in lams: T, th calculate(lam*1e-9, 5.30.02j, 0.02, 500e-6) T_list.append(T) th_list.append(th)逻辑说明t 矩阵的左上 2×2 块对应线偏振基下的琼斯矩阵。透射率取第一列x 偏振入射的模方和。旋转角用左右旋相位差的一半来算这是正入射下的标准做法。参数 500e-6 是 500 微米厚Ce_YIG 典型器件厚度在 200~1000 微米之间。3.4 结果解读与参数敏感性跑出来的透射谱在 1550 nm 附近如果出现明显干涉条纹说明厚度和折射率匹配条纹间距 Δλ ≈ λ²/(2nd)。旋转角随 ε2 线性增大但透射率随 ε2 增大而下降因为吸收项被放大。实际设计要在旋转角和插入损耗之间取折中Ce_YIG 的 ε2 通常选在 0.015~0.025 之间。提示如果旋转角算出来是负的检查 ε2 的符号和磁场方向定义是否一致。符号约定不统一是磁光计算里最常见的玄学问题。4. 透射反射联算时最容易翻车的五个地方4.1 现象透射率大于 1原因出射半空间折射率没有归一化或者 t 矩阵取了错误的块。解决确认 T |t|² * Re(n_out)/Re(n_in)正入射下 n_in n_out 1 时退化为 |t|²。4.2 现象旋转角随厚度线性增大但实验不线性原因厚度超过吸收长度后多次反射和吸收导致有效旋转饱和。解决在传输矩阵里保留后向模式不要只用前向近似。吸收长度 1/α 在 Ce_YIG 里约 1~5 mm厚度接近这个量级时必须算全矩阵。4.3 现象反射率算出来和实验差一个数量级原因界面矩阵 D 用了错误的本征模式排序导致前后向模式混淆。解决检查 D 的列顺序是否和 P 的对角元顺序一致前向两个、后向两个不能交叉。4.4 现象改变磁场方向后旋转角不变原因介电张量里 ε2 的符号没有随磁化方向翻转。解决磁化反向时 ε2 → -ε2法拉第旋转角随之反号这是法拉第效应的非互易性来源代码里要显式处理。4.5 现象波长扫描出现非物理尖峰原因k0*d 在某些波长下使矩阵接近奇异数值求逆不稳定。解决改用解线性方程组而不是显式求逆或者把波长步长减小到 0.1 nm 以下。5. 进阶用群论约化参数并做实验对标5.1 用对称性减少独立参数Ce_YIG 是立方晶系磁化沿 [111] 和沿 [100] 时介电张量的非对角元结构不同。沿 [111] 磁化时ε2 的有效值要乘一个方向因子。如果做的是磁场角度依赖实验这一步不能省否则拟合出来的 ε2 会随角度漂移看起来像材料不均匀其实是坐标没转对。def rotate_eps(eps, theta, phi): 把介电张量从磁化坐标系转到实验室坐标系 # theta, phi 是磁化方向球坐标 # 构造旋转矩阵R返回 R eps R.T ...参数 theta 是磁化与 z 轴夹角phi 是方位角。对 [111] 方向theta ≈ 54.7 度phi 45 度。旋转后非对角元不再只有 ε2 一个独立量会出现 ε_xy 和 ε_xz 同时非零。5.2 实验对标椭偏仪和磁光克尔谱怎么对椭偏仪给的是 ψ 和 Δ对应反射系数比 r_p/r_s。磁光克尔谱给的是克尔旋转角和椭圆率。把计算出的反射矩阵转成 ψ、Δ 和克尔角和实验曲线叠在一起看。如果透射谱对得上但克尔谱对不上问题多半在界面层——Ce_YIG 表面容易形成非磁性的死层厚度几纳米到几十纳米对透射影响小但对反射影响大。我一般会在模型里加一层 5~20 nm 的界面层ε2 0ε1 取体材料值然后看克尔谱能不能压下去。这个死层参数没有通用值必须用自己的样品去卡。5.3 一个具体技巧用透射极小值定位磁光共振Ce_YIG 在近红外有一个磁光共振增强区表现为 ε2 的色散峰。如果只扫透射谱这个峰被吸收背景淹没看不出来。技巧是同时算透射率和旋转角取旋转角/吸收系数的比值这个比值在共振波长附近会出现极大值。用这个比值定位共振比直接看旋转角曲线准得多因为旋转角本身也受厚度干涉调制。具体做法对每个波长算 θ_F 和 α_eff -ln(T)/d然后画 θ_F/α_eff 随波长的曲线。峰位就是共振中心。这个技巧在 Ce 浓度较低、共振不明显时尤其有用。注意α_eff 在干涉条纹存在时会有振荡画图前先对 T 做平滑或者取包络否则比值曲线全是毛刺。我自己在这个方向上踩过最深的坑是早期直接用费尔德常数公式算 Ce_YIG 的旋转角结果和实验差了近一倍后来老老实实回到 4×4 矩阵把吸收和多次反射都算进去才对上。希望帮到你。本文还有配套的精品资源点击获取
