src/gray_core.f90: refactor [wip]

src/gray_params.f90

@@ -135,14 +135,14 @@ module gray_params
 
   ! Raytracing parameters
   type raytracing_parameters
-    real(wp_)                :: dst                  ! Integration step size
-    real(wp_)                :: rwmax = 1            ! Normalized maximum radius of beam power
-    integer                  :: nray, nrayr, nrayth  ! Total/radial/angular number of rays
-    integer                  :: igrad = 0            ! Complex eikonal switch
-    integer                  :: nstep = 12000        ! Max number of integration steps
-    integer(kind(step_enum)) :: idst  = STEP_TIME    ! Choice of the integration variable
-    integer                  :: ipass = 1            ! Number of plasma passes
-    logical                  :: ipol  = .false.      ! Whether to compute wave polarisation (from chi, psi)
+    real(wp_)                :: dst                   ! Integration step size
+    real(wp_)                :: rwmax = 1             ! Normalized maximum radius of beam power
+    integer                  :: nray, nrayr, nrayth   ! Total/radial/angular number of rays
+    logical                  :: igrad = .false.       ! Complex eikonal switch
+    integer                  :: nstep = 12000         ! Max number of integration steps
+    integer(kind(step_enum)) :: idst  = STEP_TIME     ! Choice of the integration variable
+    integer                  :: ipass = 1             ! Number of plasma passes
+    logical                  :: ipol  = .false.       ! Whether to compute wave polarisation (from chi, psi)
   end type
 
   ! EC resonant heating & Current Drive parameters
@@ -341,7 +341,7 @@ contains
     call append(header, line)
 
     ! code parameters
-    write(line, '("# COD igrad idst ipass ipol:",3(1x,i4),1x,l4)') &
+    write(line, '("# COD igrad idst ipass ipol:",l4,2(1x,i4),1x,l4)') &
       params%raytracing%igrad, params%raytracing%idst, &
       params%raytracing%ipass, params%raytracing%ipol
     call append(header, line)

src/gray_plasma.f90

@@ -21,7 +21,7 @@ module gray_plasma
   end type
 
   abstract interface
-    subroutine density_sub(self, psin, dens, ddens)
+    pure subroutine density_sub(self, psin, dens, ddens)
       ! Computes the density its first derivative as a function of
       ! normalised poloidal flux.
       !
@@ -32,7 +32,7 @@ module gray_plasma
       real(wp_),              intent(out) :: dens, ddens  ! density and first derivative
     end subroutine density_sub
 
-    function temp_fun(self, psin) result(temp)
+    pure function temp_fun(self, psin) result(temp)
       ! Computes the temperature as a function of the
       ! normalised poloidal flux.
       !
@@ -43,7 +43,7 @@ module gray_plasma
       real(wp_)                          :: temp
     end function temp_fun
 
-    function zeff_fun(self, psin) result(zeff)
+    pure function zeff_fun(self, psin) result(zeff)
       ! Computes the effective charge Z_eff as a
       ! function of the normalised poloidal flux.
       import :: abstract_plasma, wp_
@@ -97,7 +97,7 @@ module gray_plasma
 
 contains
 
-  subroutine analytic_density(self, psin, dens, ddens)
+  pure subroutine analytic_density(self, psin, dens, ddens)
     ! subroutine arguments
     class(analytic_plasma), intent(in)  :: self
     real(wp_),              intent(in)  :: psin         ! normalised poloidal flux
@@ -116,7 +116,7 @@ contains
   end subroutine analytic_density
 
 
-  function analytic_temp(self, psin) result(temp)
+  pure function analytic_temp(self, psin) result(temp)
     ! function arguments
     class(analytic_plasma), intent(in) :: self
     real(wp_),              intent(in) :: psin
@@ -139,7 +139,7 @@ contains
   end function analytic_temp
 
 
-  function analytic_zeff(self, psin) result(zeff)
+  pure function analytic_zeff(self, psin) result(zeff)
     ! function arguments
     class(analytic_plasma), intent(in) :: self
     real(wp_),              intent(in) :: psin
@@ -155,7 +155,7 @@ contains
   end function analytic_zeff
 
 
-  subroutine numeric_density(self, psin, dens, ddens)
+  pure subroutine numeric_density(self, psin, dens, ddens)
     use logger, only : log_error
 
     ! subroutine arguments
@@ -221,16 +221,16 @@ contains
       end block
     end if
 
-    if (dens < 0) then
-      write (msg, '("negative density:", 2(x,a,"=",g0.3))') &
-        'ne', dens, 'ψ', psin
-      call log_error(msg, mod='coreprofiles', proc='density')
-    end if
+    ! if (dens < 0) then
+    !   write (msg, '("negative density:", 2(x,a,"=",g0.3))') &
+    !     'ne', dens, 'ψ', psin
+    !   call log_error(msg, mod='coreprofiles', proc='density')
+    ! end if
 
   end subroutine numeric_density
 
 
-  function numeric_temp(self, psin) result(temp)
+  pure function numeric_temp(self, psin) result(temp)
     ! function arguments
     class(numeric_plasma), intent(in) :: self
     real(wp_),             intent(in) :: psin
@@ -246,7 +246,7 @@ contains
   end function numeric_temp
 
 
-  function numeric_zeff(self, psin) result(zeff)
+  pure function numeric_zeff(self, psin) result(zeff)
     ! function arguments
     class(numeric_plasma), intent(in) :: self
     real(wp_),             intent(in) :: psin

src/main.f90

@@ -77,10 +77,10 @@ program main
 
   ! Do some checks on the inputs
   associate (p => params%raytracing)
-    if (p%igrad == 1 .and. p%nrayr < 5) then
-      p%igrad = 0
+    if (p%igrad .and. p%nrayr < 5) then
+      p%igrad = .false.
       call log_message(level=WARNING, mod='main', &
-                       msg='igrad = 1 but nrayr < 5: disabling beamtracing')
+                       msg='igrad = .true. but nrayr < 5: disabling beamtracing')
     end if
     if (p%nrayr == 1) p%nrayth = 1
 

src/main_convert.f90

@@ -156,8 +156,8 @@ contains
     write(u, fmt) "nrayth", params%raytracing%nrayth
     if (params%raytracing%rwmax /= defaults%raytracing%rwmax) &
       write(u, fmt) "rwmax", params%raytracing%rwmax
-    if (params%raytracing%igrad /= defaults%raytracing%igrad) &
-      write(u, fmt) "igrad", params%raytracing%igrad
+    if (params%raytracing%igrad .neqv. defaults%raytracing%igrad) &
+      write(u, fmt) "igrad", format_bool(params%raytracing%igrad)
     if (params%raytracing%ipass /= defaults%raytracing%ipass) &
       write(u, fmt) "ipass", params%raytracing%ipass
     if (params%raytracing%ipol .neqv. defaults%raytracing%ipol) &

tests/06-ITER-startup/inputs/gray.ini

@@ -1,7 +1,7 @@
 [raytracing]
 nrayr = 11
 nrayth = 16
-igrad = 1
+igrad = true
 dst = 0.1
 nstep = 8000