From b2925afbabfedc5b233f25c9b754f3a45659081b Mon Sep 17 00:00:00 2001
From: Antoine Hoffmann <antoine.hoffmann@epfl.ch>
Date: Tue, 25 Jul 2023 14:43:11 +0200
Subject: [PATCH] generic scripts to run and scan parameter files
---
wk/lin_run_script.m | 86 +++++++++++++++++++++++++++
wk/lin_scan_script.m | 134 +++++++++++++++++++++++++++++++++++++++++++
2 files changed, 220 insertions(+)
create mode 100644 wk/lin_run_script.m
create mode 100644 wk/lin_scan_script.m
diff --git a/wk/lin_run_script.m b/wk/lin_run_script.m
new file mode 100644
index 00000000..0eb275e6
--- /dev/null
+++ b/wk/lin_run_script.m
@@ -0,0 +1,86 @@
+%% 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);
+mpirun = 'mpirun';
+% mpirun = '/opt/homebrew/bin/mpirun'; % for macos
+addpath(genpath([gyacomodir,'matlab'])) % Add matlab folder
+addpath(genpath([gyacomodir,'matlab/plot'])) % Add plot folder
+addpath(genpath([gyacomodir,'matlab/compute'])) % Add compute folder
+addpath(genpath([gyacomodir,'matlab/load'])) % Add load folder
+addpath(genpath([gyacomodir,'wk/parameters'])) % Add parameters folder
+
+%% Setup run or load an executable
+RUN = 1; % To run or just to load
+default_plots_options
+EXECNAME = 'gyacomo23_sp'; % single precision
+% EXECNAME = 'gyacomo23_dp'; % double precision
+
+%% Setup parameters
+% run lin_DTT_AB_rho85_PT
+% run lin_Entropy
+% run lin_ITG
+
+%%-------------------------------------------------------------------------
+%% 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 8 ',gyacomodir,'bin/',EXECNAME,' 1 4 2 0;'];
+% RUN =['time ',mpirun,' -np 1 ',gyacomodir,'bin/',EXECNAME,' 1 1 1 0;'];
+ % RUN = ['./../../../bin/gyacomo23_sp 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 0 % 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
+
+if 1
+%% Ballooning plot
+[data.PHI, data.Ts3D] = compile_results_3D(LOCALDIR,J0,J1,'phi');
+if data.inputs.BETA > 0
+[data.PSI, data.Ts3D] = compile_results_3D(LOCALDIR,J0,J1,'psi');
+end
+options.time_2_plot = [120];
+options.kymodes = [0.25];
+options.normalized = 1;
+options.PLOT_KP = 0;
+% options.field = 'phi';
+fig = plot_ballooning(data,options);
+end
+
diff --git a/wk/lin_scan_script.m b/wk/lin_scan_script.m
new file mode 100644
index 00000000..a10d35fa
--- /dev/null
+++ b/wk/lin_scan_script.m
@@ -0,0 +1,134 @@
+%% 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);
+mpirun = 'mpirun';
+% mpirun = '/opt/homebrew/bin/mpirun'; % for macos
+addpath(genpath([gyacomodir,'matlab'])) % Add matlab folder
+addpath(genpath([gyacomodir,'matlab/plot'])) % Add plot folder
+addpath(genpath([gyacomodir,'matlab/compute'])) % Add compute folder
+addpath(genpath([gyacomodir,'matlab/load'])) % Add load folder
+addpath(genpath([gyacomodir,'wk/parameters'])) % Add parameters folder
+
+%% Setup run or load an executable
+RUN = 1; % To run or just to load
+RERUN = 0; % rerun if the data does not exist
+default_plots_options
+EXECNAME = 'gyacomo23_sp'; % single precision
+% EXECNAME = 'gyacomo23_dp'; % double precision
+
+%% Setup basic parameters
+run lin_DTT_AB_rho85_PT
+% run lin_Entropy
+% run lin_ITG
+
+%% Modify parameters
+% NZ = 1;
+NY = 2;
+DT = 2e-3;
+
+%% Scan parameters
+SIMID = [SIMID,'_scan'];
+P_a = [2 4 8];
+ky_a= 0.1:0.1:0.5;
+
+%% Scan loop
+% arrays for the result
+g_ky = zeros(numel(ky_a),numel(P_a),2);
+g_avg= g_ky*0;
+g_std= g_ky*0;
+j = 1;
+for PMAX = P_a
+ JMAX = P/2;
+ i = 1;
+ for ky = ky_a
+ LY = 2*pi/ky;
+ %% RUN
+ setup
+ % naming
+ filename = [SIMID,'/',PARAMS,'/'];
+ LOCALDIR = [gyacomodir,'results/',filename,'/'];
+ % check if data exist to run if no data
+ data_ = {};
+ try
+ data_ = compile_results_low_mem(data_,LOCALDIR,00,00);
+ Ntime = numel(data_.Ts0D);
+ catch
+ data_.outfilenames = [];
+ end
+ if RUN && (RERUN || isempty(data_.outfilenames) || Ntime < 10)
+ % system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 2 ',gyacomodir,'bin/',EXECNAME,' 1 2 1 0; cd ../../../wk'])
+ system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 4 ',gyacomodir,'bin/',EXECNAME,' 1 2 2 0; cd ../../../wk'])
+ % system(['cd ../results/',SIMID,'/',PARAMS,'/; mpirun -np 6 ',gyacomodir,'bin/',EXECNAME,' 3 2 1 0; cd ../../../wk'])
+ end
+ data_ = compile_results_low_mem(data_,LOCALDIR,00,00);
+ [data_.PHI, data_.Ts3D] = compile_results_3D(LOCALDIR,00,00,'phi');
+ if numel(data_.Ts3D)>10
+ if numel(data_.Ts3D)>5
+ % Load results after trying to run
+ filename = [SIMID,'/',PARAMS,'/'];
+ LOCALDIR = [gyacomodir,'results/',filename,'/'];
+
+ data_ = compile_results_low_mem(data_,LOCALDIR,00,00);
+ [data_.PHI, data_.Ts3D] = compile_results_3D(LOCALDIR,00,00,'phi');
+
+ % linear growth rate (adapted for 2D zpinch and fluxtube)
+ options.TRANGE = [0.5 1]*data_.Ts3D(end);
+ options.NPLOTS = 0; % 1 for only growth rate and error, 2 for omega local evolution, 3 for plot according to z
+ options.GOK = 0; %plot 0: gamma 1: gamma/k 2: gamma^2/k^3
+
+ [~,it1] = min(abs(data_.Ts3D-0.5*data_.Ts3D(end))); % start of the measurement time window
+ [~,it2] = min(abs(data_.Ts3D-1.0*data_.Ts3D(end))); % end of ...
+ field = 0;
+ field_t = 0;
+ for ik = 2:NY/2+1
+ field = squeeze(sum(abs(data_.PHI),3)); % take the sum over z
+ field_t = squeeze(field(ik,1,:)); % take the kx =0, ky = ky mode only
+ to_measure = log(field_t(it1:it2));
+ tw = double(data_.Ts3D(it1:it2));
+ % gr = polyfit(tw,to_measure,1);
+ gr = fit(tw,to_measure,'poly1');
+ err= confint(gr);
+ g_ky(i,j,ik) = gr.p1;
+ g_std(i,j,ik) = abs(err(2,1)-err(1,1))/2;
+ end
+ [gmax, ikmax] = max(g_ky(i,j,:));
+
+ msg = sprintf('gmax = %2.2f, kmax = %2.2f',gmax,data_.grids.ky(ikmax)); disp(msg);
+ end
+ end
+ i = i + 1;
+end
+j = j + 1;
+end
+
+%% take max growth rate among z coordinate
+y_ = g_ky(:,:,2);
+e_ = g_std(:,:,2);
+
+%% Save scan results (gamma)
+if(numel(ky_a)>1 && numel(P_a)>1)
+ pmin = num2str(min(P_a)); pmax = num2str(max(P_a));
+ kymin = num2str(min(ky_a)); kymax= num2str(max(ky_a));
+ filename = [num2str(NX),'x',num2str(NZ),'_ky_',kymin,'_',kymax,...
+ '_P_',pmin,'_',pmax,'_',CONAME,'_',num2str(NU),'_be_',num2str(BETA),'.mat'];
+ metadata.name = filename;
+ metadata.kymin = ky;
+ metadata.title = ['$\nu_{',CONAME,'}=$',num2str(NU),'$\kappa_T=$',num2str(K_Ti),', $\kappa_N=$',num2str(K_Ni)];
+ metadata.par = [num2str(NX),'x1x',num2str(NZ)];
+ metadata.nscan = 2;
+ metadata.s2name = '$P$';
+ metadata.s2 = P_a;
+ metadata.s1name = '$ky$';
+ metadata.s1 = ky_a;
+ metadata.dname = '$\gamma c_s/R$';
+ metadata.data = y_;
+ metadata.err = e_;
+ save([SIMDIR,filename],'-struct','metadata');
+ disp(['saved in ',SIMDIR,filename]);
+ clear metadata tosave
+end
--
GitLab