diff --git a/matlab/photomaton.m b/matlab/photomaton.m
index 6f4343921e1a6ee8b1feb4b4bea2994322674d58..d5508bfc87b4f3311e722c3f50656325b67f3e91 100644
--- a/matlab/photomaton.m
+++ b/matlab/photomaton.m
@@ -4,9 +4,9 @@ if 0
 % FIELD = ne00; FNAME = 'ne'; FIELDLTX = '$N_e^{00}$'; XX = RR; YY = ZZ;
 FIELD = phi; FNAME = 'phi'; FIELDLTX = '$\phi$'; XX = RR; YY = ZZ;
 % FIELD = dr2phi; FNAME = 'dr2phi'; XX = RR; YY = ZZ;
-tf = 800;  [~,it1] = min(abs(Ts2D-tf));
-tf = 900;  [~,it2] = min(abs(Ts2D-tf)); 
-tf = 1000; [~,it3] = min(abs(Ts2D-tf));
+tf = 100;  [~,it1] = min(abs(Ts2D-tf));
+tf = 285;  [~,it2] = min(abs(Ts2D-tf)); 
+tf = 350; [~,it3] = min(abs(Ts2D-tf));
 tf = 1100; [~,it4] = min(abs(Ts2D-tf));
 fig = figure; FIGNAME = [FNAME,'_snaps','_',PARAMS]; set(gcf, 'Position',  [100, 100, 1500, 350])
 plt = @(x) x;%./max(max(x));
diff --git a/matlab/setup.m b/matlab/setup.m
index a6b17e7b5f76e012ce3a806de627b2141c2c500b..dd1765268c914e2234c2c0027f90f1efb103f96c 100644
--- a/matlab/setup.m
+++ b/matlab/setup.m
@@ -70,7 +70,8 @@ elseif (CO == -2) % Write matrice filename for Sugama
     cmat_jmaxe = 5;
     cmat_pmaxi = 10;
     cmat_jmaxi = 5;
-    INITIAL.mat_file = ['''../../../iCa/SG_P_10_J_5_N_200_dk_0.05236_MFLR_0.h5'''];
+%     INITIAL.mat_file = ['''../../../iCa/SG_P_10_J_5_N_200_dk_0.05236_MFLR_0.h5'''];
+    INITIAL.mat_file = ['''../../../iCa/gk_sugama_P_10_J_5.h5'''];
     INITIAL.selfmat_file = ...
         ['''../../../iCa/self_Coll_GKE_0_GKI_0_ESELF_3_ISELF_3_Pmaxe_',num2str(cmat_pmaxe),...
         '_Jmaxe_',num2str(cmat_jmaxe),'_Pmaxi_',num2str(cmat_pmaxi),'_Jmaxi_',...
@@ -88,7 +89,11 @@ elseif (CO == 2) % Write matrice filename for Sugama GK
     cmat_jmaxe = 5;
     cmat_pmaxi = 10;
     cmat_jmaxi = 5;
-    INITIAL.mat_file = ['''../../../iCa/SG_P_10_J_5_N_200_dk_0.05236_MFLR_0.h5'''];
+%     INITIAL.mat_file = ['''../../../iCa/SG_P_10_J_5_N_200_dk_0.05236_MFLR_0.h5'''];
+%     INITIAL.mat_file = ['''../../../iCa/gk_sugama_P_4_J_2_N_10_kpm_1.h5'''];
+%     INITIAL.mat_file = ['''../../../iCa/gk_sugama_P_4_J_2_N_150_kpm_8.h5'''];
+    INITIAL.mat_file = ['''../../../iCa/gk_sugama_P_20_J_10_N_150_kpm_8.h5'''];
+%     INITIAL.mat_file = ['''../../../iCa/gk_sugama_P_10_J_5.h5'''];
     INITIAL.selfmat_file = ...
         ['''../../../iCa/self_Coll_GKE_1_GKI_1_ESELF_3_ISELF_3_Pmaxe_',num2str(cmat_pmaxe),...
         '_Jmaxe_',num2str(cmat_jmaxe),'_Pmaxi_',num2str(cmat_pmaxi),'_Jmaxi_',...
@@ -171,8 +176,8 @@ end
 %% Compile and WRITE input file
 INPUT = write_fort90(OUTPUTS,GRID,MODEL,INITIAL,TIME_INTEGRATION,BASIC);
 nproc = 1;
-MAKE  = 'cd ..; make; cd wk';
-system(MAKE);
+% MAKE  = 'cd ..; make; cd wk';
+% system(MAKE);
 %%
 disp(['Set up ',SIMID]);
 disp([resolution,gridname,degngrad]);
diff --git a/wk/analysis_2D.m b/wk/analysis_2D.m
index 29245fc1903201ea874a6a3eecf0b06564befd86..ee087aa66a141d78d88500acac76412b78239509 100644
--- a/wk/analysis_2D.m
+++ b/wk/analysis_2D.m
@@ -1,11 +1,13 @@
 addpath(genpath('../matlab')) % ... add
+% for ETA_ =[0.6:0.1:0.9]
 %% Load results
-outfile ='';
 if 1 % Local results
 outfile ='';
 outfile ='';
 outfile ='';
-outfile ='v2.6_P_2_J_1/200x100_L_120_P_2_J_1_eta_0.4_nu_1e-01_DGGK_CLOS_0_mu_1e-01';
+outfile ='kobayashi/100x50_L_50_P_4_J_2_eta_0.71429_nu_5e-03_DGGK_CLOS_0_mu_1e-02';
+% outfile ='kobayashi/100x50_L_50_P_4_J_2_eta_0.71429_nu_5e-03_SGGK_CLOS_0_mu_1e-02';
+% outfile ='v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
 
     BASIC.RESDIR      = ['../results/',outfile,'/'];
     BASIC.MISCDIR     = ['../results/',outfile,'/'];
@@ -13,9 +15,8 @@ end
 if 0 % Marconi results
 outfile ='';
 outfile ='';
-outfile ='';
-outfile ='';
-outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt';
+outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/kobayashi/100x50_L_50_P_6_J_3_eta_0.71429_nu_5e-03_DGGK_CLOS_0_mu_1e-02/out.txt';
+% outfile =['/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_',num2str(ETA_),'_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt'];
 % load_marconi(outfile);
 BASIC.RESDIR      = ['../',outfile(46:end-8),'/'];
 BASIC.MISCDIR     = ['/misc/HeLaZ_outputs/',outfile(46:end-8),'/'];
@@ -23,7 +24,7 @@ end
 
 %%
 % JOBNUM = 1; load_results;
-JOBNUMMAX = 20; compile_results %Compile the results from first output found to JOBNUMMAX if existing
+JOBNUMMAX = 0; compile_results %Compile the results from first output found to JOBNUMMAX if existing
 
 %% Retrieving max polynomial degree and sampling info
 Npe = numel(Pe); Nje = numel(Je); [JE,PE] = meshgrid(Je,Pe);
@@ -55,17 +56,19 @@ disp('- iFFT')
 % IFFT (Lower case = real space, upper case = frequency space)
 ne00   = zeros(Nr,Nz,Ns2D);         % Gyrocenter density
 ni00   = zeros(Nr,Nz,Ns2D);
+dzni00 = zeros(Nr,Nz,Ns2D);
 np_i   = zeros(Nr,Nz,Ns5D); % Ion particle density
 si00   = zeros(Nr,Nz,Ns5D);
 phi    = zeros(Nr,Nz,Ns2D);
 drphi  = zeros(Nr,Nz,Ns2D);
-dr2phi = zeros(Nr,Nz,Ns2D);
 dzphi  = zeros(Nr,Nz,Ns2D);
+dr2phi = zeros(Nr,Nz,Ns2D);
 
 for it = 1:numel(Ts2D)
     NE_ = Ne00(:,:,it); NI_ = Ni00(:,:,it); PH_ = PHI(:,:,it);
     ne00(:,:,it)  = real(fftshift(ifft2((NE_),Nr,Nz)));
     ni00(:,:,it)  = real(fftshift(ifft2((NI_),Nr,Nz)));
+    dzni00(:,:,it) = real(fftshift(ifft2(1i*KZ.*(NI_),Nr,Nz)));
     phi (:,:,it)  = real(fftshift(ifft2((PH_),Nr,Nz)));
     drphi(:,:,it) = real(fftshift(ifft2(1i*KR.*(PH_),Nr,Nz)));
     dr2phi(:,:,it)= real(fftshift(ifft2(-KR.^2.*(PH_),Nr,Nz)));
@@ -248,7 +251,7 @@ end
 
 if 1
 %% Space time diagramm (fig 11 Ivanov 2020)
-TAVG = 3500; % Averaging time duration
+TAVG = 500; % Averaging time duration
 %Compute steady radial transport
 tend = Ts0D(end); tstart = tend - TAVG;
 [~,its0D] = min(abs(Ts0D-tstart));
@@ -310,7 +313,8 @@ fig = figure; FIGNAME = ['space_time','_',PARAMS];set(gcf, 'Position',  [100, 10
         pclr = pcolor(TX,TY,squeeze(mean(ni00(:,:,trange).*dzphi(:,:,trange),2))'); set(pclr, 'edgecolor','none'); colorbar;
         shading interp
         colormap hot;
-        caxis([0.0,max(max(mean(ni00(:,:,its2D:ite2D).*dzphi(:,:,its2D:ite2D),2)))]);
+        caxis([0.0,0.8*max(max(mean(ni00(:,:,its2D:ite2D).*dzphi(:,:,its2D:ite2D),2)))]);
+%         caxis([0.0,0.1]);
          xticks([]); ylabel('$r/\rho_s$')
         legend('Radial part. transport $\langle n_i^{00}\partial_z\phi\rangle_z$')
         title(['$\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
@@ -325,26 +329,22 @@ fig = figure; FIGNAME = ['space_time','_',PARAMS];set(gcf, 'Position',  [100, 10
 save_figure
 end
 
-if 0
-%% Averaged phi, ZF and shear profiles
+if 1
+%% Averaged shear and Reynold stress profiles
 figure;
 plt = @(x) squeeze(mean(mean(real(x(:,:,its2D:ite2D)),2),3))/max(abs(squeeze(mean(mean(real(x(:,:,its2D:ite2D)),2),3))));
-subplot(311)
-plot(r,plt(phi)); hold on;
-xlim([-60,60]); xticks([]); yticks([]);
-subplot(312)
-plot(r,plt(drphi)); hold on;
-xlim([-60,60]);xticks([]);yticks([]);
-subplot(313)
-plot(r,plt(dr2phi)); hold on;
-xlim([-60,60]); xlabel('$r/\rho_s$');xticks([]);yticks([]);
+plot(r,plt(-dr2phi),'-k','DisplayName','Zonal shear'); hold on;
+plot(r,plt(-drphi.*dzphi),'--','DisplayName','$\Pi_\phi$'); hold on;
+plot(r,plt(-drphi.*dzni00),'--','DisplayName','$\Pi_{ni00}$'); hold on;
+plot(r,plt(-drphi.*dzphi-drphi.*dzni00),'--','DisplayName','$\Pi_\phi+\Pi_{ni00}$'); hold on;
+xlim([-60,60]); xlabel('$r/\rho_s$');
 end
 
 %%
-t0    = 500;
+t0    = 00;
 [~, it02D] = min(abs(Ts2D-t0));
 [~, it05D] = min(abs(Ts5D-t0));
-skip_ = 5; 
+skip_ = 1; 
 DELAY = 1e-3*skip_;
 FRAMES_2D = it02D:skip_:numel(Ts2D);
 FRAMES_5D = it05D:skip_:numel(Ts5D);
@@ -354,7 +354,8 @@ if 0
 GIFNAME = ['ni',sprintf('_%.2d',JOBNUM),'_',PARAMS]; INTERP = 0;
 FIELD = real(ni00); X = RR; Y = ZZ; T = Ts2D; FRAMES = FRAMES_2D;
 FIELDNAME = '$n_i$'; XNAME = '$r/\rho_s$'; YNAME = '$z/\rho_s$';
-create_gif
+% create_gif
+create_mov
 end
 if 0
 %% Density electrons
@@ -455,4 +456,6 @@ FIELD = plt(Nipj); X = sort(kz'); Y = Pi'; T = Ts5D; FRAMES = FRAMES_5D;
 FIELDNAME = 'N_i'; XNAME = '$k_z\rho_s$'; YNAME = '$P$, $\sum_r$';
 create_gif_5D
 end
-%%
\ No newline at end of file
+%%
+ZF_fourier_analysis
+% end
\ No newline at end of file
diff --git a/wk/compile_cosolver_mat.m b/wk/compile_cosolver_mat.m
index 5a3970d26613508c74ae589f45b93cbd97178b64..983d661850673082c72d87ad101db4867481438c 100644
--- a/wk/compile_cosolver_mat.m
+++ b/wk/compile_cosolver_mat.m
@@ -41,12 +41,12 @@ kperp = unique([0:dk:kperpmax,kperpmax]);
 if CO == 1; CONAME = 'FC'; end;
 if CO == 2; CONAME = 'PA'; end;
 if CO == 3; CONAME = 'SG'; end;
-matfilename = ['../iCa/',CONAME,'_P_',num2str(P),'_J_',num2str(J),...
-    '_N_',num2str(N),'_dk_',num2str(dk),'_MFLR_',num2str(M_FLR),'.h5'];
-
-n_ = 1;
+% matfilename = ['../iCa/',CONAME,'_P_',num2str(P),'_J_',num2str(J),...
+%     '_N_',num2str(N),'_dk_',num2str(dk),'_MFLR_',num2str(M_FLR),'.h5'];
+matfilename = '../iCa/gk_sugama_P_10_J_5.h5';
+n_ = 0;
 for k_ = kperp
-    disp(['-Writing matrix for kperp = ',num2str(k_)])
+    disp(['-Writing matrices for kperp = ',num2str(k_)])
     %% Script to run COSOlver in order to create needed collision matrices
     COSOlver_path = '../../Documents/MoliSolver/COSOlver/';
     COSOLVER.pmaxe = P;
@@ -68,6 +68,7 @@ for k_ = kperp
     COSOLVER.co = CO;
     
     k_string      = sprintf('%0.4f',k_);
+    n_string      = sprintf('%0.5d',n_);
     self_mat_file_name = ['../iCa/self_Coll_GKE_',num2str(COSOLVER.gk),'_GKI_',num2str(COSOLVER.gk),...
         '_ESELF_',num2str(COSOLVER.co),'_ISELF_',num2str(COSOLVER.co),...
         '_Pmaxe_',num2str(COSOLVER.pmaxe),'_Jmaxe_',num2str(COSOLVER.jmaxe),...
@@ -94,48 +95,43 @@ for k_ = kperp
     C_self = h5read(self_mat_file_name,'/Caapj/Ceepj');
     sz_ = size(C_self);
     % Write it in the compiled file
-    h5create(matfilename,['/Caapj/',k_string],sz_)
-    h5write(matfilename,['/Caapj/',k_string],C_self)
-    
+    h5create(matfilename,['/',n_string,'/Caapj/Ceepj'],sz_)
+    h5write (matfilename,['/',n_string,'/Caapj/Ceepj'],C_self)
+    h5create(matfilename,['/',n_string,'/Caapj/Ciipj'],sz_)
+    h5write (matfilename,['/',n_string,'/Caapj/Ciipj'],C_self)    
     %% Load ei matrices
     % Field
     C_eiF = h5read(ei_mat_file_name,'/Ceipj/CeipjF');
     sz_ = size(C_eiF);
-    h5create(matfilename,['/CeipjF/',k_string],sz_)
-    h5write(matfilename,['/CeipjF/',k_string],C_eiF)
+    h5create(matfilename,['/',n_string,'/Ceipj/CeipjF'],sz_)
+    h5write (matfilename,['/',n_string,'/Ceipj/CeipjF'],C_eiF)
     % Test
     C_eiT = h5read(ei_mat_file_name,'/Ceipj/CeipjT');
     sz_ = size(C_eiT);
-    h5create(matfilename,['/CeipjT/',k_string],sz_)
-    h5write(matfilename,['/CeipjT/',k_string],C_eiT)
+    h5create(matfilename,['/',n_string,'/Ceipj/CeipjT'],sz_)
+    h5write (matfilename,['/',n_string,'/Ceipj/CeipjT'],C_eiT)
     
     %% Load ie matrices
     % Field
     C_ieF = h5read(ie_mat_file_name,'/Ciepj/CiepjF');
     sz_ = size(C_ieF);
-    h5create(matfilename,['/CiepjF/',k_string],sz_)
-    h5write(matfilename,['/CiepjF/',k_string],C_ieF)
+    h5create(matfilename,['/',n_string,'/Ciepj/CiepjF'],sz_)
+    h5write (matfilename,['/',n_string,'/Ciepj/CiepjF'],C_ieF)
     % Test
     C_ieT = h5read(ie_mat_file_name,'/Ciepj/CiepjT');
     sz_ = size(C_eiT);
-    h5create(matfilename,['/CiepjT/',k_string],sz_)
-    h5write(matfilename,['/CiepjT/',k_string],C_ieT)
+    h5create(matfilename,['/',n_string,'/Ciepj/CiepjT'],sz_)
+    h5write (matfilename,['/',n_string,'/Ciepj/CiepjT'],C_ieT)
     
     %% Copy fort.90 input file and put grid params
     if(k_ == 0)
-        h5create(matfilename,'/dk',1);
-        h5write (matfilename,'/dk',dk);   
-        h5create(matfilename,'/N',1);
-        h5write (matfilename,'/N',N);
-        h5create(matfilename,'/Pmaxe',1);
-        h5write (matfilename,'/Pmaxe',P);   
-        h5create(matfilename,'/Jmaxe',1);
-        h5write (matfilename,'/Jmaxe',J);   
-        h5create(matfilename,'/Pmaxi',1);
-        h5write (matfilename,'/Pmaxi',P);  
-        h5create(matfilename,'/Jmaxi',1);
-        h5write (matfilename,'/Jmaxi',J);   
+        h5create(matfilename,'/coordkperp',numel(kperp));
+        h5write (matfilename,'/coordkperp',kperp);   
+        h5create(matfilename,'/dims_e',2);
+        h5write (matfilename,'/dims_e',[P,J]);   
+        h5create(matfilename,'/dims_i',2);
+        h5write (matfilename,'/dims_i',[P,J]); 
     end
-    
+    n_ = n_ + 1;
 end
 disp(['File saved @ :',matfilename])
\ No newline at end of file
diff --git a/wk/continue_multiple_runs_marconi.m b/wk/continue_multiple_runs_marconi.m
index 827b4713f69d41c233ad741df2c61f6f36140e62..a52c803fb40d099604e34bfb0a500a43bd9a03a6 100644
--- a/wk/continue_multiple_runs_marconi.m
+++ b/wk/continue_multiple_runs_marconi.m
@@ -1,11 +1,14 @@
 %% Paste the list of continue_run calls
-continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
+continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.7_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.8_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.9_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
 
 %% Functions to modify preexisting fort.90 input file and launch on marconi
 function [] = continue_run(outfilename)
-    EXECNAME = 'helaz_2.62';
+    EXECNAME = 'helaz_2.6';
     %% CLUSTER PARAMETERS
-    CLUSTER.PART  = 'dbg';     % dbg or prod
+    CLUSTER.PART  = 'prod';     % dbg or prod
     CLUSTER.TIME  = '24:00:00'; % allocation time hh:mm:ss
     if(strcmp(CLUSTER.PART,'dbg')); CLUSTER.TIME  = '00:30:00'; end;
     CLUSTER.MEM   = '64GB';     % Memory
@@ -44,7 +47,7 @@ function [] = continue_run(outfilename)
         line = A{33};
         line = line(end-2:end);
         if(line(1) == '='); line = line(end); end;
-        J2L = str2num(line)+1;
+        J2L = str2num(line) + 1;
     end
     % Change job 2 load in fort.90
     A{33} = ['  job2load      = ',num2str(J2L)];
diff --git a/wk/linear_study.m b/wk/linear_study.m
index 9bf350cf332112ca36160a6f11c980c0dd10f0e8..671fa975bddd50a377bd8dfd1ce34ab4c190c37d 100644
--- a/wk/linear_study.m
+++ b/wk/linear_study.m
@@ -1,6 +1,6 @@
-for NU = [1.0]
+for NU = [1e-1]
 for ETAB = [0.5]
-for CO = [-2 2]
+for CO = [1]
 %clear all;
 addpath(genpath('../matlab')) % ... add
 default_plots_options
@@ -18,8 +18,8 @@ NU_HYP  = 0.0;   % Hyperdiffusivity coefficient
 LAMBDAD = 0.0;
 NOISE0  = 1.0e-5;
 %% GRID PARAMETERS
-N       = 20;     % Frequency gridpoints (Nkr = N/2)
-L       = 120;     % Size of the squared frequency domain
+N       = 66;     % Frequency gridpoints (Nkr = N/2)
+L       = 150;     % Size of the squared frequency domain
 KREQ0   = 1;      % put kr = 0
 MU_P    = 0.0;     % Hermite  hyperdiffusivity -mu_p*(d/dvpar)^4 f
 MU_J    = 0.0;     % Laguerre hyperdiffusivity -mu_j*(d/dvperp)^4 f
@@ -33,8 +33,9 @@ SPSCP   = 0;    % Sampling per time unit for checkpoints
 RESTART = 0;      % To restart from last checkpoint
 JOB2LOAD= 00;
 %% OPTIONS
-% SIMID   = 'test';  % Name of the simulation
-SIMID   = 'v2.6_lin_analysis';  % Name of the simulation
+SIMID   = 'test';  % Name of the simulation
+% SIMID   = 'kobayashi_lin';  % Name of the simulation
+% SIMID   = 'v2.6_lin_analysis';  % Name of the simulation
 NON_LIN = 0 *(1-KREQ0);   % activate non-linearity (is cancelled if KREQ0 = 1)
 % Collision operator
 % (0 : L.Bernstein, 1 : Dougherty, 2: Sugama, 3 : Full Couloumb ; +/- for GK/DK)
@@ -64,10 +65,10 @@ MU      = NU_HYP/(HD_CO*kmax)^4 % Hyperdiffusivity coefficient
 %% PARAMETER SCANS
 if 1
 %% Parameter scan over PJ
-% PA = [2, 3, 4, 6];%, 8, 10];
-% JA = [1, 2, 2, 3];%, 4,  5];
+% PA = [2, 4, 6, 8, 10];
+% JA = [1, 2, 3, 4,  5];
 PA = [4];
-JA = [4];
+JA = [2];
 DTA= DT./sqrt(JA)/4;
 % DTA= DT;
 mup_ = MU_P;
@@ -88,8 +89,8 @@ for i = 1:Nparam
     setup
     % Run linear simulation
 %     system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 1 ./../../../bin/helaz 1 1; cd ../../../wk'])
-%     system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ./../../../bin/helaz 1 6; cd ../../../wk'])
-    system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ./../../../bin/helaz 2 3; cd ../../../wk'])
+    system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ./../../../bin/helaz 1 6; cd ../../../wk'])
+%     system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ./../../../bin/helaz 2 3; cd ../../../wk'])
 %     Load and process results
     %%
     filename = ['../results/',SIMID,'/',PARAMS,'/outputs_00.h5'];
diff --git a/wk/load_multiple_outputs_marconi.m b/wk/load_multiple_outputs_marconi.m
index 9863a0234ebdfb8c6e9d4b2c602eac5e70fbabdb..16b333b5f8297bb31bc7776e944c38368b123343 100644
--- a/wk/load_multiple_outputs_marconi.m
+++ b/wk/load_multiple_outputs_marconi.m
@@ -1,32 +1,7 @@
 addpath(genpath('../matlab')) % ... add
 %% Paste the list of simulation results to load
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-02_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-02_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-02_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-02_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-02_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-03_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-03_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-03_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-03_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-03_SGDK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-02_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-03_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-03_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-03_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-03_SGGK_CLOS_0_mu_2e-02/out.txt')
-load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.9_nu_1e-03_SGGK_CLOS_0_mu_2e-02/out.txt')
+load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.5_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.7_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.8_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
+load_marconi('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.9_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt')
diff --git a/wk/local_run.m b/wk/local_run.m
index f3de8f68e06c90c6149e6b96fffcbf2d98719cbc..8d3c7ac35a2cc058fe605cbd8fe63531bc2a65fa 100644
--- a/wk/local_run.m
+++ b/wk/local_run.m
@@ -4,19 +4,19 @@ addpath(genpath('../matlab')) % ... add
 CLUSTER.TIME  = '99:00:00'; % allocation time hh:mm:ss
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 %% PHYSICAL PARAMETERS
-NU      = 1.0e-1;   % Collision frequency
-ETAB    = 0.4;    % Magnetic gradient
+NU      = 5.0e-3;   % Collision frequency
+ETAB    = 1/1.4;    % Magnetic gradient
 ETAN    = 1.0;    % Density gradient
 NU_HYP  = 1.0;
 %% GRID PARAMETERS
-N       = 20;     % Frequency gridpoints (Nkr = N/2)
-L       = 120;     % Size of the squared frequency domain
-P       = 0;
-J       = 1;
+N       = 100;     % Frequency gridpoints (Nkr = N/2)
+L       = 50;     % Size of the squared frequency domain
+P       = 4;
+J       = 2;
 MU_P    = 0.0;     % Hermite  hyperdiffusivity -mu_p*(d/dvpar)^4 f
 MU_J    = 0.0;     % Laguerre hyperdiffusivity -mu_j*(d/dvperp)^4 f
 %% TIME PARAMETERS
-TMAX    = 1000;  % Maximal time unit
+TMAX    = 2000;  % Maximal time unit
 DT      = 1e-2;   % Time step
 SPS0D   = 1;    % Sampling per time unit for profiler
 SPS2D   = 1;      % Sampling per time unit for 2D arrays
@@ -27,10 +27,11 @@ JOB2LOAD= 0;
 %% OPTIONS AND NAMING
 % Collision operator
 % (0 : L.Bernstein, 1 : Dougherty, 2: Sugama, 3 : Full Couloumb ; +/- for GK/DK)
-CO      = 1;
+CO      = 2;
 CLOS    = 0;   % Closure model (0: =0 truncation)
 NL_CLOS = -1;   % nonlinear closure model (-2: nmax = jmax, -1: nmax = jmax-j, >=0 : nmax = NL_CLOS)
-SIMID   = 'test';  % Name of the simulation
+% SIMID   = 'test';  % Name of the simulation
+SIMID   = 'kobayashi';  % Name of the simulation
 % SIMID   = ['v2.6_P_',num2str(P),'_J_',num2str(J)];  % Name of the simulation
 NON_LIN = 1;   % activate non-linearity (is cancelled if KREQ0 = 1)
 INIT_ZF = 0; ZF_AMP = 0.0;
@@ -53,7 +54,7 @@ KR0KH   = 0; A0KH = 0; % Background phi mode to drive Ray-Tay inst.
 KREQ0   = 0;      % put kr = 0
 KPAR    = 0.0;    % Parellel wave vector component
 LAMBDAD = 0.0;
-kmax    = 2/3*N*pi/L;% Highest fourier mode
+kmax    = N*pi/L;% Highest fourier mode
 HD_CO   = 0.5;    % Hyper diffusivity cutoff ratio
 % kmaxcut = 2.5;
 MU      = NU_HYP/(HD_CO*kmax)^4 % Hyperdiffusivity coefficient
diff --git a/wk/marconi_run.m b/wk/marconi_run.m
index ad34af7e6e4c38cd31bc63b8e071bf6b1a92ddad..92dda8fe6122e835b5660acec9b2dc3a193138bb 100644
--- a/wk/marconi_run.m
+++ b/wk/marconi_run.m
@@ -1,9 +1,9 @@
 clear all;
 addpath(genpath('../matlab')) % ... add
-SUBMIT = 1; % To submit the job automatically
+SUBMIT = 0; % To submit the job automatically
 % EXECNAME = 'helaz_dbg';
   EXECNAME = 'helaz_2.62';
-for ETAB = [0.6 0.7]
+for ETAB = [0.5]
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 %% Set Up parameters
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -17,30 +17,31 @@ NP_P          = 2;          % MPI processes along p
 NP_KR         = 24;         % MPI processes along kr
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 %% PHYSICAL PARAMETERS
-NU      = 1e-1;   % Collision frequency
+NU      = 1e-3;   % Collision frequency
 % ETAB    = 0.7;   % Magnetic gradient
 NU_HYP  = 1.0;   % Hyperdiffusivity coefficient
 NL_CLOS = -1;   % nonlinear closure model (-2: nmax = jmax, -1: nmax = jmax-j, >=0 : nmax = NL_CLOS)
 % (0 : L.Bernstein, 1 : Dougherty, 2: Sugama, 3 : Full Couloumb ; +/- for GK/DK)
-CO      = -2;
+CO      = 1;
 INIT_ZF = 0; ZF_AMP = 0.0;
 %% GRID PARAMETERS
 N       = 200;    % Frequency gridpoints (Nkr = N/2)
 L       = 120;    % Size of the squared frequency domain
-P       = 06;     % Electron and Ion highest Hermite polynomial degree
-J       = 03;     % Electron and Ion highest Laguerre polynomial degree
+P       = 10;     % Electron and Ion highest Hermite polynomial degree
+J       = 05;     % Electron and Ion highest Laguerre polynomial degree
 MU_P    = 0.0;% Hermite  hyperdiffusivity -mu_p*(d/dvpar)^4 f
 MU_J    = 0.0;% Laguerre hyperdiffusivity -mu_j*(d/dvperp)^4 f
 %% TIME PARAMETERS
-TMAX    = 5000;  % Maximal time unit
-DT      = 1e-2;  % Time step
+TMAX    = 10000;  % Maximal time unit
+DT      = 1e-3;  % Time step
 SPS0D   = 1;     % Sampling per time unit for profiler
 SPS2D   = 1/4;     % Sampling per time unit for 2D arrays
 SPS5D   = 1/100;  % Sampling per time unit for 5D arrays
 SPSCP   = 0;     % Sampling per time unit for checkpoints
-RESTART = 0;     % To restart from last checkpoint
+RESTART = 1;     % To restart from last checkpoint
 JOB2LOAD= 0;
 %% Naming
+% SIMID   = 'kobayashi';  % Name of the simulation
 % SIMID   = 'test';  % Name of the simulation
 SIMID   = ['v2.6_P_',num2str(P),'_J_',num2str(J)];  % Name of the simulation
 PREFIX  =[];
diff --git a/wk/open_figure_script.m b/wk/open_figure_script.m
index 0ec2a71fd41f5bb31c0937017701871e5b7370fb..c69ae0c029dd87801dcffa4cc5ab0e2a61a10bf6 100644
--- a/wk/open_figure_script.m
+++ b/wk/open_figure_script.m
@@ -1,19 +1,23 @@
 simname_ = '';
 fname_='';
 %% Marconi output file
-fname_='';
-fname_='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt';
-% fname_='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-02_DGGK_CLOS_0_mu_2e-02/out.txt';
+% fname_='';
+% fname_='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_1e-01_DGGK_CLOS_0_mu_2e-02/out.txt';
+fname_='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.7_nu_1e-03_DGGK_CLOS_0_mu_2e-02/out.txt';
 simname_ = fname_(54:end-8);
 
 %%
-% simname_ = 'v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.7_nu_1e+00_SGGK_CLOS_0_mu_2e-02';
-% simname_ = 'v2.5_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e+00_SGGK_CLOS_0_mu_2e-02';
+simname_ = '';
+simname_ = '';
+simname_ = '';
+simname_ = '';
+simname_ = 'kobayashi/100x50_L_50_P_2_J_1_eta_0.71429_nu_5e-03_DGGK_CLOS_0_mu_1e-02';
 
 % simname_ = 'v2.6_P_2_J_1/200x100_L_120_P_2_J_1_eta_0.6_nu_1e-02_SGGK_CLOS_0_mu_1e-01';
 % simname_ = 'v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
 % simname_ = 'v2.5_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.8_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
-% simname_ = 'v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.6_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
+% simname_ = 'v2.6_P_6_J_3/200x100_L_120_P_6_J_3_eta_0.5_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
+% simname_ = 'v2.6_P_4_J_2/200x100_L_120_P_4_J_2_eta_0.6_nu_1e-01_DGGK_CLOS_0_mu_2e-02';
 
 
 %%
diff --git a/wk/transport_results_2_5.m b/wk/transport_results_2_5.m
index 46d2362717c3c1aa38acfbb4b9a398490fe553ac..df40b342799186f9ad0761c49971f19057f883f0 100644
--- a/wk/transport_results_2_5.m
+++ b/wk/transport_results_2_5.m
@@ -1,3 +1,4 @@
+default_plots_options
 %% nuDGGK = 1.0
 figure; set(gcf, 'Position',  [100, 100, 1200, 350])
 % 6,3 nuDGGK=1.0
diff --git a/wk/transport_results_2_6.m b/wk/transport_results_2_6.m
new file mode 100644
index 0000000000000000000000000000000000000000..47e78c83ad20133d8d3e16ede258488b725edb21
--- /dev/null
+++ b/wk/transport_results_2_6.m
@@ -0,0 +1,156 @@
+default_plots_options
+%% nuDGGK = 0.1
+figure; set(gcf, 'Position',  [100, 100, 1200, 350])
+% 6,3 nuDGGK=1.0
+eta   = [0.5,0.6,0.7,0.8];
+gamma = [32.6443 3.6895 0.3744 0.056];
+std_g = [6.1529 0.7986 0.0493 0.0134];
+subplot(121)
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=1.0$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(1,:)); hold on;
+subplot(122)
+
+% 10,5 nuDGGK = 1.0
+eta = [0.5 0.6,0.7,0.8];
+gamma = [32.6 3.5546 0.3917 0.0500];
+std_g = [7.7 0.5846 0.0486 0.0088];
+shear = [1.8505 0.60866 0.048249];
+std_s = [0.1599 0.00614 0.0061403];
+subplot(121)
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=1.0$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(2,:)); hold on;
+subplot(122)
+errorbar(eta(2:end),shear,std_s,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=1.0$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(2,:)); hold on;
+
+% 12,6 nuDGGK = 1.0
+eta = [0.6];
+gamma = [4.064];
+std_g = [0.7964];
+subplot(121)
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=12,06$, $\nu_{DGGK}=1.0$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(3,:)); hold on;
+subplot(122)
+
+
+% Mix len Ricci 2006
+subplot(121)
+semilogy([0.5  1.0],[10  1e-1],'--','color',[0,0,0]+0.6, 'DisplayName','Mix. Length');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$\Gamma^\infty_{part}$') 
+subplot(122)
+grate = [0.2131 0.106 0.02021];
+semilogy([0.6 0.7 0.8],grate,'--','color',[0,0,0]+0.6, 'DisplayName','Lin. Max Growth');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$V_E$''') 
+
+%% nuDGGK = 0.5
+figure; set(gcf, 'Position',  [100, 100, 1200, 350])
+% 6,3 nuDGGK=0.5
+eta = [0.5,0.6,0.7,0.8];
+gamma = [20.511 2.6292 0.2353 0.057];
+std_g = [3.67 1.2 0.055 0.004];
+shear = [nan 1.7417  0.57345 0.25155];
+std_s = [nan 0.35991 0.041 0.00913];
+subplot(121)
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=0.5$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(4,:)); hold on;
+subplot(122)
+errorbar(eta,shear,std_s,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=0.5$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(4,:)); hold on;
+
+% 10,5 nuDGGK=0.5
+eta    = [0.6,0.7,0.8 0.9];
+gamma  = [2.4949 0.23368 0.057792 0.023572];
+std_g  = [0.8696 0.085267 0.0060463 0.0046137];
+shear = [1.7094 0.55278 0.2054 0.085678];
+std_s = [0.2428 0.068223 0.012734 0.012291];
+subplot(121)
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.5$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(5,:)); hold on;
+subplot(122)
+errorbar(eta,shear,std_s,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.5$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(5,:)); hold on;
+
+% Mix len Ricci 2006
+subplot(121)
+semilogy([0.5  1.0],[10  1e-1],'--','color',[0,0,0]+0.6, 'DisplayName','Mix. Length');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$\Gamma^\infty_{part}$') 
+subplot(122)
+grate = [0.2194 0.129 0.05084 0.01346];
+semilogy([0.6 0.7 0.8 0.9],grate,'--','color',[0,0,0]+0.6, 'DisplayName','Lin. Max Growth');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$V_E$''') 
+
+%% nuDGGK = 0.1
+figure; set(gcf, 'Position',  [100, 100, 1200, 350])
+% 6,3
+eta =   [0.6        0.7     0.8];
+gamma = [0.24321    0.085   0.0367];
+std_g  = [0.295      0.05    0.0023];
+gbmax  = [1.0        0.21    0.0367];
+gbmin  = [0.02       0.04    0.0367];
+% errorbar(eta,gamma,std_,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=0.1$',...
+%     'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(6,:)); hold on;
+subplot(121)
+plot(eta,gamma,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.1$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(6,:)); hold on;
+plot(eta,gbmax,'-.','DisplayName',' ','Color',line_colors(6,:)); hold on;
+plot(eta,gbmin,'-.','DisplayName',' ','Color',line_colors(6,:)); hold on;
+
+% 10,5
+eta =   [0.5        0.6        0.7      0.8       0.9];
+gamma = [12.2       0.19761 0.088 0.04253 0.026037];
+std_g   = [4.7      0.21328 0.065 0.009 0.00118];
+gbmax  = [12.2      0.8 0.25 0.04253 0.026037];
+gbmin  = [12.2      0.02 0.03 0.04253 0.026037];
+shear = [0.65361 0.46548 0.30645 0.21123 ];
+std_s  = [0.21288 0.10233 0.02886 0.0044664 ];
+sbmax  = [1 0.7 0.30645 0.21123 ];
+sbmin  = [0.4 0.35 0.30645 0.21123 ];
+% errorbar(eta,gamma,std_,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.1$',...
+%     'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(7,:)); hold on;
+subplot(121)
+plot(eta,gamma,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.1$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(7,:)); hold on;
+plot(eta,gbmax,'-.','DisplayName',' ','Color',line_colors(7,:)); hold on;
+plot(eta,gbmin,'-.','DisplayName',' ','Color',line_colors(7,:)); hold on;
+subplot(122)
+plot(eta(2:end),shear,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=0.1$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(7,:)); hold on;
+plot(eta(2:end),sbmax,'-.','DisplayName',' ','Color',line_colors(7,:)); hold on;
+plot(eta(2:end),sbmin,'-.','DisplayName',' ','Color',line_colors(7,:)); hold on;
+
+% Mix len Ricci 2006
+subplot(121)
+semilogy([0.5  1.0],[10  1e-1],'--','color',[0,0,0]+0.6, 'DisplayName','Mix. Length');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$\Gamma^\infty_{part}$') 
+subplot(122)
+grate = [0.236 0.174 0.112 0.053];
+semilogy([0.6 0.7 0.8 0.9],grate,'--','color',[0,0,0]+0.6, 'DisplayName','Lin. Max Growth');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$V_E$''') 
+
+%% nuSGGK = 1.0
+figure
+% 6,3 nuDGGK=1.0
+eta = [0.5 0.6,0.7,0.8];
+gamma = [2.3 0.2215 0.0133 0.0032];
+std_g   = [3.1 0.22 0.0019 0.0006];
+errorbar(eta,gamma,std_g,'-','DisplayName','$P,J=06,03$, $\nu_{SGGK}=1.0$',...
+    'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(8,:)); hold on;
+
+
+% 10,5 nuDGGK = 1.0
+eta = [0.5 0.6,0.7,0.8];
+gamma = [10 0.319 0.009 0.0026];
+std_g   = [1.34 0.228 0.001 0.001];
+% errorbar(eta,gamma,err,'-.','DisplayName','$P,J=10,5$, $\nu_{SGGK}=1.0$',...
+%     'LineWidth',2.0,'MarkerSize',5, 'Color', line_colors(1,:)); hold on;
+
+% Mix len Ricci 2006
+semilogy([0.5  1.0],[10  1e-1],'--','color',[0,0,0]+0.6, 'DisplayName','Mix. Length');
+set(gca, 'YScale', 'log'); grid on; legend('show')
+xlabel('$L_n/R$'); ylabel('$\Gamma^\infty_{part}$')