FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real, | intent(in) | :: | norm | |||
real, | intent(in), | dimension (:) | :: | wgt |
subroutine check_legendre_weights (norm, wgt)
use mp, only: mp_abort
use file_utils, only: error_unit
implicit none
real, intent (in) :: norm!, eps
real, dimension (:), intent (in) :: wgt
logical :: error=.false.
integer :: n, nh
real :: s
n = size(wgt)
error = .false.
nh = (n+1)/2
! check if weights are all positive
if (any(wgt < 0.)) then
write (error_unit(),*) 'ERROR in legendre: weights got negative'
error = .true.
end if
if (n>=2) then
! check symmetry of weights
if ( any( abs(wgt(n:n+1-n/2:-1)/wgt(1:n/2) - 1.0) > epsilon(wgt) ) ) then
write (error_unit(),*) 'WARNING in legendre: symmetry of weights broken'
error = .true.
end if
! check if weights are increasing toward center
if (n>=3) then
if (any(wgt(2:nh) <= wgt(1:nh-1))) then
write (error_unit(),*) 'ERROR in legendre: weights decreasing toward center'
error = .true.
end if
end if
end if
! check if their sum is close enough to normalized value
if (debug) then
! s = sum(wgt)
s = sum(dble(wgt)) ! ignore roundoff error arising from 8-byte summation
if (abs(norm/s-1.0) > epsilon(s)) then
write (error_unit(),*) 'WARNING in legendre: weights summation incorrect:', &
& size(wgt), s/norm-1.0
! Do not stop as it is mostly a minor issue
end if
end if
if (error) call mp_abort ('STOP in check_legendre_weights')
end subroutine check_legendre_weights