rayevolution_yp

PURPOSE ^

SYNOPSIS ^

function [P0] = rayevolution_yp(wavetype,equil_fit,waves,mode,rho_S,itn)

DESCRIPTION ^

 Analysis of the ray dynamics and the power spectrum (check is the
 spectral gap is filled

 INPUT
    - wavetype: wave type (0 -> LH, 1 -> EC, 2 ->...). May be in between if different waves types are launched into the plasma
    - equil_fit: interpolated magnetic equilibrium structure
    - waves: wave structure
    - mode: 1 -> the first wave is without fluctuations, 0 -> no distinction is done between rays [1,1]
    - itn: iteration number (optional) [1,1]

 OUTPUT:
    - P0: total rf power (MW) [1,1]

 By Yves Peysson (CEA-DRFC, yves.peysson@cea.fr) and Joan Decker (CEA-DRFC, joan.decker@cea.fr)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [P0] = rayevolution_yp(wavetype,equil_fit,waves,mode,rho_S,itn)
0002 %
0003 % Analysis of the ray dynamics and the power spectrum (check is the
0004 % spectral gap is filled
0005 %
0006 % INPUT
0007 %    - wavetype: wave type (0 -> LH, 1 -> EC, 2 ->...). May be in between if different waves types are launched into the plasma
0008 %    - equil_fit: interpolated magnetic equilibrium structure
0009 %    - waves: wave structure
0010 %    - mode: 1 -> the first wave is without fluctuations, 0 -> no distinction is done between rays [1,1]
0011 %    - itn: iteration number (optional) [1,1]
0012 %
0013 % OUTPUT:
0014 %    - P0: total rf power (MW) [1,1]
0015 %
0016 % By Yves Peysson (CEA-DRFC, yves.peysson@cea.fr) and Joan Decker (CEA-DRFC, joan.decker@cea.fr)
0017 %
0018 if nargin < 2,
0019     error('Wrong number of input arguments in rayevolution_yp.m.');
0020 end
0021 %
0022 if nargin < 3,
0023     mode = 0;
0024     rho_S = linspace(0,1,20);
0025     itn = -1;
0026 end
0027 %
0028 if nargin < 4,
0029     rho_S = linspace(0,1,20);
0030     itn = -1;
0031 end
0032 %
0033 if nargin < 5,
0034     itn = -1; 
0035 end
0036 %
0037 Rp = equil_fit.Rp;
0038 %
0039 nw = length(waves);
0040 nr = ones(1,nw);
0041 %
0042 for iw = 1:nw,nr(iw) = length(waves{iw}.rays);end
0043 %
0044 nparbin = [-10:0.1:10];
0045 hr = zeros(length(rho_S)-1,length(nparbin));
0046 hr_nofluct = zeros(length(rho_S)-1,length(nparbin));
0047 rho = (rho_S(1:end-1) + rho_S(2:end))/2;
0048 P0 = 0;
0049 %
0050 %
0051 if itn > 0,
0052     figure('Name',['Ray dynamics, itn = ',int2str(itn)]),
0053     title_fig = ['Ray evolution, iteration = ',int2str(itn)];
0054 else
0055     figure('Name','Ray dynamics'),
0056     title_fig = 'Ray evolution';
0057 end
0058 %
0059 for iw = 1:nw,
0060     for ir = 1:nr(iw),
0061         if iw == 1 & mode == 1,
0062             graph1D_jd(waves{iw}.rays{ir}.ss,real(waves{iw}.rays{ir}.sNpar),0,0,'Ray length','N_{||}',[],NaN,'',[nparbin(1),nparbin(end)],'-','none','g',2);
0063             if ~wavetype == 1,
0064                 graph1D_jd(waves{iw}.rays{ir}.ss,sign(real(waves{iw}.rays{ir}.sNpar(1)))*6.5./sqrt(waves{iw}.rays{ir}.sTe),0,0,'Ray length (m)','N_{||}',title_fig,NaN,'',[nparbin(1),nparbin(end)],'-','none','r',2);        
0065             end
0066         else
0067             graph1D_jd(waves{iw}.rays{ir}.ss,real(waves{iw}.rays{ir}.sNpar),0,0,'Ray length','N_{||}',[],NaN,'',[nparbin(1),nparbin(end)],'-','none','b',0.5);
0068             if ~wavetype == 1,
0069                 graph1D_jd(waves{iw}.rays{ir}.ss,sign(real(waves{iw}.rays{ir}.sNpar(1)))*6.5./sqrt(waves{iw}.rays{ir}.sTe),0,0,'Ray length (m)','N_{||}',title_fig,NaN,'',[nparbin(1),nparbin(end)],'-','none','r',0.5);                   
0070             end
0071             %
0072             P0 = P0 + waves{iw}.rays{ir}.P0_2piRp*2*pi*Rp/1e6;
0073         end
0074         %
0075         for irho = 1:length(rho_S)-1,
0076             if iw == 1 & mode == 1,%the first ray without fluctuations
0077                 ii = find((0 <= waves{iw}.rays{ir}.srho) & (waves{iw}.rays{ir}.srho < rho_S(irho + 1)));
0078                 hr_nofluct(irho,:) = hr_nofluct(irho,:) + hist(real(waves{iw}.rays{ir}.sNpar(ii)),nparbin);
0079             else
0080                 ii = find((0 <= waves{iw}.rays{ir}.srho) & (waves{iw}.rays{ir}.srho < rho_S(irho + 1)));
0081                 hr(irho,:) = hr(irho,:) + hist(real(waves{iw}.rays{ir}.sNpar(ii)),nparbin);                
0082             end
0083         end
0084     end
0085 end
0086 %
0087 if mode == 1,
0088     for ir = 1:nr(1),
0089         graph1D_jd(waves{1}.rays{ir}.ss,real(waves{1}.rays{ir}.sNpar),0,0,'Ray length','N_{||}',[],NaN,'',[nparbin(1),nparbin(end)],'-','none','g',2);
0090         %
0091         if ~wavetype == 1,
0092             graph1D_jd(waves{1}.rays{ir}.ss,sign(real(waves{1}.rays{ir}.sNpar(1)))*6.5./sqrt(waves{1}.rays{ir}.sTe),0,0,'Ray length (m)','N_{||}',title_fig,NaN,'',[nparbin(1),nparbin(end)],'-','none','r',2);        
0093         end
0094     end
0095     %
0096     if ~wavetype == 1,
0097         legend('no fluct.','6.5/\surdT_{e}','with fluct.');
0098     else
0099         legend('with fluct.');       
0100     end
0101 else
0102     if ~wavetype == 1,   
0103         legend('with fluct.','6.5/\surdT_{e}');
0104     else
0105         legend('with fluct.');
0106     end
0107 end
0108 %
0109 if ~wavetype == 1,
0110     if itn > 0,
0111         figure('Name',['Power density, itn = ',int2str(itn)]);
0112     else
0113         figure('Name','Power density');
0114     end
0115     %
0116     for irho = 1:length(rho),
0117         %
0118         if irho <= ceil(-1+sqrt(1+length(rho)))*(ceil(-1+sqrt(1+length(rho)))+1);
0119             subplot(ceil(-1+sqrt(1+length(rho))),ceil(-1+sqrt(1+length(rho)))+1,irho),
0120             %
0121             if rem(irho,ceil(-1+sqrt(1+length(rho)))+1) == 1,
0122                 label_y = 'dP_{rf}/dN_{||}';
0123             else
0124                 label_y = '';
0125             end
0126             %
0127             if mode == 1,
0128                 graph1D_jd(nparbin,hr_nofluct(irho,:)/max(hr_nofluct(irho,:)),0,0,'N_{||}',label_y,['\rho = ',num2str(rho(irho))],NaN,[nparbin(1),nparbin(end)],[0,1],'-','none','g',1);     
0129             end
0130             %
0131             graph1D_jd(nparbin,hr(irho,:)/max(hr(irho,:)),0,0,'N_{||}',label_y,['\rho = ',num2str(rho(irho))],NaN,[nparbin(1),nparbin(end)],[0,1],'-','none','b',1);
0132             %
0133             Npar_th = 6.5./sqrt(ppval(equil_fit.Te_fit.pp_f,rho(irho)));
0134             %
0135             graph1D_jd([Npar_th,Npar_th],[0,1],0,0,'N_{||}',label_y,['\rho = ',num2str(rho(irho))],NaN,[-Npar_th,Npar_th],[0,1],'--','none','r',2);
0136             graph1D_jd([-Npar_th,-Npar_th],[0,1],0,0,'N_{||}',label_y,['\rho = ',num2str(rho(irho))],NaN,[-Npar_th,Npar_th],[0,1],'--','none','r',2);
0137             %
0138             for iw = 1:nw,
0139                 for ir = 1:nr(iw),
0140                     graph1D_jd([real(waves{iw}.rays{ir}.sNpar(1)),real(waves{iw}.rays{ir}.sNpar(1))],[0,1],0,0,'N_{||}',label_y,['\rho = ',num2str(rho(irho))],NaN,[-Npar_th,Npar_th],[0,1],'-','none','m',1);
0141                 end
0142             end
0143         end
0144     end
0145 end
0146 %
0147 
0148
Community support and wiki are available on Redmine. Last update: 18-Apr-2019.