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

Matlab scripts update

parent aa80c0e9
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 548 additions and 186 deletions
matlab/create_mov.m
+ 2
− 2
View file @ efef5834
......@@ -16,7 +16,7 @@ figure('Color','white','Position', [100, 100, 400, 400]);
pcolor(X,Y,FIELD(:,:,1)); % to set up
colormap gray
axis tight manual % this ensures that getframe() returns a consistent size
if INTERP
if 1
shading interp;
end
in = 1;
......@@ -24,7 +24,7 @@ figure('Color','white','Position', [100, 100, 400, 400]);
for n = FRAMES % loop over selected frames
scale = max(max(abs(FIELD(:,:,n))));
pclr = pcolor(X,Y,FIELD(:,:,n)/scale); % frame plot
if INTERP
if 1
shading interp;
end
set(pclr, 'edgecolor','none'); axis square;
......
matlab/load_marconi.m
+ 10
− 2
View file @ efef5834
......@@ -3,15 +3,23 @@ function [ RESDIR ] = load_marconi( outfilename )
% Detailed explanation goes here
hostfolder = outfilename(1:end-8);
hostfile = [hostfolder,'/out*'];
MISCDIR = ['/misc/HeLaZ_outputs/results/',outfilename(46:end-8),'/'];
MISCDIR = ['/misc/HeLaZ_outputs/',outfilename(46:end-8),'/'];
RESDIR = ['../',outfilename(46:end-8),'/'];
miscfolder = [MISCDIR,'.'];
localfolder = [RESDIR,'.'];
system(['mkdir -p ',miscfolder]);
disp(['mkdir -p ',miscfolder]);
resultfolder = [RESDIR,'.'];
% SCP the output file from marconi to misc folder of SPCPC
CMD = ['scp -r ahoffman@login.marconi.cineca.it:',hostfile,' ',miscfolder];
disp(CMD);
system(CMD);
% Load the fort.90 as well in misc folder
CMD = ['scp -r ahoffman@login.marconi.cineca.it:',hostfolder,'/fort.90',' ',miscfolder];
disp(CMD);
system(CMD);
% Put it also in the result directory
CMD = ['scp -r ahoffman@login.marconi.cineca.it:',hostfolder,'/fort.90',' ',resultfolder];
disp(CMD);
system(CMD);
end
matlab/photomaton.m
+ 60
− 17
View file @ efef5834
if 0
%% Photomaton : real space
% FIELD = ni00; FNAME = 'ni'; XX = RR; YY = ZZ;
% FIELD = phi; FNAME = 'phi'; XX = RR; YY = ZZ;
FIELD = dr2phi; FNAME = 'dr2phi'; XX = RR; YY = ZZ;
% FIELD = fftshift(abs(Ni00),2); FNAME = 'Fni'; XX = fftshift(KR,2); YY = fftshift(KZ,2);
% FIELD = fftshift(abs(PHI),2); FNAME = 'Fphi'; XX = fftshift(KR,2); YY = fftshift(KZ,2);
tf = 150; [~,it1] = min(abs(Ts2D-tf));
tf = 300; [~,it2] = min(abs(Ts2D-tf));
tf = 800; [~,it3] = min(abs(Ts2D-tf));
tf = 1200; [~,it4] = min(abs(Ts2D-tf));
% FIELD = ni00; FNAME = 'ni'; FIELDLTX = '$N_i^{00}$'; XX = RR; YY = ZZ;
% 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 = 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));
subplot(141)
......@@ -35,30 +34,74 @@ plt = @(x) x;%./max(max(x));
colormap gray
xlabel('$r/\rho_s$');ylabel('$z/\rho_s$'); set(gca,'ytick',[]);
title(sprintf('$t c_s/R=%.0f$',Ts2D(it4)));
suptitle(['$\',FNAME,'$, $\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
', $P=',num2str(PMAXI),'$, $J=',num2str(JMAXI),'$']);
suptitle([FIELDLTX,', $\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
', $L=',num2str(L),'$, $N=',num2str(Nr),'$, $(P,J)=(',num2str(PMAXI),',',num2str(JMAXI),')$,',...
' $\mu_{hd}=$',num2str(MU)]);
save_figure
end
if 0
%% Photomaton : real space
% FIELD = Ni00; FNAME = 'fni00'; FIELDLTX = '$\tilde N_i^{00}$'; XX = RR; YY = ZZ;
% FIELD = Ne00; FNAME = 'fne00'; FIELDLTX = '$\tilde N_e^{00}$'; XX = RR; YY = ZZ;
FIELD = PHI; FNAME = 'fphi'; FIELDLTX = '$\tilde \phi$'; XX = RR; YY = ZZ;
tf = 700; [~,it1] = min(abs(Ts2D-tf));
tf = 800; [~,it2] = min(abs(Ts2D-tf));
tf = 900; [~,it3] = min(abs(Ts2D-tf));
tf = 1500; [~,it4] = min(abs(Ts2D-tf));
XLIM = [0,1.0]; YLIM = [-.5;.5];
fig = figure; FIGNAME = [FNAME,'_snaps','_',PARAMS]; set(gcf, 'Position', [100, 100, 1500, 400])
plt = @(x) fftshift((abs(x)),2)/max(max(abs(x)));
subplot(141)
DATA = plt(FIELD(:,:,it1));
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),DATA); set(pclr, 'edgecolor','none');pbaspect([1 1 1])
colormap gray
xlabel('$k_r \rho_s$'); ylabel('$k_z \rho_s$'); xlim(XLIM); ylim(YLIM);
title(sprintf('$t c_s/R=%.0f$',Ts2D(it1)));
subplot(142)
DATA = plt(FIELD(:,:,it2));
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),DATA); set(pclr, 'edgecolor','none');pbaspect([1 1 1])
colormap gray
xlabel('$k_r \rho_s$'); ylabel('$k_z \rho_s$'); set(gca,'ytick',[]); xlim(XLIM); ylim(YLIM);
title(sprintf('$t c_s/R=%.0f$',Ts2D(it2)));
subplot(143)
DATA = plt(FIELD(:,:,it3));
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),DATA); set(pclr, 'edgecolor','none');pbaspect([1 1 1])
colormap gray
xlabel('$k_r \rho_s$'); ylabel('$k_z \rho_s$');set(gca,'ytick',[]); xlim(XLIM); ylim(YLIM);
title(sprintf('$t c_s/R=%.0f$',Ts2D(it3)));
subplot(144)
DATA = plt(FIELD(:,:,it4));
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),DATA); set(pclr, 'edgecolor','none');pbaspect([1 1 1])
colormap gray
xlabel('$k_r \rho_s$'); ylabel('$k_z \rho_s$'); set(gca,'ytick',[]); xlim(XLIM); ylim(YLIM);
title(sprintf('$t c_s/R=%.0f$',Ts2D(it4)));
suptitle([FIELDLTX,', $\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
', $L=',num2str(L),'$, $N=',num2str(Nr),'$, $(P,J)=(',num2str(PMAXI),',',num2str(JMAXI),')$,',...
' $\mu_{hd}=$',num2str(MU)]);
save_figure
end
%%
if 0
%% Show frame in kspace
tf = 1000; [~,it2] = min(abs(Ts2D-tf)); [~,it5] = min(abs(Ts5D-tf));
tf = 800; [~,it2] = min(abs(Ts2D-tf)); [~,it5] = min(abs(Ts5D-tf));
fig = figure; FIGNAME = ['krkz_',sprintf('t=%.0f',Ts2D(it2)),'_',PARAMS];set(gcf, 'Position', [100, 100, 700, 600])
CLIM = [0,1e3]
subplot(223); plt = @(x) fftshift((abs(x)),2);
CLIM = [0,1];
subplot(223); plt = @(x) fftshift((abs(x)),2)/max(max(abs(x)));
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),plt(PHI(:,:,it2))); set(pclr, 'edgecolor','none');
caxis(CLIM); colormap hot
xlabel('$k_r$'); ylabel('$k_z$'); title(sprintf('$t c_s/R=%.0f$',Ts2D(it2))); legend('$|\hat\phi|$');
subplot(222); plt = @(x) fftshift(abs(x),2);
subplot(222);
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),plt(Ni00(:,:,it2))); set(pclr, 'edgecolor','none');
caxis(CLIM); colormap hot
xlabel('$k_r$'); ylabel('$k_z$'); legend('$|\hat n_i^{00}|$');
subplot(221); plt = @(x) fftshift(abs(x),2);
subplot(221);
pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),plt(Ne00(:,:,it2))); set(pclr, 'edgecolor','none');
caxis(CLIM); colormap hot
xlabel('$k_r$'); ylabel('$k_z$'); legend('$|\hat n_e^{00}|$');
subplot(224); plt = @(x) fftshift((abs(x)),2);
subplot(224);
colorbar;
caxis(CLIM); colormap hot
% pclr = pcolor(fftshift(KR,2),fftshift(KZ,2),plt(Si00(:,:,it5))); set(pclr, 'edgecolor','none'); colorbar;
......
matlab/setup.m
+ 9
− 4
View file @ efef5834
......@@ -38,18 +38,21 @@ MODEL.lambdaD = LAMBDAD;
% if A0KH ~= 0; SIMID = [SIMID,'_Nz_',num2str(L/2/pi*KR0KH),'_A_',num2str(A0KH)]; end;
% Time integration and intialization parameters
TIME_INTEGRATION.numerical_scheme = '''RK4''';
if INIT_PHI; INITIAL.init_noisy_phi = '.true.'; else; INITIAL.init_noisy_phi = '.false.';end;
INITIAL.init_background = 0.0e-5;
if (INIT_PHI && INIT_ZF == 0); INITIAL.init_noisy_phi = '.true.'; else; INITIAL.init_noisy_phi = '.false.';end;
INITIAL.INIT_ZF = INIT_ZF;
INITIAL.init_background = (INIT_ZF>0)*ZF_AMP;
INITIAL.init_noiselvl = NOISE0;
INITIAL.iseed = 42;
INITIAL.selfmat_file = '''null''';
INITIAL.eimat_file = '''null''';
INITIAL.iemat_file = '''null''';
if (CO == -3) % Write matrice filename for Full Coulomb DK
INITIAL.mat_file = '''null''';
if (CO == -3) % Write matrice filename for Full Coulomb
cmat_pmaxe = 25;
cmat_jmaxe = 18;
cmat_pmaxi = 25;
cmat_jmaxi = 18;
INITIAL.mat_file = ['''../../../iCa/FC_P_25_J_18_N_200_dk_0.05236_MFLR_0.h5'''];
INITIAL.selfmat_file = ...
['''../../../iCa/self_Coll_GKE_0_GKI_0_ESELF_1_ISELF_1_Pmaxe_',num2str(cmat_pmaxe),...
'_Jmaxe_',num2str(cmat_jmaxe),'_Pmaxi_',num2str(cmat_pmaxi),'_Jmaxi_',...
......@@ -62,11 +65,12 @@ if (CO == -3) % Write matrice filename for Full Coulomb DK
['''../../../iCa/ie_Coll_GKE_0_GKI_0_ITEST_1_IBACK_1_Pmaxe_',num2str(cmat_pmaxe),...
'_Jmaxe_',num2str(cmat_jmaxe),'_Pmaxi_',num2str(cmat_pmaxi),'_Jmaxi_',...
num2str(cmat_jmaxi),'_pamaxx_10_tau_1.0000_mu_0.0233.h5'''];
elseif (CO == -2) % Write matrice filename for DK Sugama
elseif (CO == -2) % Write matrice filename for Sugama
cmat_pmaxe = 10;
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.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_',...
......@@ -84,6 +88,7 @@ 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.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_',...
......
matlab/write_fort90.m
+ 2
− 1
View file @ efef5834
......@@ -64,13 +64,14 @@ fprintf(fid,'/\n');
fprintf(fid,'&INITIAL_CON\n');
fprintf(fid,[' INIT_NOISY_PHI =', INITIAL.init_noisy_phi,'\n']);
fprintf(fid,[' INIT_ZF_PHI =', INITIAL.INIT_ZF_PHI,'\n']);
fprintf(fid,[' INIT_ZF =', num2str(INITIAL.INIT_ZF),'\n']);
fprintf(fid,[' init_background =', num2str(INITIAL.init_background),'\n']);
fprintf(fid,[' init_noiselvl =', num2str(INITIAL.init_noiselvl),'\n']);
fprintf(fid,[' iseed =', num2str(INITIAL.iseed),'\n']);
fprintf(fid,[' selfmat_file =', INITIAL.selfmat_file,'\n']);
fprintf(fid,[' eimat_file =', INITIAL.eimat_file,'\n']);
fprintf(fid,[' iemat_file =', INITIAL.iemat_file,'\n']);
fprintf(fid,[' mat_file =', INITIAL.mat_file,'\n']);
fprintf(fid,'/\n');
fprintf(fid,'&TIME_INTEGRATION_PAR\n');
......
matlab/write_sbash_marconi.m
+ 1
− 1
View file @ efef5834
......@@ -36,7 +36,7 @@ fprintf(fid,[...
'#SBATCH --account=FUA35_TSVVT421\n',...
'#SBATCH --partition=skl_fua_',CLUSTER.PART,'\n',...
'module load autoload hdf5 fftw\n',...
'srun --cpu-bind=cores ./../../../bin/helaz ',num2str(NP_P),' ',num2str(NP_KR)]);
'srun --cpu-bind=cores ./../../../bin/',EXECNAME,' ',num2str(NP_P),' ',num2str(NP_KR)]);
fclose(fid);
system(['cp batch_script.sh ',BASIC.RESDIR,'/.']);
......
wk/ZF_fourier_analysis.m
+ 37
− 13
View file @ efef5834
%% Zonal flow spectral analysis
fig = figure; FIGNAME = ['zonal_flow_spectral_analysis_',PARAMS];
tend = Ts0D(end); tstart = tend - TAVG;
tend = Ts0D(end); tstart = tend-TAVG ;
[~,its0D] = min(abs(Ts0D-tstart));
[~,ite0D] = min(abs(Ts0D-tend));
[~,its2D] = min(abs(Ts2D-tstart));
[~,ite2D] = min(abs(Ts2D-tend));
set(gcf, 'Position', [100, 100, 800, 400])
% Time series analysis (burst period and time frequencies spectrum)
subplot(121)
samplerate = Ts0D(2)-Ts0D(1);
Y = log(PGAMMA_RI(its0D:ite0D)*(2*pi/Nr/Nz)^2);
[n,~] = size(Y);
Yy= fft(y);
Yy= fft(Y);
Pot = Yy .* conj(Yy) / n;
Freq = (samplerate / n * (1:n))';
Pot(1) = 0;
nmax = min(200,round(n/2));
nmax = min(20,round(n/2));
[amax, itmax] = max(Pot);
plot((0:nmax-1) , Pot(1:nmax)/amax,'DisplayName','$\Gamma_r(\omega)$');hold on;
plot([itmax-1,itmax-1],[0,1],'--k', 'DisplayName',['$T_{per}\approx',num2str(round(2*pi/Freq(itmax))),'\rho_s/c_s$']);
plot([itmax-1,itmax-1],[0,1],'--k', 'DisplayName',['$T_{per}\approx',num2str(round(1/Freq(itmax))),'L_\perp/c_s$']);
legend('show'); grid on; box on; xlabel('Period number'); yticks([]);
title('ZF temporal spectrum')
title('$\Gamma_r$ temporal spectrum')
% Space analysis (spatial period of ZF)
subplot(122)
Yy = mean(mean(1i*KR.*PHI(:,:,its2D:ite2D),3),2); hold on;
[n,~] = size(Yy);
Pot = Yy .* conj(Yy) / n;
[amax, ikmax] = max(Pot);
nmax = 20;
plot(0:nmax,Pot(1:nmax+1)/amax,'DisplayName','$\tilde V_E(k_r)$');
nmax = 20; n = numel(r);
[TT,NN] = meshgrid(Ts2D(its2D:ite2D),0:n-1);
Pot = NN;
for it = 1:ite2D-its2D+1
Y = mean(real(drphi(:,:,it)),2);
Yy = fft(Y);
[n,~] = size(Yy);
Pot(:,it) = Yy .* conj(Yy) / n;
end
[amax, ikmax] = max(mean(Pot,2));
% pclr = pcolor(NN(1:nmax,:),TT(1:nmax,:),Pot(1:nmax,:)); set(pclr, 'edgecolor','none'); hold on;
plot(0:nmax,mean(Pot(1:nmax+1,:),2)/amax,'DisplayName','$\langle\partial_r\phi\rangle_z (k_r)$'); hold on;
plot([ikmax-1,ikmax-1],[0,1],'--k', 'DisplayName',['$L_z=',num2str(2*pi/kr(ikmax)),'\rho_s$']);
grid on; box on;
title('ZF spatial spectrum')
xlabel('mode number'); yticks([]); legend('show')
xlabel('radial mode number'); yticks([]); legend('show')
save_figure
%% Shear and phi amp phase space
fig = figure; FIGNAME = ['phi_shear_phase_space_',PARAMS];
set(gcf, 'Position', [100, 100, 700, 500])
t1 = Ts2D(end); t0 = t1 - TAVG;
t1 = Ts2D(end); t0 = 0;
[~,its2D] = min(abs(Ts2D-t0)); [~,ite2D] = min(abs(Ts2D-t1));
scatter(phi_maxr_avgz(its2D:ite2D),shear_maxr_avgz(its2D:ite2D),35,Ts2D(its2D:ite2D),'.',...
'DisplayName',PARAMS); cbar = colorbar;ylabel(cbar,'$t c_s/\rho_s$','Interpreter','LaTeX')
......@@ -48,6 +56,22 @@ grid on; title('ES pot. vs Shear phase space')
% plot(phi_avgr_avgz(its2D:ite2D),shear_avgr_avgz(its2D:ite2D),'-')
save_figure
if 0
%% density and phi phase space
fig = figure; FIGNAME = ['phi_ni_phase_space_',PARAMS];
set(gcf, 'Position', [100, 100, 700, 500])
t1 = Ts2D(end); t0 = 0;
[~,its2D] = min(abs(Ts2D-t0)); [~,ite2D] = min(abs(Ts2D-t1));
scatter3(max(mean(ni00(:,:,its2D:ite2D),2),[],1),phi_maxr_avgz(its2D:ite2D),shear_maxr_avgz(its2D:ite2D),35,Ts2D(its2D:ite2D),'.',...
'DisplayName',PARAMS); cbar = colorbar;ylabel(cbar,'$t c_s/\rho_s$','Interpreter','LaTeX')
hold on
xlabel('$\langle n_i^{00} \rangle_z^r$'); ylabel('$\langle \phi \rangle_z^r$'); zlabel('$\langle dV_E/dr \rangle_z^r$')
grid on; title('ES pot. vs Shear phase space')
% plot(phi_avgr_maxz(its2D:ite2D),shear_avgr_maxz(its2D:ite2D),'-')
% plot(phi_maxr_maxz(its2D:ite2D),shear_maxr_maxz(its2D:ite2D),'-')
% plot(phi_avgr_avgz(its2D:ite2D),shear_avgr_avgz(its2D:ite2D),'-')
% save_figure
end
%% Non zonal quantities
PHI_NZ = PHI;
PHI_NZ(ikmax-1:ikmax+1,:,:) = 0;
......
wk/analysis_2D.m
+ 64
− 16
View file @ efef5834
addpath(genpath('../matlab')) % ... add
%% Load results
outfile ='';
if 0 % Local 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';
BASIC.RESDIR = ['../results/',outfile,'/'];
BASIC.MISCDIR = ['../results/',outfile,'/'];
end
if 1 % Marconi results
if 0 % Marconi results
outfile ='';
outfile ='';
outfile ='';
outfile ='';
outfile ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.9_nu_5e-01_DGGK_CLOS_0_mu_2e-02/out.txt';
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';
% load_marconi(outfile);
BASIC.RESDIR = ['../',outfile(46:end-8),'/'];
BASIC.MISCDIR = ['/misc/HeLaZ_outputs/',outfile(46:end-8),'/'];
......@@ -182,7 +183,8 @@ fig = figure; FIGNAME = ['t_evolutions',sprintf('_%.2d',JOBNUM),'_',PARAMS];
set(gcf, 'Position', [100, 100, 900, 800])
subplot(111);
suptitle(['$\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
', $L=',num2str(L),'$, $N=',num2str(Nr),'$, $(P,J)=(',num2str(PMAXI),',',num2str(JMAXI),')$']);
', $L=',num2str(L),'$, $N=',num2str(Nr),'$, $(P,J)=(',num2str(PMAXI),',',num2str(JMAXI),')$,',...
' $\mu_{hd}=$',num2str(MU)]);
subplot(421);
for ip = 1:Npe
for ij = 1:Nje
......@@ -215,25 +217,38 @@ subplot(111);
legend(['Gyro. flux';'Part. flux']);
grid on; xlabel('$t c_s/R$'); ylabel('$\Gamma_{r,i}$')
% ylim([0,2.0]);
if(1)
subplot(223)
plot(kz,g_(1,:),'-','DisplayName','Linear growth rate'); hold on;
plot([max(kz)*2/3,max(kz)*2/3],[0,10],'--k', 'DisplayName','2/3 Orszag AA');
grid on; xlabel('$k_z\rho_s$'); ylabel('$\gamma R/c_s$'); legend('show');
ylim([0,max(g_(1,:))]); xlim([0,max(kz)]);
subplot(224)
shearplot = 426; phiplot = 428;
else
shearplot = 223; phiplot = 224;
end
subplot(shearplot)
clr = line_colors(min(ip,numel(line_colors(:,1))),:);
lstyle = line_styles(min(ij,numel(line_styles)));
plot(Ts2D,shear_maxr_maxz,'DisplayName','$\max_{r,z}(s)$'); hold on;
plot(Ts2D,shear_maxr_avgz,'DisplayName','$\max_{r}\langle s \rangle_z$'); hold on;
plot(Ts2D,shear_avgr_maxz,'DisplayName','$\max_{z}\langle s \rangle_r$'); hold on;
plot(Ts2D,shear_avgr_avgz,'DisplayName','$\langle s \rangle_{r,z}$'); hold on;
grid on; xlabel('$t c_s/R$'); ylabel('$shear$');
subplot(phiplot)
clr = line_colors(min(ip,numel(line_colors(:,1))),:);
lstyle = line_styles(min(ij,numel(line_styles)));
plot(Ts2D,phi_maxr_maxz,'DisplayName','$\max_{r,z}(\phi)$'); hold on;
plot(Ts2D,phi_maxr_avgz,'DisplayName','$\max_{r}\langle\phi\rangle_z$'); hold on;
plot(Ts2D,phi_avgr_maxz,'DisplayName','$\max_{z}\langle\phi\rangle_r$'); hold on;
plot(Ts2D,phi_avgr_avgz,'DisplayName','$\langle\phi\rangle_{r,z}$'); hold on;
grid on; xlabel('$t c_s/R$'); ylabel('$T_e/e$'); legend('show');
grid on; xlabel('$t c_s/R$'); ylabel('$E.S. pot$');
save_figure
end
if 1
%% Space time diagramm (fig 11 Ivanov 2020)
TAVG = 1500; % Averaging time duration
TAVG = 3500; % Averaging time duration
%Compute steady radial transport
tend = Ts0D(end); tstart = tend - TAVG;
[~,its0D] = min(abs(Ts0D-tstart));
......@@ -283,7 +298,7 @@ fig = figure; FIGNAME = ['ZF_transport_drphi','_',PARAMS];set(gcf, 'Position',
save_figure
end
if 0
if 1
%% Space time diagramms
tstart = 0; tend = Ts2D(end);
[~,itstart] = min(abs(Ts2D-tstart));
......@@ -295,24 +310,41 @@ 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,0.02]);
caxis([0.0,max(max(mean(ni00(:,:,its2D:ite2D).*dzphi(:,:,its2D:ite2D),2)))]);
xticks([]); ylabel('$r/\rho_s$')
legend('$\Gamma_r(r)=\langle n_i^{GC}\partial_z\phi\rangle_z$')
legend('Radial part. transport $\langle n_i^{00}\partial_z\phi\rangle_z$')
title(['$\nu_{',CONAME,'}=$', num2str(NU), ', $\eta_B=$',num2str(ETAB),...
', $L=',num2str(L),'$, $N=',num2str(Nr),'$, $(P,J)=(',num2str(PMAXI),',',num2str(JMAXI),')$,',...
' $\mu_{hd}=$',num2str(MU)]);
subplot(212)
pclr = pcolor(TX,TY,squeeze(mean(drphi(:,:,trange),2))'); set(pclr, 'edgecolor','none'); colorbar;
fieldmax = max(max(mean(abs(drphi(:,:,its2D:ite2D)),2)));
caxis([-fieldmax,fieldmax]);
xlabel('$t c_s/R$'), ylabel('$r/\rho_s$')
legend('$\langle \partial_r\phi\rangle_z$')
legend('Zonal flow $\langle \partial_r\phi\rangle_z$')
save_figure
end
if 0
%% Averaged phi, ZF and shear 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([]);
end
%%
t0 = 0000;
t0 = 500;
[~, it02D] = min(abs(Ts2D-t0));
[~, it05D] = min(abs(Ts5D-t0));
skip_ = 10;
skip_ = 5;
DELAY = 1e-3*skip_;
FRAMES_2D = it02D:skip_:numel(Ts2D);
FRAMES_5D = it05D:skip_:numel(Ts5D);
......@@ -362,10 +394,26 @@ FIELDNAME = '$\langle\phi\rangle_{z}$'; XNAME = '$r/\rho_s$';
create_gif_1D_phi
end
if 0
%% flow averaged on z
GIFNAME = ['zf_z0',sprintf('_%.2d',JOBNUM),'_',PARAMS]; INTERP = 0; SCALING=1;
FIELD =(squeeze(mean(real(drphi),2))); linestyle = '-.k'; FRAMES = FRAMES_2D;
X = (r); T = Ts2D; YMIN = -1.1; YMAX = 1.1; XMIN = min(r); XMAX = max(r);
FIELDNAME = '$\langle\phi\rangle_{z}$'; XNAME = '$r/\rho_s$';
create_gif_1D_phi
end
if 0
%% shear averaged on z
GIFNAME = ['shear_z0',sprintf('_%.2d',JOBNUM),'_',PARAMS]; INTERP = 0; SCALING=1;
FIELD =(squeeze(mean(real(dr2phi),2))); linestyle = '-.k'; FRAMES = FRAMES_2D;
X = (r); T = Ts2D; YMIN = -1.1; YMAX = 1.1; XMIN = min(r); XMAX = max(r);
FIELDNAME = '$\langle\phi\rangle_{z}$'; XNAME = '$r/\rho_s$';
create_gif_1D_phi
end
if 0
%% phi @ kz = 0
GIFNAME = ['phi_kz0',sprintf('_%.2d',JOBNUM),'_',PARAMS]; INTERP = 0; SCALING = 0;
FIELD =squeeze(log10(abs(PHI(no_AA,1,:)))); linestyle = '-.'; FRAMES = FRAMES_2D;
X = kr(no_AA); T = Ts2D; YMIN = -30; YMAX = 6; XMIN = min(kr); XMAX = max(kr);
FIELD =squeeze(log10(abs(PHI(:,1,:)))); linestyle = '-.'; FRAMES = FRAMES_2D;
X = kr; T = Ts2D; YMIN = -0; YMAX = 6; XMIN = min(kr); XMAX = max(kr);
FIELDNAME = '$\log_{10}|\tilde\phi(k_z=0)|$'; XNAME = '$k_r\rho_s$';
create_gif_1D
end
......
wk/compile_cosolver_mat.m 0 → 100644
+ 141
− 0
View file @ efef5834
addpath(genpath('../matlab')) % ... add
%% Grid configuration
N = 200; % Frequency gridpoints (Nkr = N/2)
L = 120; % Size of the squared frequency domain
dk = 2.*pi/L;
kmax = N/2*dk;
kr = dk*(0:N/2);
kz = dk*(0:N/2);
[KZ, KR]= meshgrid(kz,kr);
KPERP = sqrt(KR.^2 + KZ.^2);
kperp = reshape(KPERP,[1,numel(kr)^2]);
kperp = uniquetol(kperp,1e-14);
Nperp = numel(kperp);
%% Model
% ! 0 = nothing
% ! 1 = coulomb
% ! 2 = pitch-angle (with constant coll.freq.)
% ! 3 = sugama
% ! 4 = pitch-angle with veloty dependent coll. freq.
% ! 5 = improved sugama
% ! 6 = Hirschman-Sigmar Clarke
% ! 7 = Abel (only for like species)
% ! 8 = conserving momentun pitch-angle operator (with veloty dependent coll. freq.)
% ! 9 = GK coulomb polarization term
CO = 3;
GK = 1;
P = 10;
J = 5;
M_FLR= 0; %to increase the NFLR
if 0
%% plot the kperp distribution
figure
plot(kperp)
end
% %% Check if the differences btw kperp is larger than naming precision
%%
%% We compute only on a kperp grid with dk space from 0 to kperpmax
kperpmax = sqrt(2) * kmax;
kperp = unique([0:dk:kperpmax,kperpmax]);
%% Naming
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;
for k_ = kperp
disp(['-Writing matrix for kperp = ',num2str(k_)])
%% Script to run COSOlver in order to create needed collision matrices
COSOlver_path = '../../Documents/MoliSolver/COSOlver/';
COSOLVER.pmaxe = P;
COSOLVER.jmaxe = J;
COSOLVER.pmaxi = P;
COSOLVER.jmaxi = J;
COSOLVER.kperp = k_;
COSOLVER.neFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2)))+M_FLR; % rule of thumb for sum truncation
COSOLVER.niFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2)))+M_FLR;
COSOLVER.idxT4max = 40;
COSOLVER.neFLRs = 0; % ... only for GK abel
COSOLVER.npeFLR = 0; % ... only for GK abel
COSOLVER.niFLRs = 0; % ... only for GK abel
COSOLVER.npiFLR = 0; % ... only for GK abel
COSOLVER.gk = GK;
COSOLVER.co = CO;
k_string = sprintf('%0.4f',k_);
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),...
'_Pmaxi_',num2str(COSOLVER.pmaxi),'_Jmaxi_',num2str(COSOLVER.jmaxi),...
'_JE_12',...
'_NFLR_',num2str(COSOLVER.neFLR),...
'_kperp_',k_string,'.h5'];
ei_mat_file_name = ['../iCa/ei_Coll_GKE_',num2str(COSOLVER.gk),'_GKI_',num2str(COSOLVER.gk),...
'_ETEST_',num2str(COSOLVER.co),'_EBACK_',num2str(COSOLVER.co),...
'_Pmaxe_',num2str(COSOLVER.pmaxe),'_Jmaxe_',num2str(COSOLVER.jmaxe),...
'_Pmaxi_',num2str(COSOLVER.pmaxi),'_Jmaxi_',num2str(COSOLVER.jmaxi),...
'_JE_12_tau_1.0000_mu_0.0233',...
'_NFLRe_',num2str(COSOLVER.neFLR),'_NFLRi_',num2str(COSOLVER.neFLR)...
'_kperp_',k_string,'.h5'];
ie_mat_file_name = ['../iCa/ie_Coll_GKE_',num2str(COSOLVER.gk),'_GKI_',num2str(COSOLVER.gk),...
'_ITEST_',num2str(COSOLVER.co),'_IBACK_',num2str(COSOLVER.co),...
'_Pmaxe_',num2str(COSOLVER.pmaxe),'_Jmaxe_',num2str(COSOLVER.jmaxe),...
'_Pmaxi_',num2str(COSOLVER.pmaxi),'_Jmaxi_',num2str(COSOLVER.jmaxi),...
'_JE_12_tau_1.0000_mu_0.0233',...
'_NFLRe_',num2str(COSOLVER.neFLR),'_NFLRi_',num2str(COSOLVER.neFLR)...
'_kperp_',k_string,'.h5'];
%% Load self matrix
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)
%% 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)
% 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)
%% 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)
% 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)
%% 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);
end
end
disp(['File saved @ :',matfilename])
\ No newline at end of file
wk/compute_collision_mat.m
+ 21
− 25
View file @ efef5834
......@@ -11,22 +11,28 @@ KPERP = sqrt(KR.^2 + KZ.^2);
kperp = reshape(KPERP,[1,numel(kr)^2]);
kperp = uniquetol(kperp,1e-14);
Nperp = numel(kperp);
if 0
%% plot the kperp distribution
figure
plot(kperp)
end
% %% Check if the differences btw kperp is larger than naming precision
% dkperp = diff(kperp);
% warning = sum(dkperp<0.0002);
% if warning > 0
% disp('Warning : dkperp < 0.0002');
% end
%% Model
% ! 0 = nothing
% ! 1 = coulomb
% ! 2 = pitch-angle (with constant coll.freq.)
% ! 3 = sugama
% ! 4 = pitch-angle with veloty dependent coll. freq.
% ! 5 = improved sugama
% ! 6 = Hirschman-Sigmar Clarke
% ! 7 = Abel (only for like species)
% ! 8 = conserving momentun pitch-angle operator (with veloty dependent coll. freq.)
% ! 9 = GK coulomb polarization term
CO = 3;
GK = 1;
%%
%% We compute only on a kperp grid with dk space from 0 to kperpmax
kperpmax = sqrt(2) * kmax;
kperp = unique([0:dk:kperpmax,kperpmax]);
%%
%% If DK collision, compute only kperp = 0
if (GK == 0); kperp = 0; end;
%% Kperp scan
n_ = 1;
for k_ = kperp
disp(['----------Computing matrix ',num2str(n_),'/',num2str(numel(kperp)),'----------'])
......@@ -38,8 +44,8 @@ for k_ = kperp
COSOLVER.jmaxi = 05;
COSOLVER.kperp = k_;
COSOLVER.neFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2)))+2; % rule of thumb for sum truncation
COSOLVER.niFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2)))+2;
COSOLVER.neFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2))); % rule of thumb for sum truncation
COSOLVER.niFLR = min(ceil((2/3*kperpmax)^2),max(5,ceil(COSOLVER.kperp^2)));
COSOLVER.idxT4max = 40;
COSOLVER.neFLRs = 0; % ... only for GK abel
......@@ -47,18 +53,8 @@ for k_ = kperp
COSOLVER.niFLRs = 0; % ... only for GK abel
COSOLVER.npiFLR = 0; % ... only for GK abel
COSOLVER.gk = 1;
COSOLVER.co = 3;
% ! 0 = nothing
% ! 1 = coulomb
% ! 2 = pitch-angle (with constant coll.freq.)
% ! 3 = sugama
% ! 4 = pitch-angle with veloty dependent coll. freq.
% ! 5 = improved sugama
% ! 6 = Hirschman-Sigmar Clarke
% ! 7 = Abel (only for like species)
% ! 8 = conserving momentun pitch-angle operator (with veloty dependent coll. freq.)
% ! 9 = GK coulomb polarization term
COSOLVER.gk = GK;
COSOLVER.co = CO;
write_fort90_COSOlver
......
wk/continue_multiple_runs_marconi.m
+ 4
− 11
View file @ efef5834
%% Paste the list of continue_run calls
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.5_nu_1e+00_SGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_1e+00_SGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/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/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.8_nu_1e+00_SGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.6_nu_5e-01_DGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.7_nu_5e-01_DGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.8_nu_5e-01_DGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.9_nu_5e-01_DGGK_CLOS_0_mu_2e-02/out.txt')
continue_run('/marconi_scratch/userexternal/ahoffman/HeLaZ/results/v2.5_P_10_J_5/200x100_L_120_P_10_J_5_eta_0.9_nu_1e-01_DGGK_CLOS_0_mu_2e-02/out.txt')
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')
%% Functions to modify preexisting fort.90 input file and launch on marconi
function [] = continue_run(outfilename)
EXECNAME = 'helaz_2.62';
%% CLUSTER PARAMETERS
CLUSTER.PART = 'prod'; % dbg or prod
CLUSTER.PART = 'dbg'; % 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
......@@ -51,7 +44,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)];
......
wk/linear_study.m
+ 38
− 34
View file @ efef5834
for NU = [1.0]
for ETAB = [0.5]
for CO = [1,2]
for CO = [-2 2]
%clear all;
addpath(genpath('../matlab')) % ... add
default_plots_options
......@@ -14,18 +14,18 @@ TAU = 1.0; % e/i temperature ratio
% ETAB = 0.5;
ETAN = 1.0; % Density gradient
ETAT = 0.0; % Temperature gradient
NU_HYP = 1.0; % Hyperdiffusivity coefficient
NU_HYP = 0.0; % Hyperdiffusivity coefficient
LAMBDAD = 0.0;
NOISE0 = 1.0e-5;
%% GRID PARAMETERS
N = 200; % Frequency gridpoints (Nkr = N/2)
N = 20; % Frequency gridpoints (Nkr = N/2)
L = 120; % 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
%% TIME PARMETERS
TMAX = 200; % Maximal time unit
DT = 2e-2; % Time step
TMAX = 100; % Maximal time unit
DT = 1e-2; % Time step
SPS0D = 0.5; % Sampling per time unit for 2D arrays
SPS2D = 1; % Sampling per time unit for 2D arrays
SPS5D = 1/2; % Sampling per time unit for 5D arrays
......@@ -33,7 +33,8 @@ SPSCP = 0; % Sampling per time unit for checkpoints
RESTART = 0; % To restart from last checkpoint
JOB2LOAD= 00;
%% OPTIONS
SIMID = 'linear_study'; % Name of the simulation
% SIMID = 'test'; % 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)
......@@ -41,7 +42,8 @@ NON_LIN = 0 *(1-KREQ0); % activate non-linearity (is cancelled if KREQ0 = 1)
CLOS = 0; % Closure model (0: =0 truncation, 1: semi coll, 2: Copy closure J+1 = J, P+2 = P)
NL_CLOS = 0; % nonlinear closure model (0: =0 nmax = jmax, 1: nmax = jmax-j, >1 : nmax = NL_CLOS)
KERN = 0; % Kernel model (0 : GK)
INIT_PHI= 0; % Start simulation with a noisy phi and moments
INIT_PHI= 0; % Start simulation with a noisy phi
INIT_ZF = 0; % Start simulation with a noisy phi
%% OUTPUTS
W_DOUBLE = 0;
W_GAMMA = 1;
......@@ -62,11 +64,12 @@ 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, 3, 4, 6];%, 8, 10];
% JA = [1, 2, 2, 3];%, 4, 5];
PA = [4];
JA = [2];
DTA= DT./sqrt(JA);
JA = [4];
DTA= DT./sqrt(JA)/4;
% DTA= DT;
mup_ = MU_P;
muj_ = MU_J;
Nparam = numel(PA);
......@@ -80,35 +83,36 @@ for i = 1:Nparam
PMAXE = PA(i); PMAXI = PA(i);
JMAXE = JA(i); JMAXI = JA(i);
DT = DTA(i);
MU_P = mup_/PMAXE^2;
MU_J = muj_/JMAXE^3;
% MU_P = mup_/PMAXE^2;
% MU_J = muj_/JMAXE^3;
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'])
% Load and process results
% 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'];
load_results
tend = Ts2D(end); tstart = 0.4*tend;
[~,itstart] = min(abs(Ts2D-tstart));
[~,itend] = min(abs(Ts2D-tend));
for ikr = 1:N/2+1
gamma_Ni00(i,ikr) = LinearFit_s(Ts2D(itstart:itend),squeeze(abs(Ni00(ikr,1,itstart:itend))));
gamma_Ni00(i,ikr) = (LinearFit_s(Ts2D(itstart:itend)',(squeeze(abs(Ni00(ikr,1,itstart:itend))))));
Ni00_ST(i,ikr,1:numel(Ts2D)) = squeeze((Ni00(ikr,1,:)));
end
tend = Ts5D(end); tstart = 0.4*tend;
[~,itstart] = min(abs(Ts5D-tstart));
[~,itend] = min(abs(Ts5D-tend));
for ikr = 1:N/2+1
gamma_Nipj(i,ikr) = LinearFit_s(Ts5D(itstart:itend),squeeze(max(max(abs(Nipj(:,:,ikr,1,itstart:itend)),[],1),[],2)));
gamma_Nipj(i,ikr) = LinearFit_s(Ts5D(itstart:itend)',squeeze(max(max(abs(Nipj(:,:,ikr,1,itstart:itend)),[],1),[],2)));
end
gamma_Ni00(i,:) = real(gamma_Ni00(i,:) .* (gamma_Ni00(i,:)>=0.0));
gamma_Nipj(i,:) = real(gamma_Nipj(i,:) .* (gamma_Nipj(i,:)>=0.0));
% kzmax = abs(kr(ikzmax));
% Bohm_transport(i) = ETAB/ETAN*gmax/kzmax^2;
% Clean output
system(['rm -r ',BASIC.RESDIR])
system(['rm -r ',BASIC.RESDIR]);
end
if 1
......@@ -116,7 +120,7 @@ if 1
SCALE = 1;%sqrt(2);
fig = figure; FIGNAME = 'linear_study';
plt = @(x) x;
subplot(211)
% subplot(211)
for i = 1:Nparam
clr = line_colors(mod(i-1,numel(line_colors(:,1)))+1,:);
linestyle = line_styles(floor((i-1)/numel(line_colors(:,1)))+1);
......@@ -126,22 +130,22 @@ subplot(211)
'DisplayName',['$P=$',num2str(PA(i)),', $J=$',num2str(JA(i))]);
hold on;
end
grid on; xlabel('$k_z\rho_s^{R}$'); ylabel('$\gamma(N_i^{00})\rho_s/c_s$'); xlim([0.0,max(kr)]);
title(['$\eta_B=',num2str(ETAB),'$, $\nu_{',CONAME,'}=',num2str(NU),'$, ', CLOSNAME])
legend('show')
subplot(212)
for i = 1:Nparam
clr = line_colors(mod(i-1,numel(line_colors(:,1)))+1,:);
linestyle = line_styles(floor((i-1)/numel(line_colors(:,1)))+1);
plot(plt(SCALE*kr),plt(gamma_Nipj(i,:)),...
'Color',clr,...
'LineStyle',linestyle{1},...
'DisplayName',['$P=$',num2str(PA(i)),', $J=$',num2str(JA(i))]);
hold on;
end
grid on; xlabel('$k_z\rho_s^{R}$'); ylabel('$\gamma(\max_{pj}N_i^{pj})\rho_s/c_s$'); xlim([0.0,max(kr)]);
grid on; xlabel('$k_z\rho_s^{R}$'); ylabel('$\gamma(N_i^{00})L_\perp/c_s$'); xlim([0.0,max(kr)]);
title(['$\eta_B=',num2str(ETAB),'$, $\nu_{',CONAME,'}=',num2str(NU),'$, ', CLOSNAME])
legend('show')
% subplot(212)
% for i = 1:Nparam
% clr = line_colors(mod(i-1,numel(line_colors(:,1)))+1,:);
% linestyle = line_styles(floor((i-1)/numel(line_colors(:,1)))+1);
% plot(plt(SCALE*kr),plt(gamma_Nipj(i,:)),...
% 'Color',clr,...
% 'LineStyle',linestyle{1},...
% 'DisplayName',['$P=$',num2str(PA(i)),', $J=$',num2str(JA(i))]);
% hold on;
% end
% grid on; xlabel('$k_z\rho_s^{R}$'); ylabel('$\gamma(\max_{pj}N_i^{pj})\rho_s/c_s$'); xlim([0.0,max(kr)]);
% title(['$\eta_B=',num2str(ETAB),'$, $\nu_{',CONAME,'}=',num2str(NU),'$, ', CLOSNAME])
% legend('show')
saveas(fig,[SIMDIR,'gamma_Ni_vs_',param_name,'_',PARAMS,'.fig']);
saveas(fig,[SIMDIR,'gamma_Ni_vs_',param_name,'_',PARAMS,'.png']);
end
......
wk/load_multiple_outputs_marconi.m
+ 31
− 0
View file @ efef5834
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')
wk/local_run.m
+ 12
− 12
View file @ efef5834
......@@ -4,19 +4,19 @@ addpath(genpath('../matlab')) % ... add
CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% PHYSICAL PARAMETERS
NU = 0.1; % Collision frequency
ETAB = 0.6; % Magnetic gradient
NU = 1.0e-1; % Collision frequency
ETAB = 0.4; % Magnetic gradient
ETAN = 1.0; % Density gradient
NU_HYP = 0.2;
NU_HYP = 1.0;
%% GRID PARAMETERS
N = 200; % Frequency gridpoints (Nkr = N/2)
N = 20; % Frequency gridpoints (Nkr = N/2)
L = 120; % Size of the squared frequency domain
P = 2;
P = 0;
J = 1;
MU_P = 0.0/P^2; % Hermite hyperdiffusivity -mu_p*(d/dvpar)^4 f
MU_J = 0.0/J^2; % Laguerre hyperdiffusivity -mu_j*(d/dvperp)^4 f
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 = 500; % Maximal time unit
TMAX = 1000; % 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,13 +27,13 @@ JOB2LOAD= 0;
%% OPTIONS AND NAMING
% Collision operator
% (0 : L.Bernstein, 1 : Dougherty, 2: Sugama, 3 : Full Couloumb ; +/- for GK/DK)
CO = 2;
CO = 1;
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 = 'local_low_pj_analysis'; % Name of the simulation
% SIMID = 'Sugama_no_ei_ie_coll'; % Name of the simulation
SIMID = ['v2.5_P_',num2str(P),'_J_',num2str(J)]; % Name of the simulation
SIMID = 'test'; % 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;
%% OUTPUTS
W_DOUBLE = 0;
W_GAMMA = 1;
......
wk/marconi_run.m
+ 22
− 13
View file @ efef5834
clear all;
addpath(genpath('../matlab')) % ... add
SUBMIT = 1; % To submit the job automatically
% EXECNAME = 'helaz_dbg';
EXECNAME = 'helaz_2.62';
for ETAB = [0.6 0.7]
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Set Up parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
......@@ -7,37 +11,38 @@ addpath(genpath('../matlab')) % ... add
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
CLUSTER.MEM = '128GB'; % Memory
CLUSTER.JNAME = 'HeLaZ';% Job name
NP_P = 2; % MPI processes along p
NP_KR = 24; % MPI processes along kr
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% PHYSICAL PARAMETERS
NU = 1.0; % Collision frequency
ETAB = 0.6; % Magnetic gradient
NU = 1e-1; % 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 = -3;
CO = -2;
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 = 10; % Electron and Ion highest Hermite polynomial degree
J = 05; % Electron and Ion highest Laguerre polynomial degree
MU_P = 0.0/(P/2)^4;% Hermite hyperdiffusivity -mu_p*(d/dvpar)^4 f
MU_J = 0.0/(J/2)^2;% Laguerre hyperdiffusivity -mu_j*(d/dvperp)^4 f
P = 06; % Electron and Ion highest Hermite polynomial degree
J = 03; % 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
SPS0D = 1; % Sampling per time unit for profiler
SPS2D = 1/2; % Sampling per time unit for 2D arrays
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 = 1; % To restart from last checkpoint
RESTART = 0; % To restart from last checkpoint
JOB2LOAD= 0;
%% Naming
% SIMID = 'test'; % Name of the simulation
SIMID = ['v2.5_P_',num2str(P),'_J_',num2str(J)]; % Name of the simulation
SIMID = ['v2.6_P_',num2str(P),'_J_',num2str(J)]; % Name of the simulation
PREFIX =[];
% PREFIX = sprintf('%d_%d_',NP_P, NP_KR);
%% Options
......@@ -66,7 +71,7 @@ JMAXI = J; % Highest '' Laguerre ''
kmax = N*pi/L;% Highest fourier mode
% kmax = 2/3*N*pi/L;% Highest fourier mode with AA
HD_CO = 0.5; % Hyper diffusivity cutoff ratio
MU = NU_HYP/(HD_CO*kmax)^4 % Hyperdiffusivity coefficient
MU = NU_HYP/(HD_CO*kmax)^4; % Hyperdiffusivity coefficient
NOISE0 = 1.0e-5;
ETAT = 0.0; % Temperature gradient
ETAN = 1.0; % Density gradient
......@@ -82,7 +87,11 @@ CLUSTER.CPUPT = '1'; % CPU per task
setup
write_sbash_marconi
system('rm fort.90 setup_and_run.sh batch_script.sh');
disp('done');
if(mod(NP_P*NP_KR,48)~= 0)
disp('WARNING : unused cores (ntot cores must be a 48 multiple)');
end
if(SUBMIT)
system('ssh ahoffman@login.marconi.cineca.it sh HeLaZ/wk/setup_and_run.sh');
end
disp('done');
end
\ No newline at end of file
wk/open_figure_script.m
+ 11
− 4
View file @ efef5834
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';
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_ = 'v2.5_P_2_J_1/200x100_L_120_P_2_J_1_eta_0.6_nu_1e-01_DGGK_CLOS_0_mu_2e-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.6_nu_5e-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';
%%
% figname_ = '/fig/ZF_transport_drphi_';
figname_ = '/fig/ZF_transport_drphi_';
% figname_ = '/fig/space_time_';
figname_ = '/fig/phi_shear_phase_space_';
% figname_ = '/fig/phi_shear_phase_space_';
path_ = '../results/';
......
wk/transport_results_2_5.m
+ 83
− 31
View file @ efef5834
%% nuDGGK = 1.0
figure
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_ = [6.1529 0.7986 0.0493 0.0134];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=1.0$',...
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.6,0.7,0.8];
gamma = [3.5546 0.3917 0.0421];
std_ = [0.5846 0.0486 0.0203];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=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_ = [0.7964];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=12,06$, $\nu_{DGGK}=1.0$',...
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
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_ = [3.67 1.13 0.055 0.004];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=06,03$, $\nu_{DGGK}=0.5$',...
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];
gamma = [2.4949 0.2578 0.058];
std_ = [0.8696 0.1191 0.011];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=10,05$, $\nu_{DGGK}=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
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_ = [0.295 0.05 0.0023];
bmax = [1.0 0.21 0.0367];
bmin = [0.02 0.04 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,bmax,'-.','DisplayName',' ','Color',line_colors(6,:)); hold on;
plot(eta,bmin,'-.','DisplayName',' ','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.2133 0.088 0.04253 0.02435];
std_ = [4.7 0.21328 0.065 0.009 0.00386];
bmax = [12.2 0.8 0.25 0.04253 0.02435];
bmin = [12.2 0.02 0.03 0.04253 0.02435];
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,bmax,'-.','DisplayName',' ','Color',line_colors(7,:)); hold on;
plot(eta,bmin,'-.','DisplayName',' ','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_ = [3.1 0.22 0.0019 0.0006];
errorbar(eta,gamma,std_,'-','DisplayName','$P,J=06,03$, $\nu_{SGGK}=1.0$',...
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_ = [1.34 0.228 0.001 0.001];
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;
......
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