!> A module which allows multiple values of certain parameters to be !> considered within a single simulation. !> !> In its simplest form, it starts with a given value of the parameter !> in question, runs for a given number of timesteps, and then changes !> the parameter by a given increment, until a lower limit is reached. !> In a more advanced scenario, the parameter scan continues until a !> different condition is satisfied: for example, zero heat flux. !> In a third scenario, the parameter is varied using a root finding !> algorithm. In addition, the condition for changing the parameter !> may be changed from a simple number of time steps to, for example, !> reaching a saturated state. module parameter_scan use abstract_config, only: abstract_config_type, CONFIG_MAX_NAME_LEN implicit none private public :: init_parameter_scan, finish_parameter_scan, update_scan_parameter_value public :: allocate_target_arrays, deallocate_target_arrays, target_parameter_momflux_tot public :: target_parameter_switch, scan_type_switch, scan_type_none public :: target_parameter_hflux_tot, target_parameter_phi2_tot public :: scan_restarted, scanning public :: parameter_scan_config_type, set_parameter_scan_config, get_parameter_scan_config logical :: scan_restarted integer :: scan_parameter_switch integer, parameter :: scan_parameter_tprim = 1, scan_parameter_g_exb = 2 integer :: target_parameter_switch, scan_type_switch, increment_condition_switch integer, parameter :: scan_type_none = 1, & scan_type_range = 2, scan_type_target =3, & scan_type_root_finding = 4 integer, parameter :: target_parameter_hflux_tot = 1, & target_parameter_momflux_tot = 2, & target_parameter_phi2_tot = 3 integer, parameter :: increment_condition_n_timesteps = 1, & increment_condition_delta_t = 2, & increment_condition_saturated = 3 real :: par_start, par_end, par_inc integer :: nstep_init, nstep_inc, scan_spec real :: delta_t_init, delta_t_inc, target_val integer :: scan_output_file real :: current_target_value = 0.0 logical :: scanning = .false. logical :: initialized = .false. !> Used to represent the input configuration of parameter_scan type, extends(abstract_config_type) :: parameter_scan_config_type ! namelist : parameter_scan_knobs ! indexed : false !> When the increment condition is `"delta_t"`, the parameter will be changed !> every time `delta_t_inc` time has elapsed. real :: delta_t_inc = 0.0 !> When the increment condition is `"delta_t"`, the parameter will not be !> changed until `delta_t_init` time has elapsed from the beginning of the !> simulation. Note, that if the simulation is restarted, this parameter !> will measure from beginning of original simulation. real :: delta_t_init = 0.0 !> Specifies the condition for incrementing the parameter. Possible values are: !> !> - "n_timesteps": change the parameter after a given number of time steps !> - "delta_t": change the parameter after an elapsed time !> - "saturated": change the parameter after the simulation has reached a !> saturated state (determined using the target parameter) at the current !> value of the parameter character(len = 20) :: inc_con = 'delta_t' !> When the increment condition is `'n_timesteps'`, the parameter will be !> changed every `nstep_inc`. integer :: nstep_inc = 0 !> When the increment condition is 'n_timesteps' or 'saturated', the !> parameter will not be changed until nstep_init have elapsed from the !> beginning of the simulation. Note that if the simulation is restarted, !> this parameter will measure from the restart. integer :: nstep_init = 0 !> If the scan is being run in 'range' mode, specifies the value of the !> parameter that will be reached. real :: par_end = 0.0 !> If the parameter scan is being run in 'range' or 'target' modes, !> specifies the amount by which the parameter is varied at one go. real :: par_inc = 0.0 !> Specifies the starting value for the parameter scan. real :: par_start = 0.0 !> Specify the parameter to be varied. If the parameter pertains to a !> species, the scan_spec must be specified as well. character(len = 20) :: scan_par = 'tprim' !> Specify the parameter to be varied. If the parameter pertains to a !> species, the scan_spec must be specified as well. logical :: scan_restarted = .false. !> When parameter pertains to a species, specifies the index of the !> species. integer :: scan_spec = 1 !> Specifies the way that the parameter scan is conducted. Possible values are: !> !> - 'none': do not conduct a parameter scan (default) !> - 'range': vary parameter in constant increments between 2 values: !> par_start and par_end. The step size is given by par_inc. !> - 'target': start with the parameter at par_start, and then change the !> parameter by par_inc until the target parameter has reached the target !> value !> - 'root_finding': the same as target, but the increment is changed !> intelligently using a Newton-like method. character(len = 20) :: scan_type = 'none' !> If the scan is being run in 'target' or 'root_finding' mode, specifies !> the target parameter. Possible values are: !> !> - 'hflux_tot' !> - 'momflux_tot' !> - 'phi2_tot' character(len = 20) :: target_par = 'hflux_tot' !> real :: target_val = 0.0 contains procedure, public :: read => read_parameter_scan_config procedure, public :: write => write_parameter_scan_config procedure, public :: reset => reset_parameter_scan_config procedure, public :: broadcast => broadcast_parameter_scan_config procedure, public, nopass :: get_default_name => get_default_name_parameter_scan_config procedure, public, nopass :: get_default_requires_index => get_default_requires_index_parameter_scan_config end type parameter_scan_config_type type(parameter_scan_config_type) :: parameter_scan_config contains !> FIXME : Add documentation subroutine init_parameter_scan(parameter_scan_config_in) use gs2_save, only: restore_current_scan_parameter_value use init_g, only: init_init_g use file_utils, only: open_output_file use run_parameters, only: use_old_diagnostics use mp, only: proc0, mp_abort use parameter_scan_arrays, only: write_scan_parameter, current_scan_parameter_value implicit none type(parameter_scan_config_type), intent(in), optional :: parameter_scan_config_in if (initialized) return initialized = .true. call read_parameters(parameter_scan_config_in) scanning = .true. write_scan_parameter = .true. ! Handle different cases and incompatibilities select case (scan_type_switch) case (scan_type_none) scanning = .false. write_scan_parameter = .false. return case (scan_type_target) if (.not. use_old_diagnostics) then call mp_abort('Error: Target parameter scans only supported with old_diagnostics') end if case (scan_type_root_finding) call mp_abort('Error: Root finding parameter scans not supported yet') end select if (proc0) call open_output_file(scan_output_file, ".par_scan") ! scan_restarted is set manually if (scan_restarted) then call init_init_g ! To get restart file name call restore_current_scan_parameter_value(current_scan_parameter_value) end if end subroutine init_parameter_scan !> FIXME : Add documentation subroutine allocate_target_arrays(nwrite, write_fluxes) use run_parameters, only: nstep use mp, only : mp_abort use parameter_scan_arrays, only: hflux_tot, momflux_tot, phi2_tot implicit none integer, intent(in) :: nwrite logical, intent(in) :: write_fluxes allocate(hflux_tot(nstep/nwrite+1)) allocate(momflux_tot(nstep/nwrite+1)) allocate(phi2_tot(nstep/nwrite+1)) if (.not. scanning) return select case (target_parameter_switch) case (target_parameter_hflux_tot) if (.not. write_fluxes) & call mp_abort("write_fluxes must be set to true for hflux target mode") case (target_parameter_momflux_tot) if (.not. write_fluxes) & call mp_abort("write_fluxes must be set to true for momflux target mode") end select end subroutine allocate_target_arrays !> FIXME : Add documentation subroutine deallocate_target_arrays use parameter_scan_arrays, only: hflux_tot, momflux_tot, phi2_tot implicit none if (allocated(hflux_tot)) deallocate(hflux_tot) if (allocated(momflux_tot)) deallocate(momflux_tot) if (allocated(phi2_tot)) deallocate(phi2_tot) end subroutine deallocate_target_arrays !> FIXME : Add documentation subroutine finish_parameter_scan use file_utils, only: close_output_file use mp, only: proc0 implicit none call deallocate_target_arrays if (proc0) call close_output_file(scan_output_file) initialized = .false. call parameter_scan_config%reset() end subroutine finish_parameter_scan !> FIXME : Add documentation subroutine update_scan_parameter_value(istep, reset, exit) use mp, only : mp_abort, proc0 use gs2_time, only: user_time use parameter_scan_arrays, only: current_scan_parameter_value implicit none logical, intent (inout) :: exit, reset integer, intent (in) :: istep logical :: increment_condition_satisfied, target_reached increment_condition_satisfied = .false. target_reached = .false. select case (scan_type_switch) case (scan_type_none) return case (scan_type_range) if (proc0) write (scan_output_file, fmt="('time= ',e11.4,' parameter= ',e11.4)") & user_time, current_scan_parameter_value call check_increment_condition_satisfied(istep, increment_condition_satisfied) if (.not. increment_condition_satisfied) return if (par_inc == 0.0 .or. & (par_inc > 0.0 .and. current_scan_parameter_value + par_inc > par_end) .or. & (par_inc < 0.0 .and. current_scan_parameter_value + par_inc < par_end)) & then exit = .true. else call increment_scan_parameter(par_inc, reset) end if case (scan_type_target) if (proc0) write (scan_output_file, & fmt="('time= ',e11.4,' parameter= ',e11.4,' target= ',e11.4)") & user_time, current_scan_parameter_value, current_target_value call check_increment_condition_satisfied(istep, increment_condition_satisfied) if (.not. increment_condition_satisfied) return call check_target_reached(target_reached) if (target_reached) then exit = .true. else call increment_scan_parameter(par_inc, reset) end if end select if (exit .and. proc0) write (scan_output_file, *) "Parameter scan is complete...." end subroutine update_scan_parameter_value !> FIXME : Add documentation subroutine check_target_reached(target_reached) use parameter_scan_arrays, only: hflux_tot, momflux_tot, phi2_tot, nout implicit none logical, intent (out) :: target_reached logical, save :: first = .true. real, save :: last_value = 0.0 target_reached = .false. last_value = current_target_value select case (target_parameter_switch) case (target_parameter_hflux_tot) current_target_value = hflux_tot(nout) case (target_parameter_momflux_tot) current_target_value = momflux_tot(nout) case (target_parameter_phi2_tot) current_target_value = phi2_tot(nout) end select if (first) then last_value = current_target_value first = .false. return end if ! Check if last and current values bracket the target if ((last_value < target_val .and. current_target_value > target_val) .or. & (last_value > target_val .and. current_target_value < target_val)) & target_reached = .true. end subroutine check_target_reached !> FIXME : Add documentation subroutine increment_scan_parameter(increment, reset) use parameter_scan_arrays, only: current_scan_parameter_value implicit none real, intent (in) :: increment logical, intent (inout) :: reset current_scan_parameter_value = current_scan_parameter_value + increment call set_scan_parameter(reset) end subroutine increment_scan_parameter !> FIXME : Add documentation subroutine set_scan_parameter(reset) use species, only: spec use dist_fn, only: g_exb use mp, only: proc0 use parameter_scan_arrays, only: current_scan_parameter_value implicit none logical, intent (inout) :: reset select case (scan_parameter_switch) case (scan_parameter_tprim) spec(scan_spec)%tprim = current_scan_parameter_value if (proc0) write (*,*) "Set scan parameter tprim_1 to ", spec(scan_spec)%tprim reset = .true. case (scan_parameter_g_exb) g_exb = current_scan_parameter_value if (proc0) write (*,*) "Set scan parameter g_exb to ", g_exb end select end subroutine set_scan_parameter !> FIXME : Add documentation subroutine check_increment_condition_satisfied(istep, increment_condition_satisfied) use gs2_time, only: user_time use mp, only : mp_abort, proc0 implicit none logical, intent (out) :: increment_condition_satisfied integer, intent (in) :: istep integer, save :: last_timestep = 0 real, save :: last_time = 0.0 select case (increment_condition_switch) case (increment_condition_n_timesteps) last_timestep = istep ! Note this doesn't handle the case where istep - last_timestep == nstep_inc if (istep < nstep_init) then increment_condition_satisfied = .false. else if ((istep - last_timestep) > nstep_inc) then increment_condition_satisfied = .true. if (proc0) write (*,*) 'Updating parameter at ', user_time end if case (increment_condition_delta_t) last_time = user_time ! Note this doesn't handle the case where user_time - last_timestep == delta_t_inc if (user_time < delta_t_init) then increment_condition_satisfied = .false. else if ((user_time - last_time) > delta_t_inc) then increment_condition_satisfied = .true. if (proc0) write (*,*) 'Updating parameter at ', user_time end if case (increment_condition_saturated) call mp_abort("increment_condition_saturated not implemented yet!") end select end subroutine check_increment_condition_satisfied !> FIXME : Add documentation subroutine read_parameters(parameter_scan_config_in) use file_utils, only: input_unit, error_unit, input_unit_exist use text_options, only: text_option, get_option_value use parameter_scan_arrays, only: current_scan_parameter_value implicit none type(parameter_scan_config_type), intent(in), optional :: parameter_scan_config_in type (text_option), dimension (2), parameter :: scan_parameter_opts = & (/ text_option('tprim', scan_parameter_tprim), & text_option('g_exb', scan_parameter_g_exb) /) character(20) :: scan_par type (text_option), dimension (4), parameter :: scan_type_opts = & (/ text_option('none', scan_type_none), & text_option('range', scan_type_range), & text_option('target', scan_type_target), & text_option('root_finding', scan_type_root_finding) /) character(20) :: scan_type type (text_option), dimension (3), parameter :: target_parameter_opts = & (/ text_option('hflux_tot', target_parameter_hflux_tot), & text_option('momflux_tot', target_parameter_momflux_tot), & text_option('phi2_tot', target_parameter_phi2_tot) /) character(20) :: target_par type (text_option), dimension (3), parameter :: increment_condition_opts = & (/ text_option('n_timesteps', increment_condition_n_timesteps), & text_option('delta_t', increment_condition_delta_t), & text_option('saturated', increment_condition_saturated) /) character(20) :: inc_con integer :: ierr logical :: exist if (present(parameter_scan_config_in)) parameter_scan_config = parameter_scan_config_in call parameter_scan_config%init(name = 'parameter_scan_knobs', requires_index = .false.) ! Copy out internal values into module level parameters delta_t_inc = parameter_scan_config%delta_t_inc delta_t_init = parameter_scan_config%delta_t_init inc_con = parameter_scan_config%inc_con nstep_inc = parameter_scan_config%nstep_inc nstep_init = parameter_scan_config%nstep_init par_end = parameter_scan_config%par_end par_inc = parameter_scan_config%par_inc par_start = parameter_scan_config%par_start scan_par = parameter_scan_config%scan_par scan_restarted = parameter_scan_config%scan_restarted scan_spec = parameter_scan_config%scan_spec scan_type = parameter_scan_config%scan_type target_par = parameter_scan_config%target_par target_val = parameter_scan_config%target_val exist = parameter_scan_config%exist ierr = error_unit() call get_option_value & (scan_par, scan_parameter_opts, scan_parameter_switch, & ierr, "scan_par in parameter_scan_knobs",.true.) call get_option_value & (scan_type, scan_type_opts, scan_type_switch, & ierr, "scan_type in parameter_scan_knobs",.true.) call get_option_value & (target_par, target_parameter_opts, target_parameter_switch, & ierr, "target_par in parameter_scan_knobs",.true.) call get_option_value & (inc_con, increment_condition_opts, & increment_condition_switch, & ierr, "inc_con in parameter_scan_knobs",.true.) if (.not. scan_restarted) current_scan_parameter_value = par_start end subroutine read_parameters !> Set the module level config type !> Will abort if the module has already been initialised to avoid !> inconsistencies. subroutine set_parameter_scan_config(parameter_scan_config_in) use mp, only: mp_abort type(parameter_scan_config_type), intent(in), optional :: parameter_scan_config_in if (initialized) then call mp_abort("Trying to set parameter_scan_config when already initialized.", to_screen = .true.) end if if (present(parameter_scan_config_in)) then parameter_scan_config = parameter_scan_config_in end if end subroutine set_parameter_scan_config #include "parameter_scan_auto_gen.inc" end module parameter_scan