2

Can one be found programmatically? The function I'm compiling, indicesOfMin, gives the indices of the minimum value(s) of a real-valued list, the parameter, in a single traversal (which is why I'm using this instead of Ordering or Sort).

indicesOfMin = Compile[{{list, _Real, 1}},
  Module[{bag = Internal`Bag[Most[{0}]], min = 1.7976931348623157*^308, sub = 0.},
    Do[If[(sub = list[[i]] - min) < 0, min = list[[i]]; bag = Internal`Bag[Most[{0}]]; Internal`StuffBag[bag, i];,
      If[sub == 0, Internal`StuffBag[bag, i]]],
      {i, 1, Length[list], 1}];
    Internal`BagPart[bag, All]
  ]];

As you can see, I'm using C++'s positive infinity value for my first setting of min, but I don't even know if that's the language this will be compiled to.

Potential issues:

  1. I understand I can use the first value in list for min, but that takes extra code and I would like to avoid doing so because it also creates a precondition that the list is at least one element long.
  2. I don't know what language this will be compiled to in advance (the function is part of a package, and the person running it may not have a C compiler, so the maximum field value might be wrong.
VF1
  • 4,702
  • 23
  • 31
  • 1
    Compiled functions work with machine numbers, so the same limitations apply that you have in C when working with the double type. I think the minimum/maximum depends on how your CPU works, not the language. Compile can compile to a special byte code used by Mma, or it can compile to C. Both use the same type of machine number. – Szabolcs Feb 06 '13 at 18:25
  • @Szbolcs don't various languages allocate different amounts of space for doubles regardless of your CPU? – VF1 Feb 06 '13 at 18:32
  • 2
    I'd use $MaxMachineNumber but I don't know enough about the internals to vouch for it. You need to figure out a way to throw it in to the Compile to avoid MainEvaluate perhaps with Block – ssch Feb 06 '13 at 18:35
  • Typically they don't. Only a handful have custom float implementations, like Mathematica or Maple, and this is usually to support arbitrary precision. The rest use the hardware support in your CPU. What makes things more complicated is that Intel CPUs nowadays support a number of representations ... That said, I think what you are trying to do is a hack that's not wort the risk. I'd say just use a proper implementation (take the first element). I came to the same conclusion here: I won't take the shortcut. – Szabolcs Feb 06 '13 at 18:37
  • @ssch I guess that's worth posting as answer. – Szabolcs Feb 06 '13 at 18:41
  • @Szabolcs 80-bit extended precision numbers as used in the x87 are now considered obsolete (partly because the x87 is a stack based architecture and so difficult to achieve high performance with; partly because the standard double size is and always has been 64 bits and so you throw away 16 bits of precision anyway on storing the value). In fact, x87 instructions are not supported in 64-bit mode. So, really, the only representation we have is the standard 64 bits. SSE departs slightly from IEEE754 with respect to denormal handling and certain rounding modes, but it's mostly the same. – Oleksandr R. Feb 06 '13 at 19:53
  • @VF1 by the way, your value is not positive infinity; it is the maximum possible finite number. Using this value is fine but you should bear in mind that it's not strictly equivalent to an IEEE infinity. – Oleksandr R. Feb 06 '13 at 20:02
  • @Oleksandr I wasn't sure if SSE treats the numbers the same as x87 did, thanks for clearing this up. – Szabolcs Feb 06 '13 at 20:09
  • @OleksandrR. I'm afraid I'm not familiar with what an IEEE infinity is. As long as nothing can be greater than it, it should be fine. – VF1 Feb 06 '13 at 20:32

2 Answers2

10

To my surprise, I found that it's actually not difficult to produce a true floating point infinity inside compiled code. Here's how:

cf = With[{
   $MaxMachineNumber = $MaxMachineNumber,
   $MachineEpsilon = $MachineEpsilon
  },
  Compile[{{x, _Real, 0}},
   Block[{inf = (1. + $MachineEpsilon) $MaxMachineNumber},
    inf > x
   ], RuntimeOptions -> "CompareWithTolerance" -> False
  ]
 ]

We can verify that this value is larger than the largest finite number:

cf[$MaxMachineNumber]
(* -> True *)

However, while one can use it inside compiled code, it is not possible to return it as a result: this causes a numerical error that results in the code being reevalated outside of the VM.

It would also appear that the Mathematica VM is not particularly IEEE754-compliant: attempts to produce a NaN (e.g. NaN = 0. / 0.) result in a value that is both greater and less than any other, when in fact it should be neither greater nor less than anything else. This could perhaps be a valid choice but will make it very difficult to use NaNs productively in compiled code. Furthermore, both inf / inf and inf - inf should give NaN, but while the former works as it should, the latter gives an error, possibly as a result of producing a signaling NaN value (even though exclusively quiet NaNs seem to be produced otherwise).

A final important point: to have any hope of obtaining meaningful (I won't go as far as to say consistent) results with respect to NaNs, it is important to set the option RuntimeOptions -> "CompareWithTolerance" -> False as I have done in the example. Absent this, NaN compares less (but not greater) than both itself and any strictly positive value, which is obviously semantically incorrect.

Oleksandr R.
  • 23,023
  • 4
  • 87
  • 125
8

$MaxMachineNumber has the largest machine number the current computer supports, you could force it into the Compile like this:

f = With[{$MaxMachineNumber = $MaxMachineNumber},
      Compile[{{i, _Real}},
        Module[{min = $MaxMachineNumber},
          Min[i, min]
      ]]]

I don't know if this is a reliable way though.

ssch
  • 16,590
  • 2
  • 53
  • 88
  • 2
    The same trick is used in this answer http://mathematica.stackexchange.com/a/4877/745. Works without any problems. – Ajasja Feb 06 '13 at 18:46
  • Actually machine infinity is also supported in compiled code... the problem is that there's no way to pass it in or out, so that's why I chose this value instead. – Oleksandr R. Feb 06 '13 at 19:41