proc_rundke_yp

0001 function [] = proc_rundke_yp(nfig0,filename,ir_display,save_mode)
0002 %
0003 permission = test_permissions_yp;
0004 %
0005 if ~permission 
0006     disp('Please move the script to a local folder where you have write permission before to run it')
0007     return;
0008 end
0009 %
0010 syn_mode = 0;
0011 print_mode = -1;
0012 %
0013 style0 = 'none';
0014 style = '-';
0015 style2 = '--';
0016 marker = 'none';
0017 marker1 = '+';
0018 marker2 = 'o';
0019 marker3 = 's';
0020 marker4 = '.';
0021 color1 = 'k';
0022 color2 = 'r';
0023 color3 = 'b';
0024 color4 = [0,0.75,0];
0025 color5 = 'm';
0026 color6 = 'c';
0027 %
0028 width = 0.5;
0029 width2 = 2;
0030 siz = 20;
0031 %
0032 red = 0.9;
0033 lspace = 0.7;
0034 lspace2 = 0.5;
0035 bspace = 0.7;
0036 bspace2 = 0.5;
0037 %
0038 load(filename,'ZP0_10','Zcurr10','gridDKE','mksa','equilDKE','dke_out10','dkedisplay','dkeparam')
0039 %
0040 J_rf = -Zcurr10.x_0_fsav*mksa.j_ref;
0041 P_rf = ZP0_10.x_rf_fsav.'*mksa.P_ref; 
0042 %
0043 xrho = gridDKE.xrho; 
0044 %
0045 I_rf = sum(J_rf.*equilDKE.xdA_dke);
0046 p_rf = sum(P_rf.*equilDKE.xdV_2piRp_dke);
0047 Eta_rf = I_rf/p_rf;
0048 %
0049 disp(['RF current: ',num2str(I_rf),' MA'])
0050 disp(['RF power: ',num2str(p_rf),' MW'])
0051 disp(['RF CD efficiency: ',num2str(Eta_rf),' A/W'])
0052 %
0053 % figures
0054 %
0055 figure(nfig0+1),clf
0056 %
0057 xlab = 'r/a';
0058 ylab = 'P_{rf} (MW/m^3)';
0059 tit = 'Density of RF Power Absorbed';
0060 %
0061 xlim = [0,1];
0062 xtick = 0:0.2:1;
0063 ylim = NaN;
0064 %
0065 x = xrho; 
0066 y = P_rf;
0067 %
0068 graph1D_jd(x,y,0,0,xlab,ylab,tit,NaN,xlim,ylim,style,marker1,color2,width2,siz,gca,red,lspace,bspace2);
0069 set(gca,'xtick',xtick);
0070 %
0071 if save_mode >= 0,
0072      savefig_dke_yp(['P_rf'])
0073 end
0074 %
0075 figure(nfig0+2),clf
0076 %
0077 ylab = 'J (MA/m^2)';
0078 tit = 'Current Density';
0079 %
0080 y = J_rf;
0081 %
0082 graph1D_jd(x,y,0,0,xlab,ylab,tit,NaN,xlim,ylim,style,marker1,color2,width2,siz,gca,red,lspace,bspace2);
0083 set(gca,'xtick',xtick);
0084 set(gca,'ylim',[0,max(get(gca,'ylim'))])
0085 %
0086 if save_mode >= 0,
0087     savefig_dke_yp(['J_rf'])
0088 end
0089 %
0090 figure(nfig0+3),clf
0091 %
0092 ylab = '\eta (A/W)';
0093 tit = 'CD efficiency';
0094 %
0095 y = J_rf./P_rf/(2*pi*equilDKE.Rp);
0096 %
0097 graph1D_jd(x(1:15),y(1:15),0,0,xlab,ylab,tit,NaN,xlim,ylim,style,marker1,color2,width2,siz,gca,red,lspace,bspace2);
0098 set(gca,'xtick',xtick);
0099 %
0100 if save_mode >= 0,
0101     savefig_dke_yp(['eta_rf'])
0102 end
0103 %
0104 figure(nfig0+4),clf
0105 %
0106 [logXf0_cyl_ref,ppar_cyl_ref,dppar_cyl_ref,pperp_cyl_ref,dpperp_cyl_ref] = s2c_dke_yp(log(dke_out10.XXf0(:,:,ir_display)),gridDKE.pn,gridDKE.mhu,dkedisplay.dp_cyl);%Spherical to cylindrical coordinate transformation
0107 Xf0_cyl_ref = exp(logXf0_cyl_ref);%For accurate representation in figures
0108 %
0109 [logXfinit_cyl_ref,ppar_cyl_ref,dppar_cyl_ref,pperp_cyl_ref,dpperp_cyl_ref] = s2c_dke_yp(log(dke_out10.XXfinit(:,:,ir_display)),gridDKE.pn,gridDKE.mhu,dkedisplay.dp_cyl);%Spherical to cylindrical coordinate transformation
0110 Xfinit_cyl_ref = exp(logXfinit_cyl_ref);%For accurate representation in figures
0111 %
0112 xyD0_rf = ones(dke_out10.gridindex_rf.nr,1);
0113 [ZXXD_rf] = raysum_rf_jd(dkeparam,dke_out10.ZXYD_rf,dke_out10.ZXYD_rf_tp,dke_out10.ZXYF_rf_tp,xyD0_rf,dke_out10.gridindex_rf);
0114 [XDparpar_cyl_ref,ppar_cyl_ref,dppar_cyl_ref,pperp_cyl_ref,dpperp_cyl_ref] = s2c_dke_yp(ZXXD_rf.parpar_ij(:,:,ir_display),gridDKE.pn,gridDKE.mhu,dkedisplay.dp_cyl);%Spherical to cylindrical coordinate transformation
0115 %
0116 Xppar_cyl_ref = ones(length(pperp_cyl_ref),1)*ppar_cyl_ref(:)';%Grid for the distribution functions only
0117 Xpperp_cyl_ref = pperp_cyl_ref(:)*ones(1,length(ppar_cyl_ref));
0118 Xp_cyl_ref = sqrt(Xppar_cyl_ref.*Xppar_cyl_ref + Xpperp_cyl_ref.*Xpperp_cyl_ref);
0119 %
0120 Xf0_cyl_ref = Xf0_cyl_ref.*(Xf0_cyl_ref>0);%Remove negative values
0121 Xf0_cyl_ref = Xf0_cyl_ref.*(Xp_cyl_ref<max(ppar_cyl_ref));%Remove values above max(ppar)
0122 Xf0_cyl_ref(isnan(Xf0_cyl_ref)) = 0;%Remove NaN values
0123 Xf0_cyl_ref(abs(Xf0_cyl_ref) < eps) = 0;%Remove values less than computer accuracy
0124 
0125 %
0126 Xfinit_cyl_ref = Xfinit_cyl_ref.*(Xfinit_cyl_ref>0);%Remove negative values
0127 Xfinit_cyl_ref = Xfinit_cyl_ref.*(Xp_cyl_ref<max(ppar_cyl_ref));%Remove values above max(ppar)
0128 Xfinit_cyl_ref(isnan(Xfinit_cyl_ref)) = 0;%Remove NaN values
0129 %
0130 
0131 
0132 Xppar_cyl_ref(abs(Xf0_cyl_ref) < eps) = 0;
0133 Xpperp_cyl_ref(abs(Xf0_cyl_ref) < eps) = 0;
0134 
0135 
0136 ippar = find(Xppar_cyl_ref == max(max(Xppar_cyl_ref)));
0137 ipperp = find(Xpperp_cyl_ref == max(max(Xpperp_cyl_ref)));
0138 pparmaxcode = max(Xppar_cyl_ref(ippar));
0139 
0140 pperpmaxcode = max(Xpperp_cyl_ref(ipperp));
0141 
0142 pmaxcode = max([pparmaxcode,pperpmaxcode]);
0143 
0144 %
0145 Xp_cyl_ref = sqrt(Xppar_cyl_ref.*Xppar_cyl_ref + Xpperp_cyl_ref.*Xpperp_cyl_ref);
0146 %
0147 XDparpar_cyl_ref = XDparpar_cyl_ref.*(XDparpar_cyl_ref>0);
0148 XDparpar_cyl_ref = XDparpar_cyl_ref.*(Xp_cyl_ref<max(ppar_cyl_ref));
0149 XDparpar_cyl_ref(isnan(XDparpar_cyl_ref))=0;
0150 %
0151 Fpar0_cyl_ref = 2*pi*integral_dke_jd(dpperp_cyl_ref,Xpperp_cyl_ref.*Xf0_cyl_ref,1);
0152 Fparinit_cyl_ref = 2*pi*integral_dke_jd(dpperp_cyl_ref,Xpperp_cyl_ref.*Xfinit_cyl_ref,1);        
0153 %
0154 %
0155 pnmax_S = max(gridDKE.pn_S);
0156 %
0157 x1 = ppar_cyl_ref;
0158 xlab1 = 'p_{||}/p_{Te}';
0159 %xlim1 = [-pnmax_S,pnmax_S];
0160 xlim1 = [-pmaxcode,pmaxcode];
0161 %
0162 x2 = pperp_cyl_ref;
0163 xlab2 = 'p_{\perp}/p_{Te}';
0164 %xlim2 = [0,pnmax_S];
0165 xlim2 = [0,pmaxcode];
0166 %
0167 leg = ['f_{init}';'f_0     '];
0168 %
0169 y1 = Xfinit_cyl_ref;
0170 y2 = Xf0_cyl_ref;
0171 %
0172 x9 = x1;
0173 y9 = sqrt(pnmax_S^2 - x9.^2);
0174 %
0175 graph2D_jd(x1,x2,y1','','','',xlim1,xlim2,0,mksa.betath_ref,[0:0.5:pmaxcode],0,style,color3,width,siz,red,max(max(y1)));
0176 graph2D_jd(x1,x2,y2','','','',xlim1,xlim2,0,mksa.betath_ref,[0:0.5:pmaxcode],0,style,color2,width2,siz,1,max(max(y1)));
0177 graph1D_jd(x9,y9,0,0,xlab1,xlab2,['2D Distribution function f_0 at \rho = ',num2str(gridDKE.xrho(ir_display))],NaN,xlim1,xlim2,style,marker,color1,width2,siz);
0178 axis('equal');axis([xlim1,xlim2])
0179 %
0180 ripple = 0;
0181 xpndetrap = dke_out10.Zmripple.xpndetrap;
0182 xmhusupertrap = cos((90 - dke_out10.Zmripple.xthetaloss)*pi/180);
0183 pparsupertrap = xpndetrap(ir_display)*gridDKE.mhu(abs(gridDKE.mhu)<=xmhusupertrap(ir_display));%For contour plot of the supertrapped domain
0184 pperpsupertrap = xpndetrap(ir_display)*sqrt(1-gridDKE.mhu(abs(gridDKE.mhu)<=xmhusupertrap(ir_display)).^2);%For contour plot of the supertrapped domain
0185 
0186 bounce_ripple_jd(dkeparam.bounce_mode,pnmax_S,gridDKE.xmhubounce2,ir_display,[]);
0187 %
0188 if save_mode >= 0,
0189     savefig_dke_yp(['f0_2D_ir',int2str(ir_display)])
0190 end
0191 %
0192 figure(nfig0+5),clf
0193 %
0194 y1 = XDparpar_cyl_ref;
0195 cont = 10;[0:ceil(max(max(XDparpar_cyl_ref)))];
0196 %
0197 graph2D_jd(x1,x2,y1','','','',xlim1,xlim2,0,NaN,cont,0,style,NaN,width,siz,0.95);
0198 graph1D_jd(x9,y9,0,0,xlab1,xlab2,'D_{|| ||LH} [\nu_e p_{Te}^2]',NaN,xlim1,xlim2,style,marker,color1,width2,siz);
0199 axis('equal');axis([xlim1,xlim2])
0200 bounce_ripple_jd(dkeparam.bounce_mode,pnmax_S,gridDKE.xmhubounce2,ir_display,[]);
0201 colorbar
0202 %
0203 %
0204 %figure(4),
0205 %graph2D_jd(x1,x2,y1','','','',xlim1,xlim2,0,NaN,cont,0,style2,color4,width2,siz);
0206 %legend(leg)
0207 %
0208 %
0209 % ---------- parallel distribution difference-----------
0210 %
0211 if save_mode >= 0,
0212     savefig_dke_yp(['Dparpar_ir',int2str(ir_display)])
0213 end
0214 %
0215 figure(nfig0+6),clf
0216 %
0217 y1 = Fpar0_cyl_ref - Fparinit_cyl_ref;
0218 ylab = 'F_{||0} - F_{||M}';
0219 ylim = NaN;
0220 %
0221 graph1D_jd(x1,y1,0,0,xlab1,ylab,'',NaN,xlim1,ylim,style,marker,color2,width2,siz,gca,red,lspace,bspace);
0222 %
0223 % ---------- parallel distribution -----------
0224 %
0225 if save_mode >= 0,
0226     savefig_dke_yp(['F0mFM_1D_ir',int2str(ir_display)])
0227 end
0228 
0229 figure(nfig0+7),clf
0230 %
0231 y1 = Fparinit_cyl_ref;
0232 y2 = Fpar0_cyl_ref;
0233 ylab = 'F_{0||}';
0234 ylim = [eps,1];
0235 %
0236 graph1D_jd(x1,y1,0,1,'','','',NaN,xlim1,ylim,style,marker,color3,width,siz,gca,red,lspace,bspace);
0237 graph1D_jd(x1,y2,0,1,xlab1,ylab,'',leg,xlim1,ylim,style,marker,color2,width2,siz);
0238 
0239 %
0240 if save_mode >= 0,
0241     savefig_dke_yp(['F0_1D_ir',int2str(ir_display)])
0242 end
0243  
0244    
0245 
0246 
0247