scripts update

matlab/plot/plot_ballooning.m

@@ -57,7 +57,7 @@ function [FIG] = plot_ballooning(data,options)
         % normalize real and imaginary parts at chi =0
         if options.normalized
             [~,idxLFS] = min(abs(b_angle -0));
-            normalization = phib( idxLFS);
+            normalization = abs(phib( idxLFS));
         else
             normalization = 1;
         end
@@ -71,7 +71,7 @@ function [FIG] = plot_ballooning(data,options)
         legend('real','imag','norm')
         xlabel('$\chi / \pi$')
         ylabel('$\phi_B (\chi)$');
-        title(['$k_y=',sprintf('%1.1f',data.ky(iky)),...
+        title(['$k_y=',sprintf('%2.2f',data.ky(iky)),...
                ',t_{avg}\in [',sprintf('%1.1f',data.Ts3D(it0)),',',sprintf('%1.1f',data.Ts3D(it1)),']$']);
 
         if data.BETA > 0
@@ -87,13 +87,13 @@ function [FIG] = plot_ballooning(data,options)
             end
 
             subplot(numel(ikyarray),ncol,ncol*(iplot-1)+2)
-            plot(b_angle/pi, psib_real,'-b'); hold on;
-            plot(b_angle/pi, psib_imag ,'-r');
-            plot(b_angle/pi, sqrt(psib_real .^2 + psib_imag.^2),'-k');
+            plot(b_angle/pi, psib_real/ normalization,'-b'); hold on;
+            plot(b_angle/pi, psib_imag/ normalization ,'-r');
+            plot(b_angle/pi, sqrt(psib_real .^2 + psib_imag.^2)/ normalization,'-k');
             legend('real','imag','norm')
             xlabel('$\chi / \pi$')
             ylabel('$\psi_B (\chi)$');
-            title(['$k_y=',sprintf('%1.1f',data.ky(iky)),...
+            title(['$k_y=',sprintf('%2.2f',data.ky(iky)),...
                    ',t_{avg}\in [',sprintf('%1.1f',data.Ts3D(it0)),',',sprintf('%1.1f',data.Ts3D(it1)),']$']);
         end
         

matlab/plot/show_moments_spectrum.m

@@ -6,7 +6,7 @@ Ji = DATA.Ji;
 Nipj = sum(sum(sum(abs(DATA.Nipj),3),4),5);
 Nipj = squeeze(Nipj);
 
-if DATA.K_E
+if DATA.KIN_E
 Pe = DATA.Pe;
 Je = DATA.Je;
 Nepj = sum(sum(sum(abs(DATA.Nepj),3),4),5);
@@ -100,7 +100,7 @@ else
     for ip2j = 1:numel(p2ji)
         Ni_ST(ip2j,:) = Ni_ST(ip2j,:)/weights(ip2j);
     end
-    if DATA.K_E
+    if DATA.KIN_E
     % same for electrons!!
     [JJ,PP] = meshgrid(Je,Pe);
     P2Je = PP + 2*JJ;
@@ -131,7 +131,7 @@ else
     else
     plt = @(x,ip2j) x;
     end
-    if DATA.K_E
+    if DATA.KIN_E
     subplot(2,1,1)
     end
         imagesc(DATA.Ts5D,p2ji,plt(Ni_ST,1:numel(p2ji))); 
@@ -140,7 +140,7 @@ else
 %         set(pclr, 'edgecolor','none'); hold on;
     xlabel('$t$'); ylabel('$p+2j$')
     title('$\langle\sum_k |N_i^{pj}|\rangle_{p+2j=const}$')
-    if DATA.K_E
+    if DATA.KIN_E
     subplot(2,1,2)
         imagesc(DATA.Ts5D,p2je,plt(Ne_ST,1:numel(p2ji))); 
         set(gca,'YDir','normal')

matlab/process_field.m

@@ -156,6 +156,26 @@ switch OPTIONS.NAME
                 FIELD(:,:,it) = squeeze(compr(tmp));
             end                
         end
+    case '\psi'
+        NAME = 'psi';
+        if COMPDIM == 3
+            for it = 1:numel(FRAMES)
+                tmp = squeeze(compr(DATA.PSI(:,:,:,FRAMES(it))));
+                FIELD(:,:,it) = squeeze(process(tmp));
+            end
+        else
+            if REALP
+                tmp = zeros(Ny,Nx,Nz);
+            else
+                tmp = zeros(DATA.Nky,DATA.Nkx,Nz);
+            end
+            for it = 1:numel(FRAMES)
+                for iz = 1:numel(DATA.z)
+                    tmp(:,:,iz) = squeeze(process(DATA.PSI(:,:,iz,FRAMES(it))));
+                end
+                FIELD(:,:,it) = squeeze(compr(tmp));
+            end                
+        end
     case 'N_i^{00}'
         NAME = 'Ni00';
         if COMPDIM == 3

wk/analysis_HeLaZ.m

@@ -64,16 +64,17 @@ end
 if 0
 %% 2D snapshots
 % Options
-options.INTERP    = 1;
+options.INTERP    = 0;
 options.POLARPLOT = 0;
 options.AXISEQUAL = 0;
 % options.NAME      = '\phi';
+options.NAME      = '\psi';
 % options.NAME      = 'n_e';
-options.NAME      = 'N_i^{00}';
+% options.NAME      = 'N_i^{00}';
 % options.NAME      = 'T_i';
 % options.NAME      = '\Gamma_x';
 % options.NAME      = 'k^2n_e';
-options.PLAN      = 'xy';
+options.PLAN      = 'xz';
 % options.NAME      'f_i';
 % options.PLAN      = 'sx';
 options.COMP      = 'avg';
@@ -97,16 +98,16 @@ end
 
 if 0
 %% Kinetic distribution function sqrt(<f_a^2>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        = 1;
-options.T         = [500 1000];
-options.PLT_FCT   = 'contour';
+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         = [40];
+options.PLT_FCT   = 'pcolor';
 options.ONED      = 0;
 options.non_adiab = 0;
-options.SPECIE    = 'i';
+options.SPECIE    = 'e';
 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')

wk/analysis_gene.m

+helazdir = '/home/ahoffman/HeLaZ/';
+addpath(genpath([helazdir,'matlab'])) % ... add
+addpath(genpath([helazdir,'matlab/plot'])) % ... add
+addpath(genpath([helazdir,'matlab/compute'])) % ... add
+addpath(genpath([helazdir,'matlab/load'])) % ... add
+
 % folder = '/misc/gene_results/shearless_cyclone/miller_output_1.0/';
 % folder = '/misc/gene_results/shearless_cyclone/miller_output_0.8/';
 % folder = '/misc/gene_results/shearless_cyclone/s_alpha_output_1.0/';
@@ -9,7 +15,8 @@
 % folder = '/misc/gene_results/HP_fig_2a_mu_1e-2/';
 % folder = '/misc/gene_results/HP_fig_2b_mu_5e-2/';
 % folder = '/misc/gene_results/HP_fig_2c_mu_5e-2/';
-% folder = '/misc/gene_results/LD_zpinch_1.6/';
+folder = '/misc/gene_results/LD_zpinch_1.6/';
+% folder = '/misc/gene_results/ZP_HP_kn_1.6_nuv_3.2/';
 % folder = '/misc/gene_results/ZP_HP_kn_1.6_nuv_3.2/';
 % folder = '/misc/gene_results/ZP_kn_2.5_large_box/';
 % folder = '/misc/gene_results/CBC/128x64x16x24x12/';
@@ -18,7 +25,7 @@
 % folder = '/misc/gene_results/CBC/KT_5.3_128x64x16x24x12_00/';
 % folder = '/misc/gene_results/CBC/KT_4.5_128x64x16x24x12_01/';
 % folder = '/misc/gene_results/CBC/KT_13_128x64x16x24x12/';
-folder = '/misc/gene_results/CBC/KT_13_large_box_128x64x16x24x12/';
+% folder = '/misc/gene_results/CBC/KT_13_large_box_128x64x16x24x12/';
 gene_data = load_gene_data(folder);
 gene_data = invert_kxky_to_kykx_gene_results(gene_data);
 if 1

wk/header_3D_results.m

@@ -44,7 +44,8 @@ helazdir = '/home/ahoffman/HeLaZ/';
 % outfile = 'CBC/192x96x16x3x2';
 % outfile = 'CBC/128x64x16x5x3';
 % outfile = 'CBC/kT_9_128x64x16x5x3';
-outfile = 'CBC/kT_13_large_box_128x64x16x5x3';
+outfile = 'CBC/kT_4.5_128x64x16x13x7';
+% outfile = 'CBC/kT_13_large_box_128x64x16x5x3';
 JOBNUMMIN = 00; JOBNUMMAX = 06;
 
 % outfile = 'CBC/kT_scan_128x64x16x5x3';

wk/quick_run.m

@@ -32,7 +32,7 @@ PMAXE   = P;     % Hermite basis size of electrons
 JMAXE   = J;     % Laguerre "
 PMAXI   = P;     % " ions
 JMAXI   = J;     % "
-NX      = 6;    % real space x-gridpoints
+NX      = 20;    % real space x-gridpoints
 NY      = 2;     %     ''     y-gridpoints
 LX      = 2*pi/0.1;   % Size of the squared frequency domain
 LY      = 2*pi/0.25;     % Size of the squared frequency domain
@@ -48,7 +48,7 @@ SHEAR   = 0.8;    % magnetic shear (Not implemented yet)
 EPS     = 0.18;    % inverse aspect ratio
 %% TIME PARMETERS
 TMAX    = 40;  % Maximal time unit
-DT      = 5e-3;   % Time step
+DT      = 2e-3;   % Time step
 SPS0D   = 1;      % Sampling per time unit for 2D arrays
 SPS2D   = 0;      % Sampling per time unit for 2D arrays
 SPS3D   = 1;      % Sampling per time unit for 2D arrays
@@ -60,13 +60,13 @@ LINEARITY = 'linear';   % activate non-linearity (is cancelled if KXEQ0 = 1)
 % Collision operator
 % (LB:L.Bernstein, DG:Dougherty, SG:Sugama, LR: Lorentz, LD: Landau)
 CO      = 'DG';
-GKCO    = 1; % gyrokinetic operator
+GKCO    = 0; % 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= 'phi';   % Start simulation with a noisy mom00/phi/allmom
+INIT_OPT= 'mom00';   % Start simulation with a noisy mom00/phi/allmom
 %% OUTPUTS
 W_DOUBLE = 1;
 W_GAMMA  = 1; W_HF     = 1;
@@ -82,12 +82,12 @@ INIT_BLOB = 0; WIPE_TURB = 0; ACT_ON_MODES = 0;
 MU_X    = MU;     %
 MU_Y    = MU;     %
 N_HD    = 4;
-MU_Z    = 2.0;     %
+MU_Z    = 1.0;     %
 MU_P    = 0.0;     %
 MU_J    = 0.0;     %
 LAMBDAD = 0.0;
-NOISE0  = 0.0e-5; % Init noise amplitude
-BCKGD0  = 1.0;    % Init background
+NOISE0  = 1.0e-5; % Init noise amplitude
+BCKGD0  = 0.0;    % Init background
 GRADB   = 1.0;
 CURVB   = 1.0;
 %%-------------------------------------------------------------------------
@@ -123,10 +123,10 @@ 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
+if 1
 %% Ballooning plot
 options.time_2_plot = [120];
-options.kymodes     = [0.1 0.2 0.3];
+options.kymodes     = [0.25];
 options.normalized  = 1;
 % options.field       = 'phi';
 fig = plot_ballooning(data,options);