diff --git a/matlab/evolve_tracers.m b/matlab/evolve_tracers.m index a6a738c..1363616 100644 --- a/matlab/evolve_tracers.m +++ b/matlab/evolve_tracers.m @@ -1,338 +1,338 @@ % Options SHOW_FILM = 1; field2plot ='phi'; INIT = 'lin'; % lin (for a line)/ round (for a small round)/ gauss for random -U_TIME = 40; % >0 for frozen velocity at a given time, -1 for evolving field +U_TIME = 15000; % >0 for frozen velocity at a given time, -1 for evolving field Evolve_U = 1; % 0 for frozen velocity at a given time, 1 for evolving field Tfin = 100; -dt_ = 0.01; +dt_ = 0.1; Nstep = ceil(Tfin/dt_); % Init tracers Np = 20; %number of tracers % color = tcolors; color = jet(Np); tcolors = distinguishable_colors(Np); %Their colors dimmed = 0; % to dimm the colormap in the background (infty = white, 0 normal color) Na = 10/dt_; %length of trace Traj_x = zeros(Np,Nstep); Traj_y = zeros(Np,Nstep); Disp_x = zeros(Np,Nstep); Disp_y = zeros(Np,Nstep); xmax = max(data.x); xmin = min(data.x); ymax = max(data.y); ymin = min(data.y); switch INIT case 'lin' % Evenly distributed initial positions dp_ = (xmax-xmin)/(Np-1); Xp = linspace(xmin+dp_/2,xmax-dp_/2,Np); Yp = zeros(1,Np); Zp = zeros(Np); case 'round' % All particles arround a same point xc = 0; yc = 0; theta = rand(1,Np)*2*pi; r = 0.1; Xp = xc + r*cos(theta); Yp = yc + r*sin(theta); Zp = zeros(Np); case 'gauss' % normal distribution arround a point xc = 0; yc = 0; sgm = 1.0; dx = normrnd(0,sgm,[1,Np]); dx = dx - mean(dx); Xp = xc + dx; dy = normrnd(0,sgm,[1,Np]); dy = dy - mean(dy); Yp = yc + dy; Zp = zeros(Np); end % position grid and velocity field [YY_, XX_ ,ZZ_] = meshgrid(data.y,data.x,data.z); [KX,KY] = meshgrid(data.kx,data.ky); Ux = zeros(size(XX_)); Uy = zeros(size(XX_)); Uz = zeros(size(XX_)); ni = zeros(size(XX_)); [~,itu_] = min(abs(U_TIME-data.Ts3D)); % computing the frozen velocity field for iz = 1:data.Nz % Ux(:,:,iz) = real(ifft2( 1i*KY.*(data.PHI(:,:,iz,itu_)),data.Nx,data.Ny)); % Uy(:,:,iz) = real(ifft2(-1i*KX.*(data.PHI(:,:,iz,itu_)),data.Nx,data.Ny)); % ni(:,:,iz) = real(ifft2(data.DENS_I(:,:,iz,itu_),data.Nx,data.Ny)); Ux(:,:,iz) = ifourier_GENE( 1i*KY.*(data.PHI(:,:,iz,itu_)))'; Uy(:,:,iz) = ifourier_GENE(-1i*KX.*(data.PHI(:,:,iz,itu_)))'; ni(:,:,iz) = ifourier_GENE(data.DENS_I(:,:,iz,itu_))'; end %% FILM options FPS = 30; DELAY = 1/FPS; FORMAT = '.gif'; if SHOW_FILM FILENAME = [data.localdir,'tracer_evolution',FORMAT]; switch FORMAT case '.avi' vidfile = VideoWriter(FILENAME,'Uncompressed AVI'); vidfile.FrameRate = FPS; open(vidfile); end fig = figure; end % %Time loop t_ = 0; it = 1; itu_old = 0; nbytes = fprintf(2,'frame %d/%d',it,Nstep); while(t_ data.x(ixC) % right from grid point ix0 = ixC; ix1 = ixC+1; else % left ix0 = ixC-1; ix1 = ixC; end y_ = Yp(ip); [e_y,iyC] = min(abs(y_-data.y)); if e_y == 0 % on the face iy0 = iyC; iy1 = iyC; elseif y_ > data.y(iyC) % above iy0 = iyC; iy1 = iyC+1; else % under iy0 = iyC-1; iy1 = iyC; end z_ = Zp(ip,1); [e_z,izC] = min(abs(z_-data.z)); if e_z == 0 % on the face iz0 = izC; iz1 = izC; elseif z_ > data.z(izC) % before iz0 = izC; iz1 = izC+1; else % behind iz0 = izC-1; iz1 = izC; end x0 = data.x(ix0); x1 = data.x(ix1); %left right y0 = data.y(iy0); y1 = data.y(iy1); %down top z0 = data.z(iz0); z1 = data.z(iz1); %back front if(e_x > 0) ai__ = (x_ - x0)/(x1-x0); % interp coeff x else ai__ = 0; end if(e_y > 0) a_i_ = (y_ - y0)/(y0-y1); % interp coeff y else a_i_ = 0; end if(e_z > 0) a__i = (z_ - z0)/(z1-z0); % interp coeff z else a__i = 0; end % interp velocity and density %velocity values u000 = [Ux(ix0,iy0,iz0) Uy(ix0,iy0,iz0) ni(ix0,iy0,iz0)]; u001 = [Ux(ix0,iy0,iz1) Uy(ix0,iy0,iz1) ni(ix0,iy0,iz1)]; u010 = [Ux(ix0,iy1,iz0) Uy(ix0,iy1,iz0) ni(ix0,iy1,iz0)]; u011 = [Ux(ix0,iy1,iz1) Uy(ix0,iy1,iz1) ni(ix0,iy1,iz1)]; u100 = [Ux(ix1,iy0,iz0) Uy(ix1,iy0,iz0) ni(ix1,iy0,iz0)]; u101 = [Ux(ix1,iy0,iz1) Uy(ix1,iy0,iz1) ni(ix1,iy0,iz1)]; u110 = [Ux(ix1,iy1,iz0) Uy(ix1,iy1,iz0) ni(ix1,iy1,iz0)]; u111 = [Ux(ix1,iy1,iz1) Uy(ix1,iy1,iz1) ni(ix1,iy1,iz1)]; %linear interpolation linterp = @(x1,x2,a) (1-a)*x1 + a*x2; u_00 = linterp(u000,u100,ai__); u_10 = linterp(u010,u110,ai__); u__0 = linterp(u_00,u_10,a_i_); u_01 = linterp(u001,u101,ai__); u_11 = linterp(u011,u111,ai__); u__1 = linterp(u_01,u_11,a_i_); u___ = linterp(u__0,u__1,a__i); % push the particle % q = sign(-u___(3)); q = -u___(3); % q =1; x_ = x_ + dt_*u___(1)*q; y_ = y_ + dt_*u___(2)*q; % apply periodic boundary conditions if(x_ > xmax) x_ = xmin + (x_ - xmax); elseif(x_ < xmin) x_ = xmax + (x_ - xmin); end if(y_ > ymax) y_ = ymin + (y_ - ymax); elseif(y_ < ymin) y_ = ymax + (y_ - ymin); end % store data % Trajectory Traj_x(ip,it) = x_; Traj_y(ip,it) = y_; % Displacement Disp_x(ip,it) = Disp_x(ip,it) + dt_*u___(1)*q; Disp_y(ip,it) = Disp_y(ip,it) + dt_*u___(2)*q; % update position Xp(ip) = x_; Yp(ip) = y_; end %% Movie if SHOW_FILM && (~Evolve_U || (itu_old ~= itu_)) % updating the velocity field clf(fig); switch field2plot case 'phi' F2P = ifourier_GENE(data.PHI(:,:,iz,itu_))'; case 'ni' F2P = ifourier_GENE(data.DENS_I(:,:,iz,itu_))'; case 'Ni00' F2P = ifourier_GENE(data.Ni00(:,:,iz,itu_))'; case 'none' F2P = 0*XX_; end scale = max(max(abs(F2P))); % Scaling to normalize pclr = pcolor(XX_,YY_,F2P/scale); caxis([-1 1]); colormap(bluewhitered); colorbar; set(pclr, 'edgecolor','none'); hold on; caxis([-1+dimmed,1+dimmed]); shading interp for ip = 1:Np ia0 = max(1,it-Na); plot(Traj_x(ip,ia0:it),Traj_y(ip,ia0:it),'.','Color',color(ip,:)); hold on end for ip = 1:Np plot(Traj_x(ip, 1),Traj_y(ip, 1),'xk'); hold on plot(Traj_x(ip,it),Traj_y(ip,it),'ok','MarkerFaceColor',color(ip,:)); hold on end title(['$t \approx$', sprintf('%.3d',ceil(data.Ts3D(itu_)))]); axis equal xlim([xmin xmax]); ylim([ymin ymax]); drawnow % Capture the plot as an image frame = getframe(fig); switch FORMAT case '.gif' im = frame2im(frame); [imind,cm] = rgb2ind(im,32); % Write to the GIF File if it == 1 imwrite(imind,cm,FILENAME,'gif', 'Loopcount',inf); else imwrite(imind,cm,FILENAME,'gif','WriteMode','append', 'DelayTime',DELAY); end case '.avi' writeVideo(vidfile,frame); otherwise disp('Unknown format'); break end end t_ = t_ + dt_; it = it + 1; itu_old = itu_; % terminal info while nbytes > 0 fprintf('\b') nbytes = nbytes - 1; end nbytes = fprintf(2,'frame %d/%d',it,Nstep); end disp(' ') switch FORMAT case '.gif' disp(['Gif saved @ : ',FILENAME]) case '.avi' disp(['Video saved @ : ',FILENAME]) close(vidfile); end Nt = it; %% Plot trajectories and statistics xtot = Disp_x; ytot = Disp_y; % aggregate the displacements for iq = 1:Np x_ = 0; y_ = 0; for iu = 1:Nt-1 x_ = x_ + Disp_x(iq,iu); y_ = y_ + Disp_y(iq,iu); xtot(iq,iu) = x_; ytot(iq,iu) = y_; end end ma = @(x) x;%movmean(x,100); time_ = linspace(U_TIME,U_TIME+Tfin,it-1); figure; subplot(221); for ip = 1:Np % plot(ma(time_),ma(Traj_x(ip,:)),'-','Color',tcolors(ip,:)); hold on plot(ma(time_),xtot(ip,:),'-','Color',tcolors(ip,:)); hold on % plot(ma(time_),ma(Disp_x(ip,:)),'-','Color',tcolors(ip,:)); hold on end ylabel('$x_p$'); xlim(U_TIME + [0 Tfin]); subplot(222); itf = floor(Nt/2); %fit end time % x^2 displacement plot(time_,mean(xtot.^2,1),'DisplayName','$\langle x.^2\rangle_p$'); hold on fit = polyfit(time_(1:itf),mean(xtot(:,1:itf).^2,1),1); plot(time_,fit(1)*time_+fit(2),'--k'); hold on ylabel('$\langle x^2 \rangle_p$'); % % y^2 displacement % fit = polyfit(time_(1:itf),mean(ytot(:,1:itf).^2,1),1); % plot(time_,fit(1)*time_+fit(2),'--k','DisplayName',['$\alpha=',num2str(fit(1)),'$']); hold on % plot(time_,mean(ytot.^2,1),'DisplayName','$\langle y.^2\rangle_p$'); % % % r^2 displacement % fit = polyfit(time_(1:itf),mean(xtot(:,1:itf).^2+ytot(:,1:itf).^2,1),1); % plot(time_,fit(1)*time_+fit(2),'--k','DisplayName',['$\alpha=',num2str(fit(1)),'$']); hold on % plot(time_,mean(xtot.^2+ytot.^2,1),'DisplayName','$\langle r.^2\rangle_p$'); % ylabel('$\langle x^2 \rangle_p$'); % xlim(U_TIME + [0 Tfin]); subplot(223); for ip = 1:Np % plot(ma(time_),ma(Traj_y(ip,:)),'-','Color',tcolors(ip,:)); hold on plot(ma(time_),ytot(ip,:),'-','Color',tcolors(ip,:)); hold on % plot(ma(time_),ma(Disp_y(ip,:)),'-','Color',tcolors(ip,:)); hold on end xlabel('time'); ylabel('$y_p$'); xlim(U_TIME + [0 Tfin]); subplot(224); histogram(xtot(:,1),20); hold on histogram(xtot(:,end),20) xlabel('position'); ylabel('$n$'); diff --git a/testcases/miller_example/fort_00.90 b/testcases/fort.90 similarity index 93% rename from testcases/miller_example/fort_00.90 rename to testcases/fort.90 index 728a345..abce02c 100644 --- a/testcases/miller_example/fort_00.90 +++ b/testcases/fort.90 @@ -1,86 +1,86 @@ &BASIC nrun = 100000000 dt = 0.01 tmax = 100 maxruntime = 356400 / &GRID pmaxe = 2 jmaxe = 1 pmaxi = 2 jmaxi = 1 Nx = 128 Lx = 100 - Ny = 64 - Ly = 120 + Ny = 128 + Ly = 240 Nz = 16 Npol = 1 SG = .false. / &GEOMETRY geom = 'miller' q0 = 1.4 shear = 0.0 eps = 0.18 / &OUTPUT_PAR nsave_0d = 50 nsave_1d = -1 nsave_2d = Inf nsave_3d = 100 - nsave_5d = 500 + nsave_5d = 1000 write_doubleprecision = .false. write_gamma = .true. write_hf = .true. write_phi = .true. write_Na00 = .true. write_Napj = .true. write_Sapj = .false. write_dens = .true. write_temp = .true. job2load = -1 / &MODEL_PAR ! Collisionality CLOS = 0 NL_CLOS = -1 LINEARITY = 'nonlinear' KIN_E = .false. - mu_x = 0.1 - mu_y = 0.1 + mu_x = 1.0 + mu_y = 1.0 N_HD = 4 - mu_z = 0.2 + mu_z = 1.0 mu_p = 0 mu_j = 0 nu = 0.2 tau_e = 1 tau_i = 1 sigma_e = 0.023338 sigma_i = 1 q_e = -1 q_i = 1 K_Ni = 2.22 K_Ti = 6.92 K_Ne = 0 K_Te = 0 GradB = 1 CurvB = 1 lambdaD = 0 beta = 0 / &COLLISION_PAR collision_model = 'DG' gyrokin_CO = .false. interspecies = .true. mat_file = 'null' / &INITIAL_CON INIT_OPT = 'ppj' ACT_ON_MODES = 'none' init_background = 0 init_noiselvl = 1e-3 iseed = 42 / &TIME_INTEGRATION_PAR numerical_scheme = 'RK4' / diff --git a/testcases/miller_example/fort_01.90 b/testcases/miller_example/fort_01.90 deleted file mode 100644 index 9e3c980..0000000 --- a/testcases/miller_example/fort_01.90 +++ /dev/null @@ -1,86 +0,0 @@ -&BASIC - nrun = 100000000 - dt = 0.01 - tmax = 100 - maxruntime = 356400 -/ -&GRID - pmaxe = 2 - jmaxe = 1 - pmaxi = 2 - jmaxi = 1 - Nx = 128 - Lx = 100 - Ny = 64 - Ly = 120 - Nz = 16 - Npol = 1 - SG = .false. -/ -&GEOMETRY - geom = 'miller' - q0 = 1.4 - shear = 0.0 - eps = 0.18 -/ -&OUTPUT_PAR - nsave_0d = 50 - nsave_1d = -1 - nsave_2d = Inf - nsave_3d = 100 - nsave_5d = 500 - write_doubleprecision = .false. - write_gamma = .true. - write_hf = .true. - write_phi = .true. - write_Na00 = .true. - write_Napj = .true. - write_Sapj = .false. - write_dens = .true. - write_temp = .true. - job2load = 0 -/ -&MODEL_PAR - ! Collisionality - CLOS = 0 - NL_CLOS = -1 - LINEARITY = 'nonlinear' - KIN_E = .false. - mu_x = 0.1 - mu_y = 0.1 - N_HD = 4 - mu_z = 0.2 - mu_p = 0 - mu_j = 0 - nu = 0.2 - tau_e = 1 - tau_i = 1 - sigma_e = 0.023338 - sigma_i = 1 - q_e = -1 - q_i = 1 - K_Ni = 2.22 - K_Ti = 6.92 - K_Ne = 0 - K_Te = 0 - GradB = 1 - CurvB = 1 - lambdaD = 0 - beta = 0 -/ -&COLLISION_PAR - collision_model = 'DG' - gyrokin_CO = .false. - interspecies = .true. - mat_file = 'null' -/ -&INITIAL_CON - INIT_OPT = 'ppj' - ACT_ON_MODES = 'none' - init_background = 0 - init_noiselvl = 1e-3 - iseed = 42 -/ -&TIME_INTEGRATION_PAR - numerical_scheme = 'RK4' -/ diff --git a/testcases/miller_example/fort_02.90 b/testcases/miller_example/fort_02.90 deleted file mode 100644 index e0af42a..0000000 --- a/testcases/miller_example/fort_02.90 +++ /dev/null @@ -1,86 +0,0 @@ -&BASIC - nrun = 100000000 - dt = 0.0075 - tmax = 100 - maxruntime = 356400 -/ -&GRID - pmaxe = 2 - jmaxe = 1 - pmaxi = 2 - jmaxi = 1 - Nx = 128 - Lx = 100 - Ny = 64 - Ly = 120 - Nz = 16 - Npol = 1 - SG = .false. -/ -&GEOMETRY - geom = 'miller' - q0 = 1.4 - shear = 0.0 - eps = 0.18 -/ -&OUTPUT_PAR - nsave_0d = 50 - nsave_1d = -1 - nsave_2d = Inf - nsave_3d = 100 - nsave_5d = 500 - write_doubleprecision = .false. - write_gamma = .true. - write_hf = .true. - write_phi = .true. - write_Na00 = .true. - write_Napj = .true. - write_Sapj = .false. - write_dens = .true. - write_temp = .true. - job2load = 1 -/ -&MODEL_PAR - ! Collisionality - CLOS = 0 - NL_CLOS = -1 - LINEARITY = 'nonlinear' - KIN_E = .false. - mu_x = 0.1 - mu_y = 0.1 - N_HD = 4 - mu_z = 1.0 - mu_p = 0 - mu_j = 0 - nu = 0.2 - tau_e = 1 - tau_i = 1 - sigma_e = 0.023338 - sigma_i = 1 - q_e = -1 - q_i = 1 - K_Ni = 2.22 - K_Ti = 6.92 - K_Ne = 0 - K_Te = 0 - GradB = 1 - CurvB = 1 - lambdaD = 0 - beta = 0 -/ -&COLLISION_PAR - collision_model = 'DG' - gyrokin_CO = .false. - interspecies = .true. - mat_file = 'null' -/ -&INITIAL_CON - INIT_OPT = 'ppj' - ACT_ON_MODES = 'none' - init_background = 0 - init_noiselvl = 1e-3 - iseed = 42 -/ -&TIME_INTEGRATION_PAR - numerical_scheme = 'RK4' -/ diff --git a/wk/analysis_gyacomo.m b/wk/analysis_gyacomo.m index 54353a9..d521f60 100644 --- a/wk/analysis_gyacomo.m +++ b/wk/analysis_gyacomo.m @@ -1,230 +1,230 @@ %% UNCOMMENT FOR TUTORIAL % gyacomodir = pwd; gyacomodir = gyacomodir(1:end-2); % get code directory % resdir = '.'; %Name of the directory where the results are located % JOBNUMMIN = 00; JOBNUMMAX = 10; %% addpath(genpath([gyacomodir,'matlab'])) % ... add addpath(genpath([gyacomodir,'matlab/plot'])) % ... add addpath(genpath([gyacomodir,'matlab/compute'])) % ... add addpath(genpath([gyacomodir,'matlab/load'])) % ... add %% Load the results 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); 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; %% PLOTS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% default_plots_options disp('Plots') FMT = '.fig'; if 1 %% Space time diagramm (fig 11 Ivanov 2020) % data.scale = 1;%/(data.Nx*data.Ny)^2; -i_ = 3; +i_ = 1; disp([num2str(data.TJOB_SE(i_)),' ',num2str(data.TJOB_SE(i_+1))]) disp([num2str(data.NU_EVOL(i_)),' ',num2str(data.NU_EVOL(i_+1))]) options.TAVG_0 = data.TJOB_SE(i_);%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 = 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); save_figure(data,fig,'.png') end if 0 %% statistical transport averaging for i_ = 1:2:21 % i_ = 3; disp([num2str(data.TJOB_SE(i_)),' ',num2str(data.TJOB_SE(i_+1))]) disp([num2str(data.NU_EVOL(i_)),' ',num2str(data.NU_EVOL(i_+1))]) options.T = [data.TJOB_SE(i_) data.TJOB_SE(i_+1)]; options.NPLOTS = 0; fig = statistical_transport_averaging(data,options); end end if 0 %% MOVIES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Options options.INTERP = 1; options.POLARPLOT = 0; -% options.NAME = '\phi'; -options.NAME = '\omega_z'; +options.NAME = '\phi'; +% options.NAME = '\omega_z'; % options.NAME = 'N_i^{00}'; % options.NAME = 'v_y'; % options.NAME = 'n_i^{NZ}'; % options.NAME = '\Gamma_x'; % options.NAME = 'n_i'; -options.PLAN = 'kxky'; +options.PLAN = 'xy'; % options.NAME = 'f_i'; % options.PLAN = 'sx'; options.COMP = 'avg'; % options.TIME = data.Ts5D(end-30:end); % options.TIME = data.Ts3D; options.TIME = [000:0.1:7000]; data.EPS = 0.1; data.a = data.EPS * 2000; options.RESOLUTION = 256; create_film(data,options,'.gif') end if 1 %% 2D snapshots % Options options.INTERP = 1; options.POLARPLOT = 0; options.AXISEQUAL = 1; options.NAME = '\phi'; % options.NAME = '\psi'; % options.NAME = 'n_i'; % options.NAME = 'N_i^{00}'; % options.NAME = 'T_i'; % options.NAME = '\Gamma_x'; % options.NAME = 'k^2n_e'; options.PLAN = 'xy'; % options.NAME 'f_i'; % options.PLAN = 'sx'; options.COMP = 'avg'; options.TIME = [1000 1800 2500 3000 4000]; data.a = data.EPS * 2e3; fig = photomaton(data,options); % save_figure(data,fig) end if 0 %% 3D plot on the geometry options.INTERP = 0; options.NAME = '\phi'; options.PLANES = [1]; options.TIME = [30]; options.PLT_MTOPO = 1; options.PLT_FTUBE = 0; data.EPS = 0.4; data.rho_o_R = 3e-3; % Sound larmor radius over Machine size ratio fig = show_geometry(data,options); save_figure(data,fig,'.png') end if 0 %% Kinetic distribution function sqrt(xy) (GENE vsp) options.SPAR = linspace(-3,3,32)+(6/127/2); options.XPERP = linspace( 0,6,32); % options.SPAR = gene_data.vp'; % options.XPERP = gene_data.mu'; options.iz = 'avg'; options.T = [250 600]; options.PLT_FCT = 'pcolor'; options.ONED = 0; options.non_adiab = 0; options.SPECIE = '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') end if 0 %% Hermite-Laguerre spectrum % options.TIME = 'avg'; options.P2J = 0; options.ST = 1; options.PLOT_TYPE = 'space-time'; options.NORMALIZED = 0; options.JOBNUM = 0; options.TIME = [1000]; options.specie = 'i'; options.compz = 'avg'; fig = show_moments_spectrum(data,options); % fig = show_napjz(data,options); % save_figure(data,fig,'.png'); end if 0 %% Time averaged spectrum options.TIME = [2000 3000]; options.NORM =1; options.NAME = '\phi'; % options.NAME = 'N_i^{00}'; % options.NAME ='\Gamma_x'; options.PLAN = 'kxky'; options.COMPZ = 'avg'; options.OK = 0; options.COMPXY = 'avg'; % avg/sum/max/zero/ 2D plot otherwise options.COMPT = 'avg'; options.PLOT = 'semilogy'; fig = spectrum_1D(data,options); % save_figure(data,fig,'.png') end if 0 %% 1D real plot options.TIME = [50 100 200]; options.NORM = 0; options.NAME = '\phi'; % options.NAME = 'n_i'; % options.NAME ='\Gamma_x'; % options.NAME ='s_y'; options.COMPX = 'avg'; options.COMPY = 'avg'; options.COMPZ = 1; options.COMPT = 1; options.MOVMT = 1; fig = real_plot_1D(data,options); % save_figure(data,fig,'.png') end if 0 %% Mode evolution options.NORMALIZED = 0; options.K2PLOT = [0.1 0.2 0.3 0.4]; options.TIME = [00:1200]; options.NMA = 1; options.NMODES = 5; options.iz = 'avg'; fig = mode_growth_meter(data,options); save_figure(data,fig,'.png') end if 0 %% ZF caracteristics (space time diagrams) TAVG_0 = 1200; TAVG_1 = 1500; % Averaging times duration % chose your field to plot in spacetime diag (uzf,szf,Gx) fig = ZF_spacetime(data,TAVG_0,TAVG_1); save_figure(data,fig,'.png') end if 0 %% Metric infos options.SHOW_FLUXSURF = 1; options.SHOW_METRICS = 0; fig = plot_metric(data,options); end if 0 %% linear growth rate for 3D fluxtube trange = [0 100]; nplots = 1; lg = compute_fluxtube_growth_rate(data,trange,nplots); end if 0 %% linear growth rate for 3D Zpinch trange = [5 15]; options.keq0 = 1; % chose to plot planes at k=0 or max options.kxky = 1; options.kzkx = 0; options.kzky = 1; [lg, fig] = compute_3D_zpinch_growth_rate(data,trange,options); save_figure(data,fig,'.png') end diff --git a/wk/header_2DZP_results.m b/wk/header_2DZP_results.m index 10428ea..544f479 100644 --- a/wk/header_2DZP_results.m +++ b/wk/header_2DZP_results.m @@ -1,203 +1,205 @@ %% Directory of the simulation gyacomodir = pwd; gyacomodir = gyacomodir(1:end-2); % get code directory % if 1% Local results resdir =''; resdir =''; resdir =''; resdir =''; % resdir ='debug/ppj_init'; %% nu = 5e-1 % Sugama % resdir ='Hallenbert_nu_5e-01/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_5e-01_SGGK';% also in 7x4 % resdir ='Hallenbert_nu_5e-01/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_5e-01_SGGK'; % resdir ='Hallenbert_nu_5e-01/200x32_7x4_L_120_kN_1.7_kT_0.425_nu_5e-01_SGGK';% also in 7x4 % resdir ='Hallenbert_nu_5e-01/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_5e-01_SGGK'; % resdir ='Hallenbert_nu_5e-01/200x32_5x3_L_120_kN_1.9_kT_0.475_nu_5e-01_SGGK';%also in 7x4 %% nu = 1e-1 % Landau % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/150x50_5x3_L_120_kN_1.6_kT_0.4_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.7_kT_0.425_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/150x50_5x3_L_120_kN_1.8_kT_0.45_nu_1e-01_LDGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.9_kT_0.475_nu_1e-01_LDGK'; % Sugama % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_1e-01_SGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.7_kT_0.425_nu_1e-01_SGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-01_SGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.9_kT_0.475_nu_1e-01_SGGK'; % Dougherty % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_1e-01_DGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_1e-01_DGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.7_kT_0.425_nu_1e-01_DGGK'; % resdir ='Hallenbert_nu_1e-01/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-01_DGGK'; %% nu = 5e-2 % resdir ='Hallenbert_nu_5e-02/200x32_11x6_L_120_kN_1.8_kT_0.45_nu_5e-02_SGGK';%For GENE benchmark % to analyse (added HD) % resdir ='Hallenbert_nu_5e-02/200x32_11x6_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGGK'; % testing various NL closures % resdir ='Hallenbert_nu_5e-02/200x32_7x4_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGGK'; % resdir ='Hallenbert_nu_5e-02/200x32_5x3_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/200x32_11x6_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/256x64_5x3_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/256x64_11x6_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/200x32_21x3_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/200x32_17x9_L_120_kN_1.8_kT_0.45_nu_5e-02_SGDK'; % resdir ='Hallenbert_nu_5e-02/200x32_5x3_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_DGGK'; % resdir ='Hallenbert_nu_5e-02/200x32_11x6_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_DGGK'; % resdir ='Hallenbert_nu_5e-02/128x32_5x3_Lx_120_Ly_60_kN_1.8_kT_0.45_nu_5e-02_FCGK'; %% nu = 1e-2 % Landau % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_1e-02_LDGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_1e-02_LDGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.7_kT_0.425_nu_1e-02_LDGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-02_LDGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.9_kT_0.475_nu_1e-02_LDGK'; % Sugama % resdir ='kobayashi_2015_fig1/150x150_5x3_L_100_kN_1.4_nu_5e-03_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_7x4_L_120_kN_1.5_kT_0.375_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.5_kT_0.375_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_7x4_L_120_kN_1.6_kT_0.4_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.7_kT_0.425_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/300x64_5x3_L_120_kN_1.7_kT_0.425_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_11x6_L_120_kN_1.7_kT_0.425_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-02_SGGK'; % To analyse (added HD) % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.9_kT_0.475_nu_1e-02_SGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_11x6_L_120_kN_1.9_kT_0.475_nu_1e-02_SGGK'; % Dougherty % resdir ='Hallenbert_nu_1e-02/200x32_7x4_L_120_kN_1.5_kT_0.375_nu_1e-02_DGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.6_kT_0.4_nu_1e-02_DGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_8x5_L_120_kN_1.6_kT_0.4_nu_0e+00_DGGK'; % resdir ='Hallenbert_nu_1e-02/200x32_5x3_L_120_kN_1.8_kT_0.45_nu_1e-02_DGGK'; %% nu = 5e-3 % resdir ='Hallenbert_nu_5e-03/200x32_5x3_Lx_120_Ly_60_kN_1.8_eta_0.25_nuSG_5e-03_muxy_5e-2'; % resdir ='Hallenbert_nu_5e-03/200x32_5x3_Lx_120_Ly_60_kN_1.8_eta_0.25_nuSG_5e-03_mux_5e-2_muy_6e-1'; % resdir ='Hallenbert_nu_5e-03/200x32_11x6_Lx_120_Ly_60_kN_1.8_eta_0.25_nuSG_5e-03_mux_5e-2_muy_6e-1'; % resdir ='Hallenbert_nu_5e-03/200x32_11x6_Lx_120_Ly_60_kN_1.8_eta_0.25_nuSG_5e-03_muxy_5e-2'; %% nu = 0 % resdir ='Hallenbert_fig2a/200x32_21x11_Lx_120_Ly_60_kN_1.6_eta_0.4_nu_0_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_11x6_Lx_120_Ly_60_kN_1.6_eta_0.4_nu_0_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_5x3_Lx_120_Ly_60_kN_1.6_eta_0.4_nu_0_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_5x3_Lx_120_Ly_60_kN_1.6_eta_0.4_nuDGGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_11x6_Lx_120_Ly_60_kN_1.6_eta_0.4_nuDGGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_5x3_Lx_120_Ly_60_kN_1.6_eta_0.4_nuSGGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_11x6_Lx_120_Ly_60_kN_1.6_eta_0.4_nuSGGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_5x3_Lx_120_Ly_60_kN_1.6_eta_0.4_nuLDGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2a/200x32_5x3_Lx_120_Ly_60_kN_1.6_eta_0.4_nuLRGK_0.1_muxy_1e-2'; % resdir ='Hallenbert_fig2b/200x32_11x6_Lx_240_Ly_120_kN_2.5_eta_0.25_nu_0_muxy_1e-1'; % resdir ='Hallenbert_fig2b/200x32_5x3_Lx_240_Ly_120_kN_2.5_eta_0.25_nu_0_muxy_1e-1'; % resdir ='Hallenbert_fig2b/200x32_5x3_Lx_240_Ly_120_kN_2.5_eta_0.25_nuSGGK_0.1_muxy_1e-1'; % resdir ='Hallenbert_fig2b/200x32_5x3_Lx_240_Ly_120_kN_2.5_eta_0.25_nuDGGK_0.1_muxy_1e-1'; % resdir ='Hallenbert_fig2b/200x32_5x3_Lx_240_Ly_120_kN_2.5_eta_0.25_nuLDGK_0.1_muxy_1e-1'; % resdir ='Hallenbert_fig2b/200x32_5x3_Lx_240_Ly_120_kN_2.5_eta_0.25_nuLRGK_0.1_muxy_1e-1'; % resdir ='Hallenbert_fig2c/200x32_11x6_Lx_120_Ly_60_kN_2.0_eta_0.25_nu_0_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_5x3_Lx_120_Ly_60_kN_2.0_eta_0.25_nu_0_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_5x3_Lx_120_Ly_60_kN_2.0_eta_0.25_nuSGGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_11x6_Lx_120_Ly_60_kN_2.0_eta_0.25_nuSGGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_5x3_Lx_120_Ly_60_kN_2.0_eta_0.25_nuDGGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_11x6_Lx_120_Ly_60_kN_2.0_eta_0.25_nuDGGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_5x3_Lx_120_Ly_60_kN_2.0_eta_0.25_nuLDGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_5x3_Lx_120_Ly_60_kN_2.0_eta_0.25_nuLRGK_0.1_muxy_5e-2'; % resdir ='Hallenbert_fig2c/200x32_9x5_Lx_120_Ly_60_kN_2.0_eta_0.25_nuLRGK_0.1_muxy_5e-2'; %% Transport scan % resdir = 'nu_0.1_transport_scan/colless_kn_1.7_to_2.0'; % resdir = 'nu_0.1_transport_scan/colless_kn_2.1_to_2.5'; % resdir = 'nu_0.1_transport_scan/LB_kn_2.0'; % resdir = 'nu_0.1_transport_scan/DG_kn_1.8_to_2.1'; % resdir = 'nu_0.1_transport_scan/DG_kn_2.2_to_2.5'; % resdir = 'nu_0.1_transport_scan/DG_conv_kN_1.9'; % resdir = 'nu_0.1_transport_scan/SG_kn_1.7_to_2.0'; % resdir = 'nu_0.1_transport_scan/SG_10x5_conv_test'; % resdir = 'nu_0.1_transport_scan/SG_kn_2.2_to_2.5'; % resdir = 'nu_0.1_transport_scan/LD_kn_2.0_to_2.5'; % resdir = 'nu_0.1_transport_scan/LD_kn_1.7_to_2.5'; % resdir = 'nu_0.1_transport_scan/LR_kn_1.7_to_2.0'; % resdir = 'nu_0.1_transport_scan/LR_kn_2.1_to_2.5'; % resdir = 'nu_0.1_transport_scan/colless_kn_2.2_Lx1.5'; % resdir = 'nu_0.1_transport_scan/colless_kn_2.2_HD'; % resdir = 'nu_0.1_transport_scan/colless_kn_1.6_HD'; % resdir = 'nu_0.1_transport_scan/large_box_kN_2.1_nu_0.1'; % resdir = 'nu_0.1_transport_scan/large_box_kN_2.0_nu_0.1'; % resdir = 'predator_prey_nu_scan/DG_Kn_1.7_nu_0.01'; % resdir = 'ZF_damping_linear_nu_0_20x10_kn_1.6_GK/LR_4x2_nu_0.1'; % resdir = 'ZF_damping_nu_0_20x10_kn_1.6_GK/HSG_4x2_nu_0.1'; % resdir = 'ZF_damping_nu_0_5x3_kn_2.5_GK/LR_4x2_nu_0.1'; % resdir = 'hacked_sugama/hacked_B_kn_1.6_200x32_L_120x60_nu_0.1'; % resdir = 'shearless_cyclone/200x32x24_5x4_Lx_120_Ly_60_q0_1.4_e_0.18_kN_2.22_kT_6.9_nuLR_0.01_adiab_e'; % resdir = 'shearless_cyclone/no_sg_128x32x36_6x3_Lx_120_Ly_60_q0_1.4_e_0.18_kN_2.22_kT_6.9_adiab_e'; % resdir = 'shearless_cyclone/sgrid_128x64x32_4x2_Lx_100_Ly_120_q0_1.4_e_0.18_kN_2.22_kT_6.96_adiab_e'; % resdir = 'shearless_cyclone/sgrid_128x64x32_4x2_Lx_100_Ly_120_q0_1.4_e_0.18_kN_1.78_kT_5.52_adiab_e'; % resdir = 'linear_shearless_cyclone/4_2_cyclone_1.0'; % resdir = 'linear_shearless_cyclone/test_fmom'; % else% Marconi results % resdir =''; % resdir =''; % resdir ='';fd % resdir =''; % resdir ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A_new/300x300_5x3_L_120_kN_1.6667_nu_1e-01_SGGK/out.txt'; % % resdir ='/marconi_scratch/userexternal/ahoffman/HeLaZ/results/simulation_A/300x150_L_120_P_8_J_4_eta_0.6_nu_1e-01_SGGK_mu_0e+00/out.txt'; % % BASIC.RESDIR = ['../',resdir(46:end-8),'/']; % MISCDIR = ['/misc/HeLaZ_outputs/',resdir(46:end-8),'/']; % end %% ZPINCH rerun % resdir ='Zpinch_rerun/Kn_2.5_200x48x5x3'; % resdir ='Zpinch_rerun/Kn_2.5_256x128x5x3'; % resdir ='Zpinch_rerun/Kn_2.5_312x196x5x3_Lx_400_Ly_200'; % resdir ='Zpinch_rerun/Kn_2.5_256x64x5x3'; % resdir ='Zpinch_rerun/Kn_2.0_200x48x9x5_large_box'; % resdir ='Zpinch_rerun/Kn_2.0_256x64x9x5_Lx_240_Ly_120'; % resdir ='Zpinch_rerun/Kn_1.6_256x128x7x4'; % resdir ='Zpinch_rerun/Kn_1.6_200x64x11x6'; % resdir ='Zpinch_rerun/Kn_1.6_200x64x11x6_conv'; % resdir ='Zpinch_rerun/Kn_1.6_200x64x11x6_mu_0.5'; % resdir ='Zpinch_rerun/Kn_1.6_256x128x21x11'; %% nu scan % resdir = 'Zpinch_rerun/kN_2.2_coll_scan_128x48x5x3'; % resdir = 'Zpinch_rerun/Ultra_HD_312x196x5x3'; % resdir = 'Zpinch_rerun/UHD_nu_001_LDGK'; % resdir = 'Zpinch_rerun/UHD_nu_01_LDGK'; % resdir = 'Zpinch_rerun/UHD_nu_1_LDGK'; % resdir ='Zpinch_rerun/kN_1.7_SGGK_conv_200x32x7x3_nu_0.01'; % resdir ='Zpinch_rerun/kN_1.7_LDGKii_200x32x7x3_nu_scan'; % resdir ='Zpinch_rerun/nu_0.1_LDGKii_200x48x7x4_kN_scan'; % resdir ='Zpinch_rerun/nu_0.1_FCGK_200x48x5x3_kN_scan'; % resdir ='Zpinch_rerun/nu_0.1_SGGK_200x48x5x3_kN_scan'; % resdir = 'Zpinch_rerun/kN_1.7_FCGK_200x32x5x3_nu_scan'; -% resdir = 'Zpinch_rerun/kN_1.7_SGGK_200x32x7x4_nu_scan'; -% resdir = 'Zpinch_rerun/kN_1.7_SGGK_200x32x5x3_nu_scan'; -resdir = 'Zpinch_rerun/kN_1.7_SGGK_256x64x5x3_nu_scan'; -%% +% % resdir = 'Zpinch_rerun/kN_1.7_SGGK_200x32x7x4_nu_scan'; +resdir = 'Zpinch_rerun/kN_2.2_SGGK_200x32x5x3_nu_scan'; +% resdir = 'Zpinch_rerun/kN_1.7_SGGK_256x64x5x3_nu_scan'; +%% Convergence cases kN = (1.6 2.2) nu = (0.01 1.0) +resdir = 'Zpinch_rerun/convcoll_Kn_1.6_200x32x3x2_mu_0.01'; + JOBNUMMIN = 00; JOBNUMMAX = 10; resdir = ['results/',resdir]; run analysis_gyacomo \ No newline at end of file diff --git a/wk/header_3D_results.m b/wk/header_3D_results.m index 779dd11..63cb0bf 100644 --- a/wk/header_3D_results.m +++ b/wk/header_3D_results.m @@ -1,64 +1,66 @@ % Directory of the code "mypathtoHeLaZ/HeLaZ/" gyacomodir = '/home/ahoffman/gyacomo/'; % Directory of the simulation (from results) % if 1% Local results % resdir ='volcokas/64x32x16x5x3_kin_e_npol_1'; %% Dimits % resdir ='shearless_cyclone/128x64x16x5x3_Dim_90'; % resdir ='shearless_cyclone/128x64x16x9x5_Dim_scan/128x64x16x9x5_Dim_60'; % resdir ='shearless_cyclone/128x64x16x5x3_Dim_scan/128x64x16x5x3_Dim_70'; % resdir ='shearless_cyclone/64x32x16x5x3_Dim_scan/64x32x16x5x3_Dim_70'; %% AVS % resdir = 'volcokas/64x32x16x5x3_kin_e_npol_1'; % resdir = 'volcokas/64x32x16x5x3_kin_e_npol_1'; % resdir = 'shearless_cyclone/64x32x80x5x3_CBC_Npol_5_kine'; % resdir = 'shearless_cyclone/96x32x160x5x3_CBC_Npol_10_kine'; % resdir = 'shearless_cyclone/64x32x160x5x3_CBC_Npol_10_kine'; % resdir = 'shearless_cyclone/96x32x160x5x3_CBC_Npol_10_kine'; %% shearless CBC % resdir ='shearless_cyclone/64x32x16x5x3_CBC_080'; % resdir ='shearless_cyclone/64x32x16x5x3_CBC_scan/64x32x16x5x3_CBC_100'; % resdir ='shearless_cyclone/64x32x16x5x3_CBC_120'; % resdir ='shearless_cyclone/64x32x16x9x5_CBC_080'; % resdir ='shearless_cyclone/64x32x16x9x5_CBC_100'; % resdir ='shearless_cyclone/64x32x16x9x5_CBC_120'; % resdir = 'shearless_cyclone/64x32x16x5x3_CBC_CO/64x32x16x5x3_CBC_LRGK'; %% CBC % resdir = 'CBC/64x32x16x5x3'; % resdir = 'CBC/64x128x16x5x3'; % resdir = 'CBC/128x64x16x5x3'; % resdir = 'CBC/96x96x16x3x2_Nexc_6'; % resdir = 'CBC/128x96x16x3x2'; % resdir = 'CBC/192x96x16x3x2'; % resdir = 'CBC/192x96x24x13x7'; % resdir = 'CBC/kT_11_128x64x16x5x3'; % resdir = 'CBC/kT_9_256x128x16x3x2'; % resdir = 'CBC/kT_4.5_128x64x16x13x3'; % resdir = 'CBC/kT_4.5_192x96x24x13x7'; % resdir = 'CBC/kT_4.5_128x64x16x13x7'; % resdir = 'CBC/kT_4.5_128x96x24x15x5'; % resdir = 'CBC/kT_5.3_192x96x24x13x7'; % resdir = 'CBC/kT_13_large_box_128x64x16x5x3'; % resdir = 'CBC/kT_11_96x64x16x5x3_ky_0.02'; % resdir = 'CBC/kT_scan_128x64x16x5x3'; % resdir = 'CBC/kT_scan_192x96x16x3x2'; % resdir = 'CBC/kT_13_96x96x16x3x2_Nexc_6'; % resdir = 'dbg/nexc_dbg'; % resdir = 'CBC/NM_F4_kT_4.5_192x64x24x6x4'; % resdir = 'CBC_Ke_EM/192x96x24x5x3'; % resdir = 'CBC_Ke_EM/96x48x16x5x3'; % resdir = 'CBC_Ke_EM/minimal_res'; %% KBM % resdir = 'NL_KBM/192x64x24x5x3'; %% Linear CBC % resdir = 'linear_CBC/20x2x32_21x11_Lx_62.8319_Ly_31.4159_q0_1.4_e_0.18_s_0.8_kN_2.22_kT_5.3_nu_1e-02_DGDK_adiabe'; -resdir = 'testcases/miller_example'; +% resdir = 'testcases/miller_example'; +resdir = 'Miller/128x256x3x2_CBC_dt_5e-3'; +% resdir = ['results/',resdir]; JOBNUMMIN = 00; JOBNUMMAX = 10; run analysis_gyacomo diff --git a/wk/lin_EPY.m b/wk/lin_EPY.m index 8a00601..429c132 100644 --- a/wk/lin_EPY.m +++ b/wk/lin_EPY.m @@ -1,178 +1,178 @@ %% QUICK RUN SCRIPT for linear entropy mode in a Zpinch % This script create a directory in /results and run a simulation directly % from matlab framework. It is meant to run only small problems in linear % for benchmark and debugging purpose since it makes matlab "busy" % SIMID = 'lin_EPY'; % Name of the simulation RUN = 1; % To run or just to load addpath(genpath('../matlab')) % ... add default_plots_options PROGDIR = '/home/ahoffman/gyacomo/'; EXECNAME = 'gyacomo'; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Set Up parameters CLUSTER.TIME = '99:00:00'; % allocation time hh:mm:ss %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% PHYSICAL PARAMETERS NU = 0.1; % Collision frequency TAU = 1.0; % e/i temperature ratio -K_Ne = 2.0; % ele Density ''' -K_Te = 0.5; % ele Temperature ''' -K_Ni = 2.0; % ion Density gradient drive -K_Ti = 0.5; % ion Temperature ''' +K_Ne = 2.2; % ele Density ''' +K_Te = K_Ne/4; % ele Temperature ''' +K_Ni = K_Ne; % ion Density gradient drive +K_Ti = K_Ni/4; % ion Temperature ''' SIGMA_E = 0.0233380; % mass ratio sqrt(m_a/m_i) (correct = 0.0233380) KIN_E = 1; % 1: kinetic electrons, 2: adiabatic electrons BETA = 0.0; % electron plasma beta %% GRID PARAMETERS P = 4; J = P/2; PMAXE = P; % Hermite basis size of electrons JMAXE = J; % Laguerre " PMAXI = P; % " ions JMAXI = J; % " NX = 2; % real space x-gridpoints NY = 100; % '' y-gridpoints LX = 2*pi/0.8; % Size of the squared frequency domain -LY = 2*pi/0.05; % Size of the squared frequency domain +LY = 120;%2*pi/0.05; % Size of the squared frequency domain NZ = 1; % number of perpendicular planes (parallel grid) NPOL = 1; SG = 0; % Staggered z grids option %% GEOMETRY GEOMETRY= 'Z-pinch'; % Z-pinch overwrites q0, shear and eps Q0 = 1.4; % safety factor SHEAR = 0.8; % magnetic shear NEXC = 1; % To extend Lx if needed (Lx = Nexc/(kymin*shear)) EPS = 0.18; % inverse aspect ratio %% TIME PARMETERS TMAX = 50; % Maximal time unit DT = 1e-2; % Time step SPS0D = 1; % Sampling per time unit for 2D arrays SPS2D = 0; % Sampling per time unit for 2D arrays SPS3D = 2; % Sampling per time unit for 2D arrays SPS5D = 1/5; % Sampling per time unit for 5D arrays SPSCP = 0; % Sampling per time unit for checkpoints JOB2LOAD= -1; %% OPTIONS LINEARITY = 'linear'; % activate non-linearity (is cancelled if KXEQ0 = 1) % Collision operator % (LB:L.Bernstein, DG:Dougherty, SG:Sugama, LR: Lorentz, LD: Landau) CO = 'LD'; GKCO = 1; % gyrokinetic operator ABCO = 1; % interspecies collisions INIT_ZF = 0; ZF_AMP = 0.0; CLOS = 0; % Closure model (0: =0 truncation, 1: v^Nmax 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 %% OUTPUTS W_DOUBLE = 1; W_GAMMA = 1; W_HF = 1; W_PHI = 1; W_NA00 = 1; W_DENS = 1; W_TEMP = 1; W_NAPJ = 1; W_SAPJ = 0; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % unused HD_CO = 0.0; % Hyper diffusivity cutoff ratio MU = 0.0; % Hyperdiffusivity coefficient INIT_BLOB = 0; WIPE_TURB = 0; ACT_ON_MODES = 0; MU_X = MU; % MU_Y = MU; % N_HD = 4; MU_Z = 2.0; % MU_P = 0.0; % MU_J = 0.0; % LAMBDAD = 0.0; NOISE0 = 0.0e-5; % Init noise amplitude BCKGD0 = 1.0; % Init background GRADB = 1.0; CURVB = 1.0; %%------------------------------------------------------------------------- %% RUN setup % system(['rm fort*.90']); % Run linear simulation if RUN system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ',PROGDIR,'bin/',EXECNAME,' 1 6 1 0; cd ../../../wk']) end %% Load results %% filename = [SIMID,'/',PARAMS,'/']; LOCALDIR = [PROGDIR,'results/',filename,'/']; % Load outputs from jobnummin up to jobnummax JOBNUMMIN = 00; JOBNUMMAX = 00; data = compile_results(LOCALDIR,JOBNUMMIN,JOBNUMMAX); %Compile the results from first output found to JOBNUMMAX if existing %% Short analysis if 1 %% linear growth rate (adapted for 2D zpinch and fluxtube) trange = [0.5 1]*data.Ts3D(end); nplots = 1; % 1 for only growth rate and error, 2 for omega local evolution, 3 for plot according to z lg = compute_fluxtube_growth_rate(data,trange,nplots); [gmax, kmax] = max(lg.g_ky(:,end)); [gmaxok, kmaxok] = max(lg.g_ky(:,end)./lg.ky); msg = sprintf('gmax = %2.2f, kmax = %2.2f',gmax,lg.ky(kmax)); disp(msg); msg = sprintf('gmax/k = %2.2f, kmax/k = %2.2f',gmaxok,lg.ky(kmaxok)); disp(msg); end if 0 %% Ballooning plot options.time_2_plot = [120]; options.kymodes = [0.9]; options.normalized = 1; % options.field = 'phi'; fig = plot_ballooning(data,options); end if 0 %% Hermite-Laguerre spectrum % options.TIME = 'avg'; options.P2J = 1; options.ST = 1; options.PLOT_TYPE = 'space-time'; % options.PLOT_TYPE = 'Tavg-1D'; % options.PLOT_TYPE = 'Tavg-2D'; options.NORMALIZED = 0; options.JOBNUM = 0; options.TIME = [0 50]; options.specie = 'i'; options.compz = 'avg'; fig = show_moments_spectrum(data,options); % fig = show_napjz(data,options); save_figure(data,fig) end if 0 %% linear growth rate for 3D Zpinch (kz fourier transform) trange = [0.5 1]*data.Ts3D(end); options.keq0 = 1; % chose to plot planes at k=0 or max options.kxky = 1; options.kzkx = 0; options.kzky = 0; [lg, fig] = compute_3D_zpinch_growth_rate(data,trange,options); save_figure(data,fig) end if 0 %% Mode evolution options.NORMALIZED = 0; options.K2PLOT = 1; options.TIME = [0:1000]; options.NMA = 1; options.NMODES = 1; options.iz = 'avg'; fig = mode_growth_meter(data,options); save_figure(data,fig,'.png') end if 0 %% RH TEST ikx = 2; iky = 2; t0 = 0; t1 = data.Ts3D(end); [~, it0] = min(abs(t0-data.Ts3D));[~, it1] = min(abs(t1-data.Ts3D)); plt = @(x) squeeze(mean(real(x(iky,ikx,:,it0:it1)),3))./squeeze(mean(real(x(iky,ikx,:,it0)),3)); figure plot(data.Ts3D(it0:it1), plt(data.PHI)); xlabel('$t$'); ylabel('$\phi_z(t)/\phi_z(0)$') title(sprintf('$k_x=$%2.2f, $k_y=$%2.2f',data.kx(ikx),data.ky(iky))) end