diff --git a/matlab/create_gif.m b/matlab/create_gif.m index 179a71e591aaea02bcd8b0c61eed05dfb84cf383..38bdbfed693c98bf0dd0b3080109b3223bf8cbed 100644 --- a/matlab/create_gif.m +++ b/matlab/create_gif.m @@ -23,9 +23,9 @@ fig = figure; shading interp; hold on - n = 0; + in = 1; nbytes = fprintf(2,'frame %d/%d',n,numel(FIELD(1,1,:))); - for n = 1:numel(FIELD(1,1,:)) % time loop + for n = FRAMES % loop over selected frames pclr = pcolor(X,Y,FIELD(:,:,n)); shading interp; % frame plot set(pclr, 'edgecolor','none'); title([FIELDNAME,', $t \approx$', sprintf('%.3d',ceil(T(n)))]); @@ -35,7 +35,7 @@ fig = figure; im = frame2im(frame); [imind,cm] = rgb2ind(im,64); % Write to the GIF File - if n == 1 + if in == 1 imwrite(imind,cm,GIFNAME,'gif', 'Loopcount',inf); else imwrite(imind,cm,GIFNAME,'gif','WriteMode','append', 'DelayTime',DELAY); @@ -46,6 +46,7 @@ fig = figure; nbytes = nbytes - 1; end nbytes = fprintf(2,'frame %d/%d',n,numel(FIELD(1,1,:))); + in = in + 1; end disp(' ') disp(['Gif saved @ : ',GIFNAME]) diff --git a/matlab/load_2D_data.m b/matlab/load_2D_data.m index ddb42e1c7e35694fd3209bb4afcdcee065d3f16c..32c192fe3655a8cc88946a12698b402f4f8d7b9c 100644 --- a/matlab/load_2D_data.m +++ b/matlab/load_2D_data.m @@ -1,4 +1,4 @@ -function [ data, kr, kz, time ] = load_2D_data( filename, variablename ) +function [ data, kr, kz, time, dt ] = load_2D_data( filename, variablename ) %LOAD_2D_DATA load a 2D variable stored in a hdf5 result file from HeLaZ time = h5read(filename,'/data/var2d/time'); kr = h5read(filename,['/data/var2d/',variablename,'/coordkr']); diff --git a/matlab/load_5D_data.m b/matlab/load_5D_data.m index 34cb449f5b50e0fa88489eb1c47febbe14aef195..767ca2b5e43644e826f9c4f0b72fd096b77ce0c7 100644 --- a/matlab/load_5D_data.m +++ b/matlab/load_5D_data.m @@ -1,4 +1,4 @@ -function [ data, p, j, kr, kz, time ] = load_5D_data( filename, variablename ) +function [ data, p, j, kr, kz, time, dt ] = load_5D_data( filename, variablename ) %LOAD_5D_DATA load a 5D variable stored in a hdf5 result file from HeLaZ time = h5read(filename,'/data/var5d/time'); p = h5read(filename,['/data/var5d/',variablename,'/coordp']); diff --git a/src/diagnose.F90 b/src/diagnose.F90 index 0c8a7dcb2cca0f8ad9bcadc2b38922f1fa49ac82..7b217c7eb6a6ca5de38139f34858389ce0efdaa4 100644 --- a/src/diagnose.F90 +++ b/src/diagnose.F90 @@ -175,7 +175,7 @@ SUBROUTINE diagnose(kstep) ! Terminal info IF (MOD(cstep, INT(1.0/dt)) == 0) THEN - WRITE(*,"(F4.0,A,F4.0)") time,"/",tmax + WRITE(*,"(F5.0,A,F5.0)") time,"/",tmax ENDIF ! 2.1 0d history arrays diff --git a/src/model_mod.F90 b/src/model_mod.F90 index 9d6f3ff916578f8831d95a34f19b4941fdfde6f9..58604cf58b28a9922094500c838d6c54fa8732a0 100644 --- a/src/model_mod.F90 +++ b/src/model_mod.F90 @@ -6,6 +6,7 @@ MODULE model PRIVATE INTEGER, PUBLIC, PROTECTED :: CO = 0 ! Collision Operator + LOGICAL, PUBLIC, PROTECTED :: DK = 0 ! Drift kinetic model LOGICAL, PUBLIC, PROTECTED :: NON_LIN = .true. ! To turn on non linear bracket term REAL(dp), PUBLIC, PROTECTED :: mu = 0._dp ! Hyperdiffusivity coefficient (for num. stability) REAL(dp), PUBLIC, PROTECTED :: nu = 1._dp ! Collision frequency @@ -31,7 +32,7 @@ CONTAINS USE prec_const IMPLICIT NONE - NAMELIST /MODEL_PAR/ CO, NON_LIN, mu, nu, tau_e, tau_i, sigma_e, sigma_i, & + NAMELIST /MODEL_PAR/ CO, DK, NON_LIN, mu, nu, tau_e, tau_i, sigma_e, sigma_i, & q_e, q_i, eta_n, eta_T, eta_B, lambdaD READ(lu_in,model_par) diff --git a/wk/analysis.m b/wk/analysis.m index 229dbd1963c75a2511372d0cee086616720abb2c..fbcac91a340e4d3d28b79e55ae28d506507769c1 100644 --- a/wk/analysis.m +++ b/wk/analysis.m @@ -2,14 +2,14 @@ JOBNUM = 00; filename = [BASIC.SIMID,'_','%.2d.h5']; filename = sprintf(filename,JOBNUM); disp(['Analysing ',filename]) -[Nipj, p_, j_, kr, kz, Ts] = load_5D_data(filename, 'moments_i'); -Nepj = load_5D_data(filename, 'moments_e'); +[Nipj, p_, j_, kr, kz, Ts, dt] = load_5D_data(filename, 'moments_i'); +Nepj = load_5D_data(filename, 'moments_e'); Ni = squeeze(Nipj(1,1,:,:,:)); Ne = squeeze(Nepj(1,1,:,:,:)); PH = load_2D_data(filename, 'phi'); Ts = Ts'; Ns = numel(Ts); -dt = mean(diff(Ts)); +dt_samp = mean(diff(Ts)); if strcmp(OUTPUTS.write_non_lin,'.true.') Sipj = load_5D_data(filename, 'Sipj'); Sepj = load_5D_data(filename, 'Sepj'); @@ -45,36 +45,54 @@ for it = 1:numel(PH(1,1,:)) end %% Post processing -phi_ST_r = zeros(Nr,Ns); % Space-Time diagram of ES potential -ne_ST_r = zeros(Nr,Ns); % Space-Time diagram of density -ni_ST_r = zeros(Nr,Ns); % Space-Time diagram of density -phi_ST_z = zeros(Nz,Ns); % Space-Time diagram of ES potential -ne_ST_z = zeros(Nz,Ns); % Space-Time diagram of density -ni_ST_z = zeros(Nz,Ns); % Space-Time diagram of density -phi_00 = zeros(1,Ns); % Time evolution of ES potential at origin -ne_00 = zeros(1,Ns); % Time evolution of density at origin -ni_00 = zeros(1,Ns); % Time evolution of density at origin -[~,ir0] = min(abs(r)); [~,iz0] = min(abs(z)); -Ne_11 = zeros(1,Ns); % Time evolution of F density at 1,1 -Ni_11 = zeros(1,Ns); % Time evolution of F density at 1,1 -[~,ikr1] = min(abs(kr-round(1/dkr)*dkr)); [~,ikz1] = min(abs(kz-round(1/dkz)*dkz)); -Sni_norm= zeros(1,Ns); % Time evolution of the amp of density nonlin term -Sne_norm= zeros(1,Ns); % Time evolution of the amp of vorti. nonlin term -E_pot = zeros(1,Ns); % Potential energy n^2 -E_kin = zeros(1,Ns); % Kinetic energy grad(phi)^2 -ExB = zeros(1,Ns); % ExB drift intensity \propto |\grad \phi| -CFL = zeros(1,Ns); % CFL time step +Ne_ST_kr = zeros(Nkr,Ns); % Space-Time of max_kz Ne(k) +Ne_ST_kz = zeros(Nkz,Ns); % '' max_kr Ne(k) +ne_ST_r = zeros(Nr,Ns); % Space-Time of ne(z==0) +ne_ST_z = zeros(Nz,Ns); % '' r==0 +ne_00 = zeros(1,Ns); % Time evolution of ne(r,z) at origin +Ne_gm = zeros(1,Ns); % Time evolution of Ne(k) max gamma (max real) +Ni_ST_kr = zeros(Nkr,Ns); % same for ions +Ni_ST_kz = zeros(Nkz,Ns); % . +ni_ST_r = zeros(Nr,Ns); % . +ni_ST_z = zeros(Nz,Ns); % . +ni_00 = zeros(1,Ns); % . +Ni_gm = zeros(1,Ns); % . +PH_ST_kr = zeros(Nkr,Ns); % same for ES-potential +PH_ST_kz = zeros(Nkz,Ns); % . +phi_ST_r = zeros(Nr,Ns); % . +phi_ST_z = zeros(Nz,Ns); % . +phi_00 = zeros(1,Ns); % . +Sne_norm = zeros(1,Ns); % Time evolution of the amp of e nonlin term +Sni_norm = zeros(1,Ns); % Time evolution of the amp of i nonlin term +E_pot = zeros(1,Ns); % Potential energy n^2 +E_kin = zeros(1,Ns); % Kinetic energy grad(phi)^2 +ExB = zeros(1,Ns); % ExB drift intensity \propto |\grad \phi| +CFL = zeros(1,Ns); % CFL time step Ddr = 1i*KR; Ddz = 1i*KZ; lapl = Ddr.^2 + Ddz.^2; +[~,ir0] = min(abs(r)); % index of r==0 +[~,iz0] = min(abs(z)); % index of z==0 +[~,ikr1] = min(abs(kr-round(1/dkr)*dkr)); % index of kr==1 +[~,ikz1] = min(abs(kz-round(1/dkz)*dkz)); % index of kz==1 for it = 1:numel(PH(1,1,:)) NE_ = Ne(:,:,it); NN_ = Ni(:,:,it); PH_ = PH(:,:,it); - phi_ST_r(:,it) = phi(:,iz0,it); phi_ST_z(:,it) = phi(ir0,:,it); + ne_ST_r (:,it)= ne(:,iz0,it); ne_ST_z (:,it)= ne(ir0,:,it); + Ne_ST_kr (:,it)= max(real(Ne(:,:,it)),[],2); + Ne_ST_kz (:,it)= max(real(Ne(:,:,it)),[],1); + ne_00(it) = ne(ir0,iz0,it); + Ne_gm(it) = max(max(real(Ne(:,:,it)))); + ni_ST_r (:,it)= ni(:,iz0,it); ni_ST_z (:,it)= ni(ir0,:,it); - phi_00(it) = phi(ir0,iz0,it); - ne_00(it) = ne(ir0,iz0,it); - ni_00(it) = ni(ir0,iz0,it); - Ne_11(it) = Ne(ikr1,ikz1,it); - Ni_11(it) = Ni(ikr1,ikz1,it); + Ni_ST_kr (:,it)= max(real(Ni(:,:,it)),[],2); + Ni_ST_kz (:,it)= max(real(Ni(:,:,it)),[],1); + ni_00(it) = ni(ir0,iz0,it); + Ni_gm(it) = max(max(real(Ni(:,:,it)))); + + phi_ST_r(:,it) = phi(:,iz0,it); phi_ST_z(:,it) = phi(ir0,:,it); + PH_ST_kr(:,it) = max(real(PH(:,:,it)),[],2); + PH_ST_kz(:,it) = max(real(PH(:,:,it)),[],1); + phi_00(it) = phi(ir0,iz0,it); + if strcmp(OUTPUTS.write_non_lin,'.true.') Sni_norm(it) = sum(sum(abs(SNi(:,:,it)))); Sne_norm(it) = sum(sum(abs(SNe(:,:,it)))); @@ -94,7 +112,7 @@ fig = figure; FIGNAME = ['t_evolutions',sprintf('_%.2d',JOBNUM)]; semilogy(Ts,abs(ni_00),'-','DisplayName','$n_i^{00}$'); grid on; xlabel('$t$'); ylabel('$|n_a(x=0,y=0)|$'); subplot(222) - semilogy(Ts,abs(Ni_11),'-','DisplayName','$\phi$') + semilogy(Ts,abs(Ni_gm),'-','DisplayName','$\phi$') grid on; xlabel('$t$'); ylabel('$|\tilde n(k_r\approx 1,k_z\approx 1)|$'); subplot(223) semilogy(Ts,E_kin+E_pot,'-','DisplayName','$\sum|ik\tilde\phi_i|^2+\sum|\tilde n_i|^2$') @@ -110,12 +128,12 @@ end FMT = '.fig'; save_figure %% Spectra energy -fig = figure; FIGNAME = ['Energy_kin_KZ',sprintf('_%.2d',JOBNUM)]; - semilogy(kr(floor(end/2)+1:end),E_kin_KR(floor(end/2)+1:end),'o','DisplayName','$\sum_y\langle|ik\tilde\phi_i|^2\rangle_t$') - hold on; - loglog(kz(floor(end/2)+1:end),E_kin_KZ(floor(end/2)+1:end),'o','DisplayName','$\sum_x\langle|ik\tilde\phi_i|^2\rangle_t$') - grid on; xlabel('$k$'); legend('show'); -FMT = '.fig'; save_figure +% fig = figure; FIGNAME = ['Energy_kin_KZ',sprintf('_%.2d',JOBNUM)]; +% semilogy(kr(floor(end/2)+1:end),E_kin_KR(floor(end/2)+1:end),'o','DisplayName','$\sum_y\langle|ik\tilde\phi_i|^2\rangle_t$') +% hold on; +% loglog(kz(floor(end/2)+1:end),E_kin_KZ(floor(end/2)+1:end),'o','DisplayName','$\sum_x\langle|ik\tilde\phi_i|^2\rangle_t$') +% grid on; xlabel('$k$'); legend('show'); +% FMT = '.fig'; save_figure %% CFL condition fig = figure; FIGNAME = ['CFL',sprintf('_%.2d',JOBNUM)]; @@ -155,17 +173,41 @@ subplot(223)% density xlabel('$z\,(r=0)$'); ylabel('$t$'); title('$\phi$'); FMT = '.fig'; save_figure -%% phi -fig = figure; FIGNAME = ['phi_ST',sprintf('_%.2d',JOBNUM)]; - [TY,TX] = meshgrid(Ts,z); - pclr = pcolor(TX,TY,(plt(phi_ST_r))); set(pclr, 'edgecolor','none'); colorbar; - xlabel('$x\,(y=0)$'); ylabel('$t$'); title('$\phi$'); +%% Space-Time diagrams for max_kz(Real) +plt = @(x) (x); +fig = figure; FIGNAME = ['kr_space_time_diag',sprintf('_%.2d',JOBNUM)]; + [TY,TX] = meshgrid(Ts,kr); +subplot(221)% density + pclr = pcolor(TX,TY,(plt(Ne_ST_kr))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kr$'); ylabel('$t$'); title('$\max_{k_z}(\textrm{Re}(N_e^{00}))$'); +subplot(222)% density + pclr = pcolor(TX,TY,(plt(Ni_ST_kr))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kr$'); ylabel('$t$'); title('$\max_{k_z}(\textrm{Re}(N_i^{00}))$'); +subplot(223)% density + pclr = pcolor(TX,TY,(plt(PH_ST_kr))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kr$'); ylabel('$t$'); title('$\max_{k_z}(\textrm{Re}(\tilde\phi$))'); FMT = '.fig'; save_figure +%% Space-Time diagrams at max_kr(Real) +plt = @(x) (x); +fig = figure; FIGNAME = ['kz_space_time_diag',sprintf('_%.2d',JOBNUM)]; + [TY,TX] = meshgrid(Ts,kz); +subplot(221)% density + pclr = pcolor(TX,TY,(plt(Ne_ST_kz))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kz$'); ylabel('$t$'); title('$\max_{k_r}(\textrm{Re}(N_e^{00}))$'); +subplot(222)% density + pclr = pcolor(TX,TY,(plt(Ni_ST_kz))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kz$'); ylabel('$t$'); title('$\max_{k_r}(\textrm{Re}(N_i^{00}))$'); +subplot(223)% density + pclr = pcolor(TX,TY,(plt(PH_ST_kz))); set(pclr, 'edgecolor','none'); colorbar; + xlabel('$kz$'); ylabel('$t$'); title('$\max_{k_r}(\textrm{Re}(\tilde\phi))$'); +FMT = '.fig'; save_figure + + if 0 %% Show frame -it = min(1,numel(Ts)); -fig = figure; FIGNAME = ['frame',sprintf('_%.2d',SID)]; +it = min(50,numel(Ts)); +fig = figure; FIGNAME = ['frame',sprintf('_%.2d',JOBNUM)]; subplot(221); plt = @(x) fftshift((real(x))); pclr = pcolor(fftshift(KR),fftshift(KZ),plt(PH(:,:,it))); set(pclr, 'edgecolor','none'); colorbar; xlabel('$k_r$'); ylabel('$k_z$'); title(sprintf('t=%.3d',Ts(it))); legend('$\hat\phi$'); @@ -183,19 +225,24 @@ end FMT = '.fig'; save_figure end %% -DELAY = 0.1; skip_ = 2; -if 0 +DELAY = 0.07; skip_ = 1; +FRAMES = 200:skip_:numel(Ts); +if 1 %% GIFS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Density electron GIFNAME = ['ne',sprintf('_%.2d',JOBNUM)]; FIELDNAME = '$n_e^{00}$'; -FIELD = real(ne(:,:,1:skip_:end)); X = XX; Y = YY; T = Ts(1:skip_:end); +FIELD = real(ne); X = XX; Y = YY; T = Ts; create_gif %% Density ion GIFNAME = ['ni',sprintf('_%.2d',JOBNUM)]; FIELDNAME = '$n_i^{00}$'; -FIELD = real(ni(:,:,1:skip_:end)); X = XX; Y = YY; T = Ts(1:skip_:end); +FIELD = real(ni); X = XX; Y = YY; T = Ts; create_gif %% Phi GIFNAME = ['phi',sprintf('_%.2d',JOBNUM)]; FIELDNAME = '$\phi$'; -FIELD = real(phi(:,:,1:skip_:end)); X = XX; Y = YY; T = Ts(1:skip_:end); +FIELD = real(phi); X = XX; Y = YY; T = Ts; +create_gif +%% Density electron frequency +GIFNAME = ['Ni',sprintf('_%.2d',JOBNUM)]; FIELDNAME = '$N_i^{00}$'; +FIELD = real(Ni); X = KR; Y = KZ; T = Ts; create_gif end \ No newline at end of file diff --git a/wk/fort.90 b/wk/fort.90 index c11d29567d34becbe94fc861c90d030b6a6efb84..7001202ea04c23d51eba399dad2456e2563d3a7d 100644 --- a/wk/fort.90 +++ b/wk/fort.90 @@ -1,18 +1,19 @@ &BASIC nrun = 100000000 - dt = 0.001 - tmax = 200 + dt = 0.0001 + tmax = 100 RESTART = .true. / &GRID - pmaxe =1 - jmaxe = 0 - pmaxi = 1 - jmaxi = 0 + pmaxe =5 + jmaxe = 2 + pmaxi = 5 + jmaxi = 2 Nr = 64 - Lr = 20 + Lr = 10 Nz = 64 - Lz = 20 + Lz = 10 + kpar = 0 / &OUTPUT_PAR nsave_0d = 0 @@ -22,18 +23,19 @@ write_Ni00 = .true. write_moments = .true. write_phi = .true. - write_non_lin = .true. + write_non_lin = .false. write_doubleprecision = .true. - resfile0 = 'SwoK_32x64_L_20_P_1_J_0_nB_0.5_nN_1_mu_1e-04_' - rstfile0 = '../checkpoint/cp_SwoK_32x64_L_20_P_1_J_0_nB_0.5_nN_1_mu_1e-04_' + resfile0 = 'gvskr_32x64_L_10_lin_P_5_J_2_nB_0.1_nN_1_mu_1e-02_' + rstfile0 = '../checkpoint/cp_gvskr_32x64_L_10_lin_P_5_J_2_nB_0.1_nN_1_mu_1e-02_' job2load = 0 / &MODEL_PAR ! Collisionality - CO = 0 - NON_LIN = .true. - mu = 0.0001 - nu = 0 + CO = -2 + DK = .false. + NON_LIN = .false. + mu = 0.01 + nu = 0.01 tau_e = 1 tau_i = 1 sigma_e = 0.023338 @@ -42,7 +44,7 @@ q_i = 1 eta_n = 1 eta_T = 0 - eta_B = 0.5 + eta_B = 0.1 lambdaD = 0 / &INITIAL_CON @@ -50,9 +52,9 @@ initback_moments =0.0001 initnoise_moments =5e-05 iseed =42 - selfmat_file ='../iCa/self_Coll_GKE_0_GKI_0_ESELF_1_ISELF_1_Pmaxe_1_Jmaxe_0_Pmaxi_1_Jmaxi_0_pamaxx_10.h5' - eimat_file ='../iCa/ei_Coll_GKE_0_GKI_0_ETEST_1_EBACK_1_Pmaxe_1_Jmaxe_0_Pmaxi_1_Jmaxi_0_pamaxx_10_tau_1.0000_mu_0.0233.h5' - iemat_file ='../iCa/ie_Coll_GKE_0_GKI_0_ITEST_1_IBACK_1_Pmaxe_1_Jmaxe_0_Pmaxi_1_Jmaxi_0_pamaxx_10_tau_1.0000_mu_0.0233.h5' + selfmat_file ='../iCa/self_Coll_GKE_0_GKI_0_ESELF_1_ISELF_1_Pmaxe_5_Jmaxe_2_Pmaxi_5_Jmaxi_2_pamaxx_10.h5' + eimat_file ='../iCa/ei_Coll_GKE_0_GKI_0_ETEST_1_EBACK_1_Pmaxe_5_Jmaxe_2_Pmaxi_5_Jmaxi_2_pamaxx_10_tau_1.0000_mu_0.0233.h5' + iemat_file ='../iCa/ie_Coll_GKE_0_GKI_0_ITEST_1_IBACK_1_Pmaxe_5_Jmaxe_2_Pmaxi_5_Jmaxi_2_pamaxx_10_tau_1.0000_mu_0.0233.h5' / &TIME_INTEGRATION_PAR numerical_scheme='RK4' diff --git a/wk/setup.m b/wk/setup.m index c1c45262f8e083a5cede87e93f9f2af168c30b41..022ec0005632467fde3fa10a2d7da9cd5e56565c 100644 --- a/wk/setup.m +++ b/wk/setup.m @@ -5,29 +5,32 @@ default_plots_options %% Set Up parameters %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% PHYSICAL PARAMETERS -NU = 0.0; % Collision frequency +NU = 1e-2; % Collision frequency TAU = 1.0; % e/i temperature ratio -ETAB = 0.5; % Magnetic gradient +ETAB = 0.1; % Magnetic gradient ETAN = 1.0; % Density gradient ETAT = 0.0; % Temperature gradient -MU = 1e-4; % Hyper diffusivity coefficient +MU = 1e-2; % Hyper diffusivity coefficient +LAMBDAD = 0.0; %% GRID PARAMETERS -N = 64; % Frequency gridpoints -L = 20; % Size of the squared frequency domain -PMAXE = 01; % Highest electron Hermite polynomial degree -JMAXE = 00; % Highest '' Laguerre '' -PMAXI = 01; % Highest ion Hermite polynomial degree -JMAXI = 00; % Highest '' Laguerre '' +N = 64; % Frequency gridpoints (Nr = N/2) +L = 10; % Size of the squared frequency domain +PMAXE = 05; % Highest electron Hermite polynomial degree +JMAXE = 02; % Highest '' Laguerre '' +PMAXI = 05; % Highest ion Hermite polynomial degree +JMAXI = 02; % Highest '' Laguerre '' +KPAR = 0.0; % Parellel wave vector component %% TIME PARAMETERS -TMAX = 200.0; % Maximal time unit -DT = 1e-3; % Time step -SPS = 1; % Sampling per time unit +TMAX = 100.0; % Maximal time unit +DT = 1e-4; % Time step +SPS = 10; % Sampling per time unit RESTART = 1; % To restart from last checkpoint JOB2LOAD= 0; %% OPTIONS -SIMID = 'SwoK'; % Name of the simulation -NON_LIN = 1; % activate non-linearity -CO = 0; +SIMID = 'gvskr'; % Name of the simulation +NON_LIN = 0; % activate non-linearity (is cancelled if KREQ0 = 1) +CO = -2; % Collision operator (0 : L.Bernstein, -1 : Full Coulomb, -2 : Dougherty) +DK = 0; % Drift kinetic model (put every kernel to 1) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -37,6 +40,7 @@ params = ['P_',num2str(PMAXE),'_J_',num2str(JMAXE),... '_nB_',num2str(ETAB),'_nN_',num2str(ETAN),'_mu_%0.0e_']; params = sprintf(params,MU); if ~NON_LIN; params = ['lin_',params]; end; +if DK; params = ['DK_',params]; end; resolution = ['_',num2str(N/2),'x',num2str(N),'_']; gridname = ['L_',num2str(L),'_']; BASIC.SIMID = [SIMID,resolution,gridname,params]; @@ -66,8 +70,10 @@ GRID.Nr = N; % r grid resolution GRID.Lr = L; % r length GRID.Nz = N; % z '' GRID.Lz = L; % z '' +GRID.kpar = KPAR; % Model parameters MODEL.CO = CO; % Collision operator (0 : L.Bernstein, -1 : Full Coulomb, -2 : Dougherty) +if DK; MODEL.DK = '.true.'; else; MODEL.DK = '.false.';end; if NON_LIN; MODEL.NON_LIN = '.true.'; else; MODEL.NON_LIN = '.false.';end; MODEL.mu = MU; MODEL.nu = NU; % hyper diffusive coefficient nu for HW @@ -84,7 +90,7 @@ MODEL.q_i = 1.0; MODEL.eta_n = ETAN; % source term kappa for HW MODEL.eta_T = ETAT; % Temperature MODEL.eta_B = ETAB; % Magnetic -MODEL.lambdaD = 0.0; +MODEL.lambdaD = LAMBDAD; % Time integration and intialization parameters TIME_INTEGRATION.numerical_scheme = '''RK4'''; INITIAL.only_Na00 = '.false.';