/usr/share/octave/packages/interval-3.1.0/@infsup/private/double2decimal.m is in octave-interval 3.1.0-5.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 | ## Copyright 2014-2016 Oliver Heimlich
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 3 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, see <http://www.gnu.org/licenses/>.
## -*- texinfo -*-
## @documentencoding UTF-8
## @deffun double2decimal (@var{X})
##
## Convert a binary floating point number @var{X} in double precision to a
## decimal floating point number with arbitrary precision. The number must be
## a finite number and must not be NaN.
##
## Conversion is exact, because @code{rem (10, 2) == 0}.
##
## Sign: true (-) or false (+)
##
## Mantissa: Vector that holds the decimal digits after the decimal point. The
## first digit is not zero. Trailing zeroes are removed.
##
## Exponent: Integral exponent (base 10).
##
## @example
## @group
## x = double2decimal (-200)
## @result{}
## x.s = 1
## x.m = [2]
## x.e = 2
## y = double2decimal (0.125)
## @result{}
## y.s = 0
## y.m = [1 2 5]'
## y.e = 0
## z = double2decimal (0)
## @result{}
## z.s = 0
## z.m = []
## z.e = 0
## @end group
## @end example
## @seealso{str2decimal}
## @end deffun
## Author: Oliver Heimlich
## Created: 2014-09-29
function decimal = double2decimal (binary)
[decimal.s, exponent, fraction] = parsedouble (binary);
if (sum (fraction) == 0) # this number equals zero
decimal.s = false; # normalize: remove sign from -0
decimal.m = [];
decimal.e = int64 (0);
return
endif
## Remove trailing zeroes if this might reduce the number of loop cycles below
if (exponent < length (fraction))
fraction = fraction(1:find (fraction, 1, "last"));
endif
## Move the point to the end of the mantissa and interpret mantissa
## as a binary integer number that is now in front of the
## point. Convert binary integer to decimal.
exponent -= length (fraction);
decimal.m = zeros ();
for i = 1 : length(fraction)
## Multiply by 2
decimal.m .*= 2;
## Add 1 if necessary
decimal.m(end) += fraction(i);
## Carry
decimal.m = [0; rem(decimal.m, 10)] ...
+ [(decimal.m >= 10); 0];
endfor
clear fraction;
## Normalize: Remove leading zeroes (for performance reasons not in loop)
decimal.m = decimal.m(find (decimal.m ~= 0, 1, "first"):end);
assert (length (decimal.m) > 0, "number must not equal zero at this point");
decimal.e = int64 (length (decimal.m));
## Multiply decimal integer with 2 ^ exponent
while (exponent > 0)
decimal.m .*= 2;
decimal.m = [0; rem(decimal.m, 10)] ...
+ [(decimal.m >= 10); 0];
if (decimal.m(1) == 0)
decimal.m(1) = [];
else
decimal.e ++;
endif
exponent --;
endwhile
while (exponent < 0)
## Instead of division by 2 we devide by 10 and multiply by 5
decimal.e --; # cheap division by 10
decimal.m .*= 5;
decimal.m = [0; rem(decimal.m, 10)] ...
+ [floor(decimal.m ./ 10); 0];
if (decimal.m(1) == 0)
decimal.m(1) = [];
else
decimal.e ++;
endif
exponent ++;
endwhile
## Normalize mantissa: remove trailing zeroes;
decimal.m = decimal.m(1 : find (decimal.m ~= 0, 1, "last"));
endfunction
|