Currently, I am developing an octree data structure. And I need it to be fast. To this end, I would like to have a fast way to extract the $m$ binary digits "after the dot" of each entry in a matrix of machine precision numbers. This matrix is of size $n \times d$ where $d$ is the space dimension (say $d=3$) and represents a set of $n$ points (shift and rescaled to lie in the unit box) to which the octree has to adapt. The number $n$ is typically in the range of a couple of 100 thousands and a few million. So, the input parameters are, e.g.,
n = 1000000;
d = 3;
m = 32;
x = RandomReal[{0., 1. - $MachineEpsilon}, {n, d}];
For performance, it is crucial that the resulting $n \times d \times m$ tensor is a packed array of Integer
s. As you can see below, I have tried already several things. But in view of the fact that actually nothing has to be computed as the mantissas of the machine precision numbers are actually stored in binary, the timings are really disappointing.
ToPack = Developer`ToPackedArray;
cf = Compile[{{z0, _Integer}, {base, _Integer}, {n, _Integer}},
IntegerDigits[z0, base, n],
CompilationTarget -> "C",
RuntimeAttributes -> {Listable},
Parallelization -> True,
RuntimeOptions -> "Speed"
];
digits1 = ToPack[RealDigits[x, 2, m, -1][[All, All, 1]]]; //
AbsoluteTiming // First
digits2 = ToPack[IntegerDigits[IntegerPart[x 2.^m], 2, m]]; //
AbsoluteTiming // First
digits3 = Transpose[Table[Mod[IntegerPart[x 2.^k], 2], {k, 1, m}], {3, 1, 2}]; //
AbsoluteTiming // First
digits4 = Transpose[cf[IntegerPart[x 2.^m], 2, m], {1, 2, 3}]; //
AbsoluteTiming // First
digits1 == digits2 == digits3 == digits4
11.3763
3.00342
1.7729
0.516279
True
I have also tried to use BitGet
but that was a couple of times slower than RealDigits
.
In practice, I do not know in advance how many digits I need. So it would be even better to retrieve only the $k$-th digit at a time. In this case, the timings look as follows:
k = 4;
digit1 = ToPack[RealDigits[x, 2, 1, -k][[All, All, 1, -1]]]; //
AbsoluteTiming // First
digit2 = ToPack[
IntegerDigits[IntegerPart[x 2.^k], 2, 1][[All, All, -1]]]; //
AbsoluteTiming // First
digit3 = Mod[IntegerPart[x 2.^k], 2]; // AbsoluteTiming // First
digit4 = cf[IntegerPart[x 2.^k], 2, 1][[All, All, -1]]; //
AbsoluteTiming // First
digit1 == digit2 == digit3 == digit4
9.29457
2.31723
0.023894
0.238065
True
So, Mod
and IntegerPart
seem to do quite a good job, but still, there is a multiplication involved for each digit to be extracted.
Does anybody know of a fast, low-level way to retrieve the digits?
double
into anunsigned long long int
and extract the bits directly by masking, using info from en.wikipedia.org/wiki/IEEE_754 $\endgroup$uint64_t
version of the double number. Even with all safety checks off, and MTensor processing inside the C function (i.e. only one function call for a large list of input numbers) this plugin just about as fast as yourdigit3
. Not worth the trouble. In any case, as yourdigit3
takes only about 8 nanoseconds per conversion, I'm wondering what speed goal you had in mind? $\endgroup$BitAnd
andBitShiftRight
don't do what you think they would do on reals. $\endgroup$