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 => TYPE_EQUI3D, & 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 !! Flare equilibrium type, available: !! - FLARE_TYPE_EQUI3D = 'equi3d' (default) !! - 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 :: 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_geqdsk_path !! Path to eqdsk file for GEQDSK ! 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_params / & flare_equi3d_path_mgrid, flare_equi3d_path_axis, & amplitudes, spline_order, bspline_data, flare_units_system, & rhomin, rhomax 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) read(funit, nml=flare_equi3d_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) write(funit, nml=flare_equi3d_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