!
!This file contains namelist input for the GENRAY ray tracing code.
!GENRAY will read it as a file called either genray.in or genray.in,
!in that order.   However, genray.in is presently reserved for an
!upgrade in format, so it is best to use the genray.dat name.
!

!Normalization constants:
!-------------------------------------------------------------------------
!/genr/ namelist
!-------------------------------------------------------------------------
! mnemonic, is the run designator...to help keep track of runs.
!           It is used for naming the ray data output files:
!              mnemonic.txt and mnemonic.nc 
!              (Ray data o/p also depends on rayop nmlst variable.)
!           mnemonic is character*128, default="genray"
! rayop,    Specifies which of mnemonic.txt and mnemonic.nc files
!           are to be output:
!           "both", "text", "netcdf", or "none".  
!           rayop is character*8, default="both".
!           [Previous (related) iout3d nml no longer has functionality.]
! dielectric_op="enabled", adds output of the 9 complex dielectric tensor
!               elements to  .nc ray data file.
!               "disabled", omit such data from the .nc data file.
!               dielectric_op is character*8, default="disabled"
!-------------------------------------------------------------------------
!Normalization constants:
! r0x (m) characteristic length, can be used as scale factor
! b0 (tl) characteristiic magnetic field, can be used as scale factor
!-------------------------------------------------------------------------
!Parameters for output files
!--------------------------------------------------------------------------
! outdat*20     name of output file
! stat*3        status of output file
!--------------------------------------------------------------------------
! partner  = 'disabled' to use input profiles from genray.dat or genray.in
!          = 'transport'  to use plasma profile data from the file: 
!                      genray_profs_in  created by transport code;
!                  and to output power and current profiles on transport
!                      code radial grid to file: genray_profs_out
!------------------------------------------------------------------------
 &genr
 mnemonic="genray"
 rayop="both"
 dielectric_op="enabled"
 r0x=1.0d0
 b0=1.0d0
 outdat='zrn.dat'
 stat='new'
 partner="disabled"
 &end
!--------------------------------------------------------------------------
!/tokamak/
!-------------------------------------------------------------------------
!Tokamak
!--------------------------------------------------------------------------
! eqdskin=Name (character*128) of eqdsk equilibrium input file
!         "equilib.dat" (default)
!--------------------------------------------------------------------------
! Type of the (normalized) radial coordinates
! indexrho [1 - sqrt(area), 2 - sqrt(torflux), 3 - sqrt(volume), 
!           4 - sqrt(psi-psilim), 5 - (psi-psilim)],
!           6 - (r_max(psi)-rmin(psi))/(r_max(psi_lim)-rmin(psi_lim))
!--------------------------------------------------------------------------
! ipsi=1 calculation of contours psi(z,r)=const
!     =0 -read in these contours from psi.bin file
!--------------------------------------------------------------------------
! ionetwo=1-calculation power and current radial
!           profiles, to the file onetwo.bin
!           0 - no calculations)
!--------------------------------------------------------------------------
! ieffic  choice of formula for the current drive efficiency
!        =1 asymptotic simple formula (homogeneous, nonrelativistic)
!        =2 asymptotic formula (East-Karney )
!        =3 asymptotic formula (curba subroutine)
!--------------------------------------------------------------------------
! psifactr (it should be 0 < psifactr =<1,  psifcatr ~1)
!         is the parameter for the creation of the limiter points
!         using the closed flux surface:  psi(r,z)=psisep*psifactr 
!         psifactr is a parameter (it must be .le.1) to avoid
!         problems with the psi function near the separatrix.
!--------------------------------------------------------------------------
!deltripl is the relative amplitude of the ripple field at the
!         last flux surface (at rho=1)
!
!nloop    is number of toroidal field coils 
!
!i_ripple is the index to choose the ripple model
!         bripl_phi(z,r,phi)=(dF/dphi)/r=cos(N_loop*phi)*g(r,z)*N_loop/r
!         bripl_z(z,r,phi) =(dF/dz)    =sin(N_loop*phi)*(dg/dz)
!         bripl_r(z,r,phi) =(dF/dr)    =sin(N_loop*phi)*(dg/dr)
!         models for function g:	 
!         =1 the ripple model approximating the DIII-D field
!            g=beqd*reqd*deltripl*(r/rmax)**N_loop/N_loop 
!            beqd is the toroidal magnetic field at reqd (Tl)
!            reqd is the nominal major radius of the torus.
!            rmax is the max major radius at the last closed flux surface
!            r    is the major radius
!         =2 the ripple model using modified Bessel function I_0
!            g=beqd*reqd*deltripl*I_0(N_loop*rho(z,r))/(N_loop*I_(N_loop))
!--------------------------------------------------------------------------
! eqdskin=Name (character*128) of eqdsk equilibrium input file
!         "equilib.dat" (default)
!--------------------------------------------------------------------------
! NR is the number of bin boundaries in the small radius direction
!    for the calculation of the power and current drive radial profiles.
!    Power and current is tabulated at (NR-1) bin centers.
!--------------------------------------------------------------------------
 &tokamak
 eqdskin="equilib.dat"
 indexrho=2
 ipsi=1
 ionetwo=1
 ieffic=3
 psifactr=0.85d0
 deltripl=0.00d0
 nloop=24
 i_ripple=1
 NR=101
 &end
!/wave/
!-------------------------------------------------------------------------
!Waves
!-------------------------------------------------------------------------
! frqncy frequency f=w/2pi in GHz
! ioxm ( 1 - om, -1  - xm )  wave mode
! ireflm-max number of reflections =1 for EC
! jwave  (0 - LH wave, -1 AW, 1 - EC wave) wave harmonic
! -------------------------------------------------
! istart  if start point outside the plasma=1 else=2
! if istart=1 use namelist &eccone below, =2 use &grill
! if istart=3 it use &grill and the additional calculations in dinit
! to launch the ECR ray inside the plasma in the O_X mode
! conversion point (rhoconv,theta), Theta is a poloidal angle (degree)
! for mode conversion point. It is given in dinit.f 
!--------------------------------------------------------------------------
! delpwrmn - Minimum power in each ray, as a fraction of
!            starting power in the ray, after which ray is stopped.
!--------------------------------------------------------------------------
! ibw=0 it is not the direct launch of the Bernstein waves
!    =1 the direct launch of electron Bernstein wave from dhot tensor
!       The last case works only for istart=2 and grill_lh conditions 
!--------------------------------------------------------------------------
! i_vgr_ini =+1 the wave is directed into the plasma (in the initial point)
!           =-1 the wave is directed out the plasma (in the initial point
!---------------------------------------------------------------------------
! poldist_mx is the maximal poloidal distance (m) along the ray
!            default=1.d+5
!----------------------------------------------------------------------------
 &wave
 frqncy=45.0d+0
 ioxm=+1
 ireflm=1
 jwave=+1
 istart=1
 delpwrmn=1.d-16
 ibw=0
 i_vgr_ini=+1
 poldist_mx=1.d+5  
 &end
---------------------------------------------------------------------

!/scatnper/
!-------------------------------------------------------------------------
!N_perpendicular scattering
!-------------------------------------------------------------------------
! iscat it is the switch for the n_perp scattering
!       iscat=1 the scattering switched on,
!            =0 the scattering switched off
!-------------------------------------------------
! rhoscat(1:nscat_n) small radii for the scattering location
!        The parameter nscat_n should be given in the param.i file
!
! The scattering of the polar angle deltheta will be
!      deltheta=dsqrt(2.d0*scatd)*ranorm(fseed)
! scatd(0) the mean square scattering angle (radians**2)
!          for the plasma boundary reflection points
!
! scatd(1:nscat_n) the mean square scattering angles (radians**2)
!          for the interior plasma boundary points
!-------------------------------------------------------------------
 &scatnper
 iscat=0
 scatd= 0.01, 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.01, 0.01
 rhoscat=     0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.95, 0.97
 &end

!/dispers/
!Dispersion relation
!-------------------------------------------------------------------------
! ib<=nbulk cyclotron resonance sort(=1 for ecr)
! the number in (1-y(ib)) for the multiplication of the
! dispersion relation to delete the singularity
! -------------------------------------------------
! id gives form of the dispersion relation
!            =1 AN**4+BN**2+C=0
!            =2 N**2=(-B+ioxm*Sqrt(B**2-4AC))/2A;
!	     =3 Appleton-Hartree;
!            =4 electron relativistic plasma from Hermitian Mazzucato code
!            =5 electron relativistic plasma from total Mazzucato code
!            =6 hot non-relativistic plasma, Hermitian, Forest calc 
!            =7 electron relativistic plasma from Shkarofsky code
!            =8 Ono dispersion for fast waves
!            =9 hot non-relativistic plasma, full tensor [Bernstein-Friedland]
!            =10 Westerhof-Tokman dispersion with Mazzucato
!                relativistic electron dielectric tensor 
!            =11 Eric Nelson-Melby relativistic tensor
!                Dispersion function = Re(Det)
!            =12 Westerhof-Tokman dispersion with 
!                    Ram Abhay dielectric tensor
!            =13 Westerhof-Tokman dispersion with
!                 hot non-relativistic plasma, full tensor
!            =14 Abhay Ram's dielectric tensor, using the Trubnikov integral
!              **** NOTE: irkmeth=1 should be used for id=14
!---------------------------------------------------------
! For use with Abhay Ram's dispersion relation (id=14) parameters to
! control the integration routine for the Trubnikov integral.
!
! relres offers 4 choices for the resolution to use for the relativistic
! dispersion relation using the Trubnikov integral:
! relres=1 low resolution, errabs0=1.d-4, errrel0=1.d-4, navg=3, diff_err=0.1
!       =2 medium res., errabs0=1.d-5, errrel0=1.d-5, navg=12, diff_err=1.d-3
!       =3 high res., errabs0=1.d-6, errrel0=1.d-6, navg=25, diff_err=1.d-6
!       =4 user-defined res., set the following parameters manually.
! default: relres=2
!
! The Trubnikov one-dimensional (complex) integral is performed by splitting
! up a region from 0 to 1.d8 into 10^6 pieces, and each piece is integrated
! using the SLATEC adaptive quadrature routine dqag. errabs0 and errrel0 are
! the absolute and relative error tolerances passed directly to dqaq.
! Then the adjacent pieces are compared (it is an oscillatory integrand)
! and using navg number of pieces, when the average difference between them
! are less then diff_err, the integration is presumed finished (Thus it may
! finish long before the upper limit of 1.d8).
!
! errabs0 - absolute error for dqag integration routine
! errrel0 - relative error for dqag integration routine
! navg - number of adjacent integration intervals to use in comparison
! diff_err - error tolerance using navg pieces, when the average difference
!         is less than diff_err, then the integration is done.
!
! To decide when one should use the low, medium, or high resolution 
! integration, here are some suggestions based on the behavior of the
! Trubnikov integrand: The integrand converges more slowly, and hence
! the resolutions should be set higher, for low electron temperature,
! low (i.e. near zero) magnitude of n_parallel, and for low (near or
! below the fundamental cyclotron frequency) frequency. 
! Examples: n_parallel = -0.05, Te=400 eV, omega/omega_ce=0.4 to 1.2, 
!    it was necessary to use errabs0=1.d-5,errrel0=1.d-5,navg=20,diff_err=1.d-5
!    to be completely converged.  By changing Te to 4000 eV, it was sufficient
!    to use 1.d-4,1.d-4,15 and 1.d-4. 
!   An easy case: n_parallel=0.3, Te=7 keV, omega/omega_ce=2.4 to 2.7, 
!    complete convergence already at errabs0=1.d-4,errrel0=1.d-4,navg=2,
!    diff_err=0.5
!   An intermediate case: n_parallel=0.1, omega/omega_ce=1.0, Te=300 eV
!    errabs0=1.d-5,errrel0=1.d-5,navg=12,diff_err=1.d-3 was OK.
!------------------------------------------------
! For Mazzucato plasma dispersion tensor:
! iherm =1 hermitian dielectric tensor, 2-full
!------------------------------------------------------------------------
!Absorption:
!iabsorp -choice of Imag(N_perp) (N_perp is the perpendicular refractive index)
!-------------------------------------------------------------------------
! iabsorp=1 for EC waves with Mazzucato solver (ki/kr<<1)
!        =2 for LH waves
!        =3 for FW waves, Chiu et al themal corr., NF 1989, with corrections.
!        =4 for all frequencies with Forest code (ki/kr<<1)
!        =5 for EC and EBW waves from Shkarofsky code
!        =6 for EC and BW anti-hermitian part relativistic tensor+
!                         hermitian_part (Forest code)
!        =7 for EC wave case. The complex electric field calculations
!           using Cold plasma tensor + antihermitian relativistic tensor
!           ---EC relativistic absorption 
!           dielectric tensor=hermitian part(cold plasma)+
!                             anti-hermitian part(full relativistic) 
!        =8 for id=10,12,13 and Westerhof-Tokman dispersion 
!           function  [Projection method from real part of
!           eigenvalue.]  Integration is with respect to distance.
!        =9 The absorption is calculated for hot dispersion
!           using the formula from Stix book p.74 (17,18,21)
!           Im(k_perp)= 0.5*Power_abs/(P^+T^)
!           It uses hot plasma dielectric tensor.
!           It calculates hot plasma dielectric tensor reps() and 
!           electric field polarization (cex,cey,cez) using this 
!           hot plasma tensor.
!        =10 The absorption is calculated for relativistic tensor
!           (A.Ram)
!           using the formula from Stix book p.74 (17,18,21)
!           Im(k_perp)= 0.5*Power_abs/(P^+T^)
!           It uses relativistic dielectric tensor.
!           It calculates relativistic dielectric tensor reps() and 
!           electric field polarization (cex,cey,cez) using this 
!           tensor.
!        =11 The absorption Im(N_perpendicular) is calculated 
!           for the relativistic dispersion function Determinant=0.
!           The dispersion function uses the full Abhay Ram dielectric
!           tensor (with the Trubnikov integral).
!           The projection method is used to find Im(N_perpendicular).
!           That is, Im(N_perp)=-Im(D)/(dD/dRe(N_perp)), evaluated 
!           at give Real(N_perp) with Im(N_perp)=0.
!           The real part Real(n_perp) is obtained from the trajectory
!           solution for the given id value. 
!           Then complex n_perp is used with full dielectric tensor
!           to get the polarizations.  Good for arbitrary id value,
!           like other iabsorp calcs.
!------------------------------------------------------------------------
! The change of the dispersion relation and absorption
! near the gyro-frequency points
!-------------------------------------------------
! iswitch=1   To use the change of the dispersion relation and
!             absorption
!        =0   Do not use the change of the dispersion relation and
!             absorption 
!     del_y   If the difference |1-nY(jy)|<del_y 
!             (jy=1-nbulk ,n=...-2,-1,0, 1,2,3,...)
!             then switch on the new 
!             given type of the dispersion and absorption.
!   jy_d      is the type of plasma species 1<=jy<=nbulk
!   idswitch  is the type of the dispersion function near the 
!             gyro-frequency points
!             It can be equal 1,2,3,4,5,6
!   iabswitch is the type of the absorption near the gyro-frequency point  
!-----------------------------------------------------------------------
!   n_relt_harm1 is the lowest, i.e., minimum harmonic used in the
!               anti-hermitian dielectric tensor calculations.
!               It man be positive of negative.
!               Default value is +9999, in which case this input
!               is ignored.
!   n_relt_harm (.ge.1) gives the number of EC harmonics used 
!               in anti-hermitian dielectric tensor calculations
!               If n_relt_harm1=9999, then harmonics from 
!                  -n_relt_harm to +n_relt_harm are used.
!               If n_relt_harm1.ne.9999, then harmonics from
!                  n_relt_harm1 to n_relt_harm1+n_relt_harm are used.
!   It is necessary that the harmonics used in this calculation
!     be within the range of parameters [n_relt_harm1a,n_relt_harm2a]
!     set in the param.i file.s
!     These conditions are checked in the code.
!   n_relt_intgr is the number of points for integration along the
!     resonance curve.
!---------------------------------------------------------------------
!  flux=B~.B+E~.d(omega*eps_herm)/(domega).E
!   iflux=1 the flux will be calculated using the the group velocity from
!           the chosen dispersion relation (with given id) and the electric
!           field calculated for the chosen iabsorp
!   iflux=2 the flux will be calculated using V_gr for the electron cold plasma
!           dispersion and polarization (using subroutine  grpde2)  
!-----------------------------------------------------------------------
! i_im_nperp choice of the method to find Im_N_perp 
!    for hot plasma(iabsorp=4):
! i_im_nperp=1 Im_N_perp=abs(ImD_full/(dD_hermitian/dReN_perp)) 
! i_im_nperp=2 (Re_N_perp,Im_N_perp) is the complex root 
!              (of the complex dispersion relation)
!              calculated by Newton iterations with the numerical
!              derivatives (the chord method)
!------------------------------------------------------------------
! i_geom_optic sets  the form of the ray equations
!              =1  integration in time (default):
!                  ray-tracing equations right hand side=
!		   dr^/dt=-(dD/dN^)/(dD/domega)
!                  dN^/dt=(dD/dr^)/(dD/domega)
!                  In this case rside1 gives v_group
!              =2  integration is space,
!                  ray-tracing equations right hand side=
!		   dr^/dl=- ray_direction * (dD/dN^)p
!                  dN^/dl=  ray_direction * (dD/dr^)p
!                  p=1.d0/dsqrt(deru(1)**2+deru(2)**2+(r*deru(3))**2)=
!                  deru(1)=dD/dN_z,deru(2)=dD/dN_r,deru(3)=dD/dCM,
!                  N_phi=cm/r
!----------------------------------------------------------------------
! ray_direction =+1 as default 
!            or =-1
! It is a multiplier in right hand side of ray-tracing equations
! It is used for i_geom_optic=2 case
!----------------------------------------------------------------------

 &dispers
 ib=1
 id=14
 iherm=1
 iabsorp=4
 iswitch=0
 del_y=1.d-3
 jy_d=1
 idswitch=2
 iabswitch=1
 n_relt_harm=5
 n_relt_intgr=50
 iflux=1
 i_im_nperp=1
 i_geom_optic=1
 ray_direction=-1.d0   !Only for i_geom_optic=2
 &end

!/numercl/
!------------------------------------------------------------------------
!Numerical method
!-------------------------------------------------------------------------
! irkmeth (0-constant 1-variable step 5 order,2- variable step in RK 4 order)
!  irkmeth=0: Poloidal distance of output is at intervals .ge.prmt6.
!             Checks time step for passing outside plasma and reflects.
!  irkmeth=1: Only poloidal distance for control of output point (prmt6,
!             i_output has no effect). Output at distance.ge.prmt6,
!             i.e, the first code step beyond prmt6 distance.
!             No control for being outside the plasma and reducing
!             the step.  Correction method specified by icorrect is
!             operative.
!             Time or length for integration according to i_geom_optic.
!  irkmeth=2: Most developed method of ray equation integration.
!             Time or space step in the of the equations
!             (according to setting of i_geom_optic) is controlled
!             so that output is at intervals prmt6 (meters).
!             As ray approaches the plasma edge, it is reflected
!             at the last closed flux surface.
!              
! ndim1 (number of the ray tracing equations)
! isolv=1 correction,=2 expl.solution
! idif=1 analytic differentiation, =2 numerical
! -------------------------------------------------
! nrelt   Maximum number of ray elements per ray.
!         Must be .le. nrelta (a parameter)
!--------------------------------------------------------------------------
! -------------------------------------------------
! Runge-Kutta method parameters
! -------------------------------------------------
! prmt1=tau initial= prmt(1)
! prmt2=tau final=prmt(2)
! prmt3=initial tau step=prmt(3)
! prmt4=required accuracy=prmt(4)
! prmt6=hprint=prmt(6) poloidal or total distance(meters) for results output
! prmt9=ihlf=prmt(9)
!--------------------------------------------------------------------------
!icorrect= switch for Hamiltonian correction in subroutine outpt
!          [See manual].
!icorrect=0 switch off the correction
!        =1 switch on the correction
!-------------------------------------------------------------------------
! iout3d [obsolete]   ='enable'  to write output 3d.dat file
!           ='disable' do not write 3d.dat file
!--------------------------------------------------------------------------
! maxsteps_rk the maximal number of the time steps of the Runge-Kutta
!             solver (in default =10000)
!--------------------------------------------------------------------------
! i_output is used for irkmeth=2 only
! i_output=1 output is at the equal poloidal distance prmt6
!         =2 output is at the equal total distance prmt6
!--------------------------------------------------------------------------
!
! The following has been used for OXB in cases where the UH layer
!   is very close to the plasma boundary.   Then, in the vicinity of
!   the UH layer, switch the step size along the ray to a shorter
!   value.
! i_uh_switch=1    if uh=dsrt(xe+ye**2) < uh_switch then change 
!                  the output step prmt(6) to prmt6_uh_switch 
!            =0    do not change the output step prmt(6)
!
! prmt6_uh_switch  is the output step for i_uh_switch=1 case
!
! uh_switch       if uh < uh_switch then change the outout step for
!                  i_uh_switch=1 case
!----------------------------------------------------------------------
 &numercl
 irkmeth=2
 ndim1=6
 isolv=1
 idif=2
 nrelt=500
 prmt1=0.000d+00
 prmt2=9.999d+05
 prmt3=2.000d-2
 prmt4=1.000d-3
 prmt6=0.500d-3
 icorrect=0
 iout3d='enable'
 maxsteps_rk=50000
 i_uh_switch=0
 uh_switch=1.5d0 
 prmt6_uh_switch=1.d-5
 &end

!/output/
!----------------------------------------------------------------------
! iwcntr =1 genray.f will calculate the contours wb_c=n
!        =0 genray.f will not do it
! iwopen =1 mk_grapc will calculate open contours wb_c=n  (using contrb1)
!         2 mk_grapc will calculate close contours wb_c=n (using contrb2)
! iwj    =mk_grapc will calculate contours wb_cj=n, j is a kind of the plasma
!         component must be.le.nbulk, j=1 for the electron gyrofrequency
!         j.ge.2  for the ion (j kind) gyrofrequency
! itools =0 do not use mkgrtool
!        =1 to use mkgrtool
!----------------------------------------------------------------------
!
! i_plot_b =1 create figures for the magnetic field,density and temperature 
!             profiles in plot.ps file using subroutine map_b based on PGplot
! i_plot_b =0 do not write the b,n,T figures to plot.ps file 
!
! i_plot_d =1 create the dispersion function contours d(ReNperp,ImN_perp)
!             in plot.ps file using PGplot
! i_plot_d =0 do not create the dispersion function contours
!----------------------------------------------------------------------
 &output
 iwcntr=1
 iwopen=1
 iwj=1
 itools=0
 i_plot_b=0
 i_plot_d=0
 &end

!/plasma/
!-------------------------------------------------------------------------
!Plasma parameters
!-------------------------------------------------------------------------
! nbulk>=1 quantity of plasma components
!----------------------------------------------------
! izeff =0 zeff will be calculated using the given ions;
!          electron density will be calculated using ions;
!       =1 zeff and electron density are given, 
!          the ion components will be calculated;
!       =2 zeff, electron and ion (if nbulk>1) densities are given,
!          and zeff is not recalculated from the plasma components;
!       =3 Use eqdsk pres (pressure). Let temperature T_E=T_i
!          pres=dens1(k,1)*temp1(k,1)+
!          Sum(i=2,nbulk)(dens1(k,i)*temp1(k,i)
!          In this case we will calculate Zeff(rho),
!          dens_electron(rho) and T_e(rho)=T_i(rho)
!       =4 Use eqdsk pres (pressure), the given temperature
!          profiles T_i(rho) (i=1,nbulk) and the given Z_eff(rho).
!          nbulk should be .ge. 3
!          pres=dens1(k,1)*temp1(k,1)+
!          Sum(i=2,nbulk)(dens1(k,i)*temp1(k,i)
!          In this case we will calculate dense(1)(rho),
!          dense(nbulk)(rho) and dense(nbulk-1)(rho)
! -----------------------------------------------------
! idens (0 - analytic, 1 - spline) representation of
! the density, temperature and zeff radial profiles
! -----------------------------------------------------
!   temp_scale(ncompa),den_scale(ncompa) are the parameters to multiply
!   the given temperature and density profiles
! -----------------------------------------------------
! ndens is the number of points for the input radial density and 
!       temperature profiles
! -----------------------------------------------------
 &plasma
 ndens=21
 nbulk=1
 izeff=2
 idens=1
 temp_scale(1)=1.d0
 den_scale(1)=1.d0
 &end

!/species/
! plasma component charges charge(i)=mod(charge(i)/charge_electron)
! -----------------------------------------------------
! charge(1) =1 electrons
! charge(i) i=1,nbulk   charge(i+1) must be ge.charge(i)
! charge(i) i=1,nbulk
! -----------------------------------------------------
! plasma component mass dmas(i)=Mass(i)/Mass_electron
! -----------------------------------------------------
! dmas(1) 1 electrons
! dmas(i)   i=1,nbulk
! -----------------------------------------------------
 &species
 charge(1)=1.d0
 charge(2)=1.d0
 charge(3)=6.d0
 dmas(1)=1.d0
 dmas(2)=3674.d0
 dmas(3)=22044.d0 
 &end

!/varden/
! the density variation
! -----------------------------------------------------
!   var0 is an amplitude of the density variation (del_n_0)
!   (see Manual, 3.37...)
!   denm is the number of the poloidal mode in the density variation(l_theta)
!   denn is the number of the toroidal mode in the density variation(l_phi)
!   an   is the radial localization of the variation (rho_0)
!   sigman is the parameter that characterizes the radial thickness
!          of the density fluctuation    
! -----------------------------------------------------
 &varden
 var0=0d0
 denm=1.d0
 denn=15.d0
 an=0.5d0
 sigman=0.1d0 
 &end


!/denprof/
! -----------------------------------------------------
!Analytic radial profiles (idens=0).  Splines (idens=1).
!dense(i)=(dense0(i)-denseb(i))*(1-rho**rn1de(i))**rn2de(i)+denseb(i)
!-----------------------------------------------------
!        if(izeff.eq.0) then
!           zeff will be calculated using the given ions;
!	    nbulk1=nbulk
!	 else
!           =1 zeff is given, the ions component will be calculated
!            if (nbulk.eq.1) nbulk1=1
!            if (nbulk.eq.2) then
!	         nbulk1=2
!	     endif
!            if (nbulk.gt.2) nbulk1=nbulk-2
!	 endif
! -----------------------------------------------------
! dense0(i)   central density in 10**19 m**(-3) i=1,nbulk1
! -----------------------------------------------------
! denseb(i)  edge density in 10**19 m**(-3) i=1,nbulk1
! -----------------------------------------------------
! rn1de(i) i=1,nbulk1
! -----------------------------------------------------
! rn2de(i) i=1,nbulk1
! -----------------------------------------------------
 &denprof
 dense0(1)=3.575d+0
 denseb(1)=1.22d+0
 rn1de(1)=3.2d+0
 rn2de(1)=1.2d+0
 &end

!/tpoprof/
! Ratio tpop=T_perp/T_parallel
! tpop(i)=(tp0(i)-tpb(i))*(1-rho**rn1tp(i))**rn2tp(i)+tpb(i)
! -----------------------------------------------------
! tp0(i) =           central T_perp/T_parallel i=1,nbulk
! -----------------------------------------------------
! tpb(i) =ateb(i)    boundary T_perp/T_parallel i=1,nbulk
! -----------------------------------------------------
! rn1tp(i) i=1,nbulk
! -----------------------------------------------------
! rn2tp(i)  i=1,nbulk
! -----------------------------------------------------

 &tpopprof
 tp0(1)=1.0d0
 tpb(1)=1.0d0
 rn1tp(1)=2.0d0
 rn2tp(1)=1.0d0
 &end

!/vflprof/
! drift velocity || B  
! vflow(i)=(vfl0(i)-vflb(i))*(1-rho**rn1vfl(i))**rn2vfl(i)+vflb(i)
! -----------------------------------------------------
! vfl0(i)     central vflow in m/sec  i=1,nbulk
! -----------------------------------------------------
! vflb(i)     boundary vflow in m/sec i=1,nbulk
! -----------------------------------------------------
! rn1vfl(i) i=1,nbulk
! -----------------------------------------------------
! rn2vf(i)  i=1,nbulk
! -----------------------------------------------------

 &vflprof
 vfl0(1)=0.0d+0
 vflb(1)=0.0d+0
 rn1vfl(1)=2.d0
 rn2vfl(1)=1.0d0
 &end


!/zprof/
! -----------------------------------------------------
! zeff=(zeff0-zeffb)*(1-rho**rn1zeff)**rn2zeff+zeffb
! -----------------------------------------------------
! zeff0   central Z_eff
! zeffb   boundary Z_eff
! rn1zeff zeff=(zeff0-zeffb)*
! rn2zeff      (1-rho**rn1zeff)**rn2zeff+zeffb
!-----------------------------------------------------
 &zprof
 zeff0=1.50d0
 zeffb=1.50d0
 rn1zeff=2.d0
 rn2zeff=1.d0
 &end

!/tprof/
! Average temperature tempe=(T_parallel+2*T_perp)/3
! tempe(i)=(te0(i)-teb(i))*(1-rho**rn1te(i))**rn2te(i)+teb(i)
! -----------------------------------------------------
! te0(i) =at0(i)    central temperature in kev	i=1,nbulk
! -----------------------------------------------------
! teb(i) =ateb(i)    boundary temperature in kev i=1,nbulk
! -----------------------------------------------------
! rn1te(i) i=1,nbulk
! -----------------------------------------------------
! rn2te(i)  i=1,nbulk
! -----------------------------------------------------

 &tprof
 ate0(1)=1.555d0
 ateb(1)=0.063d0
 rn1te(1)=1.4d0
 rn2te(1)=1.6d0
 &end

!/grill/
!------------------LH/EBW-Starting-inside-plasma-----------------------
!  Grill conditions  for istart=2 (start point inside the plasma)
!----------------------------------------------------------------------
! i_n_poloidal =1         The input parameter is N_parallel(from grill).
!  (by default =1)        N_phi,N_theta are calculated from given N_parallel 
!                         N_rho=N_perpendicular(N_parallel) is determined 
!                         from the dispersion relation. It is directed
!                         along +,- gradient(psi) 
!
! i_n_poloidal =2         The input parameters: N_parallel(from grill)
!                         and  n_theta_pol. By default N_theta=0. 
!                         N_perpendicular(N_parallel) is determined 
!                         from the dispersion relation. 
!                         N_phi is calculated from N_parallel and N_theta
!                         N_rho is calculated form N_perpendicular, N_parallel
!                         and N_theta. 
!                         It is directed along +,- gradient(psi)
!
! i_n_poloidal=3          The given variables: N_parallel and the angle
!                         0<ksi_nperp<180 between the vector N_perpendicular 
!                         and gradient(psi). By default ksi_nperp=0.
!                         N_perpendicular(N_parallel) is determined 
!                         from the dispersion relation.
!                         N_phi,N_theta and N_rho are calculated from
!                         N_parallel,N_perpendicular and ksi_nperp.
!
! i_n_poloidal=4          The given variables:N_toroidal and
!                         N_poloidal. This case uses i_vgr_ini set in /waves/
!                         to choose the direction of the small radial N_rho
!                         component. To launch the ray inside the plasma
!                         i_vgr_ini=1 or to the plasma edge i_vgr_ini=-1 
!---------------------------------------------------------------------
! n_theta_pol            The poloidal refractive index component
!                         It is used for i_n_poloidal =2     
!                         By_default n_theta_pol=0.
!----------------------------------------------------------------------
! ksi_nperp               (degrees) the angle 0<ksi_nperp<180
!                         between the vector N_perpendicular 
!                         and gradient(psi). By default ksi_nperp=0.
!---------------------------------------------------------------------
!  Calculation of the small radius value near the plasma edge
!  where LH or FW have cutoff:  
!  i_rho_cutoff=0 (default) no calculations
!              =1 use these calculations
!--------------------------------------------------------------------
!  ngrill  is a number of the poloidal grill angles
!          It is required that ngrill.le.ngrilla, 
!          where ngrilla is parameter in param.i
!----------------------------------------------------------------------
!  igrillpw options specifying N_parallel power spectra
!           =1 power=powers/nnkpar, =2 power=sin**2x/x**2,
!           =3 power=exp-((npar-anmin)/anmax)**2    [default=1]
!----------------------------------------------------------------------
!  igrilltw specifies the form poloidal variation of power,
!                    =1 uniform over height, =2 cos**2 variation.
!----------------------------------------------------------------------
!  rhopsi0(1:ngrill) initial small radius for wave front
!                    (0<rhopsi0<1)
!  rhopsi0(i)=...    i=1,ngrill
!----------------------------------------------------------------------
!  thgrill(1:ngrill) poloidal  angle of grill, measured counter
!                    clockwise from horizontal through the
!                    magnetic axis (degrees).
!  thgrill(i)=...    i=1,ngrill (degree)         [default=0.d0]
!---------------------------------------------------------------------
!  phigrill(1;ngrill) is a toroidal grill angle of grill
!                              (degrees)
!  phigrill(i)=... i=1,ngrill (degree)         [default=0.d0]
!----------------------------------------------------------------------
! height(1:ngrill) is a poloidal length (m) of grill
!                 (giving poloidal power distribution of each grill).
! height(i)=...   i=1,ngrill                  [default=0.2d0]
!----------------------------------------------------------------------
! nthin(1:ngrill) is a number of rays near the each poloidal
!                 center, simulating a grill
! nthin(i)=...    i=1,ngrill       [default: nthin(1)=1]
!----------------------------------------------------------------------
!  anmin(1:ngrill)  position of the left bound
!                   of power spectrum P(n_parallel) (Can be neg).
!  anmin(i)=...     i=1,ngrill
!----------------------------------------------------------------------
!  anmax(1:ngrill)  position of the right bounds
!                   of power spectrum P(n_parallel) (Can be neg).
!  anmax(1)=...     i=1,ngrill
!---------------------------------------------------------------------
!  nnkpar(1:ngrill)  number of points  of power spectrum
!                    P(n_parallel)
!  nnkpar(i)=...     i=1,ngrill
!----------------------------------------------------------------------
!  powers(1:ngrill)  power in one grill (MWatts)
!  (total power of grill(in MWatts) will be powtot=sum{powers}
!  powers(i)=...     i=1,ngrill
!-------------------------------------------------------------------
!  below are for i_n_poloidal=4 case, set (N_toroidal, N_poloidal)
----------------------------------------------------------------------
!  antormin(1:ngrill)  position of the left bound
!                   of power spectrum P(n_toroidal) (Can be neg).
!  antormin(i)=...     i=1,ngrill
!----------------------------------------------------------------------
!  antormax(1:ngrill)  position of the right bounds
!                   of power spectrum P(n_toroidal) (Can be neg).
!  antormax(1)=...     i=1,ngrill
!---------------------------------------------------------------------
!  nnktor(1:ngrill)  number of points  of power spectrum
!                    P(n_toroidal)
!  nnktor(i)=...     i=1,ngrill
!----------------------------------------------------------------------
!  anpolmin(1:ngrill)  position of the left bound
!                   of power spectrum P(n_poloidal) (Can be neg).
!  anpolmin(i)=...     i=1,ngrill
!----------------------------------------------------------------------
!  anpolmax(1:ngrill)  position of the right bounds
!                   of power spectrum P(n_poloidal) (Can be neg).
!  anpolmax(1)=...     i=1,ngrill
!---------------------------------------------------------------------
!  nnkpol(1:ngrill)  number of points  of power spectrum
!                    P(n_poloidal)
!  nnkpol(i)=...     i=1,ngrill
!---------------------------------------------------------------------
!  ilaunch=1, to launch a single ray at r0launch,phi0launch,z0launch
!             in the plasma (meters and degs)
!         =0, no effect (the default)
!  This option is added for comparison with other codes.
!  r0launch is the major radius of the launch point [m]
!  z0launch is the vertical position of the launch pimt [m]
!  phi0launch is the toroidal angle of the launch point [degree] 
!---------------------------------------------------------------------
!       i_grill_pol_mesh: option specifying the poloidal mesh wtheta(j)
!                         near the central grill angle thgrill(i)
!                         =1 equispaced mesh
!                            wtheta(j)-wtheta(j-1)=zdth=Const(default)
!                         =2 poloidal mesh will be chosen to get the equal
!			     power fpwth(j) for all rays near the central
!                            grill angle fpwth(j)=1/nthini
!
!------------------------------------------------------------------------
!       i_grill_npar_ntor_npol_mesh: option specifying the refractive
!                         index meshes.
!
!                         For i_n_poloidal=1,2,3 it specifies the
!                         n_parallel mesh anzin(n) for the power
!                         spectrum pwcpl(n) n=1,...,nnkpari
!                         =1 equispaced mesh
!                            anzin(n)-anzin(n-1)=hnpar=Const (default)
!                         =2 n_parallel mesh will be chosen to get equal
!			     power pwcpl(n) for all rays in the given power
!                            spectrum  pwcpl(n)=1.d0/nnkpari
!                            pwcpl(n)=power_spectrum(anzin(n))*
!                                     delta_npar_bin(n)= 1.d0/nnkpari
!
!                         For i_n_poloidal=4 it specifies two meshes:
!                            a) n_toroidal mesh anztorin(ntor) and
!                            b) n_poloidal mesh anzpolin(npolmesh)
!                            for the power spectrum
!                            pwcpl_tp(1:nnktori,1:nnkpoli)=pwcpl_t*pwcpl_t
!                         =1 equispaced meshs (default)
!                            anztorin(ntor)- anztorin(ntor-1)=hntor=Const
!                            anzpolin(npol)- anzpolin(npol-1)=hnpol=Const
!                         =2 the meshes anztorin(1:nntori) anzpolin(1:nnkpoli)
!                            will be chosen to get the equal
!			     power pwcpl_tp(ntor,npol) for all rays in
!                            the given power spectrum
!                            pwcpl_tp(ntor,npol)=1.d0/(nnktori*nnkpoli)
!
!------------------------------------------------------------------------
 &grill
 i_n_poloidal=1
 i_rho_cutoff=0
 ngrill=1
 igrillpw=1
 rhopsi0(1)=0.95d+00,
! thgrill(1)= 15.0d+0, 
 thgrill(1)=13.8281250d0,
 phigrill(1)=0.0d+0,
 height(1)=0.10d+0,
 nthin(1)=1
 anmin(1)= -0.5d+0
 anmax(1)= -0.6d+0
 nnkpar(1)=1
 powers(1)=1.0
 ilaunch=0
 &end

!/eccone/
!--------------------D3D--EC---------------------------------------------
!     Use equilib.dat= g521022.01000, or equivalent.
!    
!
!     This namelist section specifies ECR cones  for istart=1 
!           (ray cones start outside the plasma).
!
!     Multiple source locations and launch conditions are implemented.
!     ncone=number of source cones. [default=1] [Max is parameter nconea]. 
!
!     For multiple sources (ncone.gt.1), is is necessary to set
!     ncone values for each of the the namelist variables given below:
!     powtot,zst,rst,phist,alfast,betast,alpha1,alpha2(only for raypatt
!     ="genray").  
!     The specifications of number of rays per cone do not not vary 
!     from cone to cone (i.e., na1,na2,gzone,mray(*), cr(*)) do not 
!     vary with cone number.
!    
!     powtot= total power from antenna(MW)
!
!     Two systems for specification of the ECR cone are provided,
!       chosen by raypatt:
!
!     raypatt='genray',  specify ray pattern per following
!                        genray method:	
!     zst (m)   initial z of the cone vertex
!     rst(m)    initial r of the cone vertex
!     phist(degree) initial toroidal angle phi of cone vertex,
!                   measure from x-z plane.
!     alfast(degree) toroidal angle measure from R-vector through
!                    source
!     betast(degree) poloidal angle measured from z=constant plane,
!                      positive above plane, negative below.
!     alpha1(degree) cone width, half-width at half power of the beam.
!     alpha2(degree) starting angle along cone
!     na1 number of cones (0 for central ray only)
!     na2 number of rays at cone(for na1.ge.0)
!
!     raypatt='toray',  specify ray pattern per the following
!                       toray method:  [Defn of betast changed,
!                       and there are additional namelist inputs.]
!     zst (m)   initial z of the cone vertex
!     rst(m)    initial r of the cone vertex
!     phist(degree) initial toroidal angle phi of cone vertex,
!                   measure from x-z plane.
!     alfast(degree) toroidal angle measure about z-axis
!                    from R-vector through the source.
!     betast(degree) poloidal polar angle measured from positive
!                    z-axis. [DIFFERENT FROM RAYPATT='GENRAY'!]
!     alpha1(degree) cone width, half-width at half power of the beam.
!     gzone:    if 0 then 48 ray case, as specified by mray() below
!               if 1 then there can only be 1 ray, the central ray.
!               if .gt.1 then describes number of elements in mray
!     mray(*):  if gzone .gt.0, use gaussian formulation with this number
!        of rays in corresponding annular zone, otherwise use the
!        usual 1,5,12,12,18  (48 ray) arrangement.
!        mray(1) is effectively 1.
!     cr(*):    azimuthal phase of ray pattern for each zone, in radians;
!       same size array as mray for gaussian formulation.
!      Standard setting for gzone=0 is 0.0,0.1,0.05,-0.05,0.05.
!      
!-----------------------------------------------------------------
 &eccone
! ncone=1
 ncone=2

 raypatt='genray'
! zst(1)=0.1571712110864337d0
!!rst(1)=1.573955374170274d0 
! phist(1)=4.147062324255599d0
! betast(1)=-18.81722591231102d0
! alfast(1)=199.0213150494316d0 
 zst( 1 )=0.156478331368976D+00
 rst( 1 )=0.157400845935811D+01
 phist( 1 )=0.414628673197772D+01
 betast( 1 )=-.188609626975446D+02
 alfast( 1 )=0.199026350592679D+03
!
! zst(2)=0.2571712110864337d0
! rst(2)=1.533955374170274d0 
! phist(2)=4.147062324255599d0
! betast(2)=-18.81722591231102d0
! alfast(2)=199.0213150494316d0 
 zst( 2 )=0.256798565288347D+00
 rst( 2 )=0.153399818784218D+01
 phist( 2 )=0.403642441405337D+01
 betast( 2 )=-.244236994359805D+02
 alfast( 2 )=0.200236912352825D+03
!
 alpha1(1)=4.6500d+00
 alpha2(1)=+1.500d+1
 alpha1(2)=4.6500d+00
 alpha2(2)=+1.500d+1
 na1=0
 na2=6
 powtot(1)=1.0d-1
 powtot(2)=1.0d-1

! raypatt='toray'
! zst(1)=+0.712812d0
! rst(1)=2.42352d0
! phist(1)=+0.000d0
! betast(1)=+118.521d0
! alfast(1)=+196.2d0
! alpha1(1)=1.7d0
! gzone=5
! mray=1,5,12,12,20
! cr=0.0d0,0.1d0,0.0500d0,-0.0500d0,0.0500d0
 &end

!/dentab/
!--------------------------------------------------------------------------
! density profiles (table data, case: idens=1)	dens1(ndens,nbulk)
!--------------------------------------------------------------------------
! ndensa (a parameter in param.i) is max number of points in the 
!   small radius direction.
! ncompa (a parameter in param.i) is a max number of plasma components.
! Input of profiles is set up so spline profiles can be input in tables
!   of size specified through namelist variables ndens and nbulk.
! ndens (variable)    is  number of points in small radius direction
!                     (set in namelist /plasma/). Must be .le. ndensa.   
! nbulk (variable)    is number of plasma components, must have: 
!                     nbulk.le.ncompa, and
!                     first component is for electrons
! nbulk1 is number of density components which must be specified.
! nbulk1 is calculated in dinit_mr subroutine,
! The fragment of dinit_mr is given here to understand the
! nbulk1 value
!
! The number of columns in dentab should be equal to nbulk.
! If nbulk1 < nbulk then we should put the density profiles
!   for the first nbulk1 plasma components in the table.
! The profiles for last (nbulk-nbulk1) plasma components can be arbitrary.
!
!--------------------------------------------------------------------------
!c     calculation of nbulk1
!      if(((izeff.eq.0).or.(izeff.eq.2)).or.(izeff.eq.3)) then
!c        izeff=0, zeff will be calculated using the given ions;
!c                 electron density will be calculated using ion's densities;
!c             =1  ion densities nbulk and nbulk-1 will be calculated  using
!c                 Zeff, electron density and ion's densities(i), i=2,nbulk-1;
!c        izeff=2, zeff will not coincide with the plasma components
!c             =3  it uses eqdsk pres (pressure) and ion densities_i 
!c                 for i=2,... nbulk
!c                 Let temperature T_E=T_i
!c                 pres=dens1(k,1)*temp1(k,1)+
!c                      Sum(i=2,nbulk)(dens1(k,i)*temp1(k,i)
!c                 In this case we will calculate Zeff(rho),
!c                 dens_electron(rho) and T_e(rho)=T_i(rho)
!c                  This case works for nbulk >1 only.
!c             =4  it uses eqdsk pres (pressure), zeff,ions densities
!c                 for i=2,... nbulk-2 (nbulk>2) and T_e(rho),T_i(rho)
!c                 pres=dens1(k,1)*temp1(k,1)+
!c                      Sum(i=2,nbulk)(dens1(k,i)*temp1(k,i)
!c                 In this case we will calculate dens_electron(rho) and
!c                 ion densities for i=nbulk and i=nbulk-1)
!         nbulk1=nbulk
!      else
!c        izeff=1, zeff is given, the ions component will be calculated
!         if (nbulk.eq.1) nbulk1=1
!         if (nbulk.eq.2) then
!	    write(*,*)'nbulk=2 Zeff must be equal charge(2) control it'
!	    write(*,!)'use the option izeff=0'
!	    nbulk1=2
!	    stop
!	 endif
!         if (nbulk.gt.2) nbulk1=nbulk-2
!      endif !izeff
!
!      The case nbulk=1 is used often for ECR and EBW cases.
!      In these cases only the electron component is essential.
!
!      For nbulk=1 and izeff=2 case only the electron density
!      is used in dispersion relation.In this case  Z_effective
!      is used for current drive efficiency calculations.   
!------------------------------------------------------------------------
! dens1(ndens,nbulk) (10!!13/cm!!3)
!------------------------------------------------------------------------
! If  ((izeff.eq.0).or.(izeff.eq.3)) then the electron density
! will be calculated from the charge neutrality.
! In that case we can set the arbitrary values for the electron density
! dens1(k,1), k=1:ndens and should set nbulk1-1 ion densities:
! dens1(k,i), k=1:ndens, i=2:nbulk1.
! A constant radial step is assumed, 
! The first line (k=1, i=1:nbulk1) dens1(1,i) is for rho=0
! The last line (k=ndens, i=1:nbulk1) dens1(ndens,i) is for rho=1
! The example for  izeff=0, ndens=5, nbulk=3, ncomp=4
!         electron-1 ion-2  ion-nbulk
! prof=      0.,     1.2,   1.3,       
!            0.,     2.2,   2.3,       
!            0.,     3.2,   3.3,       
!            0.,     4.2,   4.3,       
!            0.,     5.2,   5.3,       
! ------------------------------------------------------------------------
! Here array prof(ncomp,ndens) was used for convenience in namelist.
! dens1(k,i)=prof(i,k) k=1:ndens, i=1:ncomp
! ------------------------------------------------------------------------
! If (izeff.ne.0) and (izeff.ne.3) then we should set the electron density
! dens1(k,1)  and ion densities dens1(k,i) i=2:nbulk1
!------------------------------------------------------------------------
 &dentab
 prof    =   4.2000000,
             4.1400000,
             4.0500000,
             3.9500000,
             3.8900000,
             3.8100000,
             3.7400000,
             3.6700000,
             3.6100000,
             3.5500000,
             3.4900000,
             3.3500000,
             3.1900000,
             3.0100000,
             2.7600000,
             2.4500000,
             2.1400000,
             1.8500000,
             1.5100000,
             1.2500000,
             0.8500000,
 &end

!/temtab/
!--------------------------------------------------------------------------
! temperature profiles (table data, case: idens=1)	temp1(ndens,nbulk)
! Average temperature temp1=(T_parallel+2*T_perp)/3
!--------------------------------------------------------------------------
! It this namelist we must set electron temp1(ndens,1) and all ion
! species temp1(ndens,i) temperature (keV) {i=2:nbulk}
! Constant radial step is assumed.
! The first line (k=1, i=1:nbulk) temp1(1,i) is for rho=0
! The last line (k=ndens, i=1:nbulk) temp1(ndens,i) is for rho=1
! The example for   ndens=5, nbulk=3, ncomp=4
!                 
!         electron-1 ion-2  ion-nbulk 
! prof=      0.,     1.2,   1.3,       
!            0.,     2.2,   2.3,       
!            0.,     3.2,   3.3,       
!            0.,     4.2,   4.3,       
!            0.,     5.2,   5.3,
!
! In all cases temtab should have nbulk columns with the temperature
! profiles for all nbulk plasma components.    
! ------------------------------------------------------------------------
! Here array prof(ncomp,ndens) was used for convenience in namelist.
! temp1(k,i)=prof(i,k) k=1:ndens, i=1:nbulk
! ------------------------------------------------------------------------
 &temtab
 prof =        1.3000000,
               1.3000000,
               1.2700000,
               1.2400000,
               1.2200000,
               1.1300000,
               1.0500000,
               0.9700000,
               0.8700000,
               0.7900000,
               0.7300000,
               0.6800000,
               0.6100000,
               0.5300000,
               0.4500000,
               0.3700000,
               0.3000000,
               0.2300000,
               0.1600000,
               0.1000000,
               0.0400000,
 &end

!/tpoptab/
!--------------------------------------------------------------------------
! Tpop=T_perp/T_parallel profiles (table data, case: idens=1) 
|      tpop1(ndens,ncomp)
!--------------------------------------------------------------------------
! It this namelist we must set electron tpop1(ndens,1) and all ion
! species tpop1(ndens,i)  {i=2:nbulk}
! Constant radial step is assumed.
! The first line (k=1, i=1:nbulk) tpop1(1,i) is for rho=0
! The last line (k=ndens, i=1:nbulk) tpop1(ndens,i) is for rho=1
! The example for   ndens=5, nbulk=3, ncomp=4
!               
!         electron-1 ion-2  ion-nbulk 
! prof=      1.,     1.2,   1.3,       
!            1.,     1.2,   2.3,       
!            1.,     1.2,   3.3,       
!            1.,     4.2,   4.3,       
!            1.,     5.2,   5.3,       
!
! In all cases tpoptab should has nbulk columns with Tpop
! profiles for all nbulk plasma components.
! ------------------------------------------------------------------------
! Here array prof(ncomp,ndens) was used for convenience in namelist.
! tpop1(k,i)=prof(i,k) k=1:ndens, i=1:nbulk
! ------------------------------------------------------------------------
 &tpoptab
 prof=    1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
          1.0,
 &end

!/vflowtab/
!--------------------------------------------------------------------------
! vflow  profiles (table data, case: idens=1) 
|      vdflow1(ndens,ncomp) is the drift velocity || B (v/sec)
!--------------------------------------------------------------------------
! It this namelist we must set drift velocity for electrons vflow1(ndens,1)
| and all ion species vflow1(ndens,i)  {i=2:nbulk}
! Constant radial step is assumed.
! The first line (k=1, i=1:nbulk) vflow1(1,i) is for rho=0
! The last line (k=ndens, i=1:nbulk) vflow1(ndens,i) is for rho=1
! The example for   ndens=5, nbulk=3, ncompa=4
!            
!         electron-1 ion-2  ion-nbulk 
! prof=      1.,     1.2,   1.3,       
!            1.,     1.2,   2.3,       
!            1.,     1.2,   3.3,       
!            1.,     4.2,   4.3,       
!            1.,     5.2,   5.3,       
!
! In all cases vflowtab should has nbulk columns with vflow
! profiles for all nbulk plasma components.
! ------------------------------------------------------------------------
! Here array prof(ncomp,ndens) was used for convenience in namelist.
! vflow1(k,i)=prof(i,k) k=1:ndens, i=1:nbulk
! ------------------------------------------------------------------------
 &vflowtab
 prof=    0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
          0.0,
 &end


!/zeftab/
!--------------------------------------------------------------------------
! Zeff profiles (table data, case: idens=1)	zeff1(ndens)
!--------------------------------------------------------------------------
! It this namelist we must set zeff1(ndens)
! Constant radial step is assumed.
! The first value  zeff1(1) is for rho=0
! The last value zeff2(ndens) is for rho=1
! The example for  ndens=5
! zeff1= 1., 1., 1., 1. 1.
!--------------------------------------------------------------------------
 &zeftab
 zeff1=21*1.5
 &end

!/read_diskf/
!   i_diskf=0 usage analytic Maxwellian electron distribution 
!   i_diskf=1 read in the file diskf
!   i_diskf=2 read in the file netcdfnm.nc 
!   i_diskf=3 analytic calculation of the non-Maxwellian distribution
!------------------------------------------------------
!   the data for analytic non-Maxwellian electron distribution
!    
!   jx   - the number of used normalized momentum/ mesh points
!   lrz  - the number of used radial mesh points
!   iym  - the number of used pitch-angle mesh points 
!          (here the same at each radius)
!   ngen - the number of plasma species (here we use only electron specie with 
!          the number of specie k=1) 
!   jxa,iya,lrza,ngena ! the max values for jx,iym,lrz,ngen 
!   rtem0 - ratio tem0/electron_temperature(rho=0)
!           tem0 is the max energy for the momentum normalization (KeV) 
!-----tail parameters
!     f_tail=H(rho,rt1,rt2)*f_rel_Maxw(ttail),
!     H(x,x1,x2) is the box function. H=1 for x1<x<x2 otherwise H=0  
!   r1t,r2t    small normalized radii for the tail localization  
!   rtail      the relation the tail density to the total density
!   ttail      tail temperature (KeV)!tail temperature (KeV)
!-----hot parameters
!     f_hot=H(rho,rh1,rh2)*H(epar,hotmnpar,hotmxpar)*H(eper,hotmnper,hotmxper)*
!     *(p_per/mc)**hotexp*exp{-mu(thoppar)(p_par/m_ec)**2-mu(thopper)(p_per/m_ec)**2}.
!     Here mu(T)=m_e*c**2/T  
!   r1h,r2h             - small normalized radii for the hot localization
!   rhot                - the relation of hot hot density to the total density
!   thotpar,thotper     - parallel and perpendicular hot temperatures (KeV)
!   hotmnpar,hotmxpar   - the boundaries of the parallel energy box(KeV)
!                         hotmnpar < epar < hotmxpar  
!   hotmnper,hotmxper   - the boundaries of the perpendicular energy box(KeV)
!                         hotmnper < eper < hotmxper  
!   hotexp              - the degree of the perpendicular momentum: (p_per/mc)**hotexp
!-----beam parameters
!     f_beam=H(rho,rb1,rb2)*exp{-0.5*mu(tbeam)*
!              [(p_par-p_beam_par)**2+(p_per-p_beam_per)**2]/(m_e*c)**2}
!     Here
!          (p_beam /m_e*c)**2=ebeam**2/(m_e**2*c**4)-1
!           p_beam_par=p_beam*cos(thbeam)
!           p_beam_per=p_beam*sin(thbeam)
!   r1b,r2b      - small normalized radii for the beam localization
!   rbeam        - the relation of the beam density to the total density
!   ebeam        - beam energy (KeV)
!   thbeam       - beam pitch angle (0=<degree=<180) 
!   tbeam        - beam temperature (KeV)
!--------------------------------------------------------------------
 &read_diskf
 i_diskf=0
 netcdfnm='netcdfnm.nc'
 jx=200
 lrz=20
 iym=100  
 ngen=1  
 rtem0=10.d0
 r1t=1.d0
 r2t=0.d0      
 rtail=0.0d0      
 ttail=1.d0
 r1h=1.d0
 r2h=0.d0          
 rhot=2.5d-3            
 thotpar=1.d0
 thotper=1.d0     
 hotmnpar=1.d0
 hotmxpar=2.d0   
 hotmnper=1.d0
 hotmxper=2.d0    
 hotexp=0.d0    
 r1b=1.d0
 r2b=2.d0     
 rbeam=0.d0        
 ebeam=1.d0       
 thbeam=30.d0       
 tbeam=1.d0   
 &end


!--------------------------------------------------
!/emission/ 
!the data for emission calculations
!   i_emission=0 no emission
!             =1 emission calculations
!
!   0<tol_emis=<1 tolerance parameter to add the new mesh point s_n
!                if in_0(n)>tol_emis*i_0 
!   nharm1=< nharm =<nharm2 gives the used EC harmonics (Not work now)
!   nfreq=<nfreqa (see param.i) is the number of frequencies
!         nfreq=1 gives detailed plots of emission from a single ray     
!         nfreq.gt.1 gives spectra covering the specified frequency
!         range
!   if nfreq=1 code will use the frequency determined by 
!              frqncy(GHz) given in namelist wave
!   freq00=<freq01 are ratios of the minimal and maximal emission
!         frequencies to the central electron gyro-frequency
!         f_ce(rho=0)
!          
!    wallr is the wall reflection coefficient, {0=< wallr =<1}
!
!    i_rrind chooses the subroutine to calculate N_ray (default =0)
!      i_rrind=0  N_ray=N
!      i_rrind =1 from cold electron plasma using rrind   
!      i-rrind =2 from hot non-relativistic dispersion relation  
!--------------------------------------------------
!    i_r_2nd_harm=1 to calculate the major radius of the EC 2nd harmonic
!                =0 do not calculate (default =0) 
!              (to use drawemfr.in file i_r_2nd_harm=1)
!              (to use drawemf1.in file i_r_2nd_harm=0)
!---------------------------------------------------
 &emission
 i_emission=0
 tol_emis=5.0d-3
 nharm1=1
 nharm2=1
 nfreq=5
 freq00=0.585d0
 freq01=0.885d0
 wallr=0.9d0
 i_rrind=1
 i_r_2nd_harm=0
 &end

!/ox/
!------------------------------------------------------------------------
!namelist related to calculation of EC cone vertex coordinates
!for the optimal OX mode  conversion (details on this calculation
!will be given in the Genray manual)
!------------------------------------------------------------------------
! i_ox =0 /default/ do not use these calculations.
!      =1 calculations of the optimal EC cone vertex for OX conversion. 
!         The optimal direction will give optimal N_parallel
!         in OX conversion point at Xe=1.
!         In this case should have istart=1.
!         Optimal launch angles are output to ECcone_optimal.dat.
!      =2 launch the ray using EC cone vertex coordinates calculations
!         and using OX transmission procedure.
!------------------------------------------------------------------------
! theta_bot(icone) < theta_top(icone) 
! are the poloidal angle boundaries (degree) at
! Xe=1 surface. They are used in genray.f to find the optimal ray
! direction at the given EC cone vertex icone=1,...,ncone
!------------------------------------------------------------------------
! i_ox_poloidal_max is the maximal number of the poloidal angles used
!                   in the bisection method. This method calculates the
!                   the poloidal angle theta_pol of the point M at the
!                   flux surface Xe(rho=1) : M(poloidal_angle,rho).  
!                   The ray launched from M to the plasma edge will
!                   go to the EC cone vertex.
!------------------------------------------------------------------------
! eps_antenna       is the accuracy with which the ray launched from the
!                   flux surface Xe(rho)=1 at the given poloidal angle
!                   reaches the given cone vertex (i_ox=1).
!                   This accuracy is calculated using the distance of 
!                   the vertex radial coordinate from the specified:
!                   antenna position, 
!                   sqrt((rst_ox-rst)**2+(zst_ox-zst)**2)< eps_antenna
!                                                          (meters)
!------------------------------------------------------------------------
!eps_xe   The parameter which sets the vicinity
!         of the O-mode cutoff surface.
!         If xe < (1-eps_xe) then this subroutine will 
!         creates the ray jump in small radius direction
!         and find the X mode.
!         eps_xe=1.d-2 is seted as default in dinit.f
!------------------------------------------------------------------------

&ox
i_ox=2
theta_bot(1)=0.0d0
theta_top(1)=30.d0
theta_bot(2)=0.0d0
theta_top(2)=30.d0
i_ox_poloidal_max=20
eps_antenna=1.d-4
eps_xe=1.d-2
&end