gray/src/gray_params.f90

module gray_params
  use const_and_precisions, only : wp_
  implicit none
  integer, parameter :: lenfnm=256
  integer, parameter :: headw=132,headl=21

  type antctrl_type
    real(wp_) :: alpha, beta, power
    real(wp_) :: psi, chi
    integer :: iox
    integer :: ibeam
    character(len=lenfnm) :: filenm
  end type antctrl_type

  type eqparam_type
    real(wp_) :: ssplps, ssplf, factb
    integer :: sgnb, sgni, ixp
    integer :: iequil, icocos, ipsinorm, idesc, ifreefmt
    character(len=lenfnm) :: filenm
  end type eqparam_type

  type prfparam_type
    real(wp_) :: psnbnd, sspld, factne, factte
    integer :: iscal, irho !, icrho, icte, icne, iczf
    integer :: iprof
    character(len=lenfnm) :: filenm
  end type prfparam_type

  type rtrparam_type
    real(wp_) :: rwmax, dst
    integer :: nrayr, nrayth, nstep
    integer :: igrad, idst, ipass, ipol
  end type rtrparam_type

  type hcdparam_type
    integer :: iwarm, ilarm, imx, ieccd
  end type hcdparam_type

  type outparam_type
    integer :: ipec, nrho, istpr, istpl
  end type outparam_type


  integer, save :: iequil,iprof,ipol
  integer, save :: iwarm,ilarm,imx,ieccd
  integer, save :: igrad,idst,ipass
  integer, save :: istpr0,istpl0
  integer, save :: ipec,nnd

contains
  
  subroutine print_params(rtrparam,hcdparam,antctrl,eqparam,rwall, &
                         prfparam,outparam,strout)
    implicit none
! arguments
    type(rtrparam_type), intent(in) :: rtrparam
    type(hcdparam_type), intent(in) :: hcdparam
    type(antctrl_type), intent(in) :: antctrl
    type(eqparam_type), intent(in) :: eqparam
    real(wp_), intent(in) :: rwall
    type(prfparam_type), intent(in) :: prfparam
    type(outparam_type), intent(in) :: outparam
    character(len=*), dimension(:), intent(out) :: strout   ! min len=110, dimension(21)
! local variables
    character(len=8) :: rdat
    character(len=10) :: rtim
#ifndef REVISION
    character(len=*), parameter :: REVISION="unknown"
#endif

! date and time
	call date_and_time(rdat,rtim)
    write(strout(1),'("# Run date/time: ",a4,2("/",a2),1x,2(a2,":"),a6)') &
      rdat(1:4),rdat(5:6),rdat(7:8),rtim(1:2),rtim(3:4),rtim(5:10)
! SVN revision
    write(strout(2),'("# GRAY SVN revision: ",a)') REVISION
! equilibrium input data
    if (eqparam%iequil > 0) then
      write(strout(3),'("# EQL input: ",a)') trim(eqparam%filenm)
!!!!!!! missing values
      write(strout(7),'("# EQL B0 R0 aminor Rax zax:",5(1x,e12.5))') &
        0._wp_, 0._wp_, 0._wp_, 0._wp_, 0._wp_
    else
      write(strout(3),'("# EQL input: N/A (vacuum)")')
      write(strout(7),'("# EQL B0 R0 aminor Rax zax: N/A (vacuum)")')
    end if
    write(strout(4),'("# EQL iequil sgnb sgni factb:",3(1x,i4),1x,e12.5)') &
      eqparam%iequil, eqparam%sgnb, eqparam%sgni, eqparam%factb
    if (eqparam%iequil > 1) then
      write(strout(5),'("# EQL icocos ipsinorm idesc ifreefmt:",4(1x,i4))') &
        eqparam%icocos, eqparam%ipsinorm, eqparam%idesc, eqparam%ifreefmt
      write(strout(6),'("# EQL ssplps ssplf ixp:",2(1x,e12.5),1x,i4)') &
        eqparam%ssplps, eqparam%ssplf, eqparam%ixp
    else
      write(strout(5),'("# EQL icocos ipsinorm idesc ifreefmt: N/A (analytical)")')
      write(strout(6),'("# EQL ssplps ssplf ixp: N/A (analytical)")')
    end if
! profiles input data
    if (eqparam%iequil > 0) then
      write(strout(8),'("# PRF input: ",a)') trim(prfparam%filenm)
      write(strout(9),'("# PRF iprof iscal factne factte:",2(1x,i4),2(1x,e12.5))') &
        prfparam%iprof,prfparam%iscal,prfparam%factne,prfparam%factte
      if (prfparam%iprof > 0) then
        write(strout(10),'("# PRF irho psnbnd sspld:",1x,i4,2(1x,e12.5))') &
          prfparam%irho,prfparam%psnbnd,prfparam%sspld
      else
        write(strout(10),'("# PRF irho psnbnd sspld: N/A (analytical)")')
      end if
!!!!!!! missing values
      write(strout(11),'("# PRF Te0 ne0 Zeff0:",3(1x,e12.5))') &
        0._wp_, 0._wp_, 0._wp_
    else
      write(strout(8),'("# PRF input: N/A (vacuum)")')
      write(strout(9),'("# PRF iprof iscal factne factte: N/A (vacuum)")')
      write(strout(10),'("# PRF irho psnbnd sspld: N/A (vacuum)")')
      write(strout(11),'("# PRF Te0 ne0 Zeff0: N/A (vacuum)")')
    end if
! launch parameters
    write(strout(12),'("# ANT input: ",a)') trim(antctrl%filenm)
    write(strout(13),'("# ANT ibeam iox psi chi:",2(1x,i4),2(1x,e12.5))') &
      antctrl%ibeam, antctrl%iox, antctrl%psi, antctrl%chi
    write(strout(14),'("# ANT alpha beta power:",3(1x,e12.5))') &
      antctrl%alpha, antctrl%beta, antctrl%power
!!!!!!! missing values
    write(strout(15),'("# ANT x0 y0 z0:",3(1x,e12.5))') &
      0._wp_, 0._wp_, 0._wp_
!!!!!!! missing values
    write(strout(16),'("# ANT wx wy Rcix Rciy psiw psir:",6(1x,e12.5))') &
      0._wp_, 0._wp_, 0._wp_, 0._wp_, 0._wp_, 0._wp_
! wall parameters
    write(strout(17),'("# RFL rwall:",1x,e12.5)') rwall
! code parameters
    write(strout(18),'("# COD igrad idst ipass ipol:",4(1x,i4))') &
      rtrparam%igrad, rtrparam%idst, rtrparam%ipass, rtrparam%ipol
    write(strout(19),'("# COD nrayr nrayth nstep rwmax dst:",3(1x,i4),2(1x,e12.5))') &
      rtrparam%nrayr, rtrparam%nrayth, rtrparam%nstep, rtrparam%rwmax, rtrparam%dst
    write(strout(20),'("# COD iwarm ilarm imx ieccd:",4(1x,i4))') &
      hcdparam%iwarm, hcdparam%ilarm, hcdparam%imx, hcdparam%ieccd
    write(strout(21),'("# COD ipec nrho istpr istpl:",4(1x,i4))') &
      outparam%ipec, outparam%nrho, outparam%istpr, outparam%istpl
  end subroutine print_params
  
  subroutine read_params(filenm,rtrparam,hcdparam,antctrl,eqparam,rwall, &
                         prfparam,outparam,unit)
    use utils, only : get_free_unit
    implicit none
! arguments
    character(len=*), intent(in) :: filenm
    type(rtrparam_type), intent(out) :: rtrparam
    type(hcdparam_type), intent(out) :: hcdparam
    type(antctrl_type), intent(out) :: antctrl
    type(eqparam_type), intent(out) :: eqparam
    real(wp_), intent(out) :: rwall
    type(prfparam_type), intent(out) :: prfparam
    type(outparam_type), intent(out) :: outparam
    integer, intent(in), optional :: unit
! local variables
    integer :: u

    if (present(unit)) then
      u=unit
    else
      u = get_free_unit()
    end if
    open(u,file=filenm,status= 'old',action='read')

! ==========================================================================
!     nrayr     number of rays in radial direction
!     nrayth      number of rays in angular direction
!     rwmax      normalized maximum radius of beam power
!     rwmax=1 -> last ray at radius = waist
    read(u,*) rtrparam%nrayr, rtrparam%nrayth, rtrparam%rwmax
!     igrad=0    optical ray-tracing, initial conditions as for beam
!     igrad=1    quasi-optical ray-tracing
!     igrad=-1   ray-tracing, init. condit.
!                from center of mirror and with angular spread
!     ipass=1/2 1 or 2 passes into plasma
!     ipol=0 compute mode polarization at antenna, ipol=1 use polariz angles
    read(u,*) rtrparam%igrad, rtrparam%ipass, rtrparam%ipol
!     dst       integration step
!     nstep      maximum number of integration steps
!     idst=0/1/2  0 integration in s, 1 integr. in ct, 2 integr. in Sr
    read(u,*) rtrparam%dst, rtrparam%nstep, rtrparam%idst
! ==========================================================================
!     iwarm=0  :no absorption and cd
!     iwarm=1  :weakly relativistic absorption
!     iwarm=2  :relativistic absorption, n<1 asymptotic expansion
!     iwarm=3  :relativistic absorption, numerical integration
!     ilarm    :order of larmor expansion
!     imx      :max n of iterations in dispersion, imx<0 uses 1st
!               iteration in case of failure after |imx| iterations
    read(u,*) hcdparam%iwarm,hcdparam%ilarm,hcdparam%imx

!     ieccd 0/1 NO/YES  ECCD calculation ieccd>0 different CD models
    read(u,*) hcdparam%ieccd
! ==========================================================================
!     alpha0, beta0 (cartesian) launching angles
    read(u,*) antctrl%alpha, antctrl%beta
!     p0mw      injected power (MW)
    read(u,*) antctrl%power
!     abs(iox)=1/2    OM/XM
!     psipol0,chipol0 polarization angles at the antenna (if iox<0)
    read(u,*) antctrl%iox, antctrl%psi, antctrl%chi

!     ibeam=0  :read data for beam as above
!     ibeam=1  :read data from file simple astigmatic beam
!     ibeam=2  :read data from file general astigmatic beam
    read(u,*) antctrl%ibeam
    read(u,*) antctrl%filenm
! ==========================================================================
!     iequil=0  :vacuum
!     iequil=1  :analytical equilibrium
!     iequil=2  :read eqdsk
    read(u,*) eqparam%iequil
    read(u,*) eqparam%filenm
!     icocos    :index for equilibrium from COCOS - O. Sauter Feb 2012
!     ipsinorm  :0 standard EQDSK format, 1 format Portone summer 2004
    read(u,*) eqparam%icocos, eqparam%ipsinorm, eqparam%idesc, eqparam%ifreefmt
!     ixp=0,-1,+1  : no X point , bottom/up X point
!     ssplps       : spline parameter for psi interpolation
    read(u,*) eqparam%ixp, eqparam%ssplps !, eqparam%ssplf
    eqparam%ssplf=0.01_wp_
!     signum of toroidal B and I
!     factb     factor for magnetic field (only for numerical equil)
!     scaling adopted: beta=const, qpsi=const, nustar=const
    read(u,*) eqparam%sgnb, eqparam%sgni, eqparam%factb
! ==========================================================================
    read(u,*) rwall
! ==========================================================================
!     iprof=0  :analytical density and temp. profiles
!     iprof>0  :numerical density and temp. profiles
    read(u,*) prfparam%iprof, prfparam%irho  ! irho=0,1,2 -> num profiles vs rhot,rhop,psin
    read(u,*) prfparam%filenm
!     psbnd     value of psi ( > 1 ) of density boundary
    read(u,*) prfparam%psnbnd, prfparam%sspld
!    prfparam%sspld=0.001_wp_
!     iscal    ne Te scaling 0: nustar=const,  1: n_greenw=const; 2 no rescaling
!     factT factn   factor for Te&ne scaling 
    read(u,*) prfparam%factte, prfparam%factne, prfparam%iscal
! ==========================================================================
!     ipec=0/1  :pec profiles grid in psi/rhop
!     nrho       :number of grid steps for pec profiles +1
    read(u,*) outparam%ipec, outparam%nrho
!     istpr0     projection step = dsdt*istprj
!     istpl0      plot step = dsdt*istpl
    read(u,*) outparam%istpr, outparam%istpl
    close(u)
  end subroutine read_params

  subroutine set_codepar(eqparam,prfparam,outparam,rtrparam,hcdparam)
    implicit none
    type(eqparam_type), intent(in) :: eqparam
    type(prfparam_type), intent(in) :: prfparam
    type(outparam_type), intent(in) :: outparam
    type(rtrparam_type), intent(in) :: rtrparam
    type(hcdparam_type), intent(in) :: hcdparam

    iequil=eqparam%iequil
    iprof=prfparam%iprof

    ipec=outparam%ipec
    nnd=outparam%nrho
    istpr0=outparam%istpr
    istpl0=outparam%istpl

    ipol=rtrparam%ipol
    igrad=rtrparam%igrad
    idst=rtrparam%idst
    ipass=rtrparam%ipass
    if (rtrparam%nrayr<5) then
      igrad=0
      print*,'    nrayr < 5 ! => OPTICAL CASE ONLY'
      print*,' '
    end if

    iwarm=hcdparam%iwarm
    ilarm=hcdparam%ilarm
    imx=hcdparam%imx
    ieccd=hcdparam%ieccd

  end subroutine set_codepar

end module gray_params