make_wave_RFP

 script make_wave_MST_400ka

 runs C3PO based on formatted ray launching parameters and builds the WAVE structure

 by     J. Decker (DRFC/DSM/CEA) <joan.decker@cea.fr> 
 and   Y. Peysson (DRFC/DSM/CEA) <yves.peysson@cea.fr>

0001 function [] = make_wave_RFP
0002 % script make_wave_MST_400ka
0003 %
0004 % runs C3PO based on formatted ray launching parameters and builds the WAVE structure
0005 %
0006 % by     J. Decker (DRFC/DSM/CEA) <joan.decker@cea.fr>
0007 % and   Y. Peysson (DRFC/DSM/CEA) <yves.peysson@cea.fr>
0008 %
0009 close all
0010 %
0011 make_rayinit_RFP;
0012 %
0013 % Path parameters
0014 %
0015 id_dkepath = '';%For all paths used by DKE solver
0016 path_dkepath = '';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0017 %
0018 %
0019 id_wave = 'RFP';
0020 flag_analytic = 6;
0021 p_opt = -1;
0022 %
0023 % Equilibrium parameters
0024 %
0025 id_equil = 'RFP';%For plasma equilibrium
0026 path_equil = '../EQUIL/';%if nothing is specified, the working directory is first used and then MatLab is looking in all the path
0027 %
0028 % initial ray conditions
0029 %
0030 ids_rayinit = {'LH'};
0031 paths_rayinit = {''};
0032 %
0033 % C3PO computing parameters
0034 %
0035 mdce_mode_main_C3PO_jd = 0;%MatLab distributed computing environment disabled (0), enabled with the dedicated toolbox (1), enabled with a private method (2)for the function main_C3PO_jd.m (MDC toolbox must be installed for option 1)
0036 %
0037 % Display parameters
0038 %
0039 C3POdisplay_numeric.ray = 1;
0040 C3POdisplay_numeric.equilibrium = 0;
0041 C3POdisplay_numeric.p_opt = 2;%Printing or saving option of the figures
0042 %
0043 C3POdisplay_analytic.ray = 0;
0044 C3POdisplay_analytic.equilibrium = 0;
0045 C3POdisplay_analytic.p_opt = 2;%Printing or saving option of the figures
0046 %
0047 % Wave parameters
0048 %
0049 waveparam.mmode = -1;%cold plasma mode [1] : (-1) m (1) p, p is the slow mode when kperp > 0 (ex : LH slow wave)
0050 waveparam.kmode = 0;%(0:cold,1:warm,2:hot;3:weak realtivistic,4:full relativistic)
0051 %
0052 %Option parameter for FLR effects and cross-comparison between old FP code:
0053 %    - (0): all FLR effects
0054 %    - (1): small FLR effects and 1/vpar dependence
0055 %    - (2): small FLR effects and no 1/vpar dependence and old grid technique for DQL calculations (Karney, Shoucri) (see rfdiff_dke_jd)
0056 %
0057 waveparam.opt_rf = 1;
0058 %
0059 waveparam.dsmin = Inf;%minimum size for ray fragments
0060 waveparam.cut_fac = Inf;% the ray is cut when P/P_in = exp(-10*cut_fac). Use Inf if the whole ray is to be kept
0061 %
0062 n_rf_list = 0;
0063 %
0064 % -------------------------------------------------------------------------
0065 %
0066 % Global parameters for the vectorial magnetic equilibrium
0067 %
0068 % Equilibrium parameters for the ray-tracing
0069 %
0070 fitparam.mode_equil = 1;%Magnetic equilibrium grid type: (1): (psi-theta), (2): (x-y)
0071 fitparam.method = 'spline';%nearest,spline,pchip
0072 fitparam.nharm = NaN;%Number of harmonics in the magnetic equilibrium interpolation (NaN, Inf or empty, nharm = ntheta-1)
0073 fitparam.ngridresample = 1001;%Number of grid points for resampling the radial profile of magnetic equilibrium parameters
0074 %
0075 % Global parameters for the ray-tracing
0076 %
0077 rayparam.testmode = 0;
0078 rayparam.tensortype = waveparam.kmode;%(0:cold,1:warm,2:hot;3:weak relativistic,4:full relativistic)
0079 rayparam.t0 = 0;
0080 rayparam.tfinal = 1000;
0081 rayparam.dt0 = 1.e-4;
0082 rayparam.dS = 1.e-4;
0083 rayparam.tol = 1e-12;%when tolerance is increased (less accurate calculation of D=0), tfinal must be increased accordingly
0084 rayparam.kmax = 60000;
0085 rayparam.ncyclharm = 3;%number of cyclotron harmonics (just for hot and relativistic dielectric tensors)
0086 rayparam.reflection = 1;%1:Enforce wave reflection at plasma boundary, 0: the code calculates itself if the ray must leave of not the plasma
0087 rayparam.rel_opt = 1;%option for (1) relativistic or (0) non-relativistic calculations
0088 rayparam.nperp = 10000;%number of points in pperp integration for damping calculations
0089 rayparam.pperpmax = 10;%maximum value of pperp in damping calculations
0090 rayparam.tau_lim = Inf;%value of optical depth beyond which the wave is considered absorbed (usually 20. Otherwise Inf)
0091 rayparam.kextra = 1000;%number of calculations performed beyond the full linear absorption (for quasilinear calculations which may require more points)%
0092 % -------------------------------------------------------------------------
0093 %
0094 % Load structures
0095 %
0096 [dkepath,equil,rayinits] = load_structures_yp('dkepath',id_dkepath,path_dkepath,'equil',id_equil,path_equil,'rayinits',ids_rayinit,paths_rayinit);
0097 %
0098 % Vectorial description of the magnetic equilibrium
0099 %
0100 equil_fit = fitequil_yp(equil,fitparam.mode_equil,fitparam.method,fitparam.ngridresample,fitparam.nharm);%Build vectorized magnetic equilibrium structure
0101 info_dke_yp(2,['Vectorial form of the magnetic equilibrium ',equil.id,' is calculated.']);
0102 if C3POdisplay_numeric.equilibrium,testfitequil_yp(equil,equil_fit);end
0103 %
0104 % C3PO computing parameters
0105 %
0106 C3POparam.clustermode.main_C3PO_jd.scheduler.mode = mdce_mode_main_C3PO_jd;%MatLab distributed computing environment
0107 %
0108 % C3P0 ray tracing
0109 %
0110 wave_numeric = main_C3PO_jd(dkepath,[id_wave,'_numeric'],equil,'',rayinits{1},waveparam,fitparam,rayparam,C3POdisplay_numeric,C3POparam,[],[],0);
0111 %
0112 info_dke_yp(2,'Ray trajectories calculated (interpolated magnetic equilibrium)');
0113 %
0114 rayparam.tfinal = 600;
0115 wave_analytic = main_C3PO_jd(dkepath,[id_wave,'_analytic'],equil,'',rayinits{1},waveparam,fitparam,rayparam,C3POdisplay_analytic,C3POparam,[],[],flag_analytic);
0116 %
0117 info_dke_yp(2,'Ray trajectories calculated (analytic magnetic equilibrium)');
0118 %
0119 save_str = ['WAVE_',id_wave,'.mat'];
0120 save(save_str,'id_wave','wave_numeric','wave_analytic');
0121 %
0122 info_dke_yp(2,'Wave parameters saved');
0123 %
0124 % --- display results ---
0125 %
0126 rays = {wave_numeric.rays{1},wave_analytic.rays{1}};
0127 %
0128 legs = {'Numeric','Analytic'};
0129 %
0130 filename = ['Fig_',id_wave];
0131 %
0132 opt.p_opt = C3POdisplay_numeric.p_opt;
0133 opt.ntheta_fit = 65;
0134 opt.nrho_fit = 15;
0135 opt.propvar = 1;   
0136 %
0137 graph_comp_RT_jd(rays,legs,'',filename,opt)
0138 %
0139 % diary some results for CVS validation
0140 %
0141 diaryname = ['res_',id_wave,'.log'];
0142 [LUKEversion,matver] = LUKEversion_jd;
0143 %
0144 if exist(diaryname,'file'),
0145     delete(diaryname)
0146 end
0147 %
0148 diary(diaryname)
0149 %
0150 disp(['Simulation id: ',id_wave]);
0151 disp(['Simulation date: ',timeid_jd(timeid_jd(clock))]);
0152 disp(['User id: ',userid_jd]);
0153 disp(['LUKE Version: ',LUKEversion]);
0154 disp(['MATLAB Version: ',matver]);
0155 disp(['MACHINE: ',dkepath.hostname]);
0156 %
0157 disp(['----------------------'])
0158 %
0159 for iray = 1:10,
0160     disp(['Npar(numeric) = ',num2str(wave_numeric.rays{1}.sNpar(iray),'%17.15f'),' - Npar(analytic) = ',num2str(wave_analytic.rays{1}.sNpar(iray),'%17.15f')])
0161 end
0162 disp(['----------------------'])
0163 %
0164 diary off