如何在Fortran中执行整数log2()



一个显而易见的方法是:

integer function log2_i(val) result(res)
    implicit none
    integer, intent(IN) :: val
    if (val<0) then
        print *, "ERROR in log2_i(): val cannot be negative."
    else if (val==0) then
        print *, "ERROR in log2_i(): todo: return the integer equivalent of (-inf)..."
    else if (val==1) then
        res = 0
    else
        res = FLOOR( LOG(val+0.0d0)/LOG(2.0d0) )
    endif
end function log2_i

有没有更好的方法来使用Fortran的位移位运算符?

这个问题几乎相同,但使用了无符号整数。不幸的是,符号位将禁止使用相同的算法。

正如@harold所提到的,这应该不是问题:因为对数只为正数定义,所以符号位总是零(请参阅相应的维基百科文章)。因此,链接答案中的算法可以直接移植到Fortran(2008标准):

module mod_test
contains
  function ilog2_b(val ) result( res )
    integer, intent(in) :: val
    integer             :: res
    integer             :: tmp
    res = -1 
    ! Negativ values not allowed
    if ( val < 1 ) return
    tmp = val
    do while (tmp > 0)
      res = res + 1
      tmp = shiftr( tmp, 1 )
    enddo
  end function
end module
program test
  use mod_test
  print *,'Bitshifting: ', ilog2_b(12345)
  print *,'Formula:     ', floor( log(real(12345) ) / log(2.) )
end program

这是一个基于@agentp:建议的Fortran 95内在BTEST的解决方案

module mod_test
contains
  function ilog2_b2(val ) result( res )
    integer, intent(in) :: val
    integer             :: res
    integer             :: i
    res = -1 
    ! Negativ values not allowed
    if ( val < 1 ) return
    do i=bit_size(val)-1,0,-1
      if ( btest(val, i) ) then
        res = i
        return
      endif
    enddo
  end function
end module
program test
  use mod_test
  print *,'Testing bits:', ilog2_b2(123456)
  print *,'Formula:     ', floor( log(real(123456) ) / log(2.) )
end program

感谢@IanH为我指出bit_size。。。如果您的编译器支持shiftr,我会使用第一个变体。


@IanH提到了另一种使用leadz的方法,这是Fortran 2008的一个特性:

module mod_test
contains
  function ilog2_b3(val ) result( res )
    integer, intent(in) :: val
    integer             :: res
    res = -1 
    ! Negativ values not allowed
    if ( val < 1 ) return
    res = bit_size(val)-leadz(val)-1
  end function
end module
program test
  use mod_test
  print *,'Testing bits:', ilog2_b3(123456)
  print *,'Formula:     ', floor( log(real(123456) ) / log(2.) )
end program

最新更新