3

Suppose I would like to speed up the function in FindMinimum[function[variables],startingPoint] by compiling function

I've constructed a minimal working example to describe what I would like to do.

Here is the code without compile (Minimizing some spring energy between several points)

p is a list of points.

springEnergy[p_] := 
 Total@Total[(DistanceMatrix[p, 
       DistanceFunction -> EuclideanDistance] - 1)^2]

positions = CirclePoints[4];

variables = p[#] & /@ Range[Length[positions]]

startingPoint = MapThread[{#1, #2} &, {variables, positions}]

This works:

FindMinimum[
 springEnergy[variables], startingPoint

]

I'd like to do the same thing, but with a compiled springEnergy. I believe this is functionally the same as springEnergy:

fc = FunctionCompile[
  Function[Typed[pos, TypeSpecifier["NumericArray"]["Real64", 2]],
   Module[{
     adm = 
      Typed[
         KernelFunction[
          DistanceMatrix], {TypeSpecifier["NumericArray"]["Real64", 
            2]} -> TypeSpecifier["NumericArray"]["Real64", 2]][pos] -
       Typed[KernelFunction[ConstantArray],
         {"Integer32", {"Integer32", "Integer32"}} -> 
          TypeSpecifier["NumericArray"]["Real64", 2]][
        1, {Length[pos], Length[pos]}]
     },
    Total[Total[adm^2]]
    ]
   ]
  ]

fc[NumericArray[N@positions]] == springEnergy[N@positions] (True)

However,

FindMinimum[
 fc[NumericArray[variables]], startingPoint
 ]

complains because--I think--FindMinimum calls the initial point symbolically. Normally, I get around that by specifying _?NumericQ on the calling function. I tried to so something like that with a wrapper function:

wrapper[variables_NumericArray] := fc[variables]

But, that doesn't prevent the initial symbolic call by FindMinimum:

FindMinimum[
 wrapper[variables], startingPoint
 ]
Craig Carter
  • 4,416
  • 16
  • 28
  • 1
    Perhaps wrapper[variables_?MatrixQ] := fc[variables]; FindMinimum[wrapper[variables], startingPoint] – Simon Woods Mar 18 '22 at 19:49
  • 1
    FindMinimum[Inactivate@fc[variables], startingPoint] also works – Simon Woods Mar 18 '22 at 20:13
  • @SimonWoods. Excellent. Both solutions work. Bravo. Inactivate would never occurred to me. I don't quite understand why it works. Also, the two solutions you give are about the same speed. But, good news, even for this simple problem, your solution gives a factor of two speed up. If you post this as an answer--I'll accept with thanks. – Craig Carter Mar 19 '22 at 09:44

1 Answers1

2

You had the right idea, but _NumericArray will only match explicit NumericArray expressions, not arrays of numbers in general. So your wrapper function isn't doing what you intend. Since you're working in 2D a suitable pattern for the wrapper is _?MatrixQ which will match the 2D numeric array but not the 1D symbolic array.

wrapper[variables_?MatrixQ] := fc[variables]

FindMinimum[wrapper[variables], startingPoint] (* {4.34315, {p[1] -> {0.426777, -0.426777}, p[2] -> {0.426777, 0.426777}, p[3] -> {-0.426777, 0.426777}, p[4] -> {-0.426777, -0.426777}}} *)

An alternative is to use Inactivate, which prevents the initial evaluation of fc with symbolic arguments, but gets removed (by an internal Activate) for the numeric calculations.

FindMinimum[Inactivate@fc[variables], startingPoint]
(* {4.34315, {p[1] -> {0.426777, -0.426777}, p[2] -> {0.426777, 0.426777}, 
  p[3] -> {-0.426777, 0.426777}, p[4] -> {-0.426777, -0.426777}}} *)
Simon Woods
  • 84,945
  • 8
  • 175
  • 324
  • Great advice and appreciated. Accepted as answered. – Craig Carter Mar 21 '22 at 15:26
  • How did you find that Inactivate works in this case and is internally activated? It is such a helpful feature that I think it should be more visible. Is that somehow documented or did you just find it by trial and error. Would be nice to have a list of function for which that works and confirmation from WRI that this is supported... – Albert Retey Mar 27 '22 at 10:47
  • @AlbertRetey I recalled seeing Hold used in a similar way in the past (eg https://mathematica.stackexchange.com/q/43318/862) and just tried Inactivate on a hunch. You can see that it is internally activated by the errors if you wrap the whole thing in Block[{Activate}, ... ]. I assume the lack of documentation means this is accidental, or at least untested, behaviour. I've only tried it on this example so no idea if it works more widely. – Simon Woods Mar 27 '22 at 22:52
  • OK, thanks. That makes sense, although I would wish these things would be documented (=specified) a bit clearer by WRI... – Albert Retey Mar 29 '22 at 06:46