bripfield_circ3_yp

    Calculation of the toroidal magnetic field ripple by the Yushmanov
    model (divB = 0 to the first order of the magnetic ripple perturbation)
   Reviews of Plasma Physics, p.117-241, "Diffusive Transport Processes Caused by Ripple in Tokamaks" by P.N. Yushmanov, 1990.

    Input:

       - R: Major radius of the point at which magnetic field ripple components are calculated (m) [nrho,mtheta]
       - Z: Vertical position of the point at which magnetic field ripple components are calculated (m) [nrho,mtheta]
        - N: Number of toroidal magnetic field coils [1,1]
       - B0: Vacuum toroidal magnetic field on the toroidal axis, (signed value as BPHI) (T) [1,1]
       - R0: Toroidal major radius of the toroidal magnetic field windings(m) [1,1]
       - rn: Inner radius of the toroidal magnetic field windings (m) [1,1]
       - rx: Outer radius of the toroidal magnetic field windings (m) [1,1]

    Output:

       - aBx: amplitude of the horizontal magnetic ripple field along the major radius direction(T) [nrho,ntheta]
       - aBy: amplitude of the vertical magnetic ripple field along the vertical radius direction (T) [nrho,ntheta]
       - aBz: amplitude of the toroidal magnetic ripple field along the toroidal direction (T) [nrho,ntheta] 
       - deltar: radial component of the magnetic ripple perturbation [nrho,ntheta] 
       - deltat: poloidal component of the magnetic ripple perturbation [nrho,ntheta]
       - deltaz: toroidal component of the magnetic ripple perturbation [nrho,ntheta]

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

0001 function [aBx,aBy,aBz,deltar,deltat,deltaz] = bripfield_circ3_yp(R,Z,N,B0,R0,rn,rx)
0002 %
0003 %    Calculation of the toroidal magnetic field ripple by the Yushmanov
0004 %    model (divB = 0 to the first order of the magnetic ripple perturbation)
0005 %   Reviews of Plasma Physics, p.117-241, "Diffusive Transport Processes Caused by Ripple in Tokamaks" by P.N. Yushmanov, 1990.
0006 %
0007 %    Input:
0008 %
0009 %       - R: Major radius of the point at which magnetic field ripple components are calculated (m) [nrho,mtheta]
0010 %       - Z: Vertical position of the point at which magnetic field ripple components are calculated (m) [nrho,mtheta]
0011 %        - N: Number of toroidal magnetic field coils [1,1]
0012 %       - B0: Vacuum toroidal magnetic field on the toroidal axis, (signed value as BPHI) (T) [1,1]
0013 %       - R0: Toroidal major radius of the toroidal magnetic field windings(m) [1,1]
0014 %       - rn: Inner radius of the toroidal magnetic field windings (m) [1,1]
0015 %       - rx: Outer radius of the toroidal magnetic field windings (m) [1,1]
0016 %
0017 %    Output:
0018 %
0019 %       - aBx: amplitude of the horizontal magnetic ripple field along the major radius direction(T) [nrho,ntheta]
0020 %       - aBy: amplitude of the vertical magnetic ripple field along the vertical radius direction (T) [nrho,ntheta]
0021 %       - aBz: amplitude of the toroidal magnetic ripple field along the toroidal direction (T) [nrho,ntheta]
0022 %       - deltar: radial component of the magnetic ripple perturbation [nrho,ntheta]
0023 %       - deltat: poloidal component of the magnetic ripple perturbation [nrho,ntheta]
0024 %       - deltaz: toroidal component of the magnetic ripple perturbation [nrho,ntheta]
0025 %
0026 % By Yves Peysson (CEA-DRFC, yves.peysson@cea.fr) and Joan Decker (CEA-DRFC, joan.decker@cea.fr)
0027 %
0028 x = R - R0;
0029 y = Z;
0030 %
0031 r = sqrt(x.^2 + y.^2);
0032 x(x == 0) = - eps;
0033 theta = atan(y./x) + pi*(x < 0) + 2*pi*(x > 0 & y < 0);
0034 %
0035 L0 = (rn + rx).^2./8./R0;%Radial shift of nested iso-ripple surfaces (m)
0036 a = sqrt((r.*r + L0.*L0 + 2.0.*r.*L0.*cos(theta)).*(R0 - L0)./R);%iso-ripple minor radius (m)
0037 L = L0 - 0.5*a.*a./(R0 - L0);
0038 aL = sqrt((r.*r + L.*L + 2.0.*r.*L.*cos(theta)).*(R0 - L0)./R);%iso-ripple minor radius (m)
0039 %
0040 deltan = 1.0./besseli(0,0.5*N.*(rn + rx)./(R0 - L0));
0041 %
0042 b = a.*N./(R0 - L0);
0043 %
0044 deltar = deltan.*besseli(1,b).*(r + L.*cos(theta))./aL;
0045 deltat = deltan.*besseli(1,b).*L.*sin(theta)./aL;
0046 deltaz = deltan.*besseli(0,b);
0047 %
0048 aBx = B0.*R0.*(deltar.*cos(theta) - deltat.*sin(theta))./R;%Major radius direction
0049 aBy = B0.*R0.*(deltar.*sin(theta) + deltat.*cos(theta))./R;%Vertical radius direction
0050 aBz = B0.*R0.*deltaz./R;%Toroidal direction
0051 %
0052