Fetch all condition numbers to proc0 and report on worrying values as required.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(fieldmat_type), | intent(in) | :: | self |
subroutine fm_get_condition_numbers(self)
use file_utils, only: error_unit
use gs2_time, only: code_dt
use mp, only: proc0, sum_reduce
use unit_tests, only: debug_message
implicit none
type ky_condition_type
real, dimension(:), allocatable :: conditions
end type ky_condition_type
class(fieldmat_type), intent(in) :: self
integer, parameter :: verb = 1
type(ky_condition_type), dimension(:), allocatable :: ky_conditions
integer :: ik, is
character(len = 100) :: condition_number_message
real :: condition_number
allocate(ky_conditions(self%naky))
! Note that here we communciate over comm_world. It is likely that we
! could instead reduce individual `ky_conditions(:)%conditions` arrays
! amongst just the supercell heads of the ky block and then reduce
! across the ky heads to proc0. Currently we don't worry about this as
! we just do this once per initialisation, but this could become expensive
! when running at large scale?
do ik = 1, self%naky
allocate(ky_conditions(ik)%conditions(self%kyb(ik)%nsupercell))
ky_conditions(ik)%conditions = 0.0
do is = 1, self%kyb(ik)%nsupercell
if (self%kyb(ik)%supercells(is)%is_head) then
ky_conditions(ik)%conditions(is) = self%kyb(ik)%supercells(is)%condition_number
end if
end do
call sum_reduce(ky_conditions(ik)%conditions, 0)
end do
if (proc0) then
do ik = 1, self%naky
do is = 1, self%kyb(ik)%nsupercell
condition_number = ky_conditions(ik)%conditions(is)
write(condition_number_message,'("Condition number for response with ik = ",I0,", is = ",I0," and dt = ",E14.6E3," is ",E14.6E3)') ik, is, code_dt, condition_number
call debug_message(verb, "fields_local::fm_get_condition_numbers "//trim(condition_number_message))
if (condition_number > condition_number_warning_limit) then
write(error_unit(), '("Warning: ",A)') trim(condition_number_message)
end if
! In debug mode always display the condition number message to screen
if (debug) then
write(*,'(A)') trim(condition_number_message)
end if
end do
end do
end if
do ik = 1, self%naky
deallocate(ky_conditions(ik)%conditions)
end do
deallocate(ky_conditions)
end subroutine fm_get_condition_numbers