尧图网站建设 
简介一个基于980nm泵浦的掺铒光纤放大器(EDFA)MATLAB仿真项目面向光通信专业学生、研究人员以及对光放大技术感兴趣的开发者。项目以单模光纤为传输介质模拟多波长信号在980纳米泵浦源激励下的掺铒光纤放大过程980nm波长能提供良好的泵浦效率结合Runge-Kutta算法求解光传播方程可还原铒离子吸收、受激辐射、增益累积及波长转换等物理机制。资源包共10个M文件压缩包仅7KB代码轻量涵盖主程序与若干子功能模块既便于快速定位核心算法也支持修改信号波长、泵浦功率和光纤参数来扩展仿真场景。目前已有780人学习下载可作为EDFA课程设计、毕业设计或科研预研的基础框架帮助读者理解光放大器建模、多波长复用以及光纤传输损耗补偿的实际实现。1. 从“EDFA-980nm.zip”开始MATLAB 仿真掺铒光纤放大器的三件套搜到EDFA-980nm.zip这类资源包的人大多不是缺代码而是缺一条把“980nm 泵浦、单模光纤、掺铒光纤”串起来的思路。EDFA 仿真做的是同一件事给定泵浦功率、信号功率和掺铒光纤参数算出光纤长度方向上信号怎么被放大、泵浦怎么被消耗、增益最后落在多少 dB。要动手复现至少需要三部分铒离子在 980nm 泵浦下的能级跃迁模型、Giles 或速率方程描述的功率传播方程以及 MATLAB 里稳定的数值积分方法。这篇博文不依赖任何未公开的代码包直接把这三部分拆开讲清楚并给出一段能在 MATLAB 里跑通的最小脚本。从事光通信链路预算、EDFA 课题仿真或想验证掺铒光纤参数的人按这里的过程走一遍就能自己生成增益曲线。2. 980nm 泵浦掺铒光纤的物理基础与单模光纤重叠积分2.1 掺铒离子的光跃迁过程EDFA 的放大基础是铒离子 Er³⁺ 在 1550nm 附近的受激辐射。铒离子基态为 4I₁₅/₂980nm 泵浦光将电子激发到 4I₁₁/₂ 能级。这个能级寿命只有微秒量级绝大多数电子会快速无辐射跃迁到亚稳态 4I₁₃/₂。这个亚稳态寿命约 10ms远大于泵浦态寿命因此很适合积累粒子数反转。当 1550nm 信号光子进入光纤时处于亚稳态的电子受激跃迁回基态产生与入射光子同相位、同方向的新光子信号因此得到放大。仿真中真正需要的是“宏观速率”而不是单个离子行为。你不需要追踪每一个 Er³⁺ 离子只需要知道某一小段光纤内上能级粒子数占全部掺杂离子数的比例。MATLAB 模拟通常把这个比例记作 n₂取值范围 0 到 1。n₂ 越大受激辐射越强。但 n₂ 不是常数它会随泵浦功率、信号功率和光纤位置变化因此要在积分过程中每步重新计算。2.2 980nm 与 1480nm 泵浦的选择权衡掺铒光纤泵浦还有常见的 1480nm 波长。工程上为什么很多时候优先选 980nm核心原因是粒子数反转更彻底。980nm 把离子直接泵到更高的能级再掉回亚稳态几乎不直接激励亚稳态因此可以达到很高的反转率。1480nm 泵浦是直接把离子从基态泵到亚稳态量子效率高但反向受激吸收也更明显噪声系数通常要比 980nm 高 1 到 2dB。对比项980nm 泵浦1480nm 泵浦泵浦能级路径基态到 4I₁₁/₂再无辐射跃迁到亚稳态直接泵入 4I₁₃/₂ 亚稳态量子效率约 0.63约 0.94粒子数反转能力高接近完全反转中等噪声系数可做到 4dB 以下通常高 1 到 2dB残余泵浦少光纤长度利用充分残余较多适合高功率光纤激光典型应用前置放大器、低噪声放大功率放大级在做 MATLAB 仿真参数预设时980nm 泵浦对应的受激发射系数 g_p 可以近似取 0因为受激发射在泵浦波长上几乎不发生。1480nm 泵浦则不能这样简单处理仿真时要额外引入泵浦波长的发射截面。大部分教学仿真和低噪声链路设计默认都用 980nm因此本博客后面的代码也按 980nm 展开。2.3 单模光纤中的模场重叠与高斯近似计算EDFA 的增益不仅取决于掺杂浓度还取决于泵浦光和信号光与掺杂区域的匹配程度。单模光纤中光场并不是均匀分布在光纤横截面上而是呈近似高斯分布。980nm 波长较短模场直径比 1550nm 信号更小但又因为模式截止条件不同两者与纤芯掺杂区的交叠程度各不相同。这个交叠程度用重叠积分表示。仿真时如果不做二维光场计算可以用 Marcuse 经验公式估算模场半径再算重叠积分。下面的 MATLAB 代码就是一个最小实现lambda_p 0.98e-6; lambda_s 1.55e-6; % 单位 m NA 0.12; % 单模光纤数值孔径 r_core 4.0e-6; % 纤芯半径 V (lambda) 2*pi*r_core/lambda * NA; w (lambda) r_core*(0.65 1.619*V(lambda).^(-1.5) 2.879*V(lambda).^(-6)); gamma_s 1 - exp(-2*r_core^2/w(lambda_s)^2); gamma_p 1 - exp(-2*r_core^2/w(lambda_p)^2); disp([gamma_p, gamma_s]);这段代码先由 NA 和纤芯半径算出归一化频率 V(lambda)再用 Marcuse 拟合公式得到模场半径 w。gamma_p 和 gamma_s 分别是泵浦和信号与掺杂半径相同的纤芯区域的重叠比例。你会在路测中看到 gamma_p 和 gamma_s 不一致这直接影响了速率方程里的吸收系数和增益系数。如果只告诉你“掺杂浓度高增益就大”却无视重叠积分仿真结果会和实际相差不少。3. 传输方程与 Giles 模型把粒子数反转转成可积微分方程3.1 从速率方程到功率传播方程对一段长度为 dz 的光纤信号功率变化由受激辐射减去受激吸收得到。若把铒离子看作两能级系统并用 n₂ 表示上能级粒子数占比信号功率传播方程为dP_s/dz P_s × [(α_s g_s) × n₂ − α_s]其中 α_s 是信号波长上的小信号吸收系数g_s 是完全反转时的增益系数。这个公式里的 α_s 和 g_s 都是单位长度系数单位是 1/m。厂商给的参数往往以 dB/m 为单位转换成线性系数时除以 4.343因为 1dB/m 对应的线性系数约等于 1/4.343。泵浦功率的传播则有差别。980nm 泵浦不能像信号那样用“发射−吸收”对称处理因为泵浦波长上没有明显的受激发射g_p 近似为 0。泵浦传播方程为dP_p/dz −P_p × α_p × (1 − n₂)这里的物理意义很直接泵浦只被基态粒子吸收。n₂ 越大基态粒子越少泵浦吸收越弱残余泵浦越多。把这个方程和信号方程一起积分就能得到整段掺铒光纤上的功率分布。表 3-1 列出了后续仿真用到的几个关键参数。参数含义典型参考值α_s1550nm 信号小信号吸收系数4 dB/mg_s1550nm 信号完全反转时的单位增益8 dB/mα_p980nm 泵浦吸收系数6 dB/mτ亚稳态寿命10msN_t铒离子掺杂浓度1×10²⁵ /m³该表中单位长度系数只在速率方程和传播方程中使用实际积分前要除以 4.343。如果你遇到计算结果动不动就增益几十 dB 的“数字幻觉”多半是把这个换算漏了。3.2 数值求解中的本地反转因子实现在任意位置 zn₂ 由当前功率向量决定。忽略自发辐射后的准静态速率方程可以写成n₂ (W_abs) / (1/τ W_abs W_em)其中 W_abs 是受激吸收速率W_em 是受激发射速率。计算时把泵浦和信号两个波长都贡献进去。MATLAB 中可以用一个独立函数实现function n2 local_inversion(P, lambda, alpha, g_star, p) h 6.626e-34; c 2.9979e8; photon P ./ (h * c ./ lambda); rate_abs sum(alpha .* photon ./ (p.Acore * p.Nt)); rate_em sum(g_star .* photon ./ (p.Acore * p.Nt)); n2 rate_abs / (1 / p.tau rate_abs rate_em); end输入信号中 P 是泵浦和信号的功率向量lambda 是各自波长alpha 和 g_star 是对应的吸收、增益系数向量。计算过程先把功率转换成光子数率再除以纤芯有效面积和掺杂浓度得到每秒吸收或发射的平均速率。如果 rate_abs 远大于 rate_em说明该处泵浦充分n₂ 接近最大值反之信号被吸收。这个函数没有包含 ASE 的放大自发辐射项。忽略 ASE 会让增益偏高但因为 MATLAB 的 ode45 可以很快跑通第一版通常这样处理。后续要算噪声系数时再加入 ASE 波长组。3.3 边界条件与双向泵浦的处理实际 EDFA 模块中泵浦可以从信号同方向注入也可以从输出端反向注入或者两端同时注入。同向泵浦和信号都是初值直接在 z0 处给定功率对 z 积分即可。反向泵浦则不同泵浦初值在 zL信号初值在 z0形成了两点边值问题。常见做法是先假设一组反向泵浦在 z0 处的值正向积分到 zL将计算得到的泵浦功率和已知边界条件比较再用二分法或牛顿法修正初值。这个过程叫打靶法。如果仿真同时包含前向和后向 ASE 通道矩阵维度会从 2 变成更多ODE 也会明显变刚此时要改用 ode15s 而不是 ode45。4. 用 MATLAB 编写 EDFA 仿真代码并绘制增益曲线4.1 参数定义与初始化现在把上面的模型落实到完整脚本。下面这段代码保存为edfa_forward_pump_980.m独立可运行。参数用 struct 保存方便扫参时直接改字段。%% edfa_forward_pump_980.m % 980nm 同向泵浦 EDFA 仿真第一版不含 ASE h 6.626e-34; c 2.9979e8; p.L 8; % 掺铒光纤长度 m p.Acore pi * (4.1e-6)^2; % 掺杂区域有效面积 m^2 p.Nt 8e24; % 铒离子密度 m^-3 p.tau 10e-3; % 荧光寿命 s p.alpha_s 4.0 / 4.343; % 1550nm 吸收系数 1/m p.g_s 8.0 / 4.343; % 1550nm 增益系数 1/m p.alpha_p 6.0 / 4.343; % 980nm 吸收系数 1/m p.g_p 0; % 980nm 受激发射近似为 0 Pp0 200e-3; % 泵浦功率 200mW Ps0 1e-3; % 信号输入功率 1mW参数表里需要特别注意的是 A_core它并不是普通单模光纤的纤芯面积而是掺杂离子实际分布范围内的等效面积。如果使用多模光纤或双包层结构这个值还要结合模场重叠修改。Pp0 设为 200mW 是为了让泵浦在 8m 光纤后端仍有剩余避免信号在后段被重吸收。4.2 ODE 函数实现与 ode45 求解ODE 函数负责计算每一点的 dP/dz。它内部复用 3.2 节的 local_inversion 逻辑但不单独调用而是写到嵌套函数里减少参数传递。odefun (z, P) edfa_ode(z, P, p, h, c); [z, P] ode45(odefun, [0 p.L], [Pp0 Ps0]); function dP edfa_ode(z, P, p, h, c) Pp P(1); Ps P(2); lambda [980e-9, 1550e-9]; Pvec [Pp, Ps]; alpha [p.alpha_p, p.alpha_s]; g_star [p.g_p, p.g_s]; photon Pvec ./ (h * c ./ lambda); rate_abs sum(alpha .* photon ./ (p.Acore * p.Nt)); rate_em sum(g_star .* photon ./ (p.Acore * p.Nt)); n2 rate_abs / (1 / p.tau rate_abs rate_em); dP zeros(2, 1); % 泵浦只被基态离子吸收 dP(1) -Pp * p.alpha_p * (1 - n2); % 信号受激辐射减受激吸收 dP(2) Ps * ((p.alpha_s p.g_s) * n2 - p.alpha_s); end代码里 n₂ 的位置解释了“增益饱和”的来源信号功率增大后rate_em 变大n₂ 下降导致继续增长的信号增益变小。ode45 在每个积分步中自动调整步长你可以看到泵浦功率在光纤前端消耗最快越到后段变化越平缓。4.3 增益曲线与泵浦消耗的图形输出求解完成后输出信号和泵浦沿光纤的功率变化并计算总增益Pp P(:, 1); Ps P(:, 2); figure; yyaxis left; plot(z, Ps * 1e3, b-, LineWidth, 1.5); ylabel(信号功率 mW); yyaxis right; plot(z, Pp * 1e3, r--, LineWidth, 1.5); ylabel(泵浦功率 mW); xlabel(光纤长度 m); grid on; legend(信号, 泵浦, location, best); title(980nm 同向泵浦 EDFA 功率分布); GdB 10 * log10(Ps(end) / Ps0); fprintf(输入信号: %.3f mW, 输出信号: %.3f mW\n, Ps0 * 1e3, Ps(end) * 1e3); fprintf(增益: %.2f dB\n, GdB);如果仿真结果出现信号先下降再上升说明泵浦功率太小光纤前段信号被掺铒光纤吸收的速度超过了泵浦的放大速度。解决方法是提高 Pp0 或缩短光纤长度。反之如果泵浦在光纤前 2m 内就被吸干后段信号急剧下降说明掺杂浓度过高或 α_p 设置偏大。5. 参数扫参与仿真排错让 EDFA 仿真结果可信5.1 泵浦功率-光纤长度二维扫参单点仿真只能说明“某个参数组合下可跑通”不能指导设计。更常见的做法是固定输入信号扫泵浦功率和光纤长度画出增益等高线。一段可复用的扫参代码是这样写的Pp_list 50:25:300; % 泵浦功率范围 mW L_list 1:0.5:15; % 光纤长度范围 m G_matrix zeros(length(Pp_list), length(L_list)); for i 1:length(Pp_list) Pp0 Pp_list(i) * 1e-3; for j 1:length(L_list) p.L L_list(j); [~, Pout] ode45(odefun, [0 p.L], [Pp0 Ps0]); G_matrix(i, j) 10 * log10(Pout(end, 2) / Ps0); end end扫参结果通常呈现一个“山脊”泵浦不足时增益低泵浦过高后增益不再增加反而因为非线性效应导致失真。对 980nm 泵浦 EDFA 来说最优长度往往落在泵浦刚好耗尽、信号达到峰值增益的地方。把 G_matrix 用 contour 画出能直观看出可用的增益区间在哪里。5.2 数值异常与 ASE 扩展方向ODE 求解遇到 NaN 或负功率不外乎三个原因。第一参数单位不一致最常见的是把 dB/m 直接代入方程计算导致 α_s 和 g_s 被放大好几倍。第二rate_abs 出现零除比如 Nt 设为 0 或 A_core 设为空值。第三泵浦功率在积分过程中被衰减到接近 0而信号仍很大导致局部速率极不均衡ode45 无法收敛。建议在代码开头加一段参数检查确认 alpha_p、alpha_s、g_s 都大于 0Pp0 和 Ps0 大于 0光纤长度 L 为正数。剩下的差异来自 ASE。要计算噪声系数必须把 1530nm、1550nm、1570nm 等几个 ASE 波长放入功率向量每个波长再区分前向和后向功率向量维度从 2 扩展到 22×N_ase。此时方程刚性显著增加把 ode45 换成 ode15s并在自发射项中引入 2hνΔυ 的等效输入功率。最后验证仿真是否可靠可以看两个边界条件泵浦功率极低或光纤长度为 0 时增益应为负值输出信号等于输入信号经过吸收衰减泵浦功率充分高时增益应趋于由 g_s 和 α_s 决定的上限。这两点都满足再往模型里加 ASE 才不会带偏结果。本文还有配套的精品资源点击获取 基于MATLAB的980nm泵浦EDFA掺铒光纤放大器仿真与增益分析
简介一个基于980nm泵浦的掺铒光纤放大器(EDFA)MATLAB仿真项目面向光通信专业学生、研究人员以及对光放大技术感兴趣的开发者。项目以单模光纤为传输介质模拟多波长信号在980纳米泵浦源激励下的掺铒光纤放大过程980nm波长能提供良好的泵浦效率结合Runge-Kutta算法求解光传播方程可还原铒离子吸收、受激辐射、增益累积及波长转换等物理机制。资源包共10个M文件压缩包仅7KB代码轻量涵盖主程序与若干子功能模块既便于快速定位核心算法也支持修改信号波长、泵浦功率和光纤参数来扩展仿真场景。目前已有780人学习下载可作为EDFA课程设计、毕业设计或科研预研的基础框架帮助读者理解光放大器建模、多波长复用以及光纤传输损耗补偿的实际实现。1. 从“EDFA-980nm.zip”开始MATLAB 仿真掺铒光纤放大器的三件套搜到EDFA-980nm.zip这类资源包的人大多不是缺代码而是缺一条把“980nm 泵浦、单模光纤、掺铒光纤”串起来的思路。EDFA 仿真做的是同一件事给定泵浦功率、信号功率和掺铒光纤参数算出光纤长度方向上信号怎么被放大、泵浦怎么被消耗、增益最后落在多少 dB。要动手复现至少需要三部分铒离子在 980nm 泵浦下的能级跃迁模型、Giles 或速率方程描述的功率传播方程以及 MATLAB 里稳定的数值积分方法。这篇博文不依赖任何未公开的代码包直接把这三部分拆开讲清楚并给出一段能在 MATLAB 里跑通的最小脚本。从事光通信链路预算、EDFA 课题仿真或想验证掺铒光纤参数的人按这里的过程走一遍就能自己生成增益曲线。2. 980nm 泵浦掺铒光纤的物理基础与单模光纤重叠积分2.1 掺铒离子的光跃迁过程EDFA 的放大基础是铒离子 Er³⁺ 在 1550nm 附近的受激辐射。铒离子基态为 4I₁₅/₂980nm 泵浦光将电子激发到 4I₁₁/₂ 能级。这个能级寿命只有微秒量级绝大多数电子会快速无辐射跃迁到亚稳态 4I₁₃/₂。这个亚稳态寿命约 10ms远大于泵浦态寿命因此很适合积累粒子数反转。当 1550nm 信号光子进入光纤时处于亚稳态的电子受激跃迁回基态产生与入射光子同相位、同方向的新光子信号因此得到放大。仿真中真正需要的是“宏观速率”而不是单个离子行为。你不需要追踪每一个 Er³⁺ 离子只需要知道某一小段光纤内上能级粒子数占全部掺杂离子数的比例。MATLAB 模拟通常把这个比例记作 n₂取值范围 0 到 1。n₂ 越大受激辐射越强。但 n₂ 不是常数它会随泵浦功率、信号功率和光纤位置变化因此要在积分过程中每步重新计算。2.2 980nm 与 1480nm 泵浦的选择权衡掺铒光纤泵浦还有常见的 1480nm 波长。工程上为什么很多时候优先选 980nm核心原因是粒子数反转更彻底。980nm 把离子直接泵到更高的能级再掉回亚稳态几乎不直接激励亚稳态因此可以达到很高的反转率。1480nm 泵浦是直接把离子从基态泵到亚稳态量子效率高但反向受激吸收也更明显噪声系数通常要比 980nm 高 1 到 2dB。对比项980nm 泵浦1480nm 泵浦泵浦能级路径基态到 4I₁₁/₂再无辐射跃迁到亚稳态直接泵入 4I₁₃/₂ 亚稳态量子效率约 0.63约 0.94粒子数反转能力高接近完全反转中等噪声系数可做到 4dB 以下通常高 1 到 2dB残余泵浦少光纤长度利用充分残余较多适合高功率光纤激光典型应用前置放大器、低噪声放大功率放大级在做 MATLAB 仿真参数预设时980nm 泵浦对应的受激发射系数 g_p 可以近似取 0因为受激发射在泵浦波长上几乎不发生。1480nm 泵浦则不能这样简单处理仿真时要额外引入泵浦波长的发射截面。大部分教学仿真和低噪声链路设计默认都用 980nm因此本博客后面的代码也按 980nm 展开。2.3 单模光纤中的模场重叠与高斯近似计算EDFA 的增益不仅取决于掺杂浓度还取决于泵浦光和信号光与掺杂区域的匹配程度。单模光纤中光场并不是均匀分布在光纤横截面上而是呈近似高斯分布。980nm 波长较短模场直径比 1550nm 信号更小但又因为模式截止条件不同两者与纤芯掺杂区的交叠程度各不相同。这个交叠程度用重叠积分表示。仿真时如果不做二维光场计算可以用 Marcuse 经验公式估算模场半径再算重叠积分。下面的 MATLAB 代码就是一个最小实现lambda_p 0.98e-6; lambda_s 1.55e-6; % 单位 m NA 0.12; % 单模光纤数值孔径 r_core 4.0e-6; % 纤芯半径 V (lambda) 2*pi*r_core/lambda * NA; w (lambda) r_core*(0.65 1.619*V(lambda).^(-1.5) 2.879*V(lambda).^(-6)); gamma_s 1 - exp(-2*r_core^2/w(lambda_s)^2); gamma_p 1 - exp(-2*r_core^2/w(lambda_p)^2); disp([gamma_p, gamma_s]);这段代码先由 NA 和纤芯半径算出归一化频率 V(lambda)再用 Marcuse 拟合公式得到模场半径 w。gamma_p 和 gamma_s 分别是泵浦和信号与掺杂半径相同的纤芯区域的重叠比例。你会在路测中看到 gamma_p 和 gamma_s 不一致这直接影响了速率方程里的吸收系数和增益系数。如果只告诉你“掺杂浓度高增益就大”却无视重叠积分仿真结果会和实际相差不少。3. 传输方程与 Giles 模型把粒子数反转转成可积微分方程3.1 从速率方程到功率传播方程对一段长度为 dz 的光纤信号功率变化由受激辐射减去受激吸收得到。若把铒离子看作两能级系统并用 n₂ 表示上能级粒子数占比信号功率传播方程为dP_s/dz P_s × [(α_s g_s) × n₂ − α_s]其中 α_s 是信号波长上的小信号吸收系数g_s 是完全反转时的增益系数。这个公式里的 α_s 和 g_s 都是单位长度系数单位是 1/m。厂商给的参数往往以 dB/m 为单位转换成线性系数时除以 4.343因为 1dB/m 对应的线性系数约等于 1/4.343。泵浦功率的传播则有差别。980nm 泵浦不能像信号那样用“发射−吸收”对称处理因为泵浦波长上没有明显的受激发射g_p 近似为 0。泵浦传播方程为dP_p/dz −P_p × α_p × (1 − n₂)这里的物理意义很直接泵浦只被基态粒子吸收。n₂ 越大基态粒子越少泵浦吸收越弱残余泵浦越多。把这个方程和信号方程一起积分就能得到整段掺铒光纤上的功率分布。表 3-1 列出了后续仿真用到的几个关键参数。参数含义典型参考值α_s1550nm 信号小信号吸收系数4 dB/mg_s1550nm 信号完全反转时的单位增益8 dB/mα_p980nm 泵浦吸收系数6 dB/mτ亚稳态寿命10msN_t铒离子掺杂浓度1×10²⁵ /m³该表中单位长度系数只在速率方程和传播方程中使用实际积分前要除以 4.343。如果你遇到计算结果动不动就增益几十 dB 的“数字幻觉”多半是把这个换算漏了。3.2 数值求解中的本地反转因子实现在任意位置 zn₂ 由当前功率向量决定。忽略自发辐射后的准静态速率方程可以写成n₂ (W_abs) / (1/τ W_abs W_em)其中 W_abs 是受激吸收速率W_em 是受激发射速率。计算时把泵浦和信号两个波长都贡献进去。MATLAB 中可以用一个独立函数实现function n2 local_inversion(P, lambda, alpha, g_star, p) h 6.626e-34; c 2.9979e8; photon P ./ (h * c ./ lambda); rate_abs sum(alpha .* photon ./ (p.Acore * p.Nt)); rate_em sum(g_star .* photon ./ (p.Acore * p.Nt)); n2 rate_abs / (1 / p.tau rate_abs rate_em); end输入信号中 P 是泵浦和信号的功率向量lambda 是各自波长alpha 和 g_star 是对应的吸收、增益系数向量。计算过程先把功率转换成光子数率再除以纤芯有效面积和掺杂浓度得到每秒吸收或发射的平均速率。如果 rate_abs 远大于 rate_em说明该处泵浦充分n₂ 接近最大值反之信号被吸收。这个函数没有包含 ASE 的放大自发辐射项。忽略 ASE 会让增益偏高但因为 MATLAB 的 ode45 可以很快跑通第一版通常这样处理。后续要算噪声系数时再加入 ASE 波长组。3.3 边界条件与双向泵浦的处理实际 EDFA 模块中泵浦可以从信号同方向注入也可以从输出端反向注入或者两端同时注入。同向泵浦和信号都是初值直接在 z0 处给定功率对 z 积分即可。反向泵浦则不同泵浦初值在 zL信号初值在 z0形成了两点边值问题。常见做法是先假设一组反向泵浦在 z0 处的值正向积分到 zL将计算得到的泵浦功率和已知边界条件比较再用二分法或牛顿法修正初值。这个过程叫打靶法。如果仿真同时包含前向和后向 ASE 通道矩阵维度会从 2 变成更多ODE 也会明显变刚此时要改用 ode15s 而不是 ode45。4. 用 MATLAB 编写 EDFA 仿真代码并绘制增益曲线4.1 参数定义与初始化现在把上面的模型落实到完整脚本。下面这段代码保存为edfa_forward_pump_980.m独立可运行。参数用 struct 保存方便扫参时直接改字段。%% edfa_forward_pump_980.m % 980nm 同向泵浦 EDFA 仿真第一版不含 ASE h 6.626e-34; c 2.9979e8; p.L 8; % 掺铒光纤长度 m p.Acore pi * (4.1e-6)^2; % 掺杂区域有效面积 m^2 p.Nt 8e24; % 铒离子密度 m^-3 p.tau 10e-3; % 荧光寿命 s p.alpha_s 4.0 / 4.343; % 1550nm 吸收系数 1/m p.g_s 8.0 / 4.343; % 1550nm 增益系数 1/m p.alpha_p 6.0 / 4.343; % 980nm 吸收系数 1/m p.g_p 0; % 980nm 受激发射近似为 0 Pp0 200e-3; % 泵浦功率 200mW Ps0 1e-3; % 信号输入功率 1mW参数表里需要特别注意的是 A_core它并不是普通单模光纤的纤芯面积而是掺杂离子实际分布范围内的等效面积。如果使用多模光纤或双包层结构这个值还要结合模场重叠修改。Pp0 设为 200mW 是为了让泵浦在 8m 光纤后端仍有剩余避免信号在后段被重吸收。4.2 ODE 函数实现与 ode45 求解ODE 函数负责计算每一点的 dP/dz。它内部复用 3.2 节的 local_inversion 逻辑但不单独调用而是写到嵌套函数里减少参数传递。odefun (z, P) edfa_ode(z, P, p, h, c); [z, P] ode45(odefun, [0 p.L], [Pp0 Ps0]); function dP edfa_ode(z, P, p, h, c) Pp P(1); Ps P(2); lambda [980e-9, 1550e-9]; Pvec [Pp, Ps]; alpha [p.alpha_p, p.alpha_s]; g_star [p.g_p, p.g_s]; photon Pvec ./ (h * c ./ lambda); rate_abs sum(alpha .* photon ./ (p.Acore * p.Nt)); rate_em sum(g_star .* photon ./ (p.Acore * p.Nt)); n2 rate_abs / (1 / p.tau rate_abs rate_em); dP zeros(2, 1); % 泵浦只被基态离子吸收 dP(1) -Pp * p.alpha_p * (1 - n2); % 信号受激辐射减受激吸收 dP(2) Ps * ((p.alpha_s p.g_s) * n2 - p.alpha_s); end代码里 n₂ 的位置解释了“增益饱和”的来源信号功率增大后rate_em 变大n₂ 下降导致继续增长的信号增益变小。ode45 在每个积分步中自动调整步长你可以看到泵浦功率在光纤前端消耗最快越到后段变化越平缓。4.3 增益曲线与泵浦消耗的图形输出求解完成后输出信号和泵浦沿光纤的功率变化并计算总增益Pp P(:, 1); Ps P(:, 2); figure; yyaxis left; plot(z, Ps * 1e3, b-, LineWidth, 1.5); ylabel(信号功率 mW); yyaxis right; plot(z, Pp * 1e3, r--, LineWidth, 1.5); ylabel(泵浦功率 mW); xlabel(光纤长度 m); grid on; legend(信号, 泵浦, location, best); title(980nm 同向泵浦 EDFA 功率分布); GdB 10 * log10(Ps(end) / Ps0); fprintf(输入信号: %.3f mW, 输出信号: %.3f mW\n, Ps0 * 1e3, Ps(end) * 1e3); fprintf(增益: %.2f dB\n, GdB);如果仿真结果出现信号先下降再上升说明泵浦功率太小光纤前段信号被掺铒光纤吸收的速度超过了泵浦的放大速度。解决方法是提高 Pp0 或缩短光纤长度。反之如果泵浦在光纤前 2m 内就被吸干后段信号急剧下降说明掺杂浓度过高或 α_p 设置偏大。5. 参数扫参与仿真排错让 EDFA 仿真结果可信5.1 泵浦功率-光纤长度二维扫参单点仿真只能说明“某个参数组合下可跑通”不能指导设计。更常见的做法是固定输入信号扫泵浦功率和光纤长度画出增益等高线。一段可复用的扫参代码是这样写的Pp_list 50:25:300; % 泵浦功率范围 mW L_list 1:0.5:15; % 光纤长度范围 m G_matrix zeros(length(Pp_list), length(L_list)); for i 1:length(Pp_list) Pp0 Pp_list(i) * 1e-3; for j 1:length(L_list) p.L L_list(j); [~, Pout] ode45(odefun, [0 p.L], [Pp0 Ps0]); G_matrix(i, j) 10 * log10(Pout(end, 2) / Ps0); end end扫参结果通常呈现一个“山脊”泵浦不足时增益低泵浦过高后增益不再增加反而因为非线性效应导致失真。对 980nm 泵浦 EDFA 来说最优长度往往落在泵浦刚好耗尽、信号达到峰值增益的地方。把 G_matrix 用 contour 画出能直观看出可用的增益区间在哪里。5.2 数值异常与 ASE 扩展方向ODE 求解遇到 NaN 或负功率不外乎三个原因。第一参数单位不一致最常见的是把 dB/m 直接代入方程计算导致 α_s 和 g_s 被放大好几倍。第二rate_abs 出现零除比如 Nt 设为 0 或 A_core 设为空值。第三泵浦功率在积分过程中被衰减到接近 0而信号仍很大导致局部速率极不均衡ode45 无法收敛。建议在代码开头加一段参数检查确认 alpha_p、alpha_s、g_s 都大于 0Pp0 和 Ps0 大于 0光纤长度 L 为正数。剩下的差异来自 ASE。要计算噪声系数必须把 1530nm、1550nm、1570nm 等几个 ASE 波长放入功率向量每个波长再区分前向和后向功率向量维度从 2 扩展到 22×N_ase。此时方程刚性显著增加把 ode45 换成 ode15s并在自发射项中引入 2hνΔυ 的等效输入功率。最后验证仿真是否可靠可以看两个边界条件泵浦功率极低或光纤长度为 0 时增益应为负值输出信号等于输入信号经过吸收衰减泵浦功率充分高时增益应趋于由 g_s 和 α_s 决定的上限。这两点都满足再往模型里加 ASE 才不会带偏结果。本文还有配套的精品资源点击获取