- Added saddle-free Newton method implementation to continue optimiza…

…tion after particle swarm optimization finishes. Hessian is estimated using finite differences of the gradient (slow)

- Added constrained (trust region reflective) and unconstrained (quasi-newton) optimization options to optimize EDF length and power allocation
- .mat files with nonlinear coefficients for large-effective area fibers
- Added gap to capacity to optimization
- Modified function to run simulations on cluster
- New corning fiber

.gitignore

@@ -39,3 +39,4 @@ edfa/doc/header.tex
 *.xls
 *.mexw64
 !f/GN_model_coeff_spanLengthkm=50km_Df=50GHz.mat
+!f/GN_model_coeff_spanLengthkm=50km_Df=33GHz.mat

edfa/capacity_vs_pump_power.sbatch

@@ -4,8 +4,8 @@
 #SBATCH --output=capacity_vs_Ppump_%A_%a.out
 #SBATCH --error=capacity_vs_Ppump_%A_%a.err
 #SBATCH --array=1-29
-#SBATCH --time=4-0
-#SBATCH --partition normal
+#SBATCH --time=1-0
+#SBATCH --partition=normal
 #SBATCH --ntasks=1
 #SBATCH --mem-per-cpu=4000
 #SBATCH --cpus-per-task=8

edfa/capacity_vs_pump_power_slurm.m

@@ -6,7 +6,7 @@ function capacity_vs_pump_power_slurm(task)
 addpath f/
 addpath ../f/
 
-verbose = false;
+verbose = true;
 
 % Select task
 taskList = [30:5:150 200:100:500]; % Variable to be modified in different system calls
@@ -28,9 +28,10 @@ function capacity_vs_pump_power_slurm(task)
 Pump = Channels(980e-9, pumpPower, 'forward');
 
 % SMF fiber
-SMF = fiber(spanLengthKm*1e3, @(lamb) 0.18, @(lamb) 0);
+SMF = fiber(spanLengthKm*1e3, @(l) 0.165*ones(size(l)), @(l) 20.4e-6*ones(size(l)));
 SMF.gamma = 0.8e-3;
-[~, spanAttdB] = SMF.link_attenuation(Signal.wavelength);
+[~, spanAttdB] = SMF.link_attenuation(1550e-9); % same attenuation for all wavelengths
+spanAttdB = spanAttdB + 1.5;
 
 % Filename
 filename = sprintf('results/capacity_vs_pump_power_EDF=%s_pump=%dmW_%dnm_L=%d_x_%dkm.mat',...
@@ -41,18 +42,26 @@ function capacity_vs_pump_power_slurm(task)
 % Problem variables
 problem.spanAttdB = spanAttdB;
 problem.df = df;
+problem.Gap = 10^(-1/10);
 problem.Namp = Nspans;
 problem.step_approx = @(x) 0.5*(tanh(2*x) + 1); % Smoothing factor = 2
 problem.diff_step_approx = @(x) sech(2*x).^2; % first derivative (used for computing gradient)
 problem.excess_noise_correction = 1.4; % 1.2 for 980nm, 1.6 for 1480nm
-problem.SwarmSize = min(100, 20*(Signal.N+1));
+problem.SwarmSize = min(300, 20*(Signal.N+1));
 problem.nonlinearity = true;
 S = load('../f/GN_model_coeff_spanLengthkm=50km_Df=50GHz.mat');
 problem.nonlinear_coeff = S.nonlinear_coeff;
 problem.epsilon = 0.05; % From Fig. 17 of P. Poggiolini and I. Paper, “The GN Model
 % of Non-Linear Propagation in Uncompensated Coherent Optical Systems,” 
 % J. Lightw. Technol., vol. 30, no. 24, pp. 3857–3879, 2012.
 
+options.AdaptationConstant = 0.1; 
+options.FiniteDiffStepSize = 1e-6;
+options.MaxIterations = 100;
+options.AbsTol = 1e-3;
+options.MinStep = 1e-4;
+problem.saddle_free_newton.options = options;
+
 % Define power cap according to conservation of energy and add 3 dB as
 % safety margin (multiplication by 2)
 problem.Pon = 2*(1/Signal.N)*Pump.wavelength*Pump.P/(min(Signal.wavelength)*(10^(spanAttdB/10)-1));
@@ -85,11 +94,23 @@ function capacity_vs_pump_power_slurm(task)
     [nlin.E, nlin.S, nlin.exitflag, nlin.num, nlin.approx] = ...
         optimize_power_load_and_edf_length('particle swarm', E, Pump, Signal, problem, verbose);
 
-    % local optimization after particle swarm is done
-    problem.nonlinearity = true;
-    disp('== Nonlinear regime (local optimization)')
-    [nlin_local.E, nlin_local.S, nlin_local.exitflag, nlin_local.num, nlin_local.approx] = ...
-        optimize_power_load_and_edf_length('local', nlin.E, Pump, nlin.S, problem, verbose);
+    try
+        % local optimization after particle swarm is done
+        disp('== Nonlinear regime (unconstrained optimization)')
+        [nlin_unc.E, nlin_unc.S, nlin_unc.exitflag, nlin_unc.num, nlin_unc.approx] = ...
+            optimize_power_load_and_edf_length('local unconstrained', nlin.E, Pump, nlin.S, problem, verbose);
+    catch ee
+        warning(ee.message)
+    end
+
+    try
+        % Continue optimization using Saddle-free Newton method
+        disp('== Nonlinear regime (Saddle-free Newton)')
+        [nlin_sfn.E, nlin_sfn.S, nlin_sfn.exitflag, nlin_sfn.num, nlin_sfn.approx] = ...
+            optimize_power_load_and_edf_length('saddle-free newton', nlin.E, Pump, nlin.S, problem, verbose);
+    catch ee
+        warning(ee.message)
+    end
 
 catch e
     warning(e.message)