module params_equi_circular_toroidal_m !! Module for parameter reading for the circular toroidal 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 implicit none private real(FP) :: rhomin = 0.0_FP ! Minimum of rho within the simulation domain real(FP) :: rhomax = 0.3_FP ! Maximum of rho within the simulation domain real(FP) :: q_0 = 1.0_FP ! Safety factor value on axis (rho=0) real(FP) :: q_quad_param = 12.5_FP ! Quadratic parameter of the q-profile of the form ! q(rho) = q_0 + q_quad_param * rho**2 real(FP) :: hel_amp = 0.0_FP ! Amplitude of the helical perturbation integer :: hel_m = 3 ! Poloidal mode number of the helical perturbation integer :: hel_n = 2 ! Toroidal mode number of the helical perturbation real(FP) :: hel_rho = 0.2_FP ! Resonant rho position of the helical perturbation, default corresponds ! to m/n=3/2 perturbation with default q_0 and q_quad_param real(FP) :: hel_sigma = 0.01_FP ! Variance of the Gaussian radial envelope of the helical perturbation logical :: kiss_boundary_on = .false. ! True if some boundary description in Kisslinger format is applied integer, parameter :: lmax_char = 512 ! Maximum lengths of strings in parameters integer, parameter :: nmax_bnds = 50 ! Maximum number of boundary segments ! Parameters for description of boundaries integer, public :: nbnd_seg !! Number of boundary segments character(len=lmax_char), dimension(nmax_bnds), public :: bnd_paths = '' !! Paths to files containing boundary segments integer, dimension(nmax_bnds), public :: bnd_symms = 0 !! Symmetry of boundary segments !! - 0: no symmetry (default) !! - 1: tilt symmetry public :: read_params_circular_toroidal public :: write_params_circular_toroidal public :: read_circtor_bnd_params, write_circtor_bnd_params public :: get_circular_toroidal_rhomin public :: get_circular_toroidal_rhomax public :: get_circular_toroidal_q_0 public :: get_circular_toroidal_q_quad_param public :: get_circular_toroidal_hel_amp public :: get_circular_toroidal_hel_m public :: get_circular_toroidal_hel_n public :: get_circular_toroidal_hel_rho public :: get_circular_toroidal_hel_sigma public :: get_circular_toroidal_kiss_boundary_on namelist / equi_circular_toroidal_params / & rhomin, rhomax, q_0, q_quad_param, & hel_amp, hel_m, hel_n, hel_rho, hel_sigma, & kiss_boundary_on namelist / circtor_bnd_kisslinger / nbnd_seg, bnd_paths, bnd_symms contains pure real(FP) function get_circular_toroidal_rhomin() get_circular_toroidal_rhomin = rhomin end function pure real(FP) function get_circular_toroidal_rhomax() get_circular_toroidal_rhomax = rhomax end function pure real(FP) function get_circular_toroidal_q_0() get_circular_toroidal_q_0 = q_0 end function pure real(FP) function get_circular_toroidal_q_quad_param() get_circular_toroidal_q_quad_param = q_quad_param end function pure real(FP) function get_circular_toroidal_hel_amp() get_circular_toroidal_hel_amp = hel_amp end function pure integer function get_circular_toroidal_hel_m() get_circular_toroidal_hel_m = hel_m end function pure integer function get_circular_toroidal_hel_n() get_circular_toroidal_hel_n = hel_n end function pure real(FP) function get_circular_toroidal_hel_rho() get_circular_toroidal_hel_rho = hel_rho end function pure real(FP) function get_circular_toroidal_hel_sigma() get_circular_toroidal_hel_sigma = hel_sigma end function pure logical function get_circular_toroidal_kiss_boundary_on() get_circular_toroidal_kiss_boundary_on = kiss_boundary_on end function subroutine read_params_circular_toroidal(filename) !! Reads the equi_circular_toroidal namelist from the given file character(len=*), intent(in) :: filename integer :: funit, io_error open(newunit=funit, file=filename, status='old', action='read', & iostat=io_error) if (io_error /= 0) then call handle_error("Opening parameter file "//filename//" failed!", & PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__, & additional_info= error_info_t("Iostat:", [io_error])) end if read(funit, nml=equi_circular_toroidal_params, iostat=io_error) ! END and EOR are returned if namelist can not be found. We treat ! these cases as non fatal and write a notification. if (io_error == IOSTAT_END .or. io_error == IOSTAT_EOR) then write(get_stdout(), *) & "Info: Namelist equi_circular_toroidal_params not present." elseif (io_error /= 0) then call handle_error( & "Reading circular_toroidal equi parameters failed!", & PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__, & additional_info=error_info_t("Iostat:", [io_error])) end if close(funit) end subroutine subroutine write_params_circular_toroidal(filename) !! Writes the equi_circular_toroidal namelist into the given filename character(len=*), intent(in) :: filename integer :: funit, io_error open(newunit=funit, file=filename, action='write', position='append', & iostat=io_error) if (io_error /= 0) then call handle_error("Opening parameter file "//filename//" failed!", & PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__, & additional_info= error_info_t("Iostat:", [io_error])) end if write(funit, nml=equi_circular_toroidal_params, iostat=io_error) if (io_error /= 0) then call handle_error( & "Writing circular_toroidal equi parameters failed!", & PARALLAX_ERR_PARAMETERS, __LINE__, __FILE__, & additional_info=error_info_t("Iostat:", [io_error])) end if close(funit) end subroutine subroutine read_circtor_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 if(kiss_boundary_on) then read(funit, nml=circtor_bnd_kisslinger, & iostat=io_error, iomsg=io_errmsg) endif if (io_error /= 0) then call handle_error(io_errmsg, PARALLAX_ERR_PARAMETERS, & __LINE__, __FILE__) endif close(funit) end subroutine subroutine write_circtor_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 if(kiss_boundary_on) then write(funit, nml=circtor_bnd_kisslinger, & iostat=io_error, iomsg=io_errmsg) endif 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