emission.f

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!  The @header, @table_section, @table_subsection, @item and @end_table commands
!  are custom defined commands in Doxygen.in. They are defined under ALIASES.
!  For the page created here, the 80 column limit is exceeded. Arguments of
!  aliases are separated by ','. If you intended ',' to be a string you must use
!  an escaped comma '\,'.
!
!*******************************************************************************
!
!  Note separating the Doxygen comment block here so detailed decription is
!  found in the Module not the file.
!
!*******************************************************************************
 
      MODULE emission
      USE stel_kinds
      USE profiler
      USE mpi_inc
 
      IMPLICIT NONE
 
!*******************************************************************************
!  DERIVED-TYPE DECLARATIONS
!  1) emission class
!
!*******************************************************************************
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
      TYPE emission_class
         REAL (rprec)                          :: te_start = 0.0
         REAL (rprec)                          :: te_step = 0.0
 
         REAL (rprec), DIMENSION(:,:), POINTER :: emissivity => null()
      END TYPE
 
!*******************************************************************************
!  INTERFACE BLOCKS
!*******************************************************************************
      INTERFACE emission_construct
         MODULE PROCEDURE emission_construct_netcdf
      END INTERFACE
 
      CONTAINS
!*******************************************************************************
!  CONSTRUCTION SUBROUTINES
!*******************************************************************************
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
      FUNCTION emission_construct_netcdf(filename)
      USE ezcdf
      USE v3_utilities
 
      IMPLICIT NONE
 
!  Declare Arguments
      TYPE (emission_class), POINTER :: emission_construct_netcdf
      CHARACTER (len=*), INTENT(in)  :: filename
 
!  local variables
      REAL (rprec)                   :: start_time
      INTEGER                        :: iou
      INTEGER                        :: status
      INTEGER, DIMENSION(2)          :: dim_lengths
 
!  Start of executable code
      start_time = profiler_get_start_time()
 
      ALLOCATE(emission_construct_netcdf)
 
      CALL cdf_open(iou, trim(filename), 'r', status)
      CALL assert_eq(0, status, 'failed to open emission file')
 
      CALL cdf_read(iou, 'te_start', emission_construct_netcdf%te_start)
      CALL cdf_read(iou, 'te_step', emission_construct_netcdf%te_step)
 
      CALL cdf_inquire(iou, 'emissivity', dim_lengths)
      ALLOCATE(emission_construct_netcdf%emissivity(dim_lengths(1),            &
     &                                              dim_lengths(2)))
      CALL cdf_read(iou, 'emissivity',                                        &
     &              emission_construct_netcdf%emissivity)
 
      CALL cdf_close(iou)
 
      CALL profiler_set_stop_time('emission_construct_netcdf',                 &
     &                            start_time)
 
      END FUNCTION
 
!*******************************************************************************
!  DESTRUCTION SUBROUTINES
!*******************************************************************************
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
      SUBROUTINE emission_destruct(this)
 
      IMPLICIT NONE
 
!  Declare Arguments
      TYPE (emission_class), POINTER :: this
 
!  Start of executable code
      IF (ASSOCIATED(this%emissivity)) THEN
         DEALLOCATE(this%emissivity)
         this%emissivity => null()
      END IF
 
      DEALLOCATE(this)
 
      END SUBROUTINE
 
!*******************************************************************************
!  GETTER SUBROUTINES
!*******************************************************************************
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
      FUNCTION emission_get_emission(this, te, ne, index)
 
      IMPLICIT NONE
 
!  Declare Arguments
      REAL (rprec)                      :: emission_get_emission
      TYPE (emission_class), INTENT(in) :: this
      REAL (rprec), INTENT(in)          :: te
      REAL (rprec), INTENT(in)          :: ne
      INTEGER, INTENT(in)               :: index
 
!  local variables
      INTEGER                           :: i_low, i_high
      REAL (rprec)                      :: w_low, w_high
      INTEGER                           :: num_temp
      REAL (rprec)                      :: start_time
 
!  Start of executable code
      start_time = profiler_get_start_time()
 
!  Avoid a divid by zero error.
      IF (te .eq. 0.0) THEN
         emission_get_emission = 0.0
         CALL profiler_set_stop_time('emission_get_emission',                  &
     &                               start_time)
         RETURN
      END IF
 
!  Determine the high and low indices for interpolating. If the temperature is
!  beyond the first or last temperature, linearly interpolate from the last two
!  temperatures.
      num_temp = SIZE(this%emissivity, 2)
 
      IF (te .lt. this%te_start) THEN
!  Assume no emission at zero eV. Interpolate between zero and the first
!  temperature level.
         i_low = 1
         i_high = 2
 
         w_high = 0.0
         w_low = te/this%te_start
      ELSE IF (te .gt. this%te_start + (num_temp - 1)*this%te_step) THEN
         i_high = num_temp
         i_low = i_high - 1
 
         w_high = 1.0
         w_low = 0.0
      ELSE
         i_low = (te - this%te_start)/this%te_step + 1
         i_high = i_low + 1
 
         w_high = (te - this%te_start)/this%te_step + 1.0 - i_low
         w_low = 1.0 - w_high
      END IF
 
      emission_get_emission = w_high*this%emissivity(index,i_high)             &
     &                      + w_low*this%emissivity(index,i_low)
 
      emission_get_emission = ne*ne/sqrt(te)*max(emission_get_emission,        &
     &                                           0.0)
 
      CALL profiler_set_stop_time('emission_get_emission', start_time)
 
      END FUNCTION
 
      END MODULE