From 9ab19c11680429755038cb79bad6615bc194f0d1 Mon Sep 17 00:00:00 2001
From: Antoine Hoffmann <antoine.hoffmann@epfl.ch>
Date: Thu, 26 Jan 2023 13:40:30 +0100
Subject: [PATCH] scripts save
---
wk/CBC_nu_CO_scan.m | 132 +++++++++++++++++++++++++----------------
wk/analysis_gene.m | 53 ++++++++++++-----
wk/analysis_gyacomo.m | 36 +++++------
wk/header_3D_results.m | 19 +++++-
wk/lin_3Dzpinch.m | 32 +++++-----
5 files changed, 169 insertions(+), 103 deletions(-)
diff --git a/wk/CBC_nu_CO_scan.m b/wk/CBC_nu_CO_scan.m
index 84b33df4..019a372c 100644
--- a/wk/CBC_nu_CO_scan.m
+++ b/wk/CBC_nu_CO_scan.m
@@ -1,21 +1,21 @@
-gyacomodir = '/home/ahoffman/gyacomo/';
+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';
-% EXECNAME = 'gyacomo_dbg';
+% EXECNAME = 'gyacomo_debug';
EXECNAME = 'gyacomo';
CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss
%%
-% SIMID = 'linear_CBC_nu+PJ_scan_kT_6.96_SGGK'; % Name of the simulation
-SIMID = 'convergence_pITG_dbg'; % Name of the simulation
-% SIMID = 'linear_CBC_nu_scan_kT_11_ky_0.3_DGGK'; % Name of the simulation
-RERUN = 1; % If you want to rerun the sim (bypass the check of existing data)
+SIMID = 'p2_CBC_convergence_coll_PJ'; % Name of the simulation
+RERUN = 0; % If you want to rerun the sim (bypass the check of existing data)
+RUN = 1;
+% NU_a = [0.0];
+% P_a = [8];
-% NU_a = [0.0 0.01 0.02 0.05 0.1];
-NU_a = [0.0];
-P_a = [30];
-% P_a = [2 4:4:36];
+NU_a = [0.0:0.01:0.1];
+P_a = [2:2:30];
J_a = floor(P_a/2);
% collision setting
CO = 'DG';
@@ -25,18 +25,16 @@ COLL_KCUT = 1.75;
KIN_E = 0; % 1: kinetic electrons, 2: adiabatic electrons
BETA = 1e-4; % electron plasma beta
% background gradients setting
-K_Ne = 0*2.22; % ele Density '''
-K_Te = 0*6.96; % ele Temperature '''
-K_Ni = 0*2.22; % ion Density gradient drive
-K_Ti = 6.96; % ion Temperature '''
+K_N = 2.22;
+% K_T = 6.96;
+K_T = 5.3;
% Geometry
-GEOMETRY= 'miller';
-% GEOMETRY= 's-alpha';
-SHEAR = 0.0; % magnetic shear
+GEOMETRY= 's-alpha';
+SHEAR = 0.8; % magnetic shear
% time and numerical grid
-DT = 1e-2;
+DT = 1e-3;
TMAX = 50;
-kymin = 0.4;
+kymin = 0.3;
NY = 2;
% arrays for the result
g_ky = zeros(numel(NU_a),numel(P_a),NY/2+1);
@@ -58,7 +56,11 @@ for NU = NU_a
JMAXE = J; % Laguerre "
PMAXI = P; % " ions
JMAXI = J; % "
- NX = 2; % real space x-gridpoints
+ K_Ne = K_N; % ele Density '''
+ K_Te = K_T; % ele Temperature '''
+ K_Ni = K_N; % ion Density gradient drive
+ K_Ti = K_T; % ion Temperature '''
+ NX = 12; % real space x-gridpoints
LX = 2*pi/0.8; % Size of the squared frequency domain
LY = 2*pi/kymin; % Size of the squared frequency domain
NZ = 24; % number of perpendicular planes (parallel grid)
@@ -121,8 +123,8 @@ for NU = NU_a
filename = [SIMID,'/',PARAMS,'/'];
LOCALDIR = [gyacomodir,'results/',filename,'/'];
% check if data exist to run if no data
- data = compile_results(LOCALDIR,0,0); %Compile the results from first output found to JOBNUMMAX if existing
- if (RERUN || isempty(data.NU_EVOL))
+ data_ = compile_results(LOCALDIR,0,0); %Compile the results from first output found to JOBNUMMAX if existing
+ if RUN && (RERUN || isempty(data_.NU_EVOL) || numel(data_.Ts3D)<10)
system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 4 ',gyacomodir,'bin/',EXECNAME,' 1 2 2 0; cd ../../../wk'])
% system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 1 ',gyacomodir,'bin/',EXECNAME,' 1 1 1 0; cd ../../../wk'])
end
@@ -130,34 +132,25 @@ for NU = NU_a
% Load results after trying to run
filename = [SIMID,'/',PARAMS,'/'];
LOCALDIR = [gyacomodir,'results/',filename,'/'];
- data = compile_results(LOCALDIR,0,0); %Compile the results from first output found to JOBNUMMAX if existing
+ data_ = compile_results(LOCALDIR,0,0); %Compile the results from first output found to JOBNUMMAX if existing
- % linear growth rate (adapted for 2D zpinch and fluxtube)
- options.TRANGE = [0.5 1]*data.Ts3D(end);
- options.NPLOTS = 0; % 1 for only growth rate and error, 2 for omega local evolution, 3 for plot according to z
- options.GOK = 0; %plot 0: gamma 1: gamma/k 2: gamma^2/k^3
-% lg = compute_fluxtube_growth_rate(data,options);
-% [gmax, kmax] = max(lg.g_ky(:,end));
-% [gmaxok, kmaxok] = max(lg.g_ky(:,end)./lg.ky);
-% g_ky(i,j,:) = g_ky;
-%
-% g_avg(i,j,:) = lg.avg_g;
-% g_std(i,j,:) = lg.std_g;
-
- [~,it1] = min(abs(data.Ts3D-0.5*data.Ts3D(end))); % start of the measurement time window
- [~,it2] = min(abs(data.Ts3D-1.0*data.Ts3D(end))); % end of ...
+ [~,it1] = min(abs(data_.Ts3D-0.8*data_.Ts3D(end))); % start of the measurement time window
+ [~,it2] = min(abs(data_.Ts3D-1.0*data_.Ts3D(end))); % end of ...
field = 0;
field_t = 0;
- for ik = 1:NY/2+1
- field = squeeze(sum(abs(data.PHI),3)); % take the sum over z
+ for ik = 2:NY/2+1
+ field = squeeze(sum(abs(data_.PHI),3)); % take the sum over z
field_t = squeeze(field(ik,1,:)); % take the kx =0, ky = ky mode only
- to_measure = log(field_t);
- gr = polyfit(data.Ts3D(it1:it2),to_measure(it1:it2),1);
- g_ky(i,j,ik) = gr(1);
+ to_measure = log(field_t(it1:it2));
+ tw = data_.Ts3D(it1:it2);
+ gr = fit(tw,to_measure,'poly1');
+ err= confint(gr);
+ g_ky(i,j,ik) = gr.p1;
+ g_std(i,j,ik) = abs(err(2,1)-err(1,1))/2;
end
[gmax, ikmax] = max(g_ky(i,j,:));
- msg = sprintf('gmax = %2.2f, kmax = %2.2f',gmax,data.ky(ikmax)); disp(msg);
+ msg = sprintf('gmax = %2.2f, kmax = %2.2f',gmax,data_.ky(ikmax)); disp(msg);
i = i + 1;
@@ -165,6 +158,12 @@ end
j = j + 1;
end
+if 0
+%% Check time evolution
+figure;
+plot(data_.Ts3D,to_measure); hold on
+plot(data_.Ts3D(it1:it2),to_measure(it1:it2),'--');
+end
if 1
%% Study of the peak growth rate
figure
@@ -181,7 +180,7 @@ for i = 1:numel(idx_)
e_ = g_std(:,:,idx_(i));
end
-colors_ = lines(numel(NU_a));
+colors_ = jet(numel(NU_a));
subplot(121)
for i = 1:numel(NU_a)
% errorbar(P_a,y_(i,:),e_(i,:),...
@@ -194,23 +193,19 @@ for i = 1:numel(NU_a)
'color',colors_(i,:));
hold on;
end
-title(['$\kappa_T=$',num2str(K_Ti),' $k_y=k_y^{max}$']);
+title(['$\kappa_T=$',num2str(K_Ti),' $k_y=$',num2str(kymin)]);
legend('show'); xlabel('$P$, $J=P/2$'); ylabel('$\gamma$');
colors_ = jet(numel(P_a));
subplot(122)
for j = 1:numel(P_a)
-% errorbar(NU_a,y_(:,j),e_(:,j),...
-% 'LineWidth',1.2,...
-% 'DisplayName',['(',num2str(P_a(j)),',',num2str(J_a(j)),')'],...
-% 'color',colors_(j,:));
plot(NU_a,y_(:,j),'s-',...
'LineWidth',2.0,...
'DisplayName',['(',num2str(P_a(j)),',',num2str(J_a(j)),')'],...
'color',colors_(j,:));
hold on;
end
-title(['$\kappa_T=$',num2str(K_Ti),' $k_y=k_y^{max}$']);
+title(['$\kappa_T=$',num2str(K_Ti),' $k_y=$',num2str(kymin)]);
legend('show'); xlabel(['$\nu_{',CO,'}$']); ylabel('$\gamma$');
end
@@ -218,5 +213,38 @@ if 0
%% Pcolor of the peak
figure;
[XX_,YY_] = meshgrid(NU_a,P_a);
-pclr=pcolor(XX_,YY_,y_'); set(pclr,'EdgeColor','none');
-end
\ No newline at end of file
+% pclr=pcolor(XX_,YY_,y_'); set(pclr,'EdgeColor','none'); axis ij;
+pclr=imagesc_custom(XX_,YY_,y_'.*(y_>0)');
+title(['$\kappa_T$',num2str(data_.K_T),', $\kappa_N=$',num2str(K_N),', $k_y=$',num2str(kymin)]);
+xlabel('$\nu$'); ylabel('$P$, $J=P/2$');
+colormap(jet)
+clb=colorbar;
+clb.Label.String = '$\gamma c_s/R$';
+clb.Label.Interpreter = 'latex';
+clb.Label.FontSize= 18;
+end
+
+%% Save metadata
+numin = num2str(min(NU_a)); numax = num2str(max(NU_a));
+ pmin = num2str(min(P_a)); pmax = num2str(max(P_a));
+filename = [num2str(NX),'x',num2str(NZ),'_ky_',num2str(kymin),...
+ '_nu_',numin,'_',numax,'_',...
+ '_P_',pmin,'_',pmax,'_KT_',num2str(K_Ti),'.mat'];
+metadata.name = filename;
+metadata.kymin = kymin;
+metadata.title = ['$\kappa_T$',num2str(K_Ti),', $\kappa_N=$',num2str(K_N),', $k_y=$',num2str(kymin)];
+metadata.par = data_.PARAMS;
+metadata.nscan = 2;
+metadata.s1name = '$\nu$';
+metadata.s1 = NU_a;
+metadata.s2name = '$P$, $J=P/2$';
+metadata.s2 = P_a;
+metadata.dname = '$\gamma c_s/R$';
+metadata.data = y_;
+metadata.err = e_;
+metadata.input_file = h5read([data_.localdir,'/outputs_00.h5'],'/files/STDIN.00');
+metadata.date = date;
+% tosave.data = metadata;
+save([SIMDIR,filename],'-struct','metadata');
+disp(['saved in ',SIMDIR,filename]);
+clear metadata tosave
diff --git a/wk/analysis_gene.m b/wk/analysis_gene.m
index 42da732b..0c3c014d 100644
--- a/wk/analysis_gene.m
+++ b/wk/analysis_gene.m
@@ -29,7 +29,7 @@ addpath(genpath([gyacomodir,'matlab/load'])) % ... add
% folder = '/misc/gene_results/CBC/128x64x16x24x12/';
% folder = '/misc/gene_results/CBC/196x96x20x32x16_01/';
% folder = '/misc/gene_results/CBC/256x96x24x32x12/';
-% folder = '/misc/gene_results/CBC/128x64x16x6x4/';
+folder = '/misc/gene_results/CBC/128x64x16x6x4/';
% folder = '/misc/gene_results/CBC/KT_5.3_128x64x16x24x12_01/';
% folder = '/misc/gene_results/CBC/KT_4.5_128x64x16x24x12_01/';
% folder = '/misc/gene_results/CBC/KT_9_128x64x16x24x12/';
@@ -40,14 +40,29 @@ addpath(genpath([gyacomodir,'matlab/load'])) % ... add
% folder = '/misc/gene_results/CBC/KT_4.5_192x96x24x30x16_00/';
% folder = '/misc/gene_results/CBC/KT_5.0_192x96x24x30x16_00/';
% folder = '/misc/gene_results/CBC/KT_5.3_192x96x24x30x16_00/';
-% folder = '/misc/gene_results/CBC/KT_6.96_64x32x32x24x12_Nexc_5/';
-
% folder = '/misc/gene_results/CBC/196x96x20x32x16_01/';
% folder = '/misc/gene_results/linear_miller_CBC/hdf5_miller_s0_adiabe/';
-folder = '/misc/gene_results/linear_miller_CBC/hdf5_salpha_s0_adiabe/';
+% folder = '/misc/gene_results/linear_miller_CBC/hdf5_salpha_s0_adiabe/';
+
+%Paper 2
+% folder = '/misc/gene_results/CBC/KT_6.96_64x32x32x24x12_Nexc_5/';
+% folder = '/misc/gene_results/CBC/KT_6.96_128x64x24x32x16_Nexc_5_00/';
+% folder = '/misc/gene_results/CBC/KT_6.96_128x64x24x32x16_Nexc_5_01/';
+% folder = '/misc/gene_results/CBC/KT_5.3_128x64x24x32x16_Nexc_5_00/';
+% folder = '/misc/gene_results/CBC/KT_5.3_128x64x24x32x16_Nexc_5_01/';
+
gene_data = load_gene_data(folder);
+gene_data.FIGDIR = folder;
gene_data = invert_kxky_to_kykx_gene_results(gene_data);
-if 1
+gene_data.Pmaxi = gene_data.Nvp-1;
+gene_data.Jmaxi = gene_data.Nmu-1;
+gene_data.CODENAME = 'GENE';
+
+%% Dashboard (Compilation of main plots of the sim)
+dashboard(gene_data);
+
+%% Separated plot routines
+if 0
%% Space time diagramm (fig 11 Ivanov 2020)
options.TAVG_0 = 0.1*gene_data.Ts3D(end);
options.TAVG_1 = gene_data.Ts3D(end); % Averaging times duration
@@ -56,8 +71,7 @@ options.ST_FIELD = '\phi'; % chose your field to plot in spacetime diag
options.INTERP = 1;
options.NCUT = 4; % Number of cuts for averaging and error estimation
options.RESOLUTION = 256;
-gene_data.FIGDIR = folder;
-fig = plot_radial_transport_and_spacetime(gene_data,options);
+fig = plot_radial_transport_and_spacetime(gene_data,options,'GENE');
save_figure(gene_data,fig,'.png')
end
@@ -94,10 +108,10 @@ if 0
options.INTERP = 1;
options.POLARPLOT = 0;
options.NAME = '\phi';
-% options.NAME = 'v_y';
-% options.NAME = '\Gamma_x';
+% options.NAME = 'v_{Ey}';
+% options.NAME = 'G_x';
% options.NAME = 'n_i';
-options.PLAN = 'kyz';
+options.PLAN = 'xz';
% options.NAME = 'f_e';
% options.PLAN = 'sx';
options.COMP = 'avg';
@@ -135,23 +149,26 @@ end
if 0
%% Show f_i(vpar,mu)
-options.times = 70;
-options.specie = 'i';
+options.T = [0.5 1]*gene_data.Ts3D(end);
+options.SPECIES = 'i';
% options.PLT_FCT = 'contour';
% options.PLT_FCT = 'contourf';
% options.PLT_FCT = 'surf';
options.PLT_FCT = 'surfvv';
-options.folder = folder;
+options.non_adiab = 0;
+options.RMS = 1; % Root mean square i.e. sqrt(sum_k|f_k|^2) as in Gene
+options.folder = gene_data.folder;
options.iz = 'avg';
options.FIELD = '<f_>';
-options.ONED = 0;
-% options.FIELD = 'Q_es';
+options.SPAR = linspace(-3,3,32);
+options.XPERP = linspace( 0,sqrt(6),16).^2;
+options.ONED = 1;
plot_fa_gene(options);
end
if 0
%% Time averaged spectrum
-options.TIME = [4000 5000];
+options.TIME = [100 500];
options.NORM =1;
% options.NAME = '\phi';
% options.NAME = 'n_i';
@@ -171,9 +188,13 @@ if 0
options.NORMALIZED = 0;
options.K2PLOT = 1;
options.TIME = 1:700;
+options.KX_TW = [25 55]; %kx Growth rate time window
+options.KY_TW = [0 20]; %ky Growth rate time window
options.NMA = 1;
options.NMODES = 15;
options.iz = 'avg';
+options.ik = 1; % sum, max or index
+options.fftz.flag = 0;
fig = mode_growth_meter(gene_data,options);
save_figure(gene_data,fig)
end
diff --git a/wk/analysis_gyacomo.m b/wk/analysis_gyacomo.m
index 88578a58..a7f93e8a 100644
--- a/wk/analysis_gyacomo.m
+++ b/wk/analysis_gyacomo.m
@@ -13,13 +13,14 @@ LOCALDIR = [gyacomodir,resdir,'/'];
MISCDIR = ['/misc/gyacomo_outputs/',resdir,'/']; %For long term storage
system(['mkdir -p ',MISCDIR]);
system(['mkdir -p ',LOCALDIR]);
-CMD = ['rsync ', LOCALDIR,'outputs* ',MISCDIR]; disp(CMD);
+% CMD = ['rsync ', LOCALDIR,'outputs* ',MISCDIR]; disp(CMD);
% system(CMD);
% Load outputs from jobnummin up to jobnummax
data = compile_results(MISCDIR,JOBNUMMIN,JOBNUMMAX); %Compile the results from first output found to JOBNUMMAX if existing
data.localdir = LOCALDIR;
data.FIGDIR = LOCALDIR;
-
+data.folder = LOCALDIR;
+data.CODENAME = 'GYACOMO';
%% PLOTS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
default_plots_options
disp('Plots')
@@ -34,12 +35,12 @@ disp([num2str(data.NU_EVOL(i_)),' ',num2str(data.NU_EVOL(i_+1))])
options.TAVG_0 = data.TJOB_SE(i_)+600;%0.4*data.Ts3D(end);
options.TAVG_1 = data.TJOB_SE(i_+1);%0.9*data.Ts3D(end); % Averaging times duration
options.NCUT = 4; % Number of cuts for averaging and error estimation
-options.NMVA = 50; % Moving average for time traces
+options.NMVA = 1; % Moving average for time traces
% options.ST_FIELD = '\Gamma_x'; % chose your field to plot in spacetime diag (e.g \phi,v_x,G_x)
options.ST_FIELD = '\phi'; % chose your field to plot in spacetime diag (e.g \phi,v_x,G_x)
options.INTERP = 0;
options.RESOLUTION = 256;
-fig = plot_radial_transport_and_spacetime(data,options);
+fig = plot_radial_transport_and_spacetime(data,options,'GYACOMO');
% save_figure(data,fig,'.png')
end
@@ -62,21 +63,21 @@ if 0
% Options
options.INTERP = 1;
options.POLARPLOT = 0;
-% options.NAME = '\phi';
+options.NAME = '\phi';
% options.NAME = '\omega_z';
% options.NAME = 'N_i^{00}';
% options.NAME = 'v_x';
% options.NAME = 'n_i^{NZ}';
% options.NAME = '\Gamma_x';
-options.NAME = 'n_e';
+% options.NAME = 'n_i';
% options.NAME = 'n_i-n_e';
options.PLAN = 'xy';
-% options.NAME = 'f_i';
-% options.PLAN = 'sx';
+options.NAME = 'f_i';
+options.PLAN = 'sx';
options.COMP = 'avg';
% options.TIME = data.Ts5D(end-30:end);
% options.TIME = data.Ts3D;
-options.TIME = [00:1:9000];
+options.TIME = [0:10000];
data.EPS = 0.1;
data.a = data.EPS * 2000;
options.RESOLUTION = 256;
@@ -86,7 +87,7 @@ end
if 0
%% fields snapshots
% Options
-options.INTERP = 1;
+options.INTERP = 0;
options.POLARPLOT = 0;
options.AXISEQUAL = 0;
options.NORMALIZE = 0;
@@ -133,13 +134,13 @@ options.XPERP = linspace( 0,sqrt(6),16).^2;
% options.SPAR = gene_data.vp';
% options.XPERP = gene_data.mu';
options.iz = 'avg';
-options.T = [100];
+options.T = [0.5:0.1:1]*data.Ts3D(end);
% options.PLT_FCT = 'contour';
% options.PLT_FCT = 'contourf';
options.PLT_FCT = 'surfvv';
options.ONED = 0;
options.non_adiab = 0;
-options.SPECIE = 'i';
+options.SPECIES = 'i';
options.RMS = 1; % Root mean square i.e. sqrt(sum_k|f_k|^2) as in Gene
fig = plot_fa(data,options);
% save_figure(data,fig,'.png')
@@ -149,7 +150,7 @@ if 0
%% Hermite-Laguerre spectrum
% options.TIME = 'avg';
options.P2J = 0;
-options.ST = 0;
+options.ST = 1;
options.PLOT_TYPE = 'space-time';
options.NORMALIZED = 0;
options.JOBNUM = 0;
@@ -171,7 +172,7 @@ options.NAME ='\Gamma_x';
options.PLAN = 'kxky';
options.COMPZ = 'avg';
options.OK = 0;
-options.COMPXY = '2D'; % avg/sum/max/zero/ 2D plot otherwise
+options.COMPXY = 'avg';%'2D'; % avg/sum/max/zero/ 2D plot otherwise
options.COMPT = 'avg';
options.PLOT = 'semilogy';
fig = spectrum_1D(data,options);
@@ -198,13 +199,14 @@ end
if 0
%% Mode evolution
options.NORMALIZED = 1;
-options.TIME = [6000:9000];
-options.KX_TW = [6000:9000]; %kx Z modes time window
-options.KY_TW = [6000:9000]; %ky Growth rate time window
+options.TIME = [000:9000];
+options.KX_TW = [25 55]; %kx Growth rate time window
+options.KY_TW = [0 20]; %ky Growth rate time window
options.NMA = 1;
options.NMODES = 800;
options.iz = 'avg'; % avg or index
options.ik = 1; % sum, max or index
+options.fftz.flag = 0;
fig = mode_growth_meter(data,options);
save_figure(data,fig,'.png')
end
diff --git a/wk/header_3D_results.m b/wk/header_3D_results.m
index 34e192e0..9a860e05 100644
--- a/wk/header_3D_results.m
+++ b/wk/header_3D_results.m
@@ -31,7 +31,7 @@ gyacomodir = '/home/ahoffman/gyacomo/';
%% CBC
% resdir = 'CBC/64x32x16x5x3';
% resdir = 'CBC/64x128x16x5x3';
-resdir = 'CBC/old/128x64x16x5x3';
+% resdir = 'CBC/old/128x64x16x5x3';
% resdir = 'CBC/96x96x16x3x2_Nexc_6';
% resdir = 'CBC/128x96x16x3x2_Nexc_0';
% resdir = 'CBC/old/192x96x24x13x7';
@@ -78,8 +78,21 @@ resdir = 'CBC/old/128x64x16x5x3';
%% CBC Miller
% resdir = 'GCM_CBC/daint/Miller_GX_comparison';
% resdir = 'GCM_CBC/daint/Salpha_GX_comparison';
-%% RK3
-% resdir = '';
+%% Paper 2 simulations
+% resdir = 'paper_2_nonlinear/kT_6.96/CBC_7x4x128x64x24';
+% resdir = 'paper_2_nonlinear/kT_6.96/CBC_3x2x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_6.96/CBC_5x3x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_6.96/CBC_7x4x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_6.96/CBC_9x5x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_11x6x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_3x2x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_5x3x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_7x4x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_9x5x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_11x6x128x64x24_Nexc_5';
+% resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_15x8x128x64x24_Nexc_5';
+resdir = 'paper_2_nonlinear/kT_5.3/CBC_kT_5.3_17x9x128x64x24_Nexc_5';
+
resdir = ['results/',resdir];
JOBNUMMIN = 00; JOBNUMMAX = 10;
run analysis_gyacomo
diff --git a/wk/lin_3Dzpinch.m b/wk/lin_3Dzpinch.m
index 66ae45c1..95353bc9 100644
--- a/wk/lin_3Dzpinch.m
+++ b/wk/lin_3Dzpinch.m
@@ -10,7 +10,7 @@ addpath(genpath([gyacomodir,'matlab/plot'])) % ... add
addpath(genpath([gyacomodir,'matlab/compute'])) % ... add
addpath(genpath([gyacomodir,'matlab/load'])) % ... add% EXECNAME = 'gyacomo_1.0';
SIMID = 'dbg'; % Name of the simulation
-RUN = 0; % To run or just to load
+RUN = 1; % To run or just to load
default_plots_options
% EXECNAME = 'gyacomo_debug';
% EXECNAME = 'gyacomo_alpha';
@@ -20,12 +20,12 @@ EXECNAME = 'gyacomo';
CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% PHYSICAL PARAMETERS
-NU = 0.1; % Collision frequency
+NU = 0.5; % Collision frequency
TAU = 1e-2; % e/i temperature ratio
K_Ne = 0; % ele Density '''
K_Te = 0; % ele Temperature '''
K_Ni = 0; % ion Density gradient drive
-K_Ti = 70; % ion Temperature '''
+K_Ti = 200; % ion Temperature '''
SIGMA_E = 0.0233380; % mass ratio sqrt(m_a/m_i) (correct = 0.0233380)
KIN_E = 0; % 1: kinetic electrons, 2: adiabatic electrons
BETA = 0.0; % electron plasma beta
@@ -36,12 +36,12 @@ PMAXE = P; % Hermite basis size of electrons
JMAXE = J; % Laguerre "
PMAXI = P; % " ions
JMAXI = J; % "
-NX = 20; % real space x-gridpoints
-NY = 20; % '' y-gridpoints
+NX = 2; % real space x-gridpoints
+NY = 40; % '' y-gridpoints
LX = 2*pi/0.4; % Size of the squared frequency domain
-LY = 2*pi/0.4; % Size of the squared frequency domain
-NZ = 48; % number of perpendicular planes (parallel grid)
-NPOL = 1;
+LY = 2*pi/0.2; % Size of the squared frequency domain
+NZ = 64; % number of perpendicular planes (parallel grid)
+NPOL = 2;
SG = 0; % Staggered z grids option
NEXC = 1; % To extend Lx if needed (Lx = Nexc/(kymin*shear))
%% GEOMETRY
@@ -57,8 +57,8 @@ PARALLEL_BC = 'dirichlet'; %'dirichlet','periodic','shearless','disconnected'
% PARALLEL_BC = 'periodic'; %'dirichlet','periodic','shearless','disconnected'
SHIFT_Y = 0.0;
%% TIME PARMETERS
-TMAX = 50; % Maximal time unit
-DT = 5e-2; % Time step
+TMAX = 20; % Maximal time unit
+DT = 1e-2; % Time step
SPS0D = 1; % Sampling per time unit for 2D arrays
SPS2D = -1; % Sampling per time unit for 2D arrays
SPS3D = 1; % Sampling per time unit for 2D arrays
@@ -69,7 +69,7 @@ JOB2LOAD= -1;
LINEARITY = 'linear'; % activate non-linearity (is cancelled if KXEQ0 = 1)
% Collision operator
% (LB:L.Bernstein, DG:Dougherty, SG:Sugama, LR: Lorentz, LD: Landau)
-CO = 'DG';
+CO = 'SG';
GKCO = 0; % gyrokinetic operator
ABCO = 1; % interspecies collisions
INIT_ZF = 0; ZF_AMP = 0.0;
@@ -93,7 +93,7 @@ INIT_BLOB = 0; WIPE_TURB = 0; ACT_ON_MODES = 0;
MU_X = MU; %
MU_Y = MU; %
N_HD = 4;
-MU_Z = 1.0; %
+MU_Z = 0.0; %
MU_P = 0.0; %
MU_J = 0.0; %
LAMBDAD = 0.0;
@@ -101,7 +101,7 @@ NOISE0 = 1.0e-5; % Init noise amplitude
BCKGD0 = 0.0; % Init background
GRADB = 0.1;
CURVB = 0.1;
-COLL_KCUT = 1000;
+COLL_KCUT = 1.5;
%%-------------------------------------------------------------------------
%% RUN
setup
@@ -168,14 +168,14 @@ end
if 1
%% linear growth rate for 3D Zpinch (kz fourier transform)
trange = [0.5 1]*data.Ts3D(end);
-options.INTERP = 0;
+options.INTERP = 1;
options.kxky = 0;
options.kzkx = 0;
options.kzky = 1;
[lg, fig] = compute_3D_zpinch_growth_rate(data,trange,options);
end
if 0
-%% Mode evolution
+%% ES Mode evolution
options.NORMALIZED = 0;
options.K2PLOT = 1;
options.TIME = [0:1000];
@@ -184,6 +184,8 @@ options.NMA = 1;
options.NMODES = 32;
options.iz = 'avg';
options.ik = 1;
+options.fftz.flag = 1;
+options.fftz.ky = 1.5; options.fftz.kx = 0;
fig = mode_growth_meter(data,options);
save_figure(data,fig,'.png')
end
--
GitLab