continuation_scheduler.f90

Go to the documentation of this file.

 
 
module continuation_scheduler
  use num_types, only: rp
  use json_module, only: json_file
  use json_utils, only: json_get_or_default, json_get_or_lookup, json_get
  use utils, only: neko_error
  implicit none
  private
 
  public :: continuation_scheduler_t
 
  !------------------------------------------------------------------
  ! Continuation parameter
  !------------------------------------------------------------------
  type continuation_parameter_t
     character(len=:), allocatable :: name ! Parameter name
     real(rp), pointer :: target => null() ! Parameter to update
     real(rp), allocatable :: values(:) ! Values to step through
     integer :: iterations_per_value = 0 ! How many iterations per value
   contains
     procedure :: update => continuation_parameter_update
     procedure :: free => continuation_parameter_free
  end type continuation_parameter_t
 
  !------------------------------------------------------------------
  ! Continuation scheduler
  !------------------------------------------------------------------
  type continuation_scheduler_t
     type(continuation_parameter_t), allocatable :: params(:)
     integer :: default_iterations = 1
   contains
     procedure :: init => init_scheduler
     procedure :: free => free_scheduler
     procedure :: json_get_or_register
     procedure :: register_parameter
     procedure :: update
     procedure :: get_param_name
     procedure :: get_n_params
  end type continuation_scheduler_t
 
  !------------------------------------------------------------------
  ! Global continuation object
  !------------------------------------------------------------------
  type(continuation_scheduler_t), public, target :: nekotop_continuation
 
contains
 
  subroutine init_scheduler(this, json)
    class(continuation_scheduler_t), intent(inout) :: this
    type(json_file), intent(inout) :: json
    integer :: n_default
 
    call json_get_or_default(json, &
         'optimization.solver.continuation_iterations', n_default, 1)
 
    this%default_iterations = n_default
 
    ! Initialize empty parameter array
    if (allocated(this%params)) deallocate(this%params)
  end subroutine init_scheduler
 
  subroutine free_scheduler(this)
    class(continuation_scheduler_t), intent(inout) :: this
    integer :: i
 
    if (.not. allocated(this%params)) return
 
    do i = 1, size(this%params)
       call this%params(i)%free()
    end do
 
    deallocate(this%params)
  end subroutine free_scheduler
 
  subroutine json_get_or_register(this, json, name, target, value, &
       default_value)
    class(continuation_scheduler_t), intent(inout) :: this
    type(json_file), intent(inout) :: json
    character(len=*), intent(in) :: name
    real(kind=rp), target, intent(inout) :: target
    real(kind=rp), intent(out) :: value
    real(kind=rp), intent(in), optional :: default_value
    real(kind=rp), allocatable :: values(:)
    real(kind=rp) :: scalar_parameter
    logical :: found
    integer :: var_type, iter
 
    ! Inspect JSON for the parameter
    call json%info(name, found=found, var_type=var_type)
    if (.not. found) then
       if (present(default_value)) then
          value = default_value
       else
          call neko_error(trim(name)//" not found and no default provided")
       end if
    else if (var_type == 6) then ! scalar
       if (present(default_value)) then
          call json_get_or_default(json, name, scalar_parameter, default_value)
       else
          call json_get(json, name, scalar_parameter)
       end if
       value = scalar_parameter
    else if (var_type == 3) then ! array
       call json_get_or_lookup(json, name, values)
       value = values(1)
       ! Read optional iterations per parameter step
       call json_get_or_default(json, trim(name)//'_iterations', iter, &
            this%default_iterations)
       ! Register the parameter for continuation
       call this%register_parameter(name, target, values, iter)
    else
       call neko_error(trim(name)//" can only be real variable or real array!")
    end if
 
  end subroutine json_get_or_register
 
  subroutine register_parameter(this, name, target, values, iterations)
    class(continuation_scheduler_t), intent(inout) :: this
    character(len=*), intent(in) :: name
    real(rp), target, intent(inout) :: target
    real(rp), intent(in) :: values(:)
    integer, optional, intent(in) :: iterations
 
    integer :: n, old_size
    type(continuation_parameter_t), allocatable :: tmp(:)
    integer :: iter_val
 
    ! Determine iterations per value
    if (present(iterations)) then
       iter_val = iterations
    else
       iter_val = this%default_iterations
    end if
 
    ! Append new parameter
    if (.not. allocated(this%params)) then
       allocate(this%params(1))
       n = 1
    else
       old_size = size(this%params)
       allocate(tmp(old_size))
       tmp = this%params
       deallocate(this%params)
       allocate(this%params(old_size + 1))
       this%params(1:old_size) = tmp
       deallocate(tmp)
       n = old_size + 1
    end if
 
    this%params(n)%name = name
    this%params(n)%target => target
    allocate(this%params(n)%values(size(values)))
    this%params(n)%values = values
    this%params(n)%iterations_per_value = iter_val
  end subroutine register_parameter
 
  subroutine update(this, iter)
    class(continuation_scheduler_t), intent(inout) :: this
    integer, intent(in) :: iter
    integer :: i
 
    if (.not. allocated(this%params)) return
 
    do i = 1, size(this%params)
       call this%params(i)%update(iter)
    end do
  end subroutine update
 
  character(len=32) function get_param_name(this, i)
    class(continuation_scheduler_t), intent(in) :: this
    integer, intent(in) :: i
 
    get_param_name = ''
    if (.not. allocated(this%params)) return
    if (i < 1 .or. i > size(this%params)) return
    get_param_name = this%params(i)%name
  end function get_param_name
 
  function get_n_params(this) result(n)
    class(continuation_scheduler_t), intent(in) :: this
    integer :: n
    if (allocated(this%params)) then
       n = size(this%params)
    else
       n = 0
    end if
  end function get_n_params
 
  subroutine continuation_parameter_update(this, iter)
    class(continuation_parameter_t), intent(inout) :: this
    integer, intent(in) :: iter
    integer :: idx
 
    if (.not. allocated(this%values)) return
    if (this%iterations_per_value <= 0) return
 
    idx = (iter - 1) / this%iterations_per_value + 1
    idx = min(idx, size(this%values))
 
    this%target = this%values(idx)
  end subroutine continuation_parameter_update
 
  subroutine continuation_parameter_free(this)
    class(continuation_parameter_t), intent(inout) :: this
    if (allocated(this%name)) deallocate(this%name)
    if (allocated(this%values)) deallocate(this%values)
    nullify(this%target)
  end subroutine continuation_parameter_free
 
end module continuation_scheduler