params_equi_flare_m.f90 Source File

module params_equi_flare_m
    !! Module for parameters of the flare equilibrium type
    use, intrinsic :: iso_fortran_env, only: IOSTAT_END, IOSTAT_EOR
    use precision_m, only : FP
    use screen_io_m, only: get_stdout
    use error_handling_m, only: handle_error, error_info_t
    use status_codes_m, only: PARALLAX_ERR_PARAMETERS
    use FLARE_model, only : FLARE_TYPE_EQUI3D_MGRID => TYPE_EQUI3D_MGRID, &
        FLARE_TYPE_EQUI3D_COILSET => TYPE_EQUI3D_COILSET, &
        FLARE_TYPE_EQUI3D_HINT => TYPE_EQUI3D_HINT, &
        FLARE_TYPE_EQUI2D_GEQDSK => TYPE_EQUI2D_GEQDSK
    use FLARE_bspline3d, only : &
        FLARE_BSPLINE3D_MAGNETIC_FIELD => BSPLINE3D_MAGNETIC_FIELD
    implicit none

    integer, private, parameter :: lmax_char = 512
    ! Maximum lengths of strings in parameters
    integer, private, parameter :: nmax_amplitudes = 100
    ! Maximum number of amplitudes in parameters
    integer, private, parameter :: nmax_bnds = 50
    ! Maximum number of boundary segments

    ! Parameters for selection of flare equilibrium type and boundary type
    character(len=lmax_char), protected :: equi_flare_type = &
        FLARE_TYPE_EQUI3D_MGRID
    !! Flare equilibrium type, available:
    !! - FLARE_TYPE_EQUI3D_MGRID = 'equi3d_mgrid' (default)
    !! - FLARE_TYPE_EQUI3D_COILSET = 'equi3d_coilset'
    !! - FLARE_TYPE_EQUI3D_HINT = 'equi3d_hint'
    !! - FLARE_TYPE_EQUI2D_GEQDSK = 'geqdsk'
    character(len=lmax_char), protected :: boundary_flare_type = 'none'
    !! Type for boundary description, available:
    !! - 'none' (default)
    !! - 'kisslinger'
    character(len=lmax_char), protected :: rho_flare_type = 'none'
    !! Type for flux surface (rho) description, available:
    !! - 'none' (default)
    !! - 'kisslinger'

    ! Parameters for different flare equilibria
    integer, protected :: spline_order = 4
    !! Order for spline interpolation
    character(len=lmax_char), protected :: bspline_data = &
        FLARE_BSPLINE3D_MAGNETIC_FIELD
    !! Data type to be interpolated, available:
    !! - FLARE_BSPLINE3D_MAGNETIC_FIELD = 'magnetic_field' (default)
    !! - FLARE_BSPLINE3D_VECTOR_POTENTIAL = 'vector_potential'
    real(FP), protected :: scale_bt = 1.0_FP
    !! Scaling factor for toroidal field
    real(FP), protected :: scale_ip = 1.0_FP
    !! Scaling factor for current
    real(FP), dimension(nmax_amplitudes), protected :: amplitudes = 0.0_FP
    !! Amplitudes for flare_equi3d
    character(len=3), protected :: flare_units_system = 'SI'
    !! Units used internally in flare equilibria, available
    !! - 'SI' (default)
    !! - 'CGS'
    real(FP), protected :: Rmin = 0.0_FP
    !! Coilset min R value.
    real(FP), protected :: Rmax = 1.0_FP
    !! Coilset max R value
    real(FP), protected :: Zmin = -1.0_FP
    !! Coilset min Z value
    real(FP), protected :: Zmax = 1.0_FP
    !! Coilset max Z value
    real(FP), protected :: rhomin = 0.0_FP
    !! rho value of inner limiting flux surface
    real(FP), protected :: rhomax = 1.0_FP
    !! rho value of outer limiting flux surface
    character(len=lmax_char), protected :: flare_equi3d_path_mgrid
    !! Path to mgrid file for equi3d
    character(len=lmax_char), protected :: flare_equi3d_path_axis
    !! Path to magnetic axis file for equi3d
    character(len=lmax_char), protected :: flare_equi3d_path_bfield
    !! Path to .bfield coilset data for equi3d coilset
    character(len=lmax_char), protected :: flare_geqdsk_path
    !! Path to eqdsk file for GEQDSK
    character(len=lmax_char), protected :: flare_equi3d_path_hint
    !! Path to .bfield hint file for equi3d hint
    integer, protected :: group = -1
    !! Chooses hint vacuum (0) or plasma (-1) equilibrium
    real(FP), protected :: bmax = 5.0_FP
    !! Maximum magnetic field for bfield hint equilibrium

    ! Parameter for description of boundaries
    integer, protected :: nbnd_seg
    !! Number of boundary segments
    character(len=lmax_char), dimension(nmax_bnds), protected :: bnd_paths = ''
    !! Paths to files containing boundary segments
    integer, dimension(nmax_bnds), protected :: bnd_symms = 0
    !! Symmetry of boundary segments
    !! - 0: no symmetry (default)
    !! - 1: tilt symmetry

    ! Parameters for description of flux surfaces
    character(len=lmax_char), protected :: flare_rho_dirpath = ''
    !! Path to directory containing flux surface files
    character(len=lmax_char), protected :: flare_rho_prefix = ''
    !! File prefix of flux surface files
    integer, protected :: rho_symms = 0
    !! Symmetry of flux surfaces
    !! - 0: no symmetry (default)
    !! - 1: tilt symmetry

    namelist / equi_params_flare_type / &
        equi_flare_type, boundary_flare_type, rho_flare_type

    namelist / flare_equi3d_mgrid_params / &
        flare_equi3d_path_mgrid, flare_equi3d_path_axis, &
        amplitudes, spline_order, bspline_data, flare_units_system, &
        rhomin, rhomax

    namelist / flare_equi3d_coilset_params / &
        flare_equi3d_path_axis, flare_equi3d_path_bfield, rhomin, rhomax, &
        flare_units_system, Rmin, Rmax, Zmin, Zmax

    namelist / flare_equi3d_hint_params / &
        flare_equi3d_path_hint, flare_equi3d_path_axis, rhomin, rhomax, &
        group, bmax

    namelist / flare_geqdsk_params / &
        flare_geqdsk_path, spline_order, scale_bt, scale_ip, flare_units_system

    namelist / flare_bnds_kisslinger / &
        nbnd_seg, bnd_paths, bnd_symms

    namelist / flare_rho_kisslinger / &
        flare_rho_dirpath, flare_rho_prefix, rho_symms

    public :: read_equi_params_flare_type, write_equi_params_flare_type
    public :: read_flare_equi_params, write_flare_equi_params
    public :: read_flare_bnd_params, write_flare_bnd_params
    public :: read_flare_rho_params, write_flare_rho_params

contains

    subroutine read_equi_params_flare_type(filename)
        !! Reads the equi_params_flare_type namelist from the given file
        character(len=*), intent(in) :: filename
        !! Filepath to parameterfile

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        open(newunit=funit, file=filename, status='old', action='read', &
            iostat=io_error, iomsg=io_errmsg)
        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        read(funit, nml=equi_params_flare_type, &
             iostat=io_error, iomsg=io_errmsg)
        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        close(funit)

    end subroutine

    subroutine write_equi_params_flare_type(filename)
        !! Writes the equi_params_flare_type namelist
        character(len=*), intent(in), optional :: filename
        !! If present, filename, where to write,
        !! otherwise writes parameters to stdout

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        if (present(filename)) then
            open(newunit=funit, file=filename, action='write', &
                position='append', iostat=io_error, iomsg=io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                                  __LINE__, __FILE__)
            endif
        else
            funit = get_stdout()
        endif

        write(funit, nml=equi_params_flare_type, &
              iostat=io_error, iomsg=io_errmsg)

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        if (present(filename)) then
            close(funit)
        endif

    end subroutine

    subroutine read_flare_equi_params(filename)
        !! Reads parameters for a specific equilibrium type of flare
        character(len=*), intent(in) :: filename
        !! Filepath to parameterfile

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        open(newunit=funit, file=filename, status='old', action='read', &
            iostat=io_error, iomsg=io_errmsg)
        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        select case(equi_flare_type)
            case(FLARE_TYPE_EQUI3D_MGRID)
                read(funit, nml=flare_equi3d_mgrid_params, &
                     iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI3D_COILSET)
                read(funit, nml=flare_equi3d_coilset_params, &
                     iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI3D_HINT)
                read(funit, nml=flare_equi3d_hint_params, &
                     iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI2D_GEQDSK)
                read(funit, nml=flare_geqdsk_params, &
                     iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('equi_type_flare= '// equi_flare_type // &
                                  ' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        close(funit)

    end subroutine

    subroutine write_flare_equi_params(filename)
        !! Writes parameters for a specific equilibrium type of flare
        character(len=*), intent(in), optional :: filename
        !! If present, filename, where to write,
        !! otherwise writes parameters to stdout

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        if (present(filename)) then
            open(newunit=funit, file=filename, action='write', &
                position='append', iostat=io_error, iomsg=io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                                  __LINE__, __FILE__)
            endif
        else
            funit = get_stdout()
        endif

        select case(equi_flare_type)
            case(FLARE_TYPE_EQUI3D_MGRID)
                write(funit, nml=flare_equi3d_mgrid_params, &
                      iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI3D_COILSET)
                write(funit, nml=flare_equi3d_coilset_params, &
                      iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI3D_HINT)
                write(funit, nml=flare_equi3d_hint_params, &
                      iostat=io_error, iomsg=io_errmsg)
            case(FLARE_TYPE_EQUI2D_GEQDSK)
                write(funit, nml=flare_geqdsk_params, &
                      iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('equi_type_flare= '// equi_flare_type // &
                                  ' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        if (present(filename)) then
            close(funit)
        endif

    end subroutine

    subroutine read_flare_bnd_params(filename)
        !! Reads parameters for describing boundary segments
        character(len=*), intent(in) :: filename
        !! Filepath to parameterfile

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        open(newunit=funit, file=filename, status='old', action='read', &
             iostat=io_error, iomsg=io_errmsg)
        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        select case(boundary_flare_type)
            case('none')
                ! Do nothing
                io_error = 0
            case('kisslinger')
                read(funit, nml=flare_bnds_kisslinger, &
                     iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('boundary_flare_type= '// &
                                  boundary_flare_type //' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        close(funit)

    end subroutine

    subroutine write_flare_bnd_params(filename)
        !! Writes parameters for describing boundary segments
        character(len=*), intent(in), optional :: filename
        !! If present, filename, where to write,
        !! otherwise writes parameters to stdout

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        if (present(filename)) then
            open(newunit=funit, file=filename, action='write', &
                position='append', iostat=io_error, iomsg=io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                                  __LINE__, __FILE__)
            endif
        else
            funit = get_stdout()
        endif

        select case(boundary_flare_type)
            case('none')
                ! Do nothing
                io_error = 0
            case('kisslinger')
                write(funit, nml=flare_bnds_kisslinger, &
                      iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('boundary_flare_type= ' // &
                                  boundary_flare_type // ' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        if (present(filename)) then
            close(funit)
        endif

    end subroutine

    subroutine read_flare_rho_params(filename)
        !! Reads parameters for describing flux surfaces *rho)
        character(len=*), intent(in) :: filename
        !! Filepath to parameterfile

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        open(newunit=funit, file=filename, status='old', action='read', &
             iostat=io_error, iomsg=io_errmsg)
        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        select case(rho_flare_type)
            case('none')
                ! Do nothing
                io_error = 0
            case('kisslinger')
                read(funit, nml=flare_rho_kisslinger, &
                     iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('rho_flare_type= '// rho_flare_type // &
                                  ' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        close(funit)

    end subroutine

    subroutine write_flare_rho_params(filename)
        !! Writes parameters for describing rho segments
        character(len=*), intent(in), optional :: filename
        !! If present, filename, where to write,
        !! otherwise writes parameters to stdout

        integer :: funit, io_error
        character(len=256) :: io_errmsg

        if (present(filename)) then
            open(newunit=funit, file=filename, action='write', &
                position='append', iostat=io_error, iomsg=io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                                  __LINE__, __FILE__)
            endif
        else
            funit = get_stdout()
        endif

        select case(rho_flare_type)
            case('none')
                ! Do nothing
                io_error = 0
            case('kisslinger')
                write(funit, nml=flare_rho_kisslinger, &
                      iostat=io_error, iomsg=io_errmsg)
            case default
                call handle_error('rho_flare_type= '// rho_flare_type // &
                                  ' not valid', &
                                  PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__)
        end select

        if (io_error /= 0) then
            call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, &
                              __LINE__, __FILE__)
        endif

        if (present(filename)) then
            close(funit)
        endif

    end subroutine

end module