From ec833a8fcee57a0f1468b01843c1bf27318f2b7f Mon Sep 17 00:00:00 2001
From: Antoine Hoffmann <antoine.hoffmann@epfl.ch>
Date: Thu, 23 Jun 2022 16:11:39 +0200
Subject: [PATCH] scripts update
---
matlab/compile_results.m | 6 ++--
matlab/compute/compute_fluxtube_growth_rate.m | 6 ++--
wk/analysis_HeLaZ.m | 4 +--
wk/cbc_qx_results.m | 16 ++++-----
wk/header_3D_results.m | 15 +++++---
wk/heat_flux_postproc.m | 11 +++---
wk/quick_run.m | 34 +++++++++----------
7 files changed, 50 insertions(+), 42 deletions(-)
diff --git a/matlab/compile_results.m b/matlab/compile_results.m
index 66a18e6e..00abfff9 100644
--- a/matlab/compile_results.m
+++ b/matlab/compile_results.m
@@ -72,9 +72,11 @@ while(CONTINUE)
tmp = Nipj_; sz = size(tmp);
Nipj_ = zeros(cat(1,[Pi_new,Ji_new]',sz(3:end)')');
Nipj_(1:Pi_old,1:Ji_old,:,:,:,:) = tmp;
+ if(KIN_E)
tmp = Nepj_; sz = size(tmp);
-% Nepj_ = zeros(cat(1,[Pe_new,Je_new]',sz(3:end)')');
-% Nepj_(1:Pe_old,1:Je_old,:,:,:,:) = tmp;
+ Nepj_ = zeros(cat(1,[Pe_new,Je_new]',sz(3:end)')');
+ Nepj_(1:Pe_old,1:Je_old,:,:,:,:) = tmp;
+ end
end
if W_SAPJ
tmp = Sipj_; sz = size(tmp);
diff --git a/matlab/compute/compute_fluxtube_growth_rate.m b/matlab/compute/compute_fluxtube_growth_rate.m
index fad8aecf..9d93e421 100644
--- a/matlab/compute/compute_fluxtube_growth_rate.m
+++ b/matlab/compute/compute_fluxtube_growth_rate.m
@@ -41,9 +41,9 @@ linear_gr.ce = abs(ce);
linear_gr.ky = kys;
if PLOT >0
figure
- plot(linear_gr.ky,linear_gr.g_ky(:,end),'DisplayName','$Re(\omega_{k_y})$'); hold on;
- plot(linear_gr.ky,linear_gr.w_ky(:,end),'DisplayName','$Im(\omega_{k_y})$'); hold on;
- plot(linear_gr.ky,linear_gr.ce (:,end),'DisplayName','$\epsilon$'); hold on;
+ plot(linear_gr.ky,linear_gr.g_ky(:,end),'-o','DisplayName','$Re(\omega_{k_y})$'); hold on;
+ plot(linear_gr.ky,linear_gr.w_ky(:,end),'-o','DisplayName','$Im(\omega_{k_y})$'); hold on;
+ plot(linear_gr.ky,linear_gr.ce (:,end),'-o','DisplayName','$\epsilon$'); hold on;
xlim([min(linear_gr.ky) max(linear_gr.ky)]);
xlabel('$k_y$');
legend('show');
diff --git a/wk/analysis_HeLaZ.m b/wk/analysis_HeLaZ.m
index 21d2cb6e..4df222ad 100644
--- a/wk/analysis_HeLaZ.m
+++ b/wk/analysis_HeLaZ.m
@@ -50,7 +50,7 @@ options.NAME = '\phi';
% options.NAME = 'n_i^{NZ}';
% options.NAME = '\Gamma_x';
% options.NAME = 'n_i';
-options.PLAN = 'kxky';
+options.PLAN = 'xz';
% options.NAME = 'f_i';
% options.PLAN = 'sx';
options.COMP = 'avg';
@@ -78,7 +78,7 @@ options.PLAN = 'kxky';
% options.NAME = 'f_i';
% options.PLAN = 'sx';
options.COMP = 'avg';
-options.TIME = [300 400 600];
+options.TIME = [100];
data.a = data.EPS * 2e3;
fig = photomaton(data,options);
% save_figure(data,fig)
diff --git a/wk/cbc_qx_results.m b/wk/cbc_qx_results.m
index 5f815bc0..3fa9a986 100644
--- a/wk/cbc_qx_results.m
+++ b/wk/cbc_qx_results.m
@@ -23,19 +23,19 @@ xlabel('CBC drive [\%]'); ylabel('Radial Heat Flux $Q_x^\infty$');
KN = 2.22;
KT = [1.00 0.90 0.80 0.70 0.60 0.50]*6.96;
-gm42 = [32.2 16.8 0.00 0.00 1.74 0.00];
-gm42_err = [05.2 02.8 0.00 0.00 0.53 0.00];
-gm84 = [0.00 0.00 0.00 0.00 0.00 0.00];
-gm84_err = [0.00 0.00 0.00 0.00 0.00 0.00];
+gm42 = [32.2 18.8 10.5 5.89 1.74 0.00];
+gm42_err = [05.2 03.8 2.10 1.66 0.53 0.00];
+gm84 = [0.00 13.2 7.66 3.20 0.00 0.00];
+gm84_err = [0.00 2.79 1.93 0.98 0.00 0.00];
gne = [0.00 0.00 0.00 0.00 0.00 0.00];
gne_err = [0.00 0.00 0.00 0.00 0.00 0.00];
figure
-errorbar(KT./KN,gm42,gm42_err,'o-', 'LineWidth',1.5); hold on;
-errorbar(KT./KN,gm84,gm84_err,'o-', 'LineWidth',1.5);
-errorbar(KT./KN, gne, gne_err,'x-k','LineWidth',1.5);
+errorbar(KT,gm42./KT,gm42_err./KT,'o-', 'LineWidth',1.5); hold on;
+errorbar(KT,gm84./KT,gm84_err./KT,'o-', 'LineWidth',1.5);
+errorbar(KT, gne./KT, gne_err./KT,'x-k','LineWidth',1.5);
% set(gca, 'YScale', 'log')
legend('GM (4,2)','GM (8,4)','Gene')
-xlabel('$\eta=\kappa_T/\kappa_N$'); ylabel('Radial Heat Flux $Q_x^\infty$');
+xlabel('$\eta=\kappa_T/\kappa_N$'); ylabel('$\chi_i$');
diff --git a/wk/header_3D_results.m b/wk/header_3D_results.m
index 5fc6c648..965e457a 100644
--- a/wk/header_3D_results.m
+++ b/wk/header_3D_results.m
@@ -2,21 +2,26 @@
helazdir = '/home/ahoffman/HeLaZ/';
% Directory of the simulation (from results)
% if 1% Local results
-outfile ='';
+% outfile ='volcokas/64x32x16x5x3_kin_e_npol_1';
%% Dimits
% outfile ='shearless_cyclone/128x64x16x5x3_Dim_90';
-outfile ='shearless_cyclone/64x32x16x5x3_Dim_50';
+% outfile ='shearless_cyclone/128x64x16x9x5_Dim_scan/128x64x16x9x5_Dim_60';
+% outfile ='shearless_cyclone/128x64x16x5x3_Dim_scan/128x64x16x5x3_Dim_70';
+% outfile ='shearless_cyclone/64x32x16x5x3_Dim_scan/64x32x16x5x3_Dim_70';
%% AVS
% outfile = 'volcokas/64x32x16x5x3_kin_e_npol_1';
-
-%% Bechmark
+% outfile = 'volcokas/64x32x16x5x3_kin_e_npol_1';
+outfile = 'shearless_cyclone/64x32x160x5x3_CBC_Npol_10_kine';
+%% CBC
% outfile ='shearless_cyclone/64x32x16x5x3_CBC_080';
-% outfile ='shearless_cyclone/64x32x16x5x3_CBC_100';
+% outfile ='shearless_cyclone/64x32x16x5x3_CBC_scan/64x32x16x5x3_CBC_100';
% outfile ='shearless_cyclone/64x32x16x5x3_CBC_120';
% outfile ='shearless_cyclone/64x32x16x9x5_CBC_080';
% outfile ='shearless_cyclone/64x32x16x9x5_CBC_100';
% outfile ='shearless_cyclone/64x32x16x9x5_CBC_120';
+% outfile = 'shearless_cyclone/64x32x16x5x3_CBC_CO/64x32x16x5x3_CBC_LRGK';
+
run analysis_HeLaZ
diff --git a/wk/heat_flux_postproc.m b/wk/heat_flux_postproc.m
index 500ebf9a..1e587e00 100644
--- a/wk/heat_flux_postproc.m
+++ b/wk/heat_flux_postproc.m
@@ -1,7 +1,8 @@
% Hflux_x = 0;
% Hflux_x = 0 * data.Ts5D;
-filename = '/home/ahoffman/HeLaZ/results/shearless_cyclone/64x32x16x5x3_CBC_100/outputs_00.h5';
-kernel_i = h5read(filename,'/data/metric/kernel_i');
+% filename = '/home/ahoffman/HeLaZ/results/shearless_cyclone/64x32x16x5x3_CBC_100/outputs_00.h5';
+filename = '/misc/HeLaZ_outputs/results/shearless_cyclone/64x32x160x5x3_CBC_Npol_10_kine/outputs_00.h5';
+kernel_i = h5read(filename,'/data/metric/kernel_i'); sigma = 0.0738;
Jacobian = h5read(filename,'/data/metric/Jacobian');
STEPS = 1:numel(data.Ts5D);
Hflux_x = 1:numel(STEPS);
@@ -35,9 +36,9 @@ for it = STEPS
correct_op = zeros(Nky,Nkx,Nz);
for in = 1:Nj
n = in-1;
- Kn = squeeze(kernel_i(in,:,:,:,1));
- N2n = squeeze(data.Nipj(3,in ,:,:,:,it5d));
- N0n = squeeze(data.Nipj(1,in ,:,:,:,it5d));
+ Kn = squeeze(kernel_i(in,:,:,:,1))*sigma^(2*n);
+ N2n = squeeze(data.Nepj(3,in ,:,:,:,it5d));
+ N0n = squeeze(data.Nepj(1,in ,:,:,:,it5d));
sumkernmom = sumkernmom + ...
Kn.* (c2n(n) .* N2n + c0n(n) .* N0n);
diff --git a/wk/quick_run.m b/wk/quick_run.m
index 6cb54a7e..efb966ed 100644
--- a/wk/quick_run.m
+++ b/wk/quick_run.m
@@ -13,22 +13,22 @@ EXECNAME = 'helaz3';
CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% PHYSICAL PARAMETERS
-NU = 0.05; % Collision frequency
+NU = 0.1; % Collision frequency
TAU = 1.0; % e/i temperature ratio
-K_N = 0;%2.22; % Density gradient drive
-K_T = 0;%6.96; % Temperature '''
+K_N = 2.22; % Density gradient drive
+K_T = 6.96; % Temperature '''
K_E = 0.0; % Electrostat '''
SIGMA_E = 0.05196152422706632;%0.0233380; % mass ratio sqrt(m_a/m_i) (correct = 0.0233380)
KIN_E = 0; % 1: kinetic electrons, 2: adiabatic electrons
%% GRID PARAMETERS
-PMAXE = 8; % Hermite basis size of electrons
-JMAXE = 4; % Laguerre "
-PMAXI = 8; % " ions
-JMAXI = 4; % "
-NX = 2; % real space x-gridpoints
-NY = 1; % '' y-gridpoints
+PMAXE = 4; % Hermite basis size of electrons
+JMAXE = 2; % Laguerre "
+PMAXI = 4; % " ions
+JMAXI = 2; % "
+NX = 20; % real space x-gridpoints
+NY = 2; % '' y-gridpoints
LX = 100; % Size of the squared frequency domain
-LY = 60; % Size of the squared frequency domain
+LY = 20; % Size of the squared frequency domain
NZ = 32; % number of perpendicular planes (parallel grid)
NPOL = 1;
SG = 0; % Staggered z grids option
@@ -36,11 +36,11 @@ SG = 0; % Staggered z grids option
% GEOMETRY= 'Z-pinch'; % Z-pinch overwrites q0, shear and eps
GEOMETRY= 's-alpha';
Q0 = 1.4; % safety factor
-SHEAR = 0.0; % magnetic shear (Not implemented yet)
+SHEAR = 0.8; % magnetic shear (Not implemented yet)
EPS = 0.18; % inverse aspect ratio
%% TIME PARMETERS
TMAX = 50; % Maximal time unit
-DT = 2e-2; % Time step
+DT = 1e-2; % Time step
SPS0D = 20; % Sampling per time unit for 2D arrays
SPS2D = 0; % Sampling per time unit for 2D arrays
SPS3D = 20; % Sampling per time unit for 2D arrays
@@ -48,7 +48,7 @@ SPS5D = 1; % Sampling per time unit for 5D arrays
SPSCP = 0; % Sampling per time unit for checkpoints
JOB2LOAD= -1;
%% OPTIONS
-SIMID = 'Ros_Hin_test'; % Name of the simulation
+SIMID = 'linear_CBC'; % Name of the simulation
LINEARITY = 'linear'; % activate non-linearity (is cancelled if KXEQ0 = 1)
% Collision operator
% (LB:L.Bernstein, DG:Dougherty, SG:Sugama, LR: Lorentz, LD: Landau)
@@ -59,7 +59,7 @@ INIT_ZF = 0; ZF_AMP = 0.0;
CLOS = 0; % Closure model (0: =0 truncation, 1: gyrofluid closure (p+2j<=Pmax))s
NL_CLOS = 0; % nonlinear closure model (-2:nmax=jmax; -1:nmax=jmax-j; >=0:nmax=NL_CLOS)
KERN = 0; % Kernel model (0 : GK)
-INIT_OPT= 'mom00'; % Start simulation with a noisy mom00/phi/allmom
+INIT_OPT= 'phi'; % Start simulation with a noisy mom00/phi/allmom
%% OUTPUTS
W_DOUBLE = 1;
W_GAMMA = 1; W_HF = 1;
@@ -90,8 +90,8 @@ system(['rm fort*.90']);
% Run linear simulation
if RUN
% system(['cd ../results/',SIMID,'/',PARAMS,'/; time mpirun -np 4 ',HELAZDIR,'bin/',EXECNAME,' 1 4 1 0; cd ../../../wk'])
-% system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 1 ',HELAZDIR,'bin/',EXECNAME,' 1 1 1 0; cd ../../../wk'])
- system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 4 ',HELAZDIR,'bin/',EXECNAME,' 1 1 4 0; cd ../../../wk'])
+ system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 1 ',HELAZDIR,'bin/',EXECNAME,' 1 1 1 0; cd ../../../wk'])
+% system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 4 ',HELAZDIR,'bin/',EXECNAME,' 1 1 4 0; cd ../../../wk'])
% system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ',HELAZDIR,'bin/',EXECNAME,' 1 6 1 0; cd ../../../wk'])
end
@@ -104,7 +104,7 @@ JOBNUMMIN = 00; JOBNUMMAX = 00;
data = compile_results(LOCALDIR,JOBNUMMIN,JOBNUMMAX); %Compile the results from first output found to JOBNUMMAX if existing
%% Short analysis
-if 0
+if 1
%% linear growth rate (adapted for 2D zpinch and fluxtube)
trange = [0.5 1]*data.Ts3D(end);
nplots = 1;
--
GitLab