• Spacing of real numbers around zero

    From Beliavsky@21:1/5 to All on Tue Apr 19 08:56:32 2022
    Spacing(0.0_wp) is supposed to give the number adjacent to zero, but it looks like that is actually given by epsilon(x) * tiny(x) , which is smaller. The output of the program

    program main
    use iso_fortran_env, only: int32, real32
    implicit none
    integer, parameter :: ikind = int32, wp = real32
    integer (kind=ikind) :: i, nbits
    integer, parameter :: ni = 4
    real(kind=wp) :: x
    nbits = bit_size(i)
    print*,"bit_size(i) =",nbits
    print*,"epsilon(x) =",epsilon(x)
    print*,"tiny(x) =",tiny(x)
    print*,"epsilon(x) * tiny(x) =",epsilon(x)*tiny(x)
    print*,"spacing(0.0_wp) =",spacing(0.0_wp)
    print "(/,a12,a16)","i","transfer(i,x)"
    do i=0,ni
    print*,i,transfer(i,x)
    end do
    end program main

    with gfortran and ifort is

    bit_size(i) = 32
    epsilon(x) = 1.19209290E-07
    tiny(x) = 1.17549435E-38
    epsilon(x) * tiny(x) = 1.40129846E-45
    spacing(0.0_wp) = 1.17549435E-38

    i transfer(i,x)
    0 0.00000000
    1 1.40129846E-45
    2 2.80259693E-45
    3 4.20389539E-45
    4 5.60519386E-45

    and the output of the 64-bit program

    program main
    use iso_fortran_env, only: int64, real64
    implicit none
    integer, parameter :: ikind = int64, wp = real64
    integer (kind=ikind) :: i, nbits
    integer, parameter :: ni = 4
    real(kind=wp) :: x
    nbits = bit_size(i)
    print*,"bit_size(i) =",nbits
    print*,"epsilon(x) =",epsilon(x)
    print*,"tiny(x) =",tiny(x)
    print*,"epsilon(x) * tiny(x) =",epsilon(x)*tiny(x)
    print*,"spacing(0.0_wp) =",spacing(0.0_wp)
    print "(/,a21,a16)","i","transfer(i,x)"
    do i=0,ni
    print*,i,transfer(i,x)
    end do
    end program main

    is

    bit_size(i) = 64
    epsilon(x) = 2.2204460492503131E-016
    tiny(x) = 2.2250738585072014E-308
    epsilon(x) * tiny(x) = 4.9406564584124654E-324
    spacing(0.0_wp) = 2.2250738585072014E-308

    i transfer(i,x)
    0 0.0000000000000000
    1 4.9406564584124654E-324
    2 9.8813129168249309E-324
    3 1.4821969375237396E-323
    4 1.9762625833649862E-323

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From steve kargl@21:1/5 to Beliavsky on Tue Apr 19 19:10:03 2022
    Beliavsky wrote:

    Spacing(0.0_wp) is supposed to give the number adjacent to zero,
    but it looks like that is actually given by epsilon(x) * tiny(x) , which
    is smaller. The output of the program

    If your processor supports subnormal number, first nonzero positive
    closest to 0 is 2**(-p-emin).

    % cat z.f90
    program foo
    print *, 2.**(minexponent(1.) - digits(1.))
    end

    % gfcx -o z a.f90 && ./z
    1.40129846E-45

    The IEEE facilities of Fortran allow to query about subnormal numbers.

    --
    steve

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)
  • From gah4@21:1/5 to steve kargl on Tue Apr 19 12:59:30 2022
    On Tuesday, April 19, 2022 at 12:10:06 PM UTC-7, steve kargl wrote:
    Beliavsky wrote:

    Spacing(0.0_wp) is supposed to give the number adjacent to zero,
    but it looks like that is actually given by epsilon(x) * tiny(x) , which
    is smaller.

    If your processor supports subnormal number, first nonzero positive
    closest to 0 is 2**(-p-emin).

    Which I never thought was a good idea. You get the equivalent
    of about 0.04 extra exponent bits, for a lot of extra work
    and/or hardware.

    Note that the actual smallest number has zero significant bits.

    It especially complicates parallel processors, like vector processors,
    that try to do about the same amount of work for each one.

    Even more, some require software help to get the right result.

    If someone really sees the need, give it another whole exponent bit.

    --- SoupGate-Win32 v1.05
    * Origin: fsxNet Usenet Gateway (21:1/5)