(nearly) ieee 754 floating point single precisition

[penpen]

I've apllied jebs changes to the code to the opening post, added the first arithmetic operation: add.
Tested using win xp only , i hope this works on win7, too.

penpen

[carlsomo]

I have a pure batch rounding routine for decimal numbers, but must be expressed as:

set "pi=3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086"

roundfp.bat pi 100

echo %pi%

3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170680

I can post it if it would help at all?

[penpen]

Thanks for your offer.
You may post it according to which rounding operation(s) you are using, or in other words: If it is one of rounding modes 1 to 4 of the below table.
But i probably implement these rounding modes only after i finished the math std operations (+ - * /), so it may take a while before i use it.
(That should not discourage you to post it, only that it take a while before i add it to the complete float code. I am thankful for any help.)

The reason:
In true IEE 754 is a bijective mapping between number intervalls and bit coding, and not between numbers and bit coding.
This rounding is used only to secure that the mapping of the number intervalls to strings and vice versa is distinct.
This means, if you convert a (string of a) number (a), and map it to its bit coding, and then convert it back to a (string) number (b), then a and b are in the same Intervall.
This seems to be an unuseful definition, but if you do some math with that numbers it ensures that the result is always the same for all numbers within this so defined intervall.

So other modes of roundings except of the following, defined by IEE 754 should not be implemented:

Code: Select all

:: Rounding to integers using the IEEE 754 rules
:: No Rounding mode                                +2.5 +1.5 -1.5 -2.5
:: 0  to nearest, ties to even           (default) +2.0 +2.0 -2.0 -2.0
:: 1  to nearest, ties away from zero              +3.0 +2.0 -2.0 -3.0
:: 2  toward 0                       (truncation)  +2.0 +1.0 -1.0 -2.0
:: 3  toward +inf        (rounding up, or ceiling) +3.0 +2.0 -1.0 -2.0
:: 4  toward -inf        (rounding down, or floor) +2.0 +1.0 -2.0 -3.0

Actually i am using the default rounding mode (no 0), even for that cases in which other rounding modes should be used.
Strictly speaking this is false and leads to false results when doing math on it, but this rounding mode is used in round about 75% of all cases where rounding is performed.#
So the error should be not that big.

Additionally there are some additional rounding rules for conversion from the 32 bit coding to floating point string.
I called it 'rounding to nice numbers' as this is what is its result.
Actually this is implemented in a buggy way, but the error again should be very low.

penpen

[carlsomo]

Mode 1:

Code: Select all

@echo off&goto :start
:RoundFP Float_Env precision
if defined %~1 goto :end
echo.RoundFP.CMD Float_Env precision
echo.Rounds Float_Env to desired precison (ie 12.556 2, returns 12.56 in Float_Env) 
echo.Resets Float_Env to new precision if passed as an environment variable         
echo.If numeric argument passed instead then program will echo the result           
echo.Returns errorlevel=precision of Float_env ie. 12.06 returns 2, 11.9 returns 1
echo.Only decimal numbers are processed, else returns a negative value if error
:end
EndLocal&exit /b -1
:start
SetLocal EnableDelayedExpansion
if defined %~1 (call set "float=%%%~1%%"&set precision=%2) else (call :parse %*)
if %errorlevel% gtr 0 goto :RoundFP Syntax
if not defined precision goto :RoundFP Syntax
if /i %precision% lss 0 goto :RoundFP Syntax
set/a rc=10
for /f "tokens=1-2 delims=,.0123456789-" %%a in ("%float%") do set rc=%%a
if /i "%rc%" neq "10" goto :RoundFP Syntax
if /i [%float:~0,1%] equ [-] (set sign=-&set float=%float:~1%) else (set sign=)
for /f "tokens=1-3 delims=." %%a in ("%float%") do (
  set int=%%a&set dec=%%b&set "zero=%%c"
)
if defined zero goto :RoundFP Syntax
if /i [%float:~0,1%] equ [.] (
  if defined dec goto :RoundFP Syntax
  set zero=0&set dec=%int%&set "int="
)
set "int=%zero%%int%"
call :strlen %int%&set/a intlen=!errorlevel!
set int=%int:,=%
call :strlen %int%&set/a commas=%intlen%-!errorlevel!
call :strlen %dec%&set/a len_dec=!errorlevel!
set dec2=&set "dec_add="
if /i %precision% geq %len_dec% set/a precision=%len_dec%
if /i %precision% gtr 8 (
  set "dec2=.%dec:~8%"
  set "dec=%dec:~0,8%"
  set/a precision2=%precision%-8
  set "holder=!dec2!"
  call :start holder !precision2!
  set "dec2=!holder!"
  if /i "!dec2:~0,1!" neq "." (
    set "dec_add=!dec!"
    set "dec2=!dec2:~1!"
  )
)
if /i %len_dec% equ 0 goto :quit
set/a spot=%len_dec%-%precision%
if not defined dec2 if /i %spot% leq 0 goto :quit
set/a point=%spot%-1
if not defined dec2 if /i %point% gtr 0 (call set "dec=%%dec:~0,-!point!%%")
call :RemoveFrontZeros dec
set/a numzeros=%errorlevel%
call set dec=%%dec:~!numzeros!%%
call :strlen %dec%
set/a len_dec_minus_zeros=%errorlevel%
if not defined dec2 (set/a dec+=5) else (if defined dec_add set/a dec+=1)
call :strlen %dec%
set/a new_len_dec=%errorlevel%
if not defined dec2 (
  set/a factor=%precision%+1
) else (set/a factor=8)
if /i %new_len_dec% gtr %factor% set/a int+=1&set "dec=%dec:~1%"
if not defined dec2 set/a dec/=10
set/a addzeros=%numzeros%-%new_len_dec%+%len_dec_minus_zeros%
if /i %addzeros% gtr 0 call :AddZeros dec %addzeros%
if defined dec2 set "dec2=%dec2:.=%
set "dec=%dec%%dec2%                      "
set dec=%dec: =0%
call set "dec=%%dec:~0,!precision!%%"
:quit
set "dec=.%dec%"
if /i %precision% equ 0 (set dec=)
if defined zero if /i %int% equ 0 set "int="
if /i %commas% gtr 0 call :InsertCommas int
set "float=%int%%dec%"
set "float=%sign%%float%"
if not defined %~1 EndLocal&echo.%float%&exit /b %precision%
EndLocal&call set "%~1=%float%"&exit /b %precision%

:Strlen Returns length of string in errorlevel
setlocal&set "#=%*"
if not defined # exit /b 0
set/a len=0
:loop
set/a len+=1
set "#=!#:~1!"&if not defined # endlocal&exit/b %len%
goto :loop

:parse arguments
if "%~2"=="" exit/b 1
set "parse=%*"
if not defined parse exit /b 1
for /f "tokens=1-2* delims=." %%a in ("%parse%") do (
  set front=%%a&set back=%%b&set "out=%%c"
  if defined out exit /b 1
  if /i [%parse:~0,1%] equ [.] (
    if defined back exit /b 1
    set back=!front!&set "front="
  )
  set "float=!front!.!back!"
)
for /f "tokens=1* delims= " %%i in ("%float%") do (
  set float=%%i&set "precision=%%j"
)
exit /b 0

:InsertCommas integer
if [%~1]==[] exit /b 0
SetLocal&call set "number=%%%~1%%"
call :strlen %number%
set/a length=%errorlevel%
set/a numcoms=(%length%-1)/3
for /l %%a in (1,1,%numcoms%) do (
  set/a spot=0 - %%a * 3
  set/a spot-=%%a-1
  call set "front=%%number:~0,!spot!%%"
  call set "back=%%number:~!spot!%%"
  set "number=!front!,!back!"
)
EndLocal&call set "%~1=%number%"&exit /b %numcoms%

:RemoveFrontZeros int
SetLocal
call set "first=%%%~1%%"
for /l %%a in (0,1,99) do (
  if /i "!first:~0,1!" neq "0" EndLocal&exit /b %%a
  set "first=!first:~1!
)

:AddZeros int numz
SetLocal
call set "int=%%%~1%%"
for /l %%a in (1,1,%~2) do (
  set "int=0!int!"
)
EndLocal&call set "%~1=%int%"&exit /b

[Aacini]

I had to program a string to float and vice versa conversion in pure batch: additional exes/cscript/... all not allowed. .
It was a pain to me to do such tedious fiddling ("fieser Fummelkram"), and i
want to save the trouble for others that have to do this, too, so i've posted it here.

Excuse me, perhaps this is a dumb question, but: what is the purpose of this simulation?

Batch have no numeric variables: all values are stored as strings. Every time that SET /A command execute, it converts strings into numbers before perform the arithmetic operation, and then convert the numeric result back to string in order to store it in a variable. Your method requires several conversions that uses SET /A command: from a numeric string into the hex representation of the IEEE 754 standard (stored as string). To perform arithmetic operations the parts of the number must be extracted first from the hex representation, then operated with SET /A command, and then packed back into the IEEE hex representation. This is unbelievably inefficient in terms of Batch code and execution time.

If the purpose of all this stuff is to emulate the IEEE 754 format and operations (perhaps with educative objectives) then I agree that there is no other way to achieve it in Batch. However, if the purpose is to perform floating point arithmetic operations in Batch, then there are other simpler and more efficient methods. Just my thoughts...

Antonio

[penpen]

Excuse me, perhaps this is a dumb question, but: what is the purpose of this simulation?

My opinion: No question is a dumb question, there only exists dumb answers .
I had to write it for a company because i owe a friend who works there a favour. (I haven't sold my rights on it, so it is legal to post it here.)
They have a software package where an intermediate representation are iee 754 floating point numbers (between 1 and 100000).
One of its components has failed a regularly security check there, so they have disabled it.
This component forwards messages to employees basing on this number.
So they needed a tool, that reads these numbers and calls another tool to forward these messages with this number modulo 100 and recode it to float.
As this program has to pass the security check of that company and batch script is fast to check, this was their favorite solution, so the shutdown time was minimized. (I think they will reprogram it theirselfes using C++, but their security checks seems to take some more time... .)
They have also disallowed me to use other programs like .NET, CScript, ... to imitate this software part.
(Because of that i didn't need to program any math on it, so i posted it raw as is.)
penpen

[dbenham]

That is crazy.
It must be an awfully good friend to go through all that effort.

Impressive work though.


Dave Benham

[penpen]

I have added subtraction (:floatSub) and multiplication (:floatMul).

penpen

[penpen]

Someone asked me to create a batch/JScript (hybrid) version of this conversion; here it is (parseFloatHex.bat):

Code: Select all

@if (true==false) /*
@echo off

cscript //nologo //e:JScript "%~f0" %*
goto :eof
*/
@end
/*
 * This JScript provides conversion from (and to) an ieee 754 single precision
 * string to hexadecimal representation without rounding.
 * (Unoptimized.)
 *
 * ieee 754 sp -> hex notation: float2Code (float)
 * ieee 754 sp <- hex notation: code2Float (code)
 *
 * helping function sint32 -> hex notation (needed): sint32ToHex (sint32)
 *
 * testing functions: testFloatCode (hexCode, float), test () (default called)
 *
 * You may start a demo by typing "cscript floatHexToString.js" into the
 * command shell (cmd.exe).
 *
 * Author: Ulf Schneider aka penpen
 * Version: 1.0
 */
function sint32ToHex (sint32) {
   var hexString = ((sint32 < 0) ? 0xFFFFFFFF + sint32 + 1 : sint32).toString (16).toUpperCase ();
   return "0x" + "00000000".substring (hexString.length) + hexString;
}

// float in valid single precisition notation
function float2Code (float) {
   var floatString = ("" + float).replace (/NaN/gi, "6.805646932770577e+38").replace (/inf/gi, "3.402823669209385e+38");
   if (isNaN (floatString)) return null;

   var abs_float = Math.abs (parseFloat (floatString));
   var exponent = Math.floor (Math.log (abs_float) / Math.log (2));
   var mantissa = abs_float / Math.pow (2, exponent);

   if (exponent <= -127) {
      mantissa /= Math.pow (2, -127 - exponent + 1);
      exponent += (-127 - exponent);
   }

   var s = (floatString == abs_float) ? 0 : 1;
   var e = exponent;
   var f = (exponent == -127) ? mantissa * 8388608: 0x7FFFFF & (mantissa * 8388608);

   return sint32ToHex ((s << 31) | ((e + 127) << 23) | f);
}

// code in
// "0x0" - "0xFFFFFFFF"
// "0"   - "4294967296"
// 0x0 - 0xFFFFFFFF
// 0   - 4294967296
function code2Float (code) {
   var c = parseInt ("" + ((!code) ? 0 : code));
   var s = (c >> 31) & 0x00000001;
   var e = ((c >> 23) & 0x000000FF) - 127;
   var f = c & 0x7FFFFF;
   var float = null;
   
   if (e == 128) {
      float = ((s == 1) ? "-" : "") + ((f == 0) ? "inf" : "NaN");
   } else if (e == -127) {
      var exponent = -126;
      var mantissa = f / 8388608;
      var sign = (s == 1) ? -1 : 1

      float = sign * Math.pow (2, exponent) * mantissa;
   } else if ((-127 < e) && (e < 128)) {
      var exponent = e;
      var mantissa = (0x00800000 | f) / 8388608;
      var sign = (s == 1) ? -1 : 1

      float = sign * Math.pow (2, exponent) * mantissa;
   }

   return float;
}


function testFloatCode (hexCode, float) {
   var c = parseInt ("" + ((!hexCode) ? 0 : hexCode));
   var floatString = code2Float (hexCode);
   var codeString = float2Code (float);
   var integrity = (hexCode == codeString) && (float == floatString);

   WScript.Echo ("float (" + hexCode + ") = " + "float (" + c + ") = " + floatString);
   WScript.Echo ("hexCode (" + float + ") = " + codeString);
   WScript.Echo ("integrity: " + integrity );
   WScript.Echo ("");
}

function test () {
   WScript.Echo ("code2Float (\"0x40490FDB\") = " + code2Float ("0x40490FDB"));
   WScript.Echo ("code2Float (\"1078530011\") = " + code2Float ("1078530011"));
   WScript.Echo ("code2Float (0x40490FDB)   = " + code2Float (0x40490FDB));
   WScript.Echo ("code2Float (1078530011)   = " + code2Float (1078530011));
   WScript.Echo ("");

   WScript.Echo ("float2Code (\"3.1415927410125732\") = " + float2Code ("3.1415927410125732"));
   WScript.Echo ("float2Code (3.1415927410125732)   = " + float2Code (3.1415927410125732));
   WScript.Echo ("");

   testFloatCode ("0x7FFFFFFF", "NaN");
   testFloatCode ("0xFFFFFFFF", "-NaN");
   testFloatCode ("0x7F800000", "inf");
   testFloatCode ("0xFF800000", "-inf");
   testFloatCode ("0x40490FDB", "3.1415927410125732");
   testFloatCode ("0xC0490FDB", "-3.1415927410125732");
   testFloatCode ("0x00C00000", "1.7632415262334312e-38");
   testFloatCode ("0x80C00000", "-1.7632415262334312e-38");
   testFloatCode ("0x00400000", "5.877471754111438e-39");
   testFloatCode ("0x80400000", "-5.877471754111438e-39");
   testFloatCode ("0x00000001", "1.401298464324817e-45");
   testFloatCode ("0x80000001", "-1.401298464324817e-45");
}


switch (WScript.Arguments.Unnamed.Length) {
   case 1:
      if (WScript.Arguments.Unnamed.Item (0).toLowerCase () == "test") {
         test ();
         WScript.Quit (0);
      }
      break;
   case 2:
      var command = WScript.Arguments.Unnamed.Item (0).toLowerCase ();
      var value = WScript.Arguments.Unnamed.Item (1);

      var result = (command == "parsefloat") ? code2Float (value) :
         (command == "parsehex") ? float2Code (value) :
         null;

      if (result == null) {
         WScript.Quit (1);
      } else {
         WScript.Echo ("" + result);
         WScript.Quit (0);
      }
      break;
   default:
      break;
}

WScript.Echo ("usage:");
WScript.Echo ("parseFloatHex.bat [command] [value]");
WScript.Echo ("command  test      : Calls the test function (no value).");
WScript.Echo ("         parseFloat: Converts the given hex value to float.");
WScript.Echo ("         parseHex  : Converts the given float value to hex.");
WScript.Echo ("value    A hexadecimal value or a iee 754 value is expected; see command.");
WScript.Quit (1);

You may parse values using:

Code: Select all

for /F %a in ('parseFloatHex "parseFloat" "0x80000001"') do set "result=%a"
echo %result%

for /F %a in ('parseFloatHex "parseHex" "-1.401298464324817e-45"') do set "result=%a"
echo %result%

penpen

[einstein1969]

Hi ,

I will try to calulate sin(x), cos(x), tan(x), cot(x) with this routine.

I probe with the classical maclawrin series, not optimizaed for the moment.

I need only sum, substrat, and multiply.

The division is not necessary because x/K = x * 1/K (K is know and then I can precalc 1/K)

For example:

SIN(x)=x - x^3/6 + x^5/120 - x^7/5040 + x^9/362880 - x^11/39916800 ...

Precalc 1/6, 1/120, 1/5040, 1/362880, 1/39916800

Ok, the question is:

- How many digits are necessary for this numbers for maximum precision? 1/6=?????? , 1/5040=?????



EDIT: Next step is using the Horner's_method for reducing the number of multiplications:

The previus factorized is:
SIN(x) ~ x * ( 1 - x^2 * (1/6 - x^2 * (1/120 - x^2 * (1/5040 - x^2 * (1/362880 - x^2 * 1/39916800 )))))
This is more optimizable!!!

einstein1969

[einstein1969]

Hi,

I have probed this:

Code: Select all

call :floatTestMul  "1.73" "1.73"

the result is not exact:

Code: Select all

1.73 * 1.73 = 1071476900 * 1071476900 == 0x3FDD70A4 * 0x3FDD70A4 == 0x403F8BAD == 1077906349 == 2.99290013

EDIT: The Sub seems ok , the Add no:

Code: Select all

1.73 * 1.73 = 1071476900 * 1071476900 == 0x3FDD70A4 * 0x3FDD70A4 == 0x403F8BAD == 1077906349 == 2.99290013
1.73 + 1.73 = 1071476900 + 1071476900 == 0x3FDD70A4 + 0x3FDD70A4 == 0x3FBAE148 == 1069211976 == 1.46000004
1.73 - 1.73 = 1071476900 - 1071476900 == 0x3FDD70A4 - 0x3FDD70A4 == 0x00000000 == 0 == 0.0

einstein1969

[penpen]

Ok, the question is:
- How many digits are necessary for this numbers for maximum precision? 1/6=?????? , 1/5040=?????

The mantissa of these numbers is encoded using a 24 bit number (first one is implicit), so it's relative (rounding) error is:
2^(-25) == 0,0000000298023223876953125

Depending on the exponent used the precision is something between 7 and 8 decimal digits;
if using the "1.xx e+/-yy" notation then the precision is 8 decimal digits (7th digit after the '.').

BUT actually i haven't fully implemented the rounding, so it could be, that the precision is reduced by one digit:
So you better should assume the precision is 6 to 8 decimal digits.

In addition the static drift in long term computations, described by Knuth, may occure as long as the rounding is not fully implemented.

I have probed this:

Code: Select all

call :floatTestMul  "1.73" "1.73"

the result is not exact:

Code: Select all

1.73 * 1.73 = 1071476900 * 1071476900 == 0x3FDD70A4 * 0x3FDD70A4 == 0x403F8BAD == 1077906349 == 2.99290013

The precision seems to be 7 decimal digits in this case, so the rounding error produces the digits at the end ("13").

EDIT: The Sub seems ok , the Add no:

Code: Select all

1.73 * 1.73 = 1071476900 * 1071476900 == 0x3FDD70A4 * 0x3FDD70A4 == 0x403F8BAD == 1077906349 == 2.99290013
1.73 + 1.73 = 1071476900 + 1071476900 == 0x3FDD70A4 + 0x3FDD70A4 == 0x3FBAE148 == 1069211976 == 1.46000004
1.73 - 1.73 = 1071476900 - 1071476900 == 0x3FDD70A4 - 0x3FDD70A4 == 0x00000000 == 0 == 0.0

I've overseen that if you add two 24 bit numbers the result may need 25 bits... (with rounding 26 bits are possible).
I think i've corrected this error; just add the first two (if) lines (near above the last [?i've used two?] :intToHex label):

Code: Select all

   if 0x1000000 LEQ %fR% set /A "fR>>=1", "eR+=1"
   if 0x1000000 LEQ %fR% set /A "fR>>=1", "eR+=1"

   if %fR% LSS 8388608 (
:: denormalized
      set /A "result=((sR&1)<<31)|fR"
   ) else (
:: normalized
      set /A "result=((sR&1)<<31)|((eR+127)<<23)|(fR&0x7FFFFF)"
   )

   endlocal & set "%~3=%result%"
   goto :eof


:intToHex

Sidenote: I will have to revise my code, but actually i have not much time to do this, so it may take a while.

penpen

[einstein1969]

Thanks penpen.

I will wait.

einstein1969

[einstein1969]

TODO:

In the floatadd/floatmul there is at the middle this:

Code: Select all

   if defined result endlocal & set "%~3=%result%"

this break the flow and for example the %e1% is not more defined

einstein1969

[penpen]

Ups!
A " & goto :eof" is missing at the end of both lines, ...
my test cases are suboptimal, as they don't cover it.
Thanks for discovering this bug: I've changed the source above.

penpen