rundke_Zi_remote

 This function is used to test the use of the MDCE toolbox within run_lukert
 A loop is created over various collision operators and ion charges

   - If it is called with remtimout > 0, the function waits for each calculation to be completed 
   - If it is called with remtimout = 0, the function saves the job in a LUKE_RESULTS file
   - If it is called without arguments, it loads the previously saved LUKE_RESULTS file if it exists, and check the results.

 INPUTS :
       - remtimout : timeout for remote calculations, in minutes. Use (0) to return later and check results by running this script again
       - link_mode : connection mode prescription
           o -2  : disable remote profiles
           o -1  : local profiles only
           o NaN : (default) use value from remote.link_mode
           o []  : give priority to (0) ssh key pair, (1) ssh bridge, (2) url/ftp bridge
           o 0   : ssh key pair enforced
           o 1   : ssh bridge enforced
           o 2   : url/ftp bridge enforced
       - opt_save
           o 0   : no not save LUKE_RESULTS files (default)
           o 1   : save LUKE_RESULTS files in folder LUKE_RESULTS_files/
       - display_mode
           o 0   : no display
           o 1   : display general info and time evoution (default)
           o 2   : display all MDCE steps
       - plot_results : flag to plot results (default : false)

0001 function rundke_Zi_remote(remtimout,link_mode,opt_save,display_mode,plot_results)
0002 %
0003 % This function is used to test the use of the MDCE toolbox within run_lukert
0004 % A loop is created over various collision operators and ion charges
0005 %
0006 %   - If it is called with remtimout > 0, the function waits for each calculation to be completed
0007 %   - If it is called with remtimout = 0, the function saves the job in a LUKE_RESULTS file
0008 %   - If it is called without arguments, it loads the previously saved LUKE_RESULTS file if it exists, and check the results.
0009 %
0010 % INPUTS :
0011 %       - remtimout : timeout for remote calculations, in minutes. Use (0) to return later and check results by running this script again
0012 %       - link_mode : connection mode prescription
0013 %           o -2  : disable remote profiles
0014 %           o -1  : local profiles only
0015 %           o NaN : (default) use value from remote.link_mode
0016 %           o []  : give priority to (0) ssh key pair, (1) ssh bridge, (2) url/ftp bridge
0017 %           o 0   : ssh key pair enforced
0018 %           o 1   : ssh bridge enforced
0019 %           o 2   : url/ftp bridge enforced
0020 %       - opt_save
0021 %           o 0   : no not save LUKE_RESULTS files (default)
0022 %           o 1   : save LUKE_RESULTS files in folder LUKE_RESULTS_files/
0023 %       - display_mode
0024 %           o 0   : no display
0025 %           o 1   : display general info and time evoution (default)
0026 %           o 2   : display all MDCE steps
0027 %       - plot_results : flag to plot results (default : false)
0028 %
0029 if nargin < 5,
0030     plot_results = false;% flag to plot results
0031 end
0032 %
0033 if nargin < 4,
0034     display_mode = 1;
0035 end
0036 %
0037 if nargin < 3,
0038     opt_save = 0;
0039 end
0040 %
0041 if nargin < 2 || (~isempty(link_mode) && (~isnumeric(link_mode) || ~isscalar(link_mode) || ~any(link_mode == [-1:2]))),
0042     link_mode = NaN;% bridge mode : (NaN) follow remote, ([]) priority to ssh, (0) ssh only, (1) bridge only
0043 end
0044 %
0045 if nargin < 1,
0046     remtimout = 0;% timeout for remote calculations, in minutes. Use (0) to return later and check results by running this script again
0047 end
0048 %
0049 permission = test_permissions_yp;
0050 %
0051 if ~permission 
0052     disp('Please move the script to a local folder where you have write permission before to run it')
0053     return;
0054 end
0055 %
0056 rempause = 1/60;%pause between job checks (in minutes) (only for remtimout > 0)
0057 %
0058 debugmode = 0;%1;% debug mode (0) off (1) on
0059 cleanmode = 2;%0;% clean job files/folders : both on remote and local (2), only on remote (1) or not at all (0)
0060 luke_memory = 2000;% job memory in mb
0061 luke_walltime = 1;% job walltime in hours
0062 enforce = 0;% (0) : exit if remote calculation failed, (1): if remote calculation failed, it is done locally
0063 %
0064 opt.save = opt_save;% save individual results in LUKE_RESULTS files (WARNING : Use absolute simul.path to save locally with local cluster)
0065 %
0066 mdce_opt = 0;% option for testing distributed computing inside remote calculation
0067 %
0068 p_opt = -1;
0069 %
0070 localdatadir = 'LUKE_RESULTS_files_Zi_remote/';
0071 %
0072 % ***********************This part must be specified by the user, run make files if necessary) *****************************
0073 %
0074 id_simul = 'LH_karney_Zi_remote';%Simulation ID
0075 path_simul = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0076 %
0077 psin_S = [];%Normalized poloidal flux grid where calculations are performed (0 < psin_S < 1) (If one value: local calculation only, not used if empty)
0078 rho_S = [0.5];%Normalized radial flux grid where calculations are performed (0 < rho_S < 1) (If one value: local calculation only, not used if empty)
0079 %
0080 id_path = '';%For all paths used by DKE solver
0081 path_path = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0082 %
0083 id_equil = 'TScyl';%For plasma equilibrium
0084 path_equil = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0085 %
0086 id_dkeparam = 'UNIFORM10010020';%For DKE code parameters
0087 path_dkeparam = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0088 %
0089 id_display = 'NO_DISPLAY';%For output code display
0090 path_display = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0091 %
0092 %************************************************************************************************************************************
0093 %************************************************************************************************************************************
0094 %************************************************************************************************************************************
0095 %
0096 [dkepath,equil,dkeparam,dkedisplay] = load_structures_yp('dkepath',id_path,path_path,'equil',id_equil,path_equil,'dkeparam',id_dkeparam,path_dkeparam,'dkedisplay',id_display,path_display);
0097 %
0098 %************************************************************************************************************************************
0099 %
0100 wavestruct.omega_lh = [4]*2*pi*1e9; %(GHz -> rad/s). Wave frequency [1,1] Indicative, no effect in small FLR limit opt_lh > 0
0101 %Option parameter for cross-comparison between old LH code:
0102 %    - (1): 1/vpar dependence
0103 %    - (2): no 1/vpar dependence and old grid technique for Dlh calculations (Karney, Shoucri) (see rfdiff_dke_jd)
0104 wavestruct.opt_lh = 2; % [1,1]
0105 %
0106 % Choose (vparmin_lh,vparmax_lh) or (Nparmin_lh,Nparmax_lh) for square n// LH wave power spectrum,
0107 % or (Npar_lh,dNpar_lh) for Gaussian shape
0108 %
0109 wavestruct.norm_ref = 1;%Normalization procedure for the LH quasilinear diffusion coefficient and spectrum boundaries
0110 %
0111 wavestruct.yvparmin_lh = [3];%LH wave square N// Spectrum: Lower limit of the plateau (vth_ref or vth) [1,n_scenario_lh]
0112 wavestruct.yvparmax_lh = [5];%LH wave square N// Spectrum: Upper limit of the plateau (vth_ref or vth) [1,n_scenario_lh]
0113 %
0114 wavestruct.yNparmin_lh = [NaN];%LH wave square N// Spectrum: Lower limit [1,n_scenario_lh]
0115 wavestruct.yNparmax_lh = [NaN];%LH wave square N// Spectrum: Upper limit [1,n_scenario_lh]
0116 wavestruct.yNpar_lh = [NaN];%LH wave Gaussian N// Spectrum: peak [1,n_scenario_lh]
0117 wavestruct.ydNpar_lh = [NaN];%LH wave Gaussian N// Spectrum: width [1,n_scenario_lh]
0118 %
0119 %   Note: this diffusion coefficient is different from the general QL D0. It has a benchmarking purpose only
0120 wavestruct.yD0_in_c_lh = [1];%Central LH QL diffusion coefficient (nhuth_ref*pth_ref^2 or nhuth*pth^2) [1,n_scenario_lh]
0121 %
0122 wavestruct.yD0_in_lh_prof = [0];%Quasilinear diffusion coefficient radial profile: (0) uniform, (1) gaussian radial profile [1,n_scenario_lh]
0123 wavestruct.ypeak_lh = [NaN];%Radial peak position of the LH quasi-linear diffusion coefficient (r/a on midplane) [1,n_scenario_lh]
0124 wavestruct.ywidth_lh = [NaN];%Radial width of the LH quasi-linear diffusion coefficient (r/a on midplane) [1,n_scenario_lh]
0125 %
0126 wavestruct.ythetab_lh = [0]*pi/180;%(deg -> rad). Poloidal location of LH beam [0..2pi] [1,n_scenario_lh]
0127 %               (0) from local values Te and ne, (1) from central values Te0 and ne0
0128 %
0129 %************************************************************************************************************************************
0130 %
0131 if exist('dmumpsmex');dkeparam.invproc = -2;end
0132 %
0133 dkeparam.boundary_mode_f = 0;%Number of points where the Maxwellian distribution is enforced from p = 0 (p=0, free conservative mode but param_inv(1) must be less than 1e-4, otherwise 1e-3 is OK most of the time. Sensitive to the number of points in p)
0134 dkeparam.norm_mode_f = 1;%Local normalization of f0 at each iteration: (0) no, the default value when the numerical conservative scheme is correct, (1) yes
0135 dkeparam.tn = [50000,100000];% 2 time steps to tn=100000 best for asymptotic solution
0136 dkeparam.rt_mode = 1;% necessary to use run_lukert
0137 %
0138 dkeparam.np_S = 201;
0139 dkeparam.nmhu_S = 201;
0140 %
0141 dkeparam.rho_S = rho_S;% radial flux grid points enforced
0142 dkeparam.psin_S = psin_S;% radial flux grid points enforced
0143 %
0144 dkeparam.clustermode.coll_dke_jd.scheduler.mode = mdce_opt;
0145 dkeparam.clustermode.eecoll_dke_yp.scheduler.mode = mdce_opt;
0146 %
0147 betath = 0.001;
0148 %
0149 [~,~,~,~,~,~,~,mc2] = pc_dke_yp;%Universal physics constants
0150 %
0151 equil.pTe = betath^2*mc2*ones(size(equil.pTe));
0152 equil.pzTi = betath^2*mc2*ones(size(equil.pzTi));
0153 %equil.pzTi = 1e-10*ones(size(equil.pzTi));
0154 %
0155 opt.fields = {'Zcurr','ZP0'};% select fields to be returned from lukestructs. Use NaN to return all fields.
0156 %
0157 Zi_list = [1,100];%Z = 74 is tungsten fully stripped ion
0158 coll_mode_list = 0:2;
0159 %
0160 nZi = length(Zi_list);
0161 ncoll = length(coll_mode_list);
0162 %
0163 % simulation filename
0164 %
0165 filename = [path_simul,'LUKE_RESULTS_',id_equil,'_',id_simul,'.mat'];
0166 %
0167 % check if current jobs exists to load possible existing jobs
0168 %
0169 if nargin == 0 && exist(filename,'file'),% load existing simulation
0170     %
0171     load(filename,'lukestructs','remtimout','link_mode','plot_results','remnum','dist_mode','mdce_mode');
0172     %
0173 else% new simulation
0174     %
0175     % determite remote and distributed profiles
0176     %
0177    % if ~exist('dkepath','var'),
0178    %     dkepath = setpath_luke_jd('',link_mode,0);
0179    % else
0180    %     dkepath = setpath_luke_jd(dkepath,link_mode,0);
0181    % end
0182     %
0183     if ~isempty(dkepath.remote),
0184         remopt = input_dke_yp(['Do you want to transfer the distributed calculation to a remote host (y/n)'],'n',{'y','n'});
0185         disp(' ')
0186     else
0187         remopt = 'n';
0188     end
0189     %
0190     if remopt == 'y',% use of mdce_remote
0191         remnum = iselect_jd(dkepath.profilestr,'Please select the host for the distributed calculation ?',0,struct,1);
0192         disp(' ')
0193         dist_mode = NaN;
0194         mdce_mode = input_dke_yp(['Please select a mode for the distribution on host ',dkepath.profilestr{remnum},':\n(0) sequential \n(1) jobcomputing \n(2) parcomputing \n(3) remcomputing \n'],0,0:3);
0195     else% direct distribution from localhost
0196         remnum = 0;
0197         dist_mode = iselect_jd(dkepath.distprofiles,'Please select a mode for the distribution from the local host ?',0,struct,1);
0198         if dist_mode > 3,
0199             mdce_mode = 3 - dist_mode;
0200         else
0201             mdce_mode = dist_mode - 1;
0202         end
0203     end
0204     %
0205     dkepath.clustermode.run_lukert.scheduler.remnum = remnum;
0206     dkepath.clustermode.run_lukert.scheduler.mode = mdce_mode;
0207     dkepath.clustermode.run_lukert.scheduler.remtimout = remtimout;
0208     dkepath.clustermode.run_lukert.scheduler.rempause = rempause;
0209     dkepath.clustermode.run_lukert.scheduler.memory = luke_memory;
0210     dkepath.clustermode.run_lukert.scheduler.walltime = luke_walltime;
0211     dkepath.clustermode.run_lukert.scheduler.debugmode = debugmode;
0212     dkepath.clustermode.run_lukert.scheduler.display = display_mode;
0213     dkepath.clustermode.run_lukert.scheduler.cleanmode = cleanmode;
0214     dkepath.clustermode.run_lukert.scheduler.enforce = enforce;% exit if remote calculation failed
0215     if opt.save,
0216         dkepath.clustermode.run_lukert.scheduler.localdatadir = localdatadir;
0217         if exist(localdatadir,'dir') && length(dir(localdatadir)) > 2,
0218             delete([localdatadir,'*']);
0219         end
0220     end
0221     %
0222     % generate luke_inputs structure
0223     %
0224     lukestructs = cell(ncoll,nZi);
0225     %
0226     for iZi = 1:nZi,
0227         for icoll = 1:ncoll,
0228             %
0229             Zi = Zi_list(iZi);
0230             equil.zZi = [1,1,1,Zi];
0231             equil.pzni = [0;0;0;1/Zi]*equil.pne;
0232             %
0233             waves{1} = make_idealLHwave_jd(equil,wavestruct);
0234             %
0235             dkeparam.coll_mode = coll_mode_list(icoll);% Relativistic Maxwellian background
0236             %
0237             % generate initial luke_inputs structure
0238             %
0239             lukestructs{icoll,iZi}.simul.path = path_simul;
0240             lukestructs{icoll,iZi}.simul.locid = [id_simul,'_coll_',num2str(coll_mode_list(icoll)),'_Zi_',num2str(Zi)];
0241             lukestructs{icoll,iZi}.equil = equil;
0242             lukestructs{icoll,iZi}.dkeparam = dkeparam;
0243             lukestructs{icoll,iZi}.dkedisplay = dkedisplay;
0244             lukestructs{icoll,iZi}.waves = waves;
0245             lukestructs{icoll,iZi}.opt = opt;
0246             %
0247         end
0248     end
0249 end
0250 %
0251 lukestructs = run_lukert(lukestructs,dkepath);
0252 %
0253 j = NaN(ncoll,nZi);
0254 P = NaN(ncoll,nZi);
0255 %
0256 flag_running = true(ncoll,nZi);
0257 flag_error = true(ncoll,nZi);;
0258 luke_warning = cell(ncoll,nZi);
0259 %
0260 for iZi = 1:nZi,
0261     for icoll = 1:ncoll,
0262         if isfield(lukestructs{icoll,iZi},'output'),
0263             %
0264             j(icoll,iZi) = lukestructs{icoll,iZi}.output.Zcurr.x_0;
0265             P(icoll,iZi) = lukestructs{icoll,iZi}.output.ZP0.x_rf_fsav;
0266             %
0267             flag_running(icoll,iZi) = false;
0268             flag_error(icoll,iZi) = false;
0269             %
0270             if isfield(lukestructs{icoll,iZi},'warning')
0271                 luke_warning{icoll,iZi} = [' **** WARNING : ',lukestructs{icoll,iZi}.warning];
0272             end
0273         elseif isfield(lukestructs{icoll,iZi},'err'),
0274             flag_running(icoll,iZi) = false;
0275             flag_error(icoll,iZi) = true;
0276         end
0277     end
0278 end
0279 %
0280 eta = j./P;
0281 %
0282 j_nr_k = NaN(ncoll,nZi);
0283 P_nr_k = NaN(ncoll,nZi);
0284 eta_nr_k = NaN(ncoll,nZi);
0285 %
0286 j_nr_k(1,1) = 0.05759;
0287 P_nr_k(1,1) = 0.004012;
0288 eta_nr_k(1,1) = 14.35;
0289 %
0290 j_nr_k(3,1) = 0.07092;
0291 P_nr_k(3,1) = 0.004294;
0292 eta_nr_k(3,1) = 16.52;
0293 %
0294 eta_0_nr_Karney = eta_nr_k(1,1)*6./(5 + Zi_list);
0295 eta_2_nr_Karney = eta_nr_k(3,1)*6./(5 + Zi_list);
0296 %
0297 %************************************************************************************************************************************
0298 %
0299 if ~any(flag_running(:)) && plot_results,
0300     %
0301     figure(1),clf
0302     %
0303     %leg = {'Linearized','High v limit','Maxwellian NR','Maxwellian'};
0304     leg = {'Linearized','High v limit','Maxwellian'};
0305     xlim = [0,max(Zi_list)];
0306     ylim = [0,0.1];
0307     xlab = 'Z_i';
0308     ylab = 'j';
0309     tit = '';
0310     siz = 20+14i;
0311     %
0312     graph1D_jd(Zi_list,j(3,:),0,0,xlab,ylab,tit,NaN,xlim,ylim,'-','none','r',2,siz,gca,0.9,0.7,0.7);
0313     graph1D_jd(Zi_list,j(2,:),0,0,'','','',NaN,xlim,ylim,'-','none','b',2,siz,gca);
0314     graph1D_jd(Zi_list,j(1,:),0,0,'','','',leg,xlim,ylim,'-','none','g',2,siz,gca);
0315     %graph1D_jd(xlim,[j_nr_k(3,:),j_nr_k(3,:)],0,0,'','','',NaN,xlim,ylim,'--','none','r',2,siz,gca);
0316     %graph1D_jd(xlim,[j_nr_k(1,1),j_nr_k(1,1)],0,0,'','','',NaN,xlim,ylim,'--','none','m',2,siz,gca);
0317     %
0318     set(gca,'ytick',[0:0.2:1]*ylim(2))
0319     set(gca,'xtick',[0:0.2:1]*xlim(2))
0320     %
0321     figure(2),clf
0322     %
0323     ylim = [0,0.02];
0324     ylab = 'P';
0325     %
0326     graph1D_jd(Zi_list,P(3,:),0,0,xlab,ylab,tit,NaN,xlim,ylim,'-','none','r',2,siz,gca,0.9,0.7,0.7);
0327     graph1D_jd(Zi_list,P(2,:),0,0,'','','',NaN,xlim,ylim,'-','none','b',2,siz,gca);
0328     graph1D_jd(Zi_list,P(1,:),0,0,'','','',leg,xlim,ylim,'-','none','g',2,siz,gca);
0329     %graph1D_jd(xlim,[P_nr_k(3,1),P_nr_k(3,1)],0,0,'','','',NaN,xlim,ylim,'--','none','r',2,siz,gca);
0330     %graph1D_jd(xlim,[P_nr_k(1,1),P_nr_k(1,1)],0,0,'','','',NaN,xlim,ylim,'--','none','m',2,siz,gca);
0331     %
0332     set(gca,'ytick',[0:0.2:1]*ylim(2))
0333     set(gca,'xtick',[0:0.2:1]*xlim(2))
0334     %
0335     figure(3),clf
0336     %
0337     ylim = [0,20];
0338     ylab = 'j/P';
0339     %
0340     graph1D_jd(Zi_list,eta(3,:),0,0,xlab,ylab,tit,NaN,xlim,ylim,'-','none','r',2,siz,gca,0.9,0.7,0.7);
0341     graph1D_jd(Zi_list,eta(2,:),0,0,'','','',NaN,xlim,ylim,'-','none','b',2,siz,gca);
0342     graph1D_jd(Zi_list,eta(1,:),0,0,'','','',leg,xlim,ylim,'-','none','g',2,siz,gca);
0343     graph1D_jd(Zi_list,eta_2_nr_Karney,0,0,'','','',NaN,xlim,ylim,'--','none','b',0.5,siz,gca);
0344     graph1D_jd(Zi_list,eta_0_nr_Karney,0,0,'','','',NaN,xlim,ylim,'--','none','r',0.5,siz,gca);
0345     %
0346     set(gca,'ytick',[0:0.2:1]*ylim(2))
0347     set(gca,'xtick',[0:0.2:1]*xlim(2))
0348     %
0349     print_jd(p_opt,'fig_j_Zi','./figures',1)
0350     print_jd(p_opt,'fig_P_Zi','./figures',2)
0351     print_jd(p_opt,'fig_eta_Zi','./figures',3)
0352     %
0353 end
0354 %
0355 %************************************************************************************************************************************
0356 %
0357 format
0358 %
0359 if ~isempty(link_mode) && ~isnan(link_mode) && link_mode >= 0,
0360     link_modestr = ['_link_mode_',num2str(link_mode)];
0361 elseif ~isempty(link_mode) && link_mode == -1,
0362     link_modestr = '_noremote';
0363 else
0364     link_modestr = '';
0365 end
0366 %
0367 diaryname = ['res_karney_Zi_remote_',datestr(now,29),'_remnum_',num2str(remnum),'_mdce_',num2str(mdce_mode),link_modestr,'.log'];
0368 %
0369 if input_dke_yp(['Do you want to save the results in the log file ',diaryname,' ? (y/n)'],'y',{'y','n'}) == 'y',
0370     %
0371     if exist(diaryname,'file'),
0372         delete(diaryname);
0373     end
0374     %
0375     diary(diaryname);
0376     %
0377 else
0378     diaryname = '';
0379 end
0380 %
0381 disp(' ')
0382 if remnum == 0,
0383     disp(['Calculation distributed locally with the MDCE profile : ',dkepath.distprofiles{dist_mode}]);
0384 else
0385     disp(['Calculation distributed on remote host : ',dkepath.profilestr{remnum},' with MDCE mode ',num2str(mdce_mode)]);
0386 end
0387 disp(' ')
0388 %
0389 
0390 disp('----------------------')
0391 disp(['Comparison LUKE/Karney - Non relativistic case (3,5,1) - with tn = ',int2str(dkeparam.tn),', dtn = ',int2str(dkeparam.dtn),' and nit_f = ',int2str(dkeparam.nit_f)])
0392 disp('----------------------')
0393 %
0394 for iZi = [1,nZi],
0395     for icoll = 1:ncoll,
0396         casestr = ['Zi = ',num2str(Zi_list(iZi)),' ; coll mode = ',num2str(coll_mode_list(icoll))];
0397         if flag_running(icoll,iZi),
0398             disp([casestr,' : running'])
0399         elseif ~isnan(j(icoll,iZi)),
0400             disp([casestr,' : j = ',num2str(j(icoll,iZi)),'/',num2str(j_nr_k(icoll,iZi)),' ; P = ',num2str(P(icoll,iZi)),'/',num2str(P_nr_k(icoll,iZi)),' ; j/P = ',num2str(eta(icoll,iZi)),'/',num2str(eta_nr_k(icoll,iZi)),luke_warning{icoll,iZi}]) 
0401         else      
0402             disp([casestr,' : failed'])
0403         end
0404     end
0405 end
0406 %
0407 if ~isempty(diaryname),
0408     diary off
0409 end
0410 %
0411 %************************************************************************************************************************************
0412 %
0413 if remtimout == 0,
0414     if ~any(flag_running(:)) && ~any(flag_error(:)),
0415         if exist(filename)
0416             delete(filename);% jobs completed, file deleted
0417         end
0418     else
0419         if any(flag_error(:)),
0420             save(filename);% file saved for future check
0421             info_dke_yp(2,['Data saved in ',filename]);
0422         end
0423     end
0424 end