FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(gs2_program_state_type), | intent(inout) | :: | state |
subroutine measure_timestep(state)
use gs2_main, only: gs2_program_state_type
use gs2_main, only: initialize_gs2
use gs2_main, only: initialize_equations
use gs2_main, only: initialize_diagnostics
use gs2_main, only: evolve_equations
use gs2_main, only: run_eigensolver
use gs2_main, only: finalize_diagnostics
use gs2_main, only: finalize_equations
use gs2_main, only: finalize_gs2
use overrides, only: optimisations_overrides_type
use optimisation_config, only: optimisation_results_type
use mp, only: proc0, broadcast, mp_abort
implicit none
type(gs2_program_state_type), intent(inout) :: state
type(optimisations_overrides_type), &
dimension(:), allocatable :: sorted_opts_temp
type(optimisation_results_type), &
dimension(:), allocatable :: sorted_res_temp
integer :: i,n, iresult
real :: t, cost
logical :: completed_steps
completed_steps = .true.
call initialize_gs2(state)
call initialize_equations(state)
call initialize_diagnostics(state)
call evolve_equations(state, state%optim%nstep_measure)
if (state%included .and. state%istep_end .ne. state%optim%nstep_measure) then
completed_steps = .false.
write(*,*) 'istep_end', state%istep_end, state%optim%nstep_measure
end if
call finalize_diagnostics(state)
call finalize_equations(state)
call finalize_gs2(state)
if (.not. completed_steps) &
call mp_abort('Optimisation has failed because gs2 is not completing &
& the required number of steps. It may be hitting a convergence &
& criterion, or a time limit, or it may be a numerical instability. &
& Check exit_when_converged, avail_cpu_time, omegatol, omegatinst.',&
.true.)
if (state%optim%measure_all) then
t = state%timers%advance(1)/real(state%optim%nstep_measure)
else
t = state%timers%timestep(1)/real(state%optim%nstep_measure)
!t = state%timers%timestep(1)
endif
cost = t*real(state%nproc_actual)
call broadcast(t)
call broadcast(cost)
state%optim%results%nproc = state%nproc_actual
call broadcast(state%optim%results%nproc)
state%optim%results%time = t
state%optim%results%cost = cost
!state%optim%results%cost = t
if (t .lt. state%optim%results%optimal_time .or. &
state%optim%results%optimal_time .lt. 0.0) then
state%optim%results%optimal_time = t
state%optim%results%optimal = .true.
end if
if (cost .lt. state%optim%results%optimal_cost .or. &
state%optim%results%optimal_cost .lt. 0.0) then
state%optim%results%optimal_cost = cost
!if (proc0) write(*,*) 'optimal_cost', state%optim%results%optimal_cost, cost
end if
n = size(state%optim%sorted_results)
!write (*,*) 'size2', size(state%optim%sorted_results)
allocate(sorted_opts_temp(n), sorted_res_temp(n))
if (n>0) then
do i = 1,n
sorted_opts_temp(i) = state%optim%sorted_optimisations(i)
sorted_res_temp(i) = state%optim%sorted_results(i)
end do
end if
deallocate(state%optim%sorted_optimisations)
deallocate(state%optim%sorted_results)
allocate(state%optim%sorted_optimisations(n+1))
allocate(state%optim%sorted_results(n+1))
i=1
do
if (i>n) exit
if (sorted_res_temp(i)%time > t) exit
state%optim%sorted_optimisations(i) = sorted_opts_temp(i)
state%optim%sorted_results(i) = sorted_res_temp(i)
i = i+1
end do
state%optim%sorted_optimisations(i) = state%init%opt_ov
state%optim%sorted_results(i) = state%optim%results
iresult = i
!if (proc0) write(*,*) 'cost', state%optim%sorted_results(i)%cost, cost
i = i + 1
do while (i < n+2)
state%optim%sorted_optimisations(i) = sorted_opts_temp(i-1)
state%optim%sorted_results(i) = sorted_res_temp(i-1)
i = i + 1
end do
do i = 1,size(state%optim%sorted_results)
state%optim%sorted_results(i)%optimal_cost = &
state%optim%results%optimal_cost
state%optim%sorted_results(i)%optimal_time = &
state%optim%results%optimal_time
state%optim%sorted_results(i)%efficiency = &
state%optim%sorted_results(i)%optimal_cost / &
state%optim%sorted_results(i)%cost
end do
if(proc0) call write_summary(6,&
state%optim%sorted_results(iresult), &
state%optim%sorted_optimisations(iresult))
deallocate(sorted_opts_temp, sorted_res_temp)
end subroutine measure_timestep