Skip to content
Snippets Groups Projects
Commit 40660cf3 authored by Antoine Cyril David Hoffmann's avatar Antoine Cyril David Hoffmann
Browse files

scripts update

parent 6403a434
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
matlab/compile_results.m
+ 14
− 4
View file @ 40660cf3
......@@ -94,15 +94,21 @@ while(CONTINUE)
PHI_ = cat(4,PHI_,PHI);
end
if W_NA00
if KIN_E
Ne00_ = cat(4,Ne00_,Ne00);
end
Ni00_ = cat(4,Ni00_,Ni00);
Ne00_ = cat(4,Ne00_,Ne00);
end
if W_DENS
if KIN_E
DENS_E_ = cat(4,DENS_E_,DENS_E);
end
DENS_I_ = cat(4,DENS_I_,DENS_I);
end
if W_TEMP
TEMP_E_ = cat(4,TEMP_E_,TEMP_E);
if KIN_E
TEMP_E_ = cat(4,TEMP_E_,TEMP_E);
end
TEMP_I_ = cat(4,TEMP_I_,TEMP_I);
end
if W_NAPJ || W_SAPJ
......@@ -110,11 +116,15 @@ while(CONTINUE)
end
if W_NAPJ
Nipj_ = cat(6,Nipj_,Nipj);
Nepj_ = cat(6,Nepj_,Nepj);
if KIN_E
Nepj_ = cat(6,Nepj_,Nepj);
end
end
if W_SAPJ
Sipj_ = cat(6,Sipj_,Sipj);
Sepj_ = cat(6,Sepj_,Sepj);
if KIN_E
Sepj_ = cat(6,Sepj_,Sepj);
end
end
% Evolution of simulation parameters
......
matlab/load_results.m
+ 2
− 2
View file @ 40660cf3
......@@ -14,8 +14,8 @@ W_NAPJ = strcmp(h5readatt(filename,'/data/input','write_Napj' ),'y');
W_SAPJ = strcmp(h5readatt(filename,'/data/input','write_Sapj' ),'y');
W_DENS = strcmp(h5readatt(filename,'/data/input','write_dens' ),'y');
W_TEMP = strcmp(h5readatt(filename,'/data/input','write_temp' ),'y');
% KIN_E = strcmp(h5readatt(filename,'/data/input', 'KIN_E' ),'y');
KIN_E = 1;
KIN_E = strcmp(h5readatt(filename,'/data/input', 'KIN_E' ),'y');
% KIN_E = 1;
if W_GAMMA
......
matlab/post_processing.m
+ 21
− 11
View file @ 40660cf3
......@@ -52,8 +52,22 @@ dx2phi = zeros(Nx,Ny,Nz,Ns3D); % "
for it = 1:numel(Ts3D)
for iz = 1:numel(z)
NE_ = Ne00(:,:,iz,it); NI_ = Ni00(:,:,iz,it); PH_ = PHI(:,:,iz,it);
ne00 (:,:,iz,it) = real(fftshift(ifft2((NE_),Nx,Ny)));
if KIN_E
NE_ = Ne00(:,:,iz,it);
ne00 (:,:,iz,it) = real(fftshift(ifft2((NE_),Nx,Ny)));
if(W_DENS)
DENS_E_ = DENS_E(:,:,iz,it);
dens_e (:,:,iz,it) = real(fftshift(ifft2((DENS_E_),Nx,Ny)));
Z_n_e (:,:,iz,it) = real(fftshift(ifft2((DENS_E_.*(KY==0)),Nx,Ny)));
end
if(W_TEMP)
TEMP_E_ = TEMP_E(:,:,iz,it);
dyTe(:,:,iz,it) = real(fftshift(ifft2(1i*KY.*(TEMP_E_),Nx,Ny)));
temp_e (:,:,iz,it) = real(fftshift(ifft2((TEMP_E_),Nx,Ny)));
Z_T_e (:,:,iz,it) = real(fftshift(ifft2((TEMP_E_.*(KY==0)),Nx,Ny)));
end
end
NI_ = Ni00(:,:,iz,it); PH_ = PHI(:,:,iz,it);
ni00 (:,:,iz,it) = real(fftshift(ifft2((NI_),Nx,Ny)));
phi (:,:,iz,it) = real(fftshift(ifft2((PH_),Nx,Ny)));
Z_phi (:,:,iz,it) = real(fftshift(ifft2((PH_.*(KY==0)),Nx,Ny)));
......@@ -61,20 +75,15 @@ for it = 1:numel(Ts3D)
dx2phi(:,:,iz,it) = real(fftshift(ifft2(-KX.^2.*(PH_),Nx,Ny)));
dyphi (:,:,iz,it) = real(fftshift(ifft2(1i*KY.*(PH_),Nx,Ny)));
if(W_DENS)
DENS_E_ = DENS_E(:,:,iz,it); DENS_I_ = DENS_I(:,:,iz,it);
DENS_I_ = DENS_I(:,:,iz,it);
dyni (:,:,iz,it) = real(fftshift(ifft2(1i*KY.*(DENS_I_),Nx,Ny)));
dens_e (:,:,iz,it) = real(fftshift(ifft2((DENS_E_),Nx,Ny)));
dens_i (:,:,iz,it) = real(fftshift(ifft2((DENS_I_),Nx,Ny)));
Z_n_e (:,:,iz,it) = real(fftshift(ifft2((DENS_E_.*(KY==0)),Nx,Ny)));
Z_n_i (:,:,iz,it) = real(fftshift(ifft2((DENS_I_.*(KY==0)),Nx,Ny)));
end
if(W_TEMP)
TEMP_E_ = TEMP_E(:,:,iz,it); TEMP_I_ = TEMP_I(:,:,iz,it);
dyTe(:,:,iz,it) = real(fftshift(ifft2(1i*KY.*(TEMP_E_),Nx,Ny)));
TEMP_I_ = TEMP_I(:,:,iz,it);
dyTi(:,:,iz,it) = real(fftshift(ifft2(1i*KY.*(TEMP_I_),Nx,Ny)));
temp_e (:,:,iz,it) = real(fftshift(ifft2((TEMP_E_),Nx,Ny)));
temp_i (:,:,iz,it) = real(fftshift(ifft2((TEMP_I_),Nx,Ny)));
Z_T_e (:,:,iz,it) = real(fftshift(ifft2((TEMP_E_.*(KY==0)),Nx,Ny)));
Z_T_i (:,:,iz,it) = real(fftshift(ifft2((TEMP_I_.*(KY==0)),Nx,Ny)));
end
end
......@@ -119,7 +128,6 @@ phi_kp_t = zeros(numel(kp_ip),Nz,Ns3D);
%
for it = 1:numel(Ts3D) % Loop over 2D aX_XYays
for iz = 1:numel(z)
NE_ = Ne00(:,:,iz,it); NI_ = Ni00(:,:,iz,it); PH_ = PHI(:,:,iz,it);
phi_maxx_maxy(iz,it) = max( max(squeeze(phi(:,:,iz,it))));
phi_avgx_maxy(iz,it) = max(mean(squeeze(phi(:,:,iz,it))));
phi_maxx_avgy(iz,it) = mean( max(squeeze(phi(:,:,iz,it))));
......@@ -147,7 +155,9 @@ end
%
for it = 1:numel(Ts5D) % Loop over 5D aX_XYays
[~, it2D] = min(abs(Ts3D-Ts5D(it)));
Ne_norm(:,:,it)= sum(sum(abs(Nepj(:,:,:,:,it)),3),4)/Nkx/Nky;
if KIN_E
Ne_norm(:,:,it)= sum(sum(abs(Nepj(:,:,:,:,it)),3),4)/Nkx/Nky;
end
Ni_norm(:,:,it)= sum(sum(abs(Nipj(:,:,:,:,it)),3),4)/Nkx/Nky;
end
......
matlab/setup.m
+ 1
− 1
View file @ 40660cf3
......@@ -13,7 +13,7 @@ GRID.Nz = NZ; % z resolution
GRID.q0 = Q0; % q factor
GRID.shear = SHEAR; % shear
GRID.eps = EPS; % inverse aspect ratio
if SG; GRID.SG = '.true.'; else; GRID.SG = '.false.';end;
% Model parameters
MODEL.CO = CO; % Collision operator (0 : L.Bernstein, -1 : Full Coulomb, -2 : Dougherty)
MODEL.CLOS = CLOS;
......
matlab/write_fort90.m
+ 1
− 0
View file @ 40660cf3
......@@ -23,6 +23,7 @@ fprintf(fid,[' Nz = ', num2str(GRID.Nz),'\n']);
fprintf(fid,[' q0 = ', num2str(GRID.q0),'\n']);
fprintf(fid,[' shear = ', num2str(GRID.shear),'\n']);
fprintf(fid,[' eps = ', num2str(GRID.eps),'\n']);
fprintf(fid,[' SG = ', GRID.SG,'\n']);
fprintf(fid,'/\n');
fprintf(fid,'&OUTPUT_PAR\n');
......
wk/analysis_3D.m
+ 21
− 19
View file @ 40660cf3
......@@ -4,8 +4,8 @@ outfile ='';
%% Directory of the simulation
if 1% Local results
outfile ='';
outfile ='';
outfile ='fluxtube_salphaB_s0/50x100x20_5x3_L_300_q0_2.7_e_0.18_kN_2.22_kT_6_nu_1e-01_DGGK_lin';
outfile ='fluxtube_salphaB_s0/50x100x20_5x3_L_300_q0_2.7_e_0.18_kN_2.22_kT_6_nu_1e-01_DGGK_adiabe';
% outfile ='fluxtube_salphaB_s0/50x100x20_5x3_L_300_q0_2.7_e_0.18_kN_2.22_kT_6_nu_1e-01_DGGK_lin';
% outfile ='fluxtube_salphaB_s0/64x64x16_5x3_L_200_q0_2.7_e_0.18_kN_2.22_kT_6_nu_1e-01_DGGK';
% outfile ='simulation_A/1024x256_3x2_L_120_kN_1.6667_nu_1e-01_DGGK';
% outfile ='Linear_Device/64x64x20_5x2_Lx_20_Ly_150_q0_25_kN_0.24_kT_0.03_nu_1e-02_DGGK';
......@@ -15,16 +15,15 @@ outfile ='fluxtube_salphaB_s0/50x100x20_5x3_L_300_q0_2.7_e_0.18_kN_2.22_kT_6_nu_
CMD = ['cp ', BASIC.RESDIR,'outputs* ',BASIC.MISCDIR]; disp(CMD);
system(CMD);
else% Marconi results
outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A/300x300_L_120_P_8_J_4_eta_0.6_nu_1e-01_PAGK_mu_0e+00/out.txt';
% outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_B/300x150_L_120_P_8_J_4_eta_0.6_nu_5e-01_SGGK_mu_0e+00/out.txt';
% outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A/500x500_L_120_P_4_J_2_eta_0.6_nu_1e-01_DGGK_mu_0e+00/out.txt';
outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A/300x200_L_200_P_8_J_4_eta_0.6_nu_1e-01_PAGK_mu_0e+00/out.txt';
% outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A/300x300_L_120_P_8_J_4_eta_0.6_nu_1e-01_PAGK_mu_0e+00/out.txt';
BASIC.RESDIR = ['../',outfile(46:end-8),'/'];
BASIC.MISCDIR = ['/misc/HeLaZ_outputs/',outfile(46:end-8),'/'];
end
%% Load the results
% Load outputs from jobnummin up to jobnummax
JOBNUMMIN = 00; JOBNUMMAX = 20;
JOBNUMMIN = 00; JOBNUMMAX = 00;
% JOBNUMMIN = 07; JOBNUMMAX = 20; % For CO damping sim A
compile_results %Compile the results from first output found to JOBNUMMAX if existing
......@@ -46,7 +45,7 @@ end
if 0
%% MOVIES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Options
t0 =000; iz = 2; ix = 1; iy = 1;
t0 =000; iz = 1; ix = 1; iy = 1;
skip_ =1; FPS = 30; DELAY = 1/FPS;
[~, it03D] = min(abs(Ts3D-t0)); FRAMES_3D = it03D:skip_:numel(Ts3D);
[~, it05D] = min(abs(Ts5D-t0)); FRAMES_5D = it05D:skip_:numel(Ts5D);
......@@ -64,12 +63,12 @@ FIELD = phi; NAME = 'phi'; FIELDNAME = '\phi';
% FIELD = Gamma_x; NAME = 'Gamma_x'; FIELDNAME = '\Gamma_x';
% Sliced
plt = @(x) real(x(ix, :, :,:)); X = Y_YZ; Y = Z_YZ; XNAME = 'y'; YNAME = 'z';
% plt = @(x) real(x(ix, :, :,:)); X = Y_YZ; Y = Z_YZ; XNAME = 'y'; YNAME = 'z';
% plt = @(x) real(x( :,iy, :,:)); X = X_XZ; Y = Z_XZ; XNAME = 'x'; YNAME = 'z';
% plt = @(x) real(x( :, :,iz,:)); X = X_XY; Y = Y_XY; XNAME = 'x'; YNAME = 'y';
plt = @(x) real(x( :, :,iz,:)); X = X_XY; Y = Y_XY; XNAME = 'x'; YNAME = 'y';
% K-space
% FIELD = PHI; NAME = 'PHI'; FIELDNAME = '\tilde \phi';
% FIELD = PHI.*(KY~=0); NAME = 'NZPHI'; FIELDNAME = '\tilde \phi_{k_y\neq0}';
% FIELD = Ne00; NAME = 'Ne00'; FIELDNAME = 'N_e^{00}';
% FIELD = Ni00; NAME = 'Ni00'; FIELDNAME = 'N_i^{00}';
% plt = @(x) fftshift((abs(x( :, :,1,:))),2); X = fftshift(KX,2); Y = fftshift(KY,2); XNAME = 'k_x'; YNAME = 'k_y';
......@@ -145,13 +144,14 @@ if 0
% Chose the field to plot
% FIELD = Ni00; FNAME = 'Ni00'; FIELDLTX = 'N_i^{00}';
% FIELD = Ne00; FNAME = 'Ne00'; FIELDLTX = 'N_e^{00}'
FIELD = PHI; FNAME = 'PHI'; FIELDLTX = '\tilde\phi';
FIELD = Ne00; FNAME = 'Ne00'; FIELDLTX = 'N_e^{00}'
% FIELD = PHI; FNAME = 'PHI'; FIELDLTX = '\tilde\phi';
% FIELD_ = fft2(Gamma_x); FIELD = FIELD_(1:Nx/2+1,:,:,:); FNAME = 'Gamma_x'; FIELDLTX = '\tilde\Gamma_x';
% FIELD_ = fft2(dens_e); FIELD = FIELD_(1:Nx/2+1,:,:,:); FNAME = 'FFT_Dens_e'; FIELDLTX = '\tilde n_e';
% Chose when to plot it
tf = [0 15 27 28 30];
% tf = [0 15 27 28 30];
tf = 200:200:1200;
% tf = 8000;
% Planar plot: choose a plane to plot at x0/y0/z0 coordinates
......@@ -178,7 +178,7 @@ fig = figure; FIGNAME = [FNAME,TNAME,'_snaps','_',PARAMS]; set(gcf, 'Position',
subplot(1,numel(tf),i_)
DATA = plt_2(squeeze(plt(FIELD,it)));
pclr = pcolor(plt_x(X),plt_y(Y),plt_z(DATA)); set(pclr, 'edgecolor','none');pbaspect([0.5 1 1])
% caxis([0 1]*5e3);
caxis([0 1]*1e3);
% caxis([-1 1]*5);
xlabel(latexize(XNAME)); ylabel(latexize(YNAME));
if(i_ > 1); set(gca,'ytick',[]); end;
......@@ -217,7 +217,8 @@ end
if 0
%% zonal vs nonzonal energies for phi(t)
it0 = 01; itend = Ns3D;
t0 = 200; [~, it0] = min(abs(Ts3D-t0));
itend = Ns3D;
trange = it0:itend;
pltx = @(x) x;%-x(1);
plty = @(x) x./max(squeeze(x));
......@@ -237,9 +238,10 @@ subplot(121)
xlabel('$t c_s/R$'); grid on;% xlim([0 500]);
subplot(122)
plot(plty(Ephi_Z(trange)),plty(Ephi_NZ_kgt0(trange)));
% plot(plty(Ephi_Z(trange)),plty(Ephi_NZ_kgt0(trange)));
plot3(plty(Ephi_Z(trange)),plty(Ephi_NZ_kgt0(trange)),Ts3D(trange));
title('Phase space'); legend(CONAME)
xlabel('$E_Z$'); ylabel('$E_{NZ}$'); grid on;% xlim([0 500]);
xlabel('$E_Z$'); ylabel('$E_{NZ}$'); zlabel('time'); grid on;% xlim([0 500]);
end
if 0
......@@ -281,7 +283,7 @@ FIELD = Ne00.*conj(Ne00); FNAME = 'Ne002'; FIELDLTX = '|N_e^{00}|^2'
% FIELD_ = fft2(dens_e); FIELD = FIELD_(1:Nx/2+1,:,:,:); FNAME = 'FFT_Dens_e'; FIELDLTX = '\tilde n_e';
% Chose when to plot it
tf = 1500:200:2500;
tf = 200:200:1200;
% tf = 8000;
% Sliced
......@@ -308,7 +310,7 @@ end
if 0
%% |phi_k|^2 spectra (Kobayashi 2015 fig 3)
% tstart = 0.8*Ts3D(end); tend = Ts3D(end); % Time window
tstart = 1000; tend = 4000;
tstart = 0000; tend = 1000;
% tstart = 10000; tend = 12000;
% Chose the field to plot
% FIELD = Ni00; FNAME = 'Ni00'; FIELDLTX = 'N_i^{00}';
......
wk/benchmark_HeLaZ_molix.m
+ 2
− 1
View file @ 40660cf3
......@@ -17,7 +17,8 @@ system(['cd ',SIMDIR,';',' ./../../../bin/helaz 0; cd ../../../wk'])
% Compare the results with molix at a given time
time_2_plot = 5.0;
[Y_,X_]=molix_plot_phi([SIMDIR,'molix_phi.h5'],time_2_plot);
[ PHI, Ts3D, dt3D] = load_3D_data([SIMDIR,'outputs_00.h5'], 'phi');
filename = [SIMDIR,'/outputs_00.h5'];
[ PHI, Ts3D, dt3D] = load_3D_data(filename, 'phi');
[Pe, Je, Pi, Ji, kx, ky, z] = load_grid_data(filename);
plot_phi_ballooning; hold on
......
wk/benchmark_linear_1D_entropy_mode.m
+ 6
− 2
View file @ 40660cf3
......@@ -5,8 +5,12 @@ SIMDIR = '../results/benchmarks/1D_linear_entropy_mode/';
cd ..
system('make');
cd wk
% Run HeLaZ
system(['cd ',SIMDIR,';',' ./../../../bin/helaz 0; cd ../../../wk'])
% Run HeLaZ in sequential
system(['cd ',SIMDIR,';',' ./../../../bin/helaz 0; cd ../../../wk']);
% Run HeLaZ in parallel (discrepencies can occur at marginal growth rate)
% since the random seed will not be the same)
% system(['cd ',SIMDIR,';',' mpirun -np 6 ./../../../bin/helaz 1 6 0; cd ../../../wk'])
% system(['cd ',SIMDIR,';',' mpirun -np 6 ./../../../bin/helaz 2 3 0; cd ../../../wk'])
%compute growth rate
%%
filename = [SIMDIR,'/outputs_00.h5'];
......
wk/local_run.m
+ 10
− 9
View file @ 40660cf3
......@@ -6,11 +6,11 @@ CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss
%% PHYSICAL PARAMETERS
NU = 0.1; % Collision frequency
K_N = 2.22; % Density gradient drive
K_T = 6.0; % Temperature '''
K_T = 8.0; % Temperature '''
K_E = 0.00; % Electrostat gradient
SIGMA_E = 0.05196; % mass ratio sqrt(m_a/m_i) (correct = 0.0233380)
NU_HYP = 0.0;
KIN_E = 1; % Kinetic (1) or adiabatic (2) electron model
KIN_E = 0; % Kinetic (1) or adiabatic (2) electron model
%% GRID PARAMETERS
NX = 50; % Spatial radial resolution ( = 2x radial modes)
LX = 300; % Radial window size
......@@ -23,11 +23,12 @@ J = 2;
Q0 = 2.7; % safety factor
SHEAR = 0.0; % magnetic shear
EPS = 0.18; % inverse aspect ratio
GRADB = 1.0; % Magnetic gradient
CURVB = 1.0; % Magnetic curvature
GRADB = 1.0; % Magnetic gradient
CURVB = 1.0; % Magnetic curvature
SG = 1; % Staggered z grids option
%% TIME PARAMETERS
TMAX = 10; % Maximal time unit
DT = 5e-3; % Time step
TMAX = 50; % Maximal time unit
DT = 2e-2; % Time step
SPS0D = 1; % Sampling per time unit for profiler
SPS2D = 1; % Sampling per time unit for 2D arrays
SPS3D = 5; % Sampling per time unit for 3D arrays
......@@ -39,13 +40,13 @@ JOB2LOAD= -1;
% (0 : L.Bernstein, 1 : Dougherty, 2: Sugama, 3 : Pitch angle ; 4 : Coulomb; +/- for GK/DK)
CO = 1;
CLOS = 0; % Closure model (0: =0 truncation)
NL_CLOS = 0; % nonlinear closure model (-2: nmax = jmax, -1: nmax = jmax-j, >=0 : nmax = NL_CLOS)
NL_CLOS = -1; % nonlinear closure model (-2: nmax = jmax, -1: nmax = jmax-j, >=0 : nmax = NL_CLOS)
SIMID = 'fluxtube_salphaB_s0'; % Name of the simulation
% SIMID = 'simulation_A'; % Name of the simulation
% SIMID = ['v3.0_P_',num2str(P),'_J_',num2str(J)]; % Name of the simulation
NON_LIN = 0; % activate non-linearity (is cancelled if KXEQ0 = 1)
NON_LIN = 1; % activate non-linearity (is cancelled if KXEQ0 = 1)
% INIT options
INIT_PHI= 1; % Start simulation with a noisy phi (0= noisy moments 00)
INIT_PHI= 0; % Start simulation with a noisy phi (0= noisy moments 00)
INIT_ZF = 0; ZF_AMP = 0.0;
INIT_BLOB = 0; WIPE_TURB = 0; WIPE_ZF = 0;
%% OUTPUTS
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment