Skip to content
Snippets Groups Projects
Commit 9eb97515 authored by Antoine Cyril David Hoffmann's avatar Antoine Cyril David Hoffmann :seedling:
Browse files

updating and cleaning testcases

parent d4b2bb4e
No related branches found
No related tags found
No related merge requests found
......@@ -17,17 +17,17 @@
Nexc = 1
/
&GEOMETRY
geom = 'Z-pinch'
q0 = 0.0
shear = 0.0
eps = 0.0
geom = 'z-pinch'
q0 = 1.4
shear = 0.8
eps = 0.13
kappa = 1.0
s_kappa= 0.0
delta = 0.0
s_delta= 0.0
zeta = 0.0
s_zeta = 0.0
parallel_bc = 'shearless'
parallel_bc = 'dirichlet'
shift_y= 0.0
/
&OUTPUT_PAR
......@@ -59,11 +59,11 @@
nu = 0.1
beta = 0.0
ADIAB_E = .f.
tau_e = 1.0
tau_i = 1.0
/
&CLOSURE_PAR
!hierarchy_closure='truncation'
hierarchy_closure='max_degree'
hierarchy_closure='truncation'
!hierarchy_closure='max_degree'
dmax = 2
nonlinear_closure='anti_laguerre_aliasing' !(truncation,full_sum,anti_laguerre_aliasing)
nmax = -1
......
&BASIC
nrun = 99999999
dt = 0.01
tmax = 500
maxruntime = 72000
job2load = 0
/
&GRID
pmax = 2
jmax = 1
Nx = 64
Lx = 200
Ny = 48
Ly = 120
Nz = 1
SG = .f.
Nexc = 1
/
&GEOMETRY
geom = 'Z-pinch'
q0 = 0.0
shear = 0.0
eps = 0.0
kappa = 1.0
s_kappa= 0.0
delta = 0.0
s_delta= 0.0
zeta = 0.0
s_zeta = 0.0
parallel_bc = 'shearless'
shift_y= 0.0
/
&OUTPUT_PAR
dtsave_0d = 1
dtsave_1d = -1
dtsave_2d = -1
dtsave_3d = 1
dtsave_5d = 10
write_doubleprecision = .f.
write_gamma = .t.
write_hf = .t.
write_phi = .t.
write_Na00 = .t.
write_Napj = .t.
write_dens = .t.
write_fvel = .t.
write_temp = .t.
/
&MODEL_PAR
LINEARITY = 'nonlinear'
Na = 2 ! number of species
mu_x = 1.0
mu_y = 1.0
N_HD = 4
mu_z = 0.0
HYP_V = 'hypcoll'
mu_p = 0.0
mu_j = 0.0
nu = 0.1
beta = 0.0
ADIAB_E = .f.
tau_e = 1.0
/
&CLOSURE_PAR
!hierarchy_closure='truncation'
hierarchy_closure='max_degree'
dmax = 2
nonlinear_closure='anti_laguerre_aliasing' !(truncation,full_sum,anti_laguerre_aliasing)
nmax = -1
/
&SPECIES
! ions
name_ = 'ions'
tau_ = 1.0
sigma_= 1.0
q_ = 1.0
k_N_ = 1.6
k_T_ = 0.4
/
&SPECIES
! electrons
name_ = 'electrons'
tau_ = 1.0
sigma_= 0.023338
q_ =-1.0
k_N_ = 1.6
k_T_ = 0.4
/
&COLLISION_PAR
collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
GK_CO = .t.
INTERSPECIES = .true.
mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
/
&INITIAL_CON
INIT_OPT = 'blob' !(phi,blob)
ACT_ON_MODES = 'donothing'
init_background = 0.0
init_noiselvl = 0.005
iseed = 42
/
&TIME_INTEGRATION_PAR
numerical_scheme = 'RK4'
/
&BASIC
nrun = 99999999
dt = 0.01
tmax = 700
maxruntime = 72000
job2load = 1
/
&GRID
pmax = 2
jmax = 1
Nx = 64
Lx = 200
Ny = 48
Ly = 120
Nz = 1
SG = .f.
Nexc = 1
/
&GEOMETRY
geom = 'Z-pinch'
q0 = 0.0
shear = 0.0
eps = 0.0
kappa = 1.0
s_kappa= 0.0
delta = 0.0
s_delta= 0.0
zeta = 0.0
s_zeta = 0.0
parallel_bc = 'shearless'
shift_y= 0.0
/
&OUTPUT_PAR
dtsave_0d = 1
dtsave_1d = -1
dtsave_2d = -1
dtsave_3d = 1
dtsave_5d = 10
write_doubleprecision = .f.
write_gamma = .t.
write_hf = .t.
write_phi = .t.
write_Na00 = .t.
write_Napj = .t.
write_dens = .t.
write_fvel = .t.
write_temp = .t.
/
&MODEL_PAR
LINEARITY = 'nonlinear'
Na = 2 ! number of species
mu_x = 1.0
mu_y = 1.0
N_HD = 4
mu_z = 0.0
HYP_V = 'hypcoll'
mu_p = 0.0
mu_j = 0.0
nu = 0.1
beta = 0.0
ADIAB_E = .f.
tau_e = 1.0
/
&CLOSURE_PAR
hierarchy_closure='truncation'
!hierarchy_closure='max_degree'
dmax = 2
nonlinear_closure='anti_laguerre_aliasing' !(truncation,full_sum,anti_laguerre_aliasing)
nmax = -1
/
&SPECIES
! ions
name_ = 'ions'
tau_ = 1.0
sigma_= 1.0
q_ = 1.0
k_N_ = 1.6
k_T_ = 0.4
/
&SPECIES
! electrons
name_ = 'electrons'
tau_ = 1.0
sigma_= 0.023338
q_ =-1.0
k_N_ = 1.6
k_T_ = 0.4
/
&COLLISION_PAR
collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
GK_CO = .t.
INTERSPECIES = .true.
mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
/
&INITIAL_CON
INIT_OPT = 'blob' !(phi,blob)
ACT_ON_MODES = 'donothing'
init_background = 0.0
init_noiselvl = 0.005
iseed = 42
/
&TIME_INTEGRATION_PAR
numerical_scheme = 'RK4'
/
../../bin/gyacomo22
\ No newline at end of file
../../bin/gyacomo23_sp
\ No newline at end of file
%% QUICK RUN SCRIPT
% This script creates a directory in /results and runs a simulation directly
% from the Matlab framework. It is meant to run only small problems in linear
% for benchmarking and debugging purposes since it makes Matlab "busy".
%% Set up the paths for the necessary Matlab modules
gyacomodir = pwd;
gyacomodir = gyacomodir(1:end-2);
addpath(genpath([gyacomodir,'matlab'])) % Add matlab module
addpath(genpath([gyacomodir,'matlab/plot'])) % Add plot module
addpath(genpath([gyacomodir,'matlab/compute'])) % Add compute module
addpath(genpath([gyacomodir,'matlab/load'])) % Add load module
%% Set simulation parameters
SIMID = 'ETG_adiab_i'; % Name of the simulation
RUN = 1; % To run or just to load
default_plots_options
% EXECNAME = 'gyacomo23_sp'; % single precision
% EXECNAME = 'gyacomo23_dp'; % double precision
EXECNAME = 'gyacomo23_debug'; % single precision
%% Set up physical parameters
CLUSTER.TIME = '99:00:00'; % Allocation time hh:mm:ss
NU = 0.0; % Collision frequency
TAU = 1; % e/i temperature ratio
K_Ne = 0; % ele Density '''
K_Te = 10; % ele Temperature '''
K_Ni = 0; % ion Density gradient drive
K_Ti = 0; % ion Temperature '''
SIGMA_E = 0.0233380; % mass ratio sqrt(m_a/m_i) (correct = 0.0233380)
NA = 1; % number of kinetic species
ADIAB_E = 0; % adiabatic electron model
ADIAB_I = 1; % adiabatic ion model
BETA = 1e-5; % electron plasma beta
MHD_PD = 0; % MHD pressure drift
%% Set up grid parameters
P = 2;
J = P/2;%P/2;
PMAX = P; % Hermite basis size
JMAX = J; % Laguerre basis size
NX = 2; % real space x-gridpoints
NY = 16; % real space y-gridpoints
LX = 2*pi/0.3; % Size of the squared frequency domain in x direction
LY = 2*pi/4.0; % Size of the squared frequency domain in y direction
NZ = 16; % number of perpendicular planes (parallel grid)
SG = 0; % Staggered z grids option
NEXC = 1; % To extend Lx if needed (Lx = Nexc/(kymin*shear))
%% GEOMETRY
% GEOMETRY= 's-alpha';
% GEOMETRY= 'miller';
GEOMETRY= 'z-pinch';
EPS = 1.0; % inverse aspect ratio
Q0 = 1.0; % safety factor
SHEAR = 0.0; % magnetic shear
KAPPA = 1.0; % elongation
S_KAPPA = 0.0;
DELTA = 0.0; % triangularity
S_DELTA = 0.0;
ZETA = 0.0; % squareness
S_ZETA = 0.0;
PARALLEL_BC = 'dirichlet'; % Boundary condition for parallel direction ('dirichlet','periodic','shearless','disconnected')
SHIFT_Y = 0.0; % Shift in the periodic BC in z
NPOL = 1; % Number of poloidal turns
PB_PHASE= 0;
%% TIME PARAMETERS
TMAX = 10; % Maximal time unit
DT = 1e-3; % Time step
DTSAVE0D = 1; % Sampling per time unit for 0D arrays
DTSAVE2D = -1; % Sampling per time unit for 2D arrays
DTSAVE3D = 3; % Sampling per time unit for 3D arrays
DTSAVE5D = 100; % Sampling per time unit for 5D arrays
JOB2LOAD = -1; % Start a new simulation serie
%% OPTIONS
LINEARITY = 'linear'; % activate non-linearity (is cancelled if KXEQ0 = 1)
CO = 'DG'; % Collision operator (LB:L.Bernstein, DG:Dougherty, SG:Sugama, LR: Lorentz, LD: Landau)
GKCO = 1; % Gyrokinetic operator
ABCO = 1; % INTERSPECIES collisions
INIT_ZF = 0; % Initialize zero-field quantities
HRCY_CLOS = 'truncation'; % Closure model for higher order moments
DMAX = -1;
NLIN_CLOS = 'truncation'; % Nonlinear closure model for higher order moments
NMAX = 0;
KERN = 0; % Kernel model (0 : GK)
INIT_OPT = 'phi'; % Start simulation with a noisy mom00/phi/allmom
NUMERICAL_SCHEME = 'RK4'; % Numerical integration scheme (RK2,SSPx_RK2,RK3,SSP_RK3,SSPx_RK3,IMEX_SSP2,ARK2,RK4,DOPRI5)
%% OUTPUTS
W_DOUBLE = 1; % Output flag for double moments
W_GAMMA = 1; % Output flag for gamma (Gyrokinetic Energy)
W_HF = 1; % Output flag for high-frequency potential energy
W_PHI = 1; % Output flag for potential
W_NA00 = 1; % Output flag for nalpha00 (density of species alpha)
W_DENS = 1; % Output flag for total density
W_TEMP = 1; % Output flag for temperature
W_NAPJ = 1; % Output flag for nalphaparallel (parallel momentum of species alpha)
W_SAPJ = 0; % Output flag for saparallel (parallel current of species alpha)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% UNUSED PARAMETERS
% These parameters are usually not to play with in linear runs
MU = 0.1; % Hyperdiffusivity coefficient
MU_X = MU; % Hyperdiffusivity coefficient in x direction
MU_Y = MU; % Hyperdiffusivity coefficient in y direction
N_HD = 4; % Degree of spatial-hyperdiffusivity
MU_Z = 10.0; % Hyperdiffusivity coefficient in z direction
HYP_V = 'hypcoll'; % Kinetic-hyperdiffusivity model
MU_P = 0.0; % Hyperdiffusivity coefficient for Hermite
MU_J = 0.0; % Hyperdiffusivity coefficient for Laguerre
LAMBDAD = 0.0; % Lambda Debye
NOISE0 = 0.0e-5; % Initial noise amplitude
BCKGD0 = 1.0e-5; % Initial background
k_gB = 1.0; % Magnetic gradient strength
k_cB = 1.0; % Magnetic curvature strength
COLL_KCUT = 1; % Cutoff for collision operator
%%-------------------------------------------------------------------------
%% RUN
setup
% system(['rm fort*.90']);
% Run linear simulation
if RUN
MVIN =['cd ../results/',SIMID,'/',PARAMS,'/;'];
% RUN =['time mpirun -np 2 ',gyacomodir,'bin/',EXECNAME,' 1 2 1 0;'];
% RUN =['time mpirun -np 4 ',gyacomodir,'bin/',EXECNAME,' 1 2 2 0;'];
% RUN =['time mpirun -np 6 ',gyacomodir,'bin/',EXECNAME,' 1 6 1 0;'];
RUN =['time mpirun -np 1 ',gyacomodir,'bin/',EXECNAME,' 1 1 1 0;'];
MVOUT='cd ../../../wk;';
system([MVIN,RUN,MVOUT]);
end
%% Analysis
% load
filename = [SIMID,'/',PARAMS,'/']; % Create the filename based on SIMID and PARAMS
LOCALDIR = [gyacomodir,'results/',filename,'/']; % Create the local directory path based on gyacomodir, results directory, and filename
FIGDIR = LOCALDIR; % Set FIGDIR to the same path as LOCALDIR
% Load outputs from jobnummin up to jobnummax
J0 = 0; J1 = 0;
data = {}; % Initialize data as an empty cell array
% load grids, inputs, and time traces
data = compile_results_low_mem(data,LOCALDIR,J0,J1);
if 1 % Activate or not
%% plot mode evolution and growth rates
% Load phi
[data.PHI, data.Ts3D] = compile_results_3D(LOCALDIR,J0,J1,'phi');
options.NORMALIZED = 0;
options.TIME = data.Ts3D;
% Time window to measure the growth of kx/ky modes
options.KX_TW = [0.2 1]*data.Ts3D(end);
options.KY_TW = [0.2 1]*data.Ts3D(end);
options.NMA = 1; % Set NMA option to 1
options.NMODES = 999; % Set how much modes we study
options.iz = 'avg'; % Compressing z
options.ik = 1; %
options.fftz.flag = 0; % Set fftz.flag option to 0
fig = mode_growth_meter(data,options); % Call the function mode_growth_meter with data and options as input arguments, and store the result in fig
end
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment