%% 光纤SBS效应仿真平台 - 完整代码 (含能量守恒分析)
clear; clc;
% ========== 光纤参数 ==========
fiber.length = 20; % 光纤长度(m)
fiber.alpha = 0.2e-3; % 损耗系数(dB/m) [^1]
fiber.Aeff = 80e-12; % 有效模场面积(m²)
fiber.g_B = 5e-11; % 布里渊增益系数(m/W)
fiber.delta_nu_B = 50e6; % 布里渊带宽(Hz)
fiber.n2 = 2.3e-20; % 非线性折射率系数(m²/W)
% ========== 双音系统参数 ==========
f0 = 193.1e12; % 中心频率(Hz)
delta_f = [0.1e9, 1e9]; % 双音频差(Hz)
P_input = [20, 20]; % 输入功率(W) 每个分量
% ========== 物理常数 ==========
c = 3e8; % 光速(m/s)
h = 6.626e-34; % 普朗克常数(J·s)
% ========== 主仿真流程 ==========
% 1. 计算SBS阈值
[P_th, L_eff] = calculate_SBS_threshold(fiber, delta_f);
% 2. 生成双音信号
[t, signal] = generate_two_tone_signal(f0, delta_f, P_input);
% 3. 分布式功率演化
[z, P_pump, P_stokes] = simulate_power_evolution(fiber, P_input, delta_f);
% 4. 能量守恒分析
analyze_energy_conservation(z, P_pump, P_stokes, fiber);
% 5. 可视化结果
plot_results(z, P_pump, P_stokes, P_th, delta_f);
% ========== 函数定义 ==========
function [P_th, L_eff] = calculate_SBS_threshold(fiber, delta_f)
% 等效线宽: Δν_eff = |f1 - f2|
delta_nu_eff = abs(delta_f(1) - delta_f(2));
% 将dB/m转换为线性衰减系数: α_lin = α_dB * ln(10)/10
alpha_lin = fiber.alpha * log(10)/10;
% 有效长度: L_eff = [1 - exp(-αL)]/α
L_eff = (1 - exp(-alpha_lin * fiber.length)) / alpha_lin;
% SBS阈值公式: P_th = 21 * A_eff * Δν_eff / (g_B * L_eff * Δν_B)
P_th = 21 * fiber.Aeff * delta_nu_eff / (fiber.g_B * L_eff * fiber.delta_nu_B);
fprintf('SBS阈值 = %.2f W (双音频差 %.2f GHz)\n', P_th, delta_nu_eff/1e9);
end
function [t, signal] = generate_two_tone_signal(f0, delta_f, P_input)
fs = 2 * max(f0 + delta_f); % 采样频率
t = 0:1/fs:1e-9; % 时间向量(1ns)
% 双音信号: s(t) = √P1·cos(2πf1t) + √P2·cos(2πf2t)
signal = sqrt(P_input(1)) * cos(2*pi*(f0 + delta_f(1)) * t) + ...
sqrt(P_input(2)) * cos(2*pi*(f0 + delta_f(2)) * t);
end
function [z, P_pump, P_stokes] = simulate_power_evolution(fiber, P_input, delta_f)
dz = 1; % 空间步长(m)
z = 0:dz:fiber.length; % 位置向量
P_total = sum(P_input); % 总输入功率
% 线性衰减系数转换
alpha_lin = fiber.alpha * log(10)/10;
% 泵浦光功率分布 (仅考虑损耗)
P_pump = P_total * exp(-alpha_lin * z);
% 反向Stokes光功率计算
P_stokes = zeros(size(z));
P_stokes(end) = 1e-3; % 末端初始噪声(1nW)
for i = length(z)-1:-1:1
% SBS耦合方程: dPs/dz = g_B·Pp·Ps - α·Ps
gain_term = fiber.g_B * P_pump(i) * P_stokes(i+1) * dz;
loss_term = alpha_lin * P_stokes(i+1) * dz;
P_stokes(i) = P_stokes(i+1) + gain_term - loss_term;
end
end
function analyze_energy_conservation(z, P_pump, P_stokes, fiber)
% 计算总输入能量
input_energy = trapz(z, P_pump(1) * exp(-fiber.alpha*z)) + ...
trapz(z, P_stokes(end) * exp(-fiber.alpha*(fiber.length-z)));
% 计算输出能量分布
output_energy = struct();
output_energy.pump_end = P_pump(end); % 泵浦光末端输出
output_energy.stokes_start = P_stokes(1); % Stokes光始端输出
output_energy.loss = input_energy - (output_energy.pump_end + output_energy.stokes_start);
fprintf('\n===== 能量守恒分析 =====\n');
fprintf('总输入能量: %.4f J/m\n', input_energy);
fprintf('泵浦光末端输出: %.4f J/m\n', output_energy.pump_end);
fprintf('Stokes光始端输出: %.4f J/m\n', output_energy.stokes_start);
fprintf('光纤损耗: %.4f J/m\n', output_energy.loss);
fprintf('能量守恒误差: %.2f%%\n', ...
100*abs(input_energy - (output_energy.pump_end + output_energy.stokes_start + output_energy.loss))/input_energy);
end
function plot_results(z, P_pump, P_stokes, P_th, delta_f)
% 功率分布图
figure('Position', [100, 100, 800, 600])
subplot(2,1,1)
plot(z, P_pump, 'b-', 'LineWidth', 2)
hold on
plot(z, P_stokes, 'r--', 'LineWidth', 2)
yline(P_th, 'k:', 'SBS阈值', 'LineWidth', 1.5, 'FontSize', 12)
xlabel('光纤位置 (m)')
ylabel('功率 (W)')
title('泵浦光与Stokes光功率演化')
legend('泵浦光', '反向Stokes光', 'Location', 'best')
grid on
% 能量转移分析
subplot(2,1,2)
energy_transfer = P_pump(1) - P_pump - (P_stokes - P_stokes(end));
plot(z, energy_transfer, 'm-', 'LineWidth', 2)
xlabel('光纤位置 (m)')
ylabel('能量转移量 (W)')
title(sprintf('SBS能量转移 (双音频差 %.2f GHz)', abs(diff(delta_f))/1e9))
grid on
end
请在这个代码的基础上进行修改改进