rfdiff_dke_fw_jd

PURPOSE ^

SYNOPSIS ^

function [ZXYD_rf,ZXYF_rf,ZXYD_rf_tp,ZXYF_rf_tp,gridindex_rf] = rfdiff_dke_fw_jd(dkeparam,display,equilDKE,gridDKE,Zmomcoef,Zbouncecoef,mksa,radialDKE,waves);%RF waves flux coefficients

DESCRIPTION ^

 full wave diffusion tensor calculation

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [ZXYD_rf,ZXYF_rf,ZXYD_rf_tp,ZXYF_rf_tp,gridindex_rf] = rfdiff_dke_fw_jd(dkeparam,display,equilDKE,gridDKE,Zmomcoef,Zbouncecoef,mksa,radialDKE,waves);%RF waves flux coefficients
0002 %
0003 % full wave diffusion tensor calculation
0004 %
0005 method = 'spline';
0006 %
0007 xrho = equilDKE.xrho;
0008 xrho_S = radialDKE.xrho_S_dke;
0009 xdrho = diff(radialDKE.xrho_S_dke);
0010 %
0011 % grids for interpolation
0012 %
0013 xrho_I = [0,xrho,1];
0014 xrhoP_I = [0,equilDKE.xrhoP,1];
0015 xrhoT_I = [0,equilDKE.xrhoT,1];
0016 xrhoV_I = [0,equilDKE.xrhoV,1];
0017 %
0018 XXmhup = gridDKE.XXmhup;
0019 XXmhu = gridDKE.XXmhu;
0020 XXmhum = gridDKE.XXmhum;
0021 %
0022 XX1mmhu2p = gridDKE.XX1mmhu2p;
0023 XX1mmhu2 = gridDKE.XX1mmhu2;
0024 XX1mmhu2m = gridDKE.XX1mmhu2m;
0025 %
0026 XXpn2p = gridDKE.XXpn2p;
0027 XXpn2 = gridDKE.XXpn2;
0028 XXpn2m = gridDKE.XXpn2m;
0029 %
0030 XXgammap = Zmomcoef.XXgammap;
0031 XXgamma = Zmomcoef.XXgamma;
0032 XXgammam = Zmomcoef.XXgammam;
0033 %
0034 XXheavisidep = Zbouncecoef.XXheaviside_p;
0035 XXheaviside = Zbouncecoef.XXheaviside;
0036 XXheavisidem = Zbouncecoef.XXheaviside_m;
0037 %
0038 XXv2p = XXpn2p./XXgammap./XXgammap;
0039 XXv2 = XXpn2./XXgamma./XXgamma;
0040 XXv2m = XXpn2m./XXgammam./XXgammam;
0041 %
0042 % parallel velocity
0043 %
0044 XXvpar.ipj = XXmhu.*gridDKE.XXpnp./XXgammap;
0045 XXvpar.ij = XXmhu.*gridDKE.XXpn./XXgamma;
0046 XXvpar.imj = XXmhu.*gridDKE.XXpnm./XXgammam;
0047 XXvpar.ijp = XXmhup.*gridDKE.XXpn./XXgamma;
0048 XXvpar.ijm = XXmhum.*gridDKE.XXpn./XXgamma;
0049 %
0050 dvpar = mean([mean(diff(XXvpar.ij(:,:,1))),mean(diff(XXvpar.ij(:,:,1).'))]);
0051 %
0052 SXX = size(XXmhu);
0053 %
0054 npn = SXX(1);
0055 nmhu = SXX(2);
0056 nr = SXX(3);
0057 %
0058 XXyn = zeros(SXX);%Cerenkov resonance only
0059 %
0060 XXDld.ipj = zeros(SXX);
0061 XXDld.ij = zeros(SXX);
0062 XXDld.imj = zeros(SXX);
0063 XXDld.ijp = zeros(SXX);
0064 XXDld.ijm = zeros(SXX);
0065 %
0066 XXDmx.ipj = zeros(SXX);
0067 XXDmx.ij = zeros(SXX);
0068 XXDmx.imj = zeros(SXX);
0069 XXDmx.ijp = zeros(SXX);
0070 XXDmx.ijm = zeros(SXX);
0071 %
0072 XXDtt.ipj = zeros(SXX);
0073 XXDtt.ij = zeros(SXX);
0074 XXDtt.imj = zeros(SXX);
0075 XXDtt.ijp = zeros(SXX);
0076 XXDtt.ijm = zeros(SXX);
0077 %
0078 % calculation of the DQL components on the DKE grid
0079 % each contribution on the wave.Xvpar grid is accounted for.
0080 %
0081 for iw = 1:length(waves);
0082     %
0083     wave = waves{iw};
0084     %
0085     % normalization
0086     %
0087     % The parallel quasilinear coefficient D// is considered here, with
0088     % D// = D_LD + vperp^2 D_MX +  vperp^4 D_TT.
0089     %
0090     % We assume the following normalization.
0091     %
0092     % - xvpar and vperp are normalized to c
0093     % - xDql_LD, xDql_MX and xDql_TT are  normalized to (mc)^2
0094     %
0095     Xvpar = wave.Xvpar;
0096     %Xvpar_S = wave.Xvpar_S;
0097     XDql_LD = wave.XDql_LD;
0098     XDql_MX = wave.XDql_MX;
0099     XDql_TT = wave.XDql_TT;
0100     %rho = wave.rho;
0101     rho_S = wave.rho_S;
0102     rhotype = wave.rhotype;
0103     %
0104     [nv,nr_fw] = size(Xvpar);
0105     %
0106     % In LUKE, the velocities are normalized to vT = betath_ref*c
0107     % and the diffusion coefficients to (nhu_ref pT^2)
0108     %
0109     Xvpar = Xvpar/mksa.betath_ref;
0110     XDql_LD = XDql_LD/(mksa.nhu_ref*mksa.betath_ref^2);
0111     XDql_MX = XDql_MX/(mksa.nhu_ref*mksa.betath_ref^2)*mksa.betath_ref^2;
0112     XDql_TT = XDql_TT/(mksa.nhu_ref*mksa.betath_ref^2)*mksa.betath_ref^4;
0113     %
0114 %     Xdvpar = [Xvpar(2,:) - Xvpar(1,:);...
0115 %               (Xvpar(3:end,:) - Xvpar(1:end-2,:))/2;...
0116 %               Xvpar(end,:) - Xvpar(end-1,:)];
0117     %
0118     if rhotype == 'g',
0119         %wrho = rho;
0120         wrho_S = rho_S;
0121     elseif rhotype == 'p',
0122         %wrho = interp1(xrhoP_I,xrho_I,rho,method);
0123         wrho_S = interp1(xrhoP_I,xrho_I,rho_S,method);
0124     elseif rhotype == 't',
0125         %wrho = interp1(xrhoT_I,xrho_I,rho,method);
0126         wrho_S = interp1(xrhoT_I,xrho_I,rho_S,method);
0127     elseif rhotype == 'v',
0128         %wrho = interp1(xrhoV_I,xrho_I,rho,method);
0129         wrho_S = interp1(xrhoV_I,xrho_I,rho_S,method);
0130     else
0131         error('rho type not recognized')
0132     end
0133     %
0134     % There are several possible methods :
0135     % - smooth XDql (in both (or either) v// and (or) rho)
0136     % - take every contribution from XDql and map it on (rho,p,mhu) space
0137     % - a combination of the two (smooth in v// for each rho contribution...)
0138     %
0139     for ir_fw = 1:nr_fw,
0140         %
0141         disp(['Contruction of DQL from FW data : ir_rw = ',num2str(ir_fw),'/',num2str(nr_fw)]);
0142         %
0143         xmask = max([zeros(1,nr);(min([repmat(wrho_S(ir_fw + 1),[1,nr]);xrho_S(2:end)]) - max([repmat(wrho_S(ir_fw),[1,nr]);xrho_S(1:end-1)]))])./xdrho;
0144         XXmask = repmat(reshape(xmask,[1,1,nr]),[npn,nmhu,1]);
0145         %
0146         yvpar = Xvpar(:,ir_fw);
0147         yDql_LD = smooth_jd(yvpar,XDql_LD(:,ir_fw),dvpar);
0148         yDql_MX = smooth_jd(yvpar,XDql_MX(:,ir_fw),dvpar);
0149         yDql_TT = smooth_jd(yvpar,XDql_TT(:,ir_fw),dvpar);
0150         %
0151         XXDld.ipj = XXDld.ipj + XXmask.*interp1(yvpar,yDql_LD,XXvpar.ipj,'linear',0);
0152         XXDld.ij = XXDld.ipj + XXmask.*interp1(yvpar,yDql_LD,XXvpar.ij,'linear',0);
0153         XXDld.imj = XXDld.ipj + XXmask.*interp1(yvpar,yDql_LD,XXvpar.imj,'linear',0);
0154         XXDld.ijp = XXDld.ipj + XXmask.*interp1(yvpar,yDql_LD,XXvpar.ijp,'linear',0);
0155         XXDld.ijm = XXDld.ipj + XXmask.*interp1(yvpar,yDql_LD,XXvpar.ijm,'linear',0);
0156         %
0157         XXDmx.ipj = XXDmx.ipj + XXmask.*interp1(yvpar,yDql_MX,XXvpar.ipj,'linear',0);
0158         XXDmx.ij = XXDmx.ipj + XXmask.*interp1(yvpar,yDql_MX,XXvpar.ij,'linear',0);
0159         XXDmx.imj = XXDmx.ipj + XXmask.*interp1(yvpar,yDql_MX,XXvpar.imj,'linear',0);
0160         XXDmx.ijp = XXDmx.ipj + XXmask.*interp1(yvpar,yDql_MX,XXvpar.ijp,'linear',0);
0161         XXDmx.ijm = XXDmx.ipj + XXmask.*interp1(yvpar,yDql_MX,XXvpar.ijm,'linear',0);
0162         %
0163         XXDtt.ipj = XXDtt.ipj + XXmask.*interp1(yvpar,yDql_TT,XXvpar.ipj,'linear',0);
0164         XXDtt.ij = XXDtt.ipj + XXmask.*interp1(yvpar,yDql_TT,XXvpar.ij,'linear',0);
0165         XXDtt.imj = XXDtt.ipj + XXmask.*interp1(yvpar,yDql_TT,XXvpar.imj,'linear',0);
0166         XXDtt.ijp = XXDtt.ipj + XXmask.*interp1(yvpar,yDql_TT,XXvpar.ijp,'linear',0);
0167         XXDtt.ijm = XXDtt.ipj + XXmask.*interp1(yvpar,yDql_TT,XXvpar.ijm,'linear',0);
0168         %
0169     end
0170 end
0171 %
0172 XXDrf.ipj = XXmhu.*XXmhu.*(XXDld.ipj + XXv2p.*XX1mmhu2.*XXDmx.ipj + XXv2p.*XXv2p.*XX1mmhu2.*XX1mmhu2.*XXDtt.ipj)./XX1mmhu2;
0173 XXDrf.ij = XXmhu.*XXmhu.*(XXDld.ij + XXv2.*XX1mmhu2.*XXDmx.ij + XXv2.*XXv2.*XX1mmhu2.*XX1mmhu2.*XXDtt.ij)./XX1mmhu2;
0174 XXDrf.imj = XXmhu.*XXmhu.*(XXDld.imj + XXv2m.*XX1mmhu2.*XXDmx.imj + XXv2m.*XXv2m.*XX1mmhu2.*XX1mmhu2.*XXDtt.imj)./XX1mmhu2;
0175 XXDrf.ijp = XXmhup.*XXmhup.*(XXDld.ijp + XXv2.*XX1mmhu2p.*XXDmx.ijp + XXv2.*XXv2.*XX1mmhu2p.*XX1mmhu2p.*XXDtt.ijp)./XX1mmhu2p;
0176 XXDrf.ijm = XXmhum.*XXmhum.*(XXDld.ijm + XXv2.*XX1mmhu2m.*XXDmx.ijm + XXv2.*XXv2.*XX1mmhu2m.*XX1mmhu2m.*XXDtt.ijm)./XX1mmhu2m;
0177 %
0178 XXDrf_mhuflip.ipj = flipdim(XXDrf.ipj,2);
0179 XXDrf_mhuflip.ij = flipdim(XXDrf.ij,2);
0180 XXDrf_mhuflip.imj = flipdim(XXDrf.imj,2);
0181 XXDrf_mhuflip.ijp = flipdim(XXDrf.ijm,2);
0182 XXDrf_mhuflip.ijm = flipdim(XXDrf.ijp,2);
0183 %
0184 XXDrf.ipj(~XXheaviside) = (XXDrf.ipj(~XXheaviside) + XXDrf_mhuflip.ipj(~XXheaviside))/2;
0185 XXDrf.ij(~XXheaviside) = (XXDrf.ij(~XXheaviside) + XXDrf_mhuflip.ij(~XXheaviside))/2;
0186 XXDrf.imj(~XXheaviside) = (XXDrf.imj(~XXheaviside) + XXDrf_mhuflip.imj(~XXheaviside))/2;
0187 XXDrf.ijp(~XXheavisidep) = (XXDrf.ijp(~XXheavisidep) + XXDrf_mhuflip.ijp(~XXheavisidep))/2;
0188 XXDrf.ijm(~XXheavisidem) = (XXDrf.ijm(~XXheavisidem) + XXDrf_mhuflip.ijm(~XXheavisidem))/2;
0189 %
0190 % QL diffusion tensor elements
0191 %
0192 ZXXD_rf.pp_ipj = XX1mmhu2.*XXDrf.ipj;   
0193 ZXXD_rf.pp_imj = XX1mmhu2.*XXDrf.imj;
0194 ZXXD_rf.pm_ipj = - sqrt(XX1mmhu2).*(XX1mmhu2 - XXyn).*XXDrf.ipj./XXmhu;
0195 ZXXD_rf.pm_imj = - sqrt(XX1mmhu2).*(XX1mmhu2 - XXyn).*XXDrf.imj./XXmhu;
0196 ZXXD_rf.pr_ipj = zeros(size(XXmhu));
0197 ZXXD_rf.pr_imj = zeros(size(XXmhu));
0198 ZXXF_rf.p_ipj = zeros(size(XXmhu));
0199 ZXXF_rf.p_imj = zeros(size(XXmhu));
0200 %
0201 ZXXD_rf.mm_ijp = (XX1mmhu2p - XXyn).^2.*XXDrf.ijp./XXmhup./XXmhup;
0202 ZXXD_rf.mm_ijm = (XX1mmhu2m - XXyn).^2.*XXDrf.ijm./XXmhum./XXmhum;
0203 ZXXD_rf.mp_ijp = - sqrt(XX1mmhu2p).*(XX1mmhu2p - XXyn).*XXDrf.ijp./XXmhup;
0204 ZXXD_rf.mp_ijm = - sqrt(XX1mmhu2m).*(XX1mmhu2m - XXyn).*XXDrf.ijm./XXmhum;
0205 ZXXD_rf.mr_ijp = zeros(size(XXmhu));
0206 ZXXD_rf.mr_ijm = zeros(size(XXmhu));
0207 ZXXF_rf.m_ijp = zeros(size(XXmhu));
0208 ZXXF_rf.m_ijm = zeros(size(XXmhu));
0209 %
0210 ZXXD_rf.rr_lp = zeros(size(XXmhu));
0211 ZXXD_rf.rr_lm = zeros(size(XXmhu));
0212 ZXXD_rf.rp_lp = zeros(size(XXmhu));
0213 ZXXD_rf.rp_lm = zeros(size(XXmhu));
0214 ZXXD_rf.rm_lp = zeros(size(XXmhu));
0215 ZXXD_rf.rm_lm = zeros(size(XXmhu));
0216 ZXXF_rf.r_lp = zeros(size(XXmhu));
0217 ZXXF_rf.r_lm = zeros(size(XXmhu));
0218 %
0219 ZXXD_rf.pp_ij = XX1mmhu2.*XXDrf.ij; 
0220 ZXXD_rf.pm_ij = - sqrt(XX1mmhu2).*(XX1mmhu2 - XXyn).*XXDrf.ij./XXmhu;
0221 ZXXD_rf.mp_ij = - sqrt(XX1mmhu2).*(XX1mmhu2 - XXyn).*XXDrf.ij./XXmhu;
0222 %
0223 ZXXD_rf.parpar_ij = XX1mmhu2.*(1 - XXyn).^2.*XXDrf.ij./XXmhu./XXmhu; 
0224 ZXXD_rf.ij = XXDrf.ij; 
0225 %
0226 ZXYD_rf.pp_ipj = ZXXD_rf.pp_ipj(:);
0227 ZXYD_rf.pp_imj = ZXXD_rf.pp_imj(:);
0228 ZXYD_rf.pm_ipj = ZXXD_rf.pm_ipj(:);
0229 ZXYD_rf.pm_imj = ZXXD_rf.pm_imj(:);
0230 ZXYD_rf.pr_ipj = ZXXD_rf.pr_ipj(:);
0231 ZXYD_rf.pr_imj = ZXXD_rf.pr_imj(:);
0232 ZXYF_rf.p_ipj = ZXXF_rf.p_ipj(:);
0233 ZXYF_rf.p_imj = ZXXF_rf.p_imj(:);
0234 %
0235 ZXYD_rf.mm_ijp = ZXXD_rf.mm_ijp(:);
0236 ZXYD_rf.mm_ijm = ZXXD_rf.mm_ijm(:);
0237 ZXYD_rf.mp_ijp = ZXXD_rf.mp_ijp(:);
0238 ZXYD_rf.mp_ijm = ZXXD_rf.mp_ijm(:);
0239 ZXYD_rf.mr_ijp = ZXXD_rf.mr_ijp(:);
0240 ZXYD_rf.mr_ijm = ZXXD_rf.mr_ijm(:);
0241 ZXYF_rf.m_ijp = ZXXF_rf.m_ijp(:);
0242 ZXYF_rf.m_ijm = ZXXF_rf.m_ijm(:);
0243 %
0244 %
0245 ZXYD_rf.rr_lp = ZXXD_rf.rr_lp(:);
0246 ZXYD_rf.rr_lm = ZXXD_rf.rr_lm(:);
0247 ZXYD_rf.rp_lp = ZXXD_rf.rp_lp(:);
0248 ZXYD_rf.rp_lm = ZXXD_rf.rp_lm(:);
0249 ZXYD_rf.rm_lp = ZXXD_rf.rm_lp(:);
0250 ZXYD_rf.rm_lm = ZXXD_rf.rm_lm(:);
0251 ZXYF_rf.r_lp = ZXXF_rf.r_lp(:);
0252 ZXYF_rf.r_lm = ZXXF_rf.r_lm(:);
0253 %
0254 ZXYD_rf.pp_ij = ZXXD_rf.pp_ij(:);
0255 ZXYD_rf.pm_ij = ZXXD_rf.pm_ij(:);
0256 ZXYD_rf.mp_ij = ZXXD_rf.mp_ij(:);
0257 %
0258 ZXYD_rf.parpar_ij = ZXXD_rf.parpar_ij(:);
0259 ZXYD_rf.ij = ZXXD_rf.ij(:);
0260 %
0261 ZXYD_rf_tp = NaN;
0262 ZXYF_rf_tp = NaN;
0263 %
0264 XX_index = find(ones(SXX));
0265 %
0266 gridindex_rf.mask_r = ceil(XX_index/(npn*nmhu));
0267 gridindex_rf.mask_m = ceil((XX_index - (gridindex_rf.mask_r - 1)*npn*nmhu)/npn);
0268 gridindex_rf.mask_p = XX_index - ((gridindex_rf.mask_r - 1)*nmhu + (gridindex_rf.mask_m - 1))*npn;
0269 gridindex_rf.mask_y = ones(prod(SXX),1);
0270 gridindex_rf.mask_n = ones(prod(SXX),1);
0271 %
0272 gridindex_rf.npn = npn;
0273 gridindex_rf.nmhu = nmhu;
0274 gridindex_rf.nr = nr;
0275 gridindex_rf.nr_dke = nr;
0276 gridindex_rf.ny = 1;
0277 gridindex_rf.nn_rf = 1;
0278 
0279
Community support and wiki are available on Redmine. Last update: 18-Apr-2019.