Introspection into floats/NV

[LanX]

Hi

I was surprised to find out that 27**(1/3) == int 27**(1/3) and tried to introspect what the double float in the NV contains. I spend hours fiddling with printf , pack , unpack and (worst) the bit-order on intel processors to create this

(compare https://bartaz.github.io/ieee754-visualization/ )

$ perl
printf "%5s %.13a\n      %s\n", $_, eval($_), unpack("B*",pack("F>",ev
+al $_)) for qw(1 -1 2 -2 1/3 -1/3)
__END__
    1 0x1.0000000000000p+0
      0011111111110000000000000000000000000000000000000000000000000000
   -1 -0x1.0000000000000p+0
      1011111111110000000000000000000000000000000000000000000000000000
    2 0x1.0000000000000p+1
      0100000000000000000000000000000000000000000000000000000000000000
   -2 -0x1.0000000000000p+1
      1100000000000000000000000000000000000000000000000000000000000000
  1/3 0x1.5555555555555p-2
      0011111111010101010101010101010101010101010101010101010101010101
 -1/3 -0x1.5555555555555p-2
      1011111111010101010101010101010101010101010101010101010101010101

2 questions:

• I wasn't able yet to make unpack "B1B11B52" work to distinguish between sign, exp and significant (mantissa), because of some byte alignment problems bits were missing at the end
• Does this work on other builds and platforms too? ubuntu 64bit perl intel processor

  Platform:
    osname=linux
    osvers=6.1.0
    archname=x86_64-linux-gnu-thread-multi
  Compiler:
    cc='x86_64-linux-gnu-gcc'

update

tested on my ARM mobile which has another endianess, and got the same result

[ikegami]

I've used this:

sub double_to_bin {
   my $n = shift;

   my ( $s, $e, $m ) =
      unpack "a1 a11 a52",
         unpack "B64",
            pack "d>", $n;

   $e = oct( "0b$e" );
   if ( $e == 0x7FF ) {
      if ( $m =~ /^0+\z/ ) {
         return $s ? "-inf" : "+inf" ;
      } else {
         return "nan";
      }
   }

   my $n_bin = ( $e ? "1" : "" ) . $m;

   $e -= 1023;
   if ( $e < 0 ) {
      $n_bin = "0." . ( "0" x ( -$e - 1 ) ) . $n_bin;
   } else {
      my $z = ( $e + 1 ) - length( $n_bin );
      if ( $z >= 0 ) {
         $n_bin .= "0" x $z;
      } else {
         substr( $n_bin, $e+1, 0, "." );
      }
   }

   $n_bin =~ s/^0+(?=[01])//;
   if ( $n_bin =~ /\./ ) {
      #$n_bin =~ s/0+\z//;
      $n_bin =~ s/\.\z//;
   }

   $n_bin = ( $s ? "-" : "" ) . $n_bin;

   return $n_bin;
}

Handles subnormals, infinities and NaNs. [Updated to handle the latter two.]

---

To print the exact value in decimal, you need up to 1074 decimal places (%.1074f) or 751 in scientific notation (%.751g) for IEEE doubles.

$ perl -e'printf "%.751g\n", unpack "d", pack "Q", 1'
4.94065645841246544176568792868221372365059802614324764425585682500675
+507270208751865299836361635992379796564695445717730926656710355939796
+398774796010781878126300713190311404527845817167848982103688718636056
+998730723050006387409153564984387312473397273169615140031715385398074
+126238565591171026658556686768187039560310624931945271591492455329305
+456544401127480129709999541931989409080416563324524757147869014726780
+159355238611550134803526493472019379026810710749170333222684475333572
+083243193609238289345836806010601150616980975307834227731832924790498
+252473077637592724787465608477820373446969953364701797267771758512566
+055119913150489110145103786273816725095583738973359899366480994116420
+5702637090279242767544565229087538682506419718265533447265625e-324

[syphilis]

To print the exact value in decimal, you need up to 1074 decimal places (%.1074f) or 751 in scientific notation (%.751g) for IEEE doubles.

https://www.exploringbinary.com/maximum-number-of-decimal-digits-in-binary-floating-point-numbers explains that the maximum number of decimal places required for a double precision value is 767.

UPDATE:

D:\>perl -le "$x = (2 ** -1021) + (2 ** -1074); $s = sprintf '%.767g',
+ $x; print $s; print length $s;"
4.4501...5625e-308
773

Cheers,
Rob

[LanX]

Thanks, the second unpack is a good workaround! :)

printf "%a" should hopefully handle all edge cases.

[ikegami]

I wrote this before %a was introduced. But %a produces scientific notation, whereas the above doesn't.

[karlgoethebier]

For curiosity and as i have a new license i asked ChatGPT:

#!/usr/bin/env perl                                                   
+                                                                     
+                                                                     
+                                                                for m
+y $val (qw(1 -1 2 -2 1/3 -1/3)) {
    my ($s, $e, $m) = double_to_bin_parts(eval $val);
    printf "%5s: sign=%s exp=%s mantissa=%s\n", $val, $s, $e, $m;
  }

for my $expr (qw(1 -1 2 -2 1/3 -1/3)) {
    my $val = eval $expr;
    my ($s, $e, $m) = double_to_bin_parts($val);
    printf "%5s %.13a\n      %s %s %s\n", $expr, $val, $s, $e, $m;
}

sub double_to_bin_parts {
    my $n = shift;
    my $bits = unpack "B64", pack "d>", $n;
    my $sign     = substr $bits, 0, 1;
    my $exp      = substr $bits, 1, 11;
    my $mantissa = substr $bits, 12, 52;
    return ($sign, $exp, $mantissa);
}

__END__

1: sign=0 exp=01111111111 mantissa=00000000000000000000000000000000000
+00000000000000000
-1: sign=1 exp=01111111111 mantissa=0000000000000000000000000000000000
+000000000000000000
 2: sign=0 exp=10000000000 mantissa=0000000000000000000000000000000000
+000000000000000000
-2: sign=1 exp=10000000000 mantissa=0000000000000000000000000000000000
+000000000000000000
1/3: sign=0 exp=01111111101 mantissa=010101010101010101010101010101010
+1010101010101010101
-1/3: sign=1 exp=01111111101 mantissa=01010101010101010101010101010101
+01010101010101010101
 1 0x1.0000000000000p+0
   0 01111111111 0000000000000000000000000000000000000000000000000000
-1 -0x1.0000000000000p+0
   1 01111111111 0000000000000000000000000000000000000000000000000000
 2 0x1.0000000000000p+1
   0 10000000000 0000000000000000000000000000000000000000000000000000
-2 -0x1.0000000000000p+1
   1 10000000000 0000000000000000000000000000000000000000000000000000
1/3 0x1.5555555555555p-2
   0 01111111101 0101010101010101010101010101010101010101010101010101
-1/3 -0x1.5555555555555p-2
   1 01111111101 0101010101010101010101010101010101010101010101010101

[LanX]

What was your prompt?

This unpack "B64", pack "d>", $n; looks like taken straight from ikegamis post or probably the docs.

Update

See also https://stackoverflow.com/questions/7669010/binary-representation-of-float-to-decimal-conversions-in-perl

[karlgoethebier]

I “discussed” your one liner and ikegamis sub and asked a question on pack/unpack similar to yours. Also impressive: A version of ikegamis sub with POD. The oracle liked the sub. A while ago i asked something about FFI::Platypus with Rust - which can be very tricky. The oracle provided a solution that i never would have found in the (good) docs.