350 lines
13 KiB
Matlab
350 lines
13 KiB
Matlab
%% GPU vs CPU Comparison Test for DP_Fiber
|
|
% This script runs the fiber simulation with and without GPU acceleration
|
|
% and compares the numerical results.
|
|
|
|
% clear; clc;
|
|
|
|
%% Setup (same as gpu_processing_dpfiber.m but simplified)
|
|
s.wavelengthplan = calcWavelengthPlan(4, 400e9, 1310);
|
|
link_length = 10;
|
|
s.pmd = 0.1;
|
|
s.gamma = 0.0023;
|
|
|
|
s.M = 4;
|
|
fsym = 112e9;
|
|
fdac = 2*fsym;
|
|
fadc = 120000000000;
|
|
s.random_key = 1;
|
|
|
|
% Laser / Modulator
|
|
vbias_rel = 0.5;
|
|
u_pi = 4.6;
|
|
vbias = -vbias_rel*u_pi;
|
|
laser_linewidth = 0e6;
|
|
|
|
% DB Stuff
|
|
duob_mode = db_mode.no_db;
|
|
|
|
rcalpha = 0.05;
|
|
Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rc","pulselength",16,"alpha",rcalpha);
|
|
|
|
s.chirpalpha = 0;
|
|
s.p_launch = 3;
|
|
s.p = "co";
|
|
|
|
N = numel(s.wavelengthplan);
|
|
|
|
switch s.p
|
|
case "co"
|
|
pol_rot = 100.*ones(1,N);
|
|
d_local = 0;
|
|
case "pair"
|
|
pol_rot = repmat([100,100,0,0],1,N/4);
|
|
d_local = 0;
|
|
case "alt"
|
|
pol_rot = repmat([100,0,100,0],1,N/4);
|
|
d_local = 0;
|
|
case "seg"
|
|
pol_rot = 100.*ones(1,N);
|
|
d_local = 3;
|
|
otherwise
|
|
error('Unknown fwm_mitigation_technique: %s', string(s.p));
|
|
end
|
|
|
|
f_plan = physconst('lightspeed')./(s.wavelengthplan.*1e-9);
|
|
margin = 5e12;
|
|
f_span = (max(f_plan)+margin)-(min(f_plan)-margin);
|
|
f_nyq = f_span/2;
|
|
|
|
kover = 4;
|
|
upsample_required = f_nyq./(fdac*kover/2);
|
|
upsample_pow = 2^nextpow2(upsample_required);
|
|
|
|
s.f_opt = fdac*kover*upsample_pow;
|
|
s.f_opt_nyq = s.f_opt/2;
|
|
|
|
%% TX per channel
|
|
for l = 1:N
|
|
[Digi_sig,Symbols{l},Tx_bits{l}] = PAMsource( ...
|
|
"fsym",fsym,"M",s.M,"order",15,"useprbs",0, ...
|
|
"fs_out",fdac, ...
|
|
"applyclipping",0,"clipfactor",1.5, ...
|
|
"applypulseform",1,"pulseformer",Pform, ...
|
|
"randkey",s.random_key+l, ...
|
|
"mrds_code",0,"mrds_blocklength",512,"duobinary_mode",duob_mode ...
|
|
).process();
|
|
|
|
Lp_awg = Filter('filtdegree',3,"f_cutoff",56e9,"fs",fdac*kover, ...
|
|
"filterType",filtertypes.gaussian,"active",true);
|
|
|
|
El_sig = AWG("fdac",fdac,"f_cutoff",fsym,"lpf_active",1,"kover",kover, ...
|
|
"bit_resolution",6,"upsampling_method","samplehold","precomp_sinc_rolloff",0, ...
|
|
"H_lpf",Lp_awg,"dac_max",0.6,"dac_min",-0.6).process(Digi_sig);
|
|
|
|
clear Digi_sig
|
|
|
|
El_sig = El_sig.normalize("mode","oneone");
|
|
scaling = 0.6*(u_pi/2-abs(vbias-u_pi/2));
|
|
El_sig = El_sig .* scaling;
|
|
|
|
Eml_out = EML("mode",eml_mode.im_cosinus,"power",3,"fsimu",El_sig.fs, ...
|
|
"lambda",s.wavelengthplan(l),"bias",vbias,"u_pi",u_pi, ...
|
|
"linewidth",laser_linewidth,"randomkey",s.random_key+l,"alpha",s.chirpalpha).process(El_sig);
|
|
|
|
clear El_sig
|
|
|
|
signal_cell{l} = Polarization_Controller("mode","rot_power","desired_power",pol_rot(l)).process(Eml_out);
|
|
|
|
clear Eml_out Lp_awg
|
|
end
|
|
|
|
disp('Signal generated for all channels.');
|
|
|
|
%% WDM mux + launch
|
|
Opt_sig_wdm = Optical_Multiplex("fs_in",fdac*kover,"fs_out",upsample_pow*fdac*kover, ...
|
|
"lambda_center",1310,"random_key",0,"filtype",1,"B",120e9).process(signal_cell);
|
|
|
|
Opt_sig_wdm = Amplifier("amp_mode","ideal_no_noise","gain_mode","output_power", ...
|
|
"amplification_db",s.p_launch+10*log10(N)).process(Opt_sig_wdm);
|
|
|
|
% Save input for both runs
|
|
Opt_sig_input = Opt_sig_wdm;
|
|
|
|
segment_length = 1;
|
|
nSegments = link_length/segment_length;
|
|
if abs(nSegments - round(nSegments)) > 1e-12
|
|
error('fiber_length_km=%g must be an integer multiple of segment_length=%g km.', link_length, segment_length);
|
|
end
|
|
nSegments = round(nSegments);
|
|
|
|
zdw = 1310;
|
|
randomize_D = true;
|
|
|
|
if nSegments > 0
|
|
Dvec = getDispersionVector(nSegments, d_local, zdw, randomize_D, s.random_key);
|
|
else
|
|
Dvec = [];
|
|
end
|
|
|
|
%% Run WITHOUT GPU
|
|
fprintf('\n========== Running WITHOUT GPU (CPU) ==========\n');
|
|
Opt_sig_cpu = Opt_sig_input;
|
|
|
|
tic;
|
|
for seg = 1:nSegments
|
|
fprintf('CPU Segment %d/%d \n',seg, nSegments);
|
|
|
|
Opt_sig_cpu = DP_Fiber("L",segment_length,"D",Dvec(seg),"Dpmd",s.pmd,"Ds",0.07, ...
|
|
"beat_len",10,"corr_len",100,"dz",1,"manakov",0, ...
|
|
"gamma",s.gamma,"lambda",zdw,"n_waveplates",10,"SS_dphimax",0.01, ...
|
|
"SS_dzmax",50,"SS_dzmin",10,"X_alpha",0.3,"X_beta",0,"rng",1,"useGPU",false).process(Opt_sig_cpu);
|
|
end
|
|
time_cpu = toc;
|
|
fprintf('CPU Time: %.3f seconds\n', time_cpu);
|
|
|
|
%% Run WITH GPU (Double Precision)
|
|
fprintf('\n========== Running WITH GPU (Double Precision) ==========\n');
|
|
Opt_sig_gpu_double = Opt_sig_input;
|
|
|
|
tic;
|
|
for seg = 1:nSegments
|
|
fprintf('GPU-Double Segment %d/%d \n',seg, nSegments);
|
|
|
|
Opt_sig_gpu_double = DP_Fiber("L",segment_length,"D",Dvec(seg),"Dpmd",s.pmd,"Ds",0.07, ...
|
|
"beat_len",10,"corr_len",100,"dz",1,"manakov",0, ...
|
|
"gamma",s.gamma,"lambda",zdw,"n_waveplates",10,"SS_dphimax",0.01, ...
|
|
"SS_dzmax",50,"SS_dzmin",10,"X_alpha",0.3,"X_beta",0,"rng",1,"useGPU",true,"useSingle",false).process(Opt_sig_gpu_double);
|
|
end
|
|
time_gpu_double = toc;
|
|
fprintf('GPU Double Time: %.3f seconds\n', time_gpu_double);
|
|
|
|
%% Run WITH GPU (Single Precision)
|
|
fprintf('\n========== Running WITH GPU (Single Precision) ==========\n');
|
|
Opt_sig_gpu_single = Opt_sig_input;
|
|
|
|
tic;
|
|
for seg = 1:nSegments
|
|
fprintf('GPU-Single Segment %d/%d \n',seg, nSegments);
|
|
|
|
Opt_sig_gpu_single = DP_Fiber("L",segment_length,"D",Dvec(seg),"Dpmd",s.pmd,"Ds",0.07, ...
|
|
"beat_len",10,"corr_len",100,"dz",1,"manakov",0, ...
|
|
"gamma",s.gamma,"lambda",zdw,"n_waveplates",10,"SS_dphimax",0.01, ...
|
|
"SS_dzmax",50,"SS_dzmin",10,"X_alpha",0.3,"X_beta",0,"rng",1,"useGPU",true,"useSingle",true).process(Opt_sig_gpu_single);
|
|
end
|
|
time_gpu_single = toc;
|
|
fprintf('GPU Single Time: %.3f seconds\n', time_gpu_single);
|
|
|
|
%% Compare Results
|
|
fprintf('\n========== Numerical Comparison ==========\n');
|
|
|
|
sig_cpu = Opt_sig_cpu.signal;
|
|
sig_gpu_double = Opt_sig_gpu_double.signal;
|
|
sig_gpu_single = Opt_sig_gpu_single.signal;
|
|
|
|
% Check dimensions
|
|
fprintf('CPU signal size: [%d x %d]\n', size(sig_cpu,1), size(sig_cpu,2));
|
|
fprintf('GPU Double signal size: [%d x %d]\n', size(sig_gpu_double,1), size(sig_gpu_double,2));
|
|
fprintf('GPU Single signal size: [%d x %d]\n', size(sig_gpu_single,1), size(sig_gpu_single,2));
|
|
|
|
% CPU vs GPU Double
|
|
fprintf('\n--- CPU vs GPU Double ---\n');
|
|
diff_cpu_double = abs(sig_cpu - sig_gpu_double);
|
|
max_diff_cpu_double = max(diff_cpu_double(:));
|
|
mean_diff_cpu_double = mean(diff_cpu_double(:));
|
|
rel_diff_cpu_double = max_diff_cpu_double / max(abs(sig_cpu(:)));
|
|
fprintf('Max absolute difference: %.6e\n', max_diff_cpu_double);
|
|
fprintf('Mean absolute difference: %.6e\n', mean_diff_cpu_double);
|
|
fprintf('Max relative difference: %.6e\n', rel_diff_cpu_double);
|
|
|
|
% CPU vs GPU Single
|
|
fprintf('\n--- CPU vs GPU Single ---\n');
|
|
diff_cpu_single = abs(sig_cpu - sig_gpu_single);
|
|
max_diff_cpu_single = max(diff_cpu_single(:));
|
|
mean_diff_cpu_single = mean(diff_cpu_single(:));
|
|
rel_diff_cpu_single = max_diff_cpu_single / max(abs(sig_cpu(:)));
|
|
fprintf('Max absolute difference: %.6e\n', max_diff_cpu_single);
|
|
fprintf('Mean absolute difference: %.6e\n', mean_diff_cpu_single);
|
|
fprintf('Max relative difference: %.6e\n', rel_diff_cpu_single);
|
|
|
|
% GPU Double vs GPU Single
|
|
fprintf('\n--- GPU Double vs GPU Single ---\n');
|
|
diff_double_single = abs(sig_gpu_double - sig_gpu_single);
|
|
max_diff_double_single = max(diff_double_single(:));
|
|
mean_diff_double_single = mean(diff_double_single(:));
|
|
rel_diff_double_single = max_diff_double_single / max(abs(sig_gpu_double(:)));
|
|
fprintf('Max absolute difference: %.6e\n', max_diff_double_single);
|
|
fprintf('Mean absolute difference: %.6e\n', mean_diff_double_single);
|
|
fprintf('Max relative difference: %.6e\n', rel_diff_double_single);
|
|
|
|
% Check tolerances
|
|
fprintf('\n--- Tolerance Check ---\n');
|
|
tol_double = 1e-10;
|
|
tol_single = 1e-5; % Single precision has ~7 significant digits
|
|
|
|
if max_diff_cpu_double < tol_double
|
|
fprintf('✓ CPU vs GPU Double: EQUIVALENT (diff < %.0e)\n', tol_double);
|
|
else
|
|
fprintf('✗ CPU vs GPU Double: DIFFER beyond tolerance (%.0e)\n', tol_double);
|
|
end
|
|
|
|
if max_diff_cpu_single < tol_single
|
|
fprintf('✓ CPU vs GPU Single: ACCEPTABLE (diff < %.0e)\n', tol_single);
|
|
else
|
|
fprintf('⚠ CPU vs GPU Single: Precision loss detected (diff = %.2e, tol = %.0e)\n', max_diff_cpu_single, tol_single);
|
|
end
|
|
|
|
% Performance comparison
|
|
fprintf('\n========== Performance Summary ==========\n');
|
|
fprintf('CPU Time: %.3f s\n', time_cpu);
|
|
fprintf('GPU Double Time: %.3f s\n', time_gpu_double);
|
|
fprintf('GPU Single Time: %.3f s\n', time_gpu_single);
|
|
fprintf('\n');
|
|
fprintf('Speedup (GPU Double vs CPU): %.2fx\n', time_cpu/time_gpu_double);
|
|
fprintf('Speedup (GPU Single vs CPU): %.2fx\n', time_cpu/time_gpu_single);
|
|
fprintf('Speedup (GPU Single vs GPU Double): %.2fx\n', time_gpu_double/time_gpu_single);
|
|
|
|
%% ========== BER Comparison ==========
|
|
% Process each fiber output through simplified receiver to check if
|
|
% single-precision affects actual BER performance
|
|
|
|
fprintf('\n========== BER Comparison ==========\n');
|
|
fprintf('Processing signals through receiver chain...\n');
|
|
|
|
% Receiver parameters
|
|
rop = -7; % Received optical power [dBm]
|
|
len_tr = 4096; % Training length
|
|
mu_dc = 0.005;
|
|
mu_ffe = [0.0001 0.0008 0.001];
|
|
mu_dfe = 0.0004;
|
|
|
|
% Helper function to process through receiver and get BER
|
|
function ber = process_receiver(Opt_sig_fib, l, Symbols, Tx_bits, ...
|
|
fdac, kover, upsample_pow, fsym, fadc, rop, s, len_tr, mu_dc, mu_ffe, mu_dfe, duob_mode)
|
|
|
|
% Demux single channel
|
|
Opt_sig_demux = Optical_Demultiplex("attenuation",0,"B",200e9,"filtype",1, ...
|
|
"fs_out",fdac*kover,"fs_in",fdac*kover*upsample_pow,"lambda_center",1310).process(Opt_sig_fib);
|
|
|
|
% ROP amplifier
|
|
Opt_sig_rx = Amplifier("amp_mode","ideal_no_noise","gain_mode","output_power", ...
|
|
"amplification_db",rop).process(Opt_sig_demux{l});
|
|
|
|
% Photodiode
|
|
PD_sig = Photodiode("fsimu",fdac*kover,"dark_current",2e-08,"responsivity",1,"temperature",20, ...
|
|
"nep",1.8e-11,"randomkey",s.random_key+l).process(Opt_sig_rx);
|
|
|
|
% Low-pass filter
|
|
rx_bwl = 100e9;
|
|
PD_sig = Filter('filtdegree',4,"f_cutoff",rx_bwl,"fs",fdac*kover, ...
|
|
"filterType",filtertypes.butterworth,"active",true).process(PD_sig);
|
|
|
|
% Scope
|
|
Lp_scpe = Filter('filtdegree',4,"f_cutoff",80e9,"fs",fadc,"filterType",filtertypes.butterworth,"active",true);
|
|
Scpe_sig = Scope("fsimu",fdac*kover,"fadc",fadc, ...
|
|
"delay",0,"fixed_delay",0,"filtertype",filtertypes.butterworth, ...
|
|
"samplingdelay",0,"rand_samplingdelay",0,"freq_offset",0,"samp_jitter",0, ...
|
|
"adcresolution",8,"quantbuffer",0.1,'block_dc',1,'lpf_active',0,'H_lpf',Lp_scpe).process(PD_sig);
|
|
|
|
% Resample to 2 sps
|
|
Scpe_sig_2sps = Scpe_sig.resample("fs_out",2*fsym);
|
|
|
|
% Time sync
|
|
[~, Scpe_cell, ~, ~] = Scpe_sig_2sps.tsynch("reference", Symbols{l}, "fs_ref", fsym, "debug_plots", 0);
|
|
Rx_sig = Scpe_cell{1};
|
|
Rx_sig = Rx_sig.normalize("mode","rms");
|
|
|
|
% FFE Equalizer
|
|
ffe_order = [50, 0, 0];
|
|
eq_ffe = EQ("Ne",ffe_order,"Nb",[0,0,0],"training_length",len_tr,"training_loops",5,"dd_loops",5, ...
|
|
"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",0);
|
|
|
|
ffe_results = ffe(eq_ffe,s.M,Rx_sig,Symbols{l},Tx_bits{l}, ...
|
|
"precode_mode",duob_mode, ...
|
|
'showAnalysis',0, ...
|
|
"postFFE",[], ...
|
|
"eth_style_symbol_mapping",0);
|
|
|
|
ber = ffe_results.metrics.BER;
|
|
end
|
|
|
|
% Process each mode for channel 1
|
|
l = 4; % Use first channel for comparison
|
|
|
|
fprintf('Processing CPU result...\n');
|
|
ber_cpu = process_receiver(Opt_sig_cpu, l, Symbols, Tx_bits, ...
|
|
fdac, kover, upsample_pow, fsym, fadc, rop, s, len_tr, mu_dc, mu_ffe, mu_dfe, duob_mode);
|
|
|
|
fprintf('Processing GPU Double result...\n');
|
|
ber_gpu_double = process_receiver(Opt_sig_gpu_double, l, Symbols, Tx_bits, ...
|
|
fdac, kover, upsample_pow, fsym, fadc, rop, s, len_tr, mu_dc, mu_ffe, mu_dfe, duob_mode);
|
|
|
|
fprintf('Processing GPU Single result...\n');
|
|
ber_gpu_single = process_receiver(Opt_sig_gpu_single, l, Symbols, Tx_bits, ...
|
|
fdac, kover, upsample_pow, fsym, fadc, rop, s, len_tr, mu_dc, mu_ffe, mu_dfe, duob_mode);
|
|
|
|
% Display BER results
|
|
fprintf('\n========== BER Results (Channel %d, ROP = %d dBm) ==========\n', l, rop);
|
|
fprintf('CPU: BER = %.4e\n', ber_cpu);
|
|
fprintf('GPU Double: BER = %.4e\n', ber_gpu_double);
|
|
fprintf('GPU Single: BER = %.4e\n', ber_gpu_single);
|
|
|
|
fprintf('\n--- BER Comparison ---\n');
|
|
if ber_cpu == 0 && ber_gpu_double == 0 && ber_gpu_single == 0
|
|
fprintf('✓ All BERs are zero (no errors detected)\n');
|
|
else
|
|
ber_diff_double = abs(ber_cpu - ber_gpu_double);
|
|
ber_diff_single = abs(ber_cpu - ber_gpu_single);
|
|
fprintf('|BER_cpu - BER_gpu_double| = %.4e\n', ber_diff_double);
|
|
fprintf('|BER_cpu - BER_gpu_single| = %.4e\n', ber_diff_single);
|
|
|
|
if ber_diff_double < 1e-6 && ber_diff_single < 1e-6
|
|
fprintf('✓ BER differences are negligible\n');
|
|
elseif ber_diff_single > ber_diff_double * 10
|
|
fprintf('⚠ Single precision shows measurable BER impact\n');
|
|
else
|
|
fprintf('✓ BER differences within acceptable range\n');
|
|
end
|
|
end
|
|
|
|
fprintf('\n========== Test Complete ==========\n');
|