From d8a3f777e174dd025920917d03e905e7ae9b2258 Mon Sep 17 00:00:00 2001
From: Antoine Hoffmann <antoine.hoffmann@epfl.ch>
Date: Tue, 25 Jul 2023 14:43:24 +0200
Subject: [PATCH] save

---
 wk/fast_analysis.m                            |   8 +-
 wk/load_metadata_scan.m                       | 174 ++++++++++++++++++
 wk/paper_2_scripts_and_results/fort_00.90     |  45 +++--
 .../load_metadata_scan.m                      | 170 +++++++++++++++++
 .../nonlin_kT_scan_analysis.m                 |   8 +-
 .../nonlin_nu_scan_analysis.m                 |   6 +-
 6 files changed, 385 insertions(+), 26 deletions(-)
 create mode 100644 wk/load_metadata_scan.m
 create mode 100644 wk/paper_2_scripts_and_results/load_metadata_scan.m

diff --git a/wk/fast_analysis.m b/wk/fast_analysis.m
index 02f1f918..a30b1677 100644
--- a/wk/fast_analysis.m
+++ b/wk/fast_analysis.m
@@ -71,7 +71,7 @@ PARTITION  = '/misc/gyacomo23_outputs/';
 % resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/9x5x192x96x24/nu_0.005';
 % resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24/nu_0.005';
 % resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24_large_box/nu_0.005';
-resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24_ckcut_0.8/nu_0.005';
+% resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24_ckcut_0.8/nu_0.005';
 %% kT eff study
 % resdir = 'paper_2_GYAC23/kT_eff_study/1x3x128x64x24_kT_3.0/Lx120';
 % resdir = 'paper_2_GYAC23/kT_eff_study/1x3x128x64x24_kT_3.0/Lx240';
@@ -82,7 +82,7 @@ resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24_ckcut
 % resdir = 'paper_2_GYAC23/kT_eff_study/7x3x128x64x24_kT_3.6/dmax';
 
 %% dev
-% PARTITION='';
+PARTITION='';
 % resdir = '/home/ahoffman/gyacomo/testcases/zpinch_3D';
 % resdir = '/home/ahoffman/gyacomo/results/dev/3D_kine_zpinch_test';
 %% CBC benchmark
@@ -93,7 +93,7 @@ resdir = 'paper_2_GYAC23/collision_study/nuLDGK_scan_kT_5.3/17x9x128x64x24_ckcut
 % resdir = '/misc/gyacomo23_outputs/paper_2_GYAC23/kT_scan_nu_1e-3/31x16x128x64x24/kT_7.0';
 
 %% Paper 3
-% resdir = '/misc/gyacomo23_outputs/paper_3_GYAC23/nonlin_KBM';
+resdir = '/misc/gyacomo23_outputs/paper_3_GYAC23/CBC_kine';
  %%
 J0 = 00; J1 = 20;
 
@@ -134,7 +134,7 @@ options.NAME      = '\phi';
 % options.NAME      = 'Q_x';
 % options.NAME      = 'n_i';
 % options.NAME      = 'n_i-n_e';
-options.PLAN      = 'xz';
+options.PLAN      = 'yz';
 % options.NAME      = 'f_i';
 % options.PLAN      = 'sx';
 options.COMP      = 'avg';
diff --git a/wk/load_metadata_scan.m b/wk/load_metadata_scan.m
new file mode 100644
index 00000000..eaec3e66
--- /dev/null
+++ b/wk/load_metadata_scan.m
@@ -0,0 +1,174 @@
+% Metadata path
+gyacomodir  = pwd;
+gyacomodir = gyacomodir(1:end-2);
+addpath(genpath([gyacomodir,'matlab'])) % ... add
+addpath(genpath([gyacomodir,'matlab/plot'])) % ... add
+addpath(genpath([gyacomodir,'matlab/compute'])) % ... add
+addpath(genpath([gyacomodir,'matlab/load'])) % ... add% EXECNAME = 'gyacomo_1.0';
+
+% datafname = 'lin_KBM/12x24_ky_0.05_0.75_P_2_16_DGGK_0.01_be_0.03.mat';
+% datafname = 'lin_KBM/12x24_ky_0.05_0.75_P_2_16_DGGK_0.01_be_0.03.mat';
+% datafname = 'lin_AUG_scan/12x24_ky_0.05_0.75_P_2_16_DGGK_0.01_be_0.000152.mat';
+% datafname = 'lin_Entropy_scan/2x1_ky_0.05_0.75_P_2_8_DGDK_0_be_0.mat';
+datafname = 'lin_DTT_AB_rho85_PT_scan/16x24_ky_0.1_1.5_P_2_8_DGGK_0.05_be_0.0034.mat';
+%% Chose if we filter gamma>0.05
+FILTERGAMMA = 0;
+
+%% Load data
+fname = ['../results/',datafname];
+d = load(fname);
+if FILTERGAMMA
+    d.data = d.data.*(d.data>0.025);
+    d.err  = d.err.*(d.data>0.025);
+end
+if 0
+%% Pcolor of the peak
+figure;
+% [XX_,YY_] = meshgrid(d.s1,d.s2);
+[XX_,YY_] = meshgrid(1:numel(d.s1),1:numel(d.s2));
+pclr=imagesc_custom(XX_,YY_,d.data'.*(d.data>0)');
+% pclr=contourf(1:numel(d.s1),1:numel(d.s2),d.data'.*(d.data>0)');
+% pclr=surf(1:numel(d.s1),1:numel(d.s2),d.data'.*(d.data>0)');
+title(d.title);
+xlabel(d.s1name); ylabel(d.s2name);
+set(gca,'XTick',1:numel(d.s1),'XTicklabel',d.s1)
+set(gca,'YTick',1:numel(d.s2),'YTicklabel',d.s2)
+colormap(jet)
+colormap(bluewhitered)
+clb=colorbar; 
+clb.Label.String = '$\gamma c_s/R$';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+if 1
+%% Scan along first dimension
+figure
+colors_ = jet(numel(d.s2));
+for i = 1:numel(d.s2)
+    % plot(d.s1,d.data(:,i),'s-',...
+    plot(d.s1(d.data(:,i)>0),d.data((d.data(:,i)>0),i),'s-',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s2name,'=',num2str(d.s2(i))],...
+        'color',colors_(i,:)); 
+    % errorbar(d.s1,d.data(:,i),d.err(:,i),'s-',...
+    % 'LineWidth',2.0,...
+    % 'DisplayName',[d.s2name,'=',num2str(d.s2(i))],...
+    % 'color',colors_(i,:)); 
+    hold on;
+end
+xlabel(d.s1name); ylabel(d.dname);title(d.title);
+xlim([d.s1(1) d.s1(end)]);
+colormap(colors_);
+clb = colorbar;
+% caxis([d.s2(1)-0.5,d.s2(end)+0.5]);
+clim([1 numel(d.s2)+1]);
+clb.Ticks=linspace(d.s2(1),d.s2(end),numel(d.s2));
+clb.Ticks    =1.5:numel(d.s2)+1.5;
+clb.TickLabels=d.s2;
+clb.Label.String = d.s2name;
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+if 0
+%% Scan along second dimension
+figure
+colors_ = jet(numel(d.s1));
+for i = 1:numel(d.s1)
+    plot(d.s2,d.data(i,:),'s-',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+        'color',colors_(i,:)); 
+%     errorbar(d.s2,d.data(i,:),d.err(i,:),'s-',...
+%         'LineWidth',2.0,...
+%         'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+%         'color',colors_(i,:)); 
+    hold on;
+end
+xlabel(d.s2name); ylabel(d.dname);title(d.title);
+xlim([d.s2(1) d.s2(end)]);
+colormap(jet(numel(d.s1)));
+clb = colorbar;
+caxis([d.s1(1)-0.5,d.s1(end)+0.5]);
+clb.Ticks=linspace(d.s1(1),d.s1(end),numel(d.s1));
+clb.YTick=d.s1;
+clb.Label.String = d.s1name;
+clb.TickLabelInterpreter = 'latex';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+
+if 0
+%% Convergence analysis
+figure
+% target_ = 0.25*(d.data(end,end)+d.data(end-i_,end)+d.data(end,end-i_)+d.data(end-i_,end-i_));
+colors_ = jet(numel(d.s1));
+for i = 1:numel(d.s1)
+%     target_ = d.data(i,end);
+    target_ = 2.79666916212537142172e-01; % Value for nuDGDK = 0.05, kT=6.96, (40,20), Nkx=8
+    if target_ > 0
+    eps_    = abs(target_ - d.data(i,1:end-1))/abs(target_);
+    semilogy(d.s2(1:end-1),eps_,'-s',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+        'color',colors_(i,:));
+    hold on;
+    end
+end
+xlabel(d.s2name); ylabel('$\epsilon_r$');title(d.title);
+xlim([d.s2(1) d.s2(end)]);
+colormap(colors_);
+clb = colorbar;
+caxis([d.s1(1)-0.5,d.s1(end)+0.5]);
+clb.Ticks=linspace(d.s1(1),d.s1(end),numel(d.s1));
+clb.YTick=d.s1;
+clb.Label.String = d.s1name;
+clb.TickLabelInterpreter = 'latex';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+grid on;
+end
+if 0
+%% Pcolor of the error
+figure;  i_ = 0;
+% target_ = 2.72724991618068013377e-01; % Value for nuDGDK = 1.0, kT=6.96, (40,20), Nkx=8
+% target_ = 2.79666916212537142172e-01; % Value for nuDGDK = 0.05, kT=6.96, (40,20), Nkx=8
+% target_ = 2.71455e-01; % Value for nuDGDK = 0.001, kT=6.96, (40,20), Nkx=8
+% target_ = 1.39427e-01; % Value for nuDGDK = 0.001, kT=5.3, (30,16), Nkx=8
+target_ = 2.72510405826983714839e-01  % Value for nuDGDK = 0.001, kT=6.96, (50,25), Nkx=8
+% target_ = 2.73048910051283844069e-01; % Value for nuDGDK = 0.0, kT=6.96, (40,20), Nkx=8
+% target_ = 0.25*(d.data(end,end)+d.data(end-i_,end)+d.data(end,end-i_)+d.data(end-i_,end-i_));
+% eps_    = log(abs(target_ - d.data)/abs(target_));
+eps_    = max(-10,log(abs(target_ - d.data)/abs(target_)));
+sign_   = 1;%sign(d.data - target_);
+eps_ = d.data;
+for i = 1:numel(d.s1)
+    for j = 1:numel(d.s2)
+        % target_ = d.data(i,end);
+        % target_ = d.data(i,end);
+        % target_ = d.data(end,j);
+        eps_(i,j) = log(abs(target_ - d.data(i,j))/target);
+        if target_ > 0
+    %     eps_(i,:) = max(-12,log(abs(target_ - d.data(i,1:end))/abs(target_)));
+    %     eps_(i,:) = log(abs(target_ - d.data(i,1:end))/abs(target_));
+    %     eps_(i,:) = min(100,100*abs(target_ - d.data(i,1:end))/abs(target_));
+        else
+        end
+    end
+end
+[XX_,YY_] = meshgrid(d.s1,d.s2);
+[XX_,YY_] = meshgrid(1:numel(d.s1),1:numel(d.s2));
+pclr=imagesc_custom(XX_,YY_,eps_'.*(d.data>0)'.*sign_');
+% pclr=contourf(1:numel(d.s1),1:numel(d.s2),eps_'.*(d.data>0)'.*sign_',5);
+% pclr=surf(1:numel(d.s1),1:numel(d.s2),eps_'.*(d.data>0)'.*sign_');
+title(d.title);
+xlabel(d.s1name); ylabel(d.s2name);
+set(gca,'XTick',1:numel(d.s1),'XTicklabel',d.s1)
+set(gca,'YTick',1:numel(d.s2),'YTicklabel',d.s2)
+% colormap(jet)
+colormap(bluewhitered)
+% caxis([-10, 0]);
+clb=colorbar; 
+clb.Label.String = '$\log(\epsilon_r)$';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
\ No newline at end of file
diff --git a/wk/paper_2_scripts_and_results/fort_00.90 b/wk/paper_2_scripts_and_results/fort_00.90
index c5c6dcd5..fe18d5ae 100644
--- a/wk/paper_2_scripts_and_results/fort_00.90
+++ b/wk/paper_2_scripts_and_results/fort_00.90
@@ -1,17 +1,17 @@
 &BASIC
   nrun       = 100000000
-  dt         = 0.0017678
-  tmax       = 50
+  dt         = 0.0022361
+  tmax       = 30
   maxruntime = 356400
   job2load   = -1
 /
 &GRID
-  pmax  = 16
-  jmax  = 8
+  pmax  = 10
+  jmax  = 5
   Nx     = 8
   Lx     = 7.854
   Ny     = 2
-  Ly     = 6.2832
+  Ly     = 20.944
   Nz     = 24
   SG     = .false.
   Nexc   = 1
@@ -22,11 +22,15 @@
   shear  = 0.8
   eps    = 0.18
   kappa  = 1
+  s_kappa  = 0
   delta  = 0
+  s_delta  = 0
   zeta   = 0
+  s_zeta   = 0
   parallel_bc = 'dirichlet'
   shift_y = 0
   Npol    = 1
+  PB_PHASE= .false.
 /
 &OUTPUT_PAR
   dtsave_0d = 1
@@ -46,7 +50,7 @@
 &MODEL_PAR
 LINEARITY = 'linear'
 RM_LD_T_EQ= .false.
-  Na      = 1
+  Na      = 2
   mu_x    = 0
   mu_y    = 0
   N_HD    = 4
@@ -54,12 +58,15 @@ RM_LD_T_EQ= .false.
   HYP_V   = 'none'
   mu_p    = 0
   mu_j    = 0
-  nu      = 0.01
+  nu      = 0.001
   k_gB    = 1
   k_cB    = 1
   lambdaD = 0
-  beta    = 0.0001
-  ADIAB_E = .true.
+  beta    = 0
+  ADIAB_E = .false.
+  ADIAB_I = .false.
+  tau_i   = 1
+  MHD_PD  = .false.
 /
 &CLOSURE_PAR
   hierarchy_closure='truncation'
@@ -68,18 +75,26 @@ RM_LD_T_EQ= .false.
   nmax             =0
 /
 &SPECIES
-  name_  = ions 
+  name_  = 'ions' 
   tau_   = 1
   sigma_ = 1
   q_     = 1
   K_N_   = 2.22
-  K_T_   = 5.3
+  K_T_   = 6.96
+/
+&SPECIES
+  name_  = 'electrons' 
+  tau_   = 1
+  sigma_ = 0.023338
+  q_     = -1
+  K_N_   = 2.22
+  K_T_   = 6.96
 /
 &COLLISION_PAR
-  collision_model = 'SG'
-  GK_CO      = .true.
+  collision_model = 'DG'
+  GK_CO      = .false.
   INTERSPECIES    = .true.
-  mat_file        = '/home/ahoffman/gyacomo/wk/paper_2_scripts_and_resuliCa/gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
+  mat_file        = '/home/ahoffman/gyacomo/wk/paper_2_scripts_and_resuliCa/null'
   collision_kcut  = 1.75
 /
 &INITIAL_CON
@@ -90,4 +105,4 @@ RM_LD_T_EQ= .false.
 /
 &TIME_INTEGRATION_PAR
   numerical_scheme = 'RK4'
-/
\ No newline at end of file
+/
diff --git a/wk/paper_2_scripts_and_results/load_metadata_scan.m b/wk/paper_2_scripts_and_results/load_metadata_scan.m
new file mode 100644
index 00000000..0f0293ee
--- /dev/null
+++ b/wk/paper_2_scripts_and_results/load_metadata_scan.m
@@ -0,0 +1,170 @@
+% Metadata path
+gyacomodir  = pwd;
+gyacomodir = gyacomodir(1:end-2);
+addpath(genpath([gyacomodir,'matlab'])) % ... add
+addpath(genpath([gyacomodir,'matlab/plot'])) % ... add
+addpath(genpath([gyacomodir,'matlab/compute'])) % ... add
+addpath(genpath([gyacomodir,'matlab/load'])) % ... add% EXECNAME = 'gyacomo_1.0';
+
+datafname = 'p2_linear_new/8x24_ky_0.3_kT_3_6.96_P_2_30_DGDK_0.001.mat';
+%% Chose if we filter gamma>0.05
+FILTERGAMMA = 0;
+
+%% Load data
+fname = ['../results/',datafname];
+d = load(fname);
+if FILTERGAMMA
+    d.data = d.data.*(d.data>0.025);
+    d.err  = d.err.*(d.data>0.025);
+end
+if 1
+%% Pcolor of the peak
+figure;
+% [XX_,YY_] = meshgrid(d.s1,d.s2);
+[XX_,YY_] = meshgrid(1:numel(d.s1),1:numel(d.s2));
+% pclr=imagesc_custom(XX_,YY_,d.data'.*(d.data>0)');
+pclr=contourf(1:numel(d.s1),1:numel(d.s2),d.data'.*(d.data>0)');
+% pclr=surf(1:numel(d.s1),1:numel(d.s2),d.data'.*(d.data>0)');
+title(d.title);
+xlabel(d.s1name); ylabel(d.s2name);
+set(gca,'XTick',1:numel(d.s1),'XTicklabel',d.s1)
+set(gca,'YTick',1:numel(d.s2),'YTicklabel',d.s2)
+colormap(jet)
+colormap(bluewhitered)
+clb=colorbar; 
+clb.Label.String = '$\gamma c_s/R$';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+if 0
+%% Scan along first dimension
+figure
+colors_ = jet(numel(d.s2));
+for i = 1:numel(d.s2)
+    % plot(d.s1,d.data(:,i),'s-',...
+    plot(d.s1(d.data(:,i)>0),d.data((d.data(:,i)>0),i),'s-',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s2name,'=',num2str(d.s2(i))],...
+        'color',colors_(i,:)); 
+    % errorbar(d.s1,d.data(:,i),d.err(:,i),'s-',...
+    % 'LineWidth',2.0,...
+    % 'DisplayName',[d.s2name,'=',num2str(d.s2(i))],...
+    % 'color',colors_(i,:)); 
+    hold on;
+end
+xlabel(d.s1name); ylabel(d.dname);title(d.title);
+xlim([d.s1(1) d.s1(end)]);
+colormap(colors_);
+clb = colorbar;
+% caxis([d.s2(1)-0.5,d.s2(end)+0.5]);
+clim([1 numel(d.s2)+1]);
+clb.Ticks=linspace(d.s2(1),d.s2(end),numel(d.s2));
+clb.Ticks    =1.5:numel(d.s2)+1.5;
+clb.TickLabels=d.s2;
+clb.Label.String = d.s2name;
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+if 0
+%% Scan along second dimension
+figure
+colors_ = jet(numel(d.s1));
+for i = 1:numel(d.s1)
+    plot(d.s2,d.data(i,:),'s-',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+        'color',colors_(i,:)); 
+%     errorbar(d.s2,d.data(i,:),d.err(i,:),'s-',...
+%         'LineWidth',2.0,...
+%         'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+%         'color',colors_(i,:)); 
+    hold on;
+end
+xlabel(d.s2name); ylabel(d.dname);title(d.title);
+xlim([d.s2(1) d.s2(end)]);
+colormap(jet(numel(d.s1)));
+clb = colorbar;
+caxis([d.s1(1)-0.5,d.s1(end)+0.5]);
+clb.Ticks=linspace(d.s1(1),d.s1(end),numel(d.s1));
+clb.YTick=d.s1;
+clb.Label.String = d.s1name;
+clb.TickLabelInterpreter = 'latex';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+
+if 0
+%% Convergence analysis
+figure
+% target_ = 0.25*(d.data(end,end)+d.data(end-i_,end)+d.data(end,end-i_)+d.data(end-i_,end-i_));
+colors_ = jet(numel(d.s1));
+for i = 1:numel(d.s1)
+%     target_ = d.data(i,end);
+    target_ = 2.79666916212537142172e-01; % Value for nuDGDK = 0.05, kT=6.96, (40,20), Nkx=8
+    if target_ > 0
+    eps_    = abs(target_ - d.data(i,1:end-1))/abs(target_);
+    semilogy(d.s2(1:end-1),eps_,'-s',...
+        'LineWidth',2.0,...
+        'DisplayName',[d.s1name,'=',num2str(d.s1(i))],...
+        'color',colors_(i,:));
+    hold on;
+    end
+end
+xlabel(d.s2name); ylabel('$\epsilon_r$');title(d.title);
+xlim([d.s2(1) d.s2(end)]);
+colormap(colors_);
+clb = colorbar;
+caxis([d.s1(1)-0.5,d.s1(end)+0.5]);
+clb.Ticks=linspace(d.s1(1),d.s1(end),numel(d.s1));
+clb.YTick=d.s1;
+clb.Label.String = d.s1name;
+clb.TickLabelInterpreter = 'latex';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+grid on;
+end
+if 0
+%% Pcolor of the error
+figure;  i_ = 0;
+% target_ = 2.72724991618068013377e-01; % Value for nuDGDK = 1.0, kT=6.96, (40,20), Nkx=8
+% target_ = 2.79666916212537142172e-01; % Value for nuDGDK = 0.05, kT=6.96, (40,20), Nkx=8
+% target_ = 2.71455e-01; % Value for nuDGDK = 0.001, kT=6.96, (40,20), Nkx=8
+% target_ = 1.39427e-01; % Value for nuDGDK = 0.001, kT=5.3, (30,16), Nkx=8
+target_ = 2.72510405826983714839e-01  % Value for nuDGDK = 0.001, kT=6.96, (50,25), Nkx=8
+% target_ = 2.73048910051283844069e-01; % Value for nuDGDK = 0.0, kT=6.96, (40,20), Nkx=8
+% target_ = 0.25*(d.data(end,end)+d.data(end-i_,end)+d.data(end,end-i_)+d.data(end-i_,end-i_));
+% eps_    = log(abs(target_ - d.data)/abs(target_));
+eps_    = max(-10,log(abs(target_ - d.data)/abs(target_)));
+sign_   = 1;%sign(d.data - target_);
+eps_ = d.data;
+for i = 1:numel(d.s1)
+    for j = 1:numel(d.s2)
+        % target_ = d.data(i,end);
+        % target_ = d.data(i,end);
+        % target_ = d.data(end,j);
+        eps_(i,j) = log(abs(target_ - d.data(i,j))/target);
+        if target_ > 0
+    %     eps_(i,:) = max(-12,log(abs(target_ - d.data(i,1:end))/abs(target_)));
+    %     eps_(i,:) = log(abs(target_ - d.data(i,1:end))/abs(target_));
+    %     eps_(i,:) = min(100,100*abs(target_ - d.data(i,1:end))/abs(target_));
+        else
+        end
+    end
+end
+[XX_,YY_] = meshgrid(d.s1,d.s2);
+[XX_,YY_] = meshgrid(1:numel(d.s1),1:numel(d.s2));
+pclr=imagesc_custom(XX_,YY_,eps_'.*(d.data>0)'.*sign_');
+% pclr=contourf(1:numel(d.s1),1:numel(d.s2),eps_'.*(d.data>0)'.*sign_',5);
+% pclr=surf(1:numel(d.s1),1:numel(d.s2),eps_'.*(d.data>0)'.*sign_');
+title(d.title);
+xlabel(d.s1name); ylabel(d.s2name);
+set(gca,'XTick',1:numel(d.s1),'XTicklabel',d.s1)
+set(gca,'YTick',1:numel(d.s2),'YTicklabel',d.s2)
+% colormap(jet)
+colormap(bluewhitered)
+% caxis([-10, 0]);
+clb=colorbar; 
+clb.Label.String = '$\log(\epsilon_r)$';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
\ No newline at end of file
diff --git a/wk/paper_2_scripts_and_results/nonlin_kT_scan_analysis.m b/wk/paper_2_scripts_and_results/nonlin_kT_scan_analysis.m
index d4b8e36c..9248201a 100644
--- a/wk/paper_2_scripts_and_results/nonlin_kT_scan_analysis.m
+++ b/wk/paper_2_scripts_and_results/nonlin_kT_scan_analysis.m
@@ -1,10 +1,10 @@
 kN=2.22;
 figure
-ERRBAR = 0; LOGSCALE = 1;
-% nustr = '1e-3'; mrkstyl='o';
+ERRBAR = 0; LOGSCALE = 0;
+nustr = '1e-3'; mrkstyl='o';
 % nustr = '1e-2'; mrkstyl='^';
-nustr = '5e-2'; mrkstyl='s';
-GENE = 1;
+% nustr = '5e-2'; mrkstyl='s';
+GENE = 0;
 if ~GENE
     rootdir = ['/misc/gyacomo23_outputs/paper_2_GYAC23/kT_scan_nu_',nustr];
     % mrkstyl='v';
diff --git a/wk/paper_2_scripts_and_results/nonlin_nu_scan_analysis.m b/wk/paper_2_scripts_and_results/nonlin_nu_scan_analysis.m
index a3223ace..fd1555a9 100644
--- a/wk/paper_2_scripts_and_results/nonlin_nu_scan_analysis.m
+++ b/wk/paper_2_scripts_and_results/nonlin_nu_scan_analysis.m
@@ -1,6 +1,6 @@
 kN=2.22;
 figure
-ERRBAR = 0; LOGSCALE = 0; AU = 0;
+ERRBAR = 1; LOGSCALE = 0; AU = 0;
 resstr={};
 msz = 10; lwt = 2.0;
 % CO = 'DGGK'; mrkstyl='d';
@@ -68,7 +68,7 @@ for j = 1:numel(directories)
         subdir = subdirectories{i};
         data    = compile_results_low_mem(data,subdir,00,20);
         try
-            Trange  = data.Ts0D(end)*[0.5 1.0];
+            Trange  = data.Ts0D(end)*[0.5 0.75];
         catch % if data does not exist put 0 everywhere
             data.Ts0D = 0;
             data.HFLUX_X = 0;
@@ -79,7 +79,7 @@ for j = 1:numel(directories)
             data.inputs.K_N  = kN;
             data.inputs.NU   = nus(i);
         end
-            Trange  = data.Ts0D(end)*[0.5 1.0];
+            Trange  = data.Ts0D(end)*[0.25 0.75];
             % Trange  = [200 400];
         %
         [~,it0] = min(abs(Trange(1)  -data.Ts0D)); 
-- 
GitLab