5
$\begingroup$

I'm minimizing a large function. Evaluating the function is very slow symbolically but using Compile I can get ~100x speed-up on this step. However, despite this, the NMinimize process is only marginally faster for large equations, and very much slower for small ones.

Is there a way to pass NMinimize the information it would normally get from a symbolic function but in compiled form? An example in this would be passing the function and its Jacobian to a root finder, avoiding the need to use finite differences.

Here is the example I'm working on. Any comments/pedantry as to how I've written this would also be appreciated; I still have no idea about the best ways to do things, and have trouble forcing/maintaining consistency.

centralangle[{q1_, p1_}, {q2_, p2_}] := 
  ArcCos[Sin[q1] Sin[q2] + Cos[q1] Cos[q2] Cos[p1 - p2]];
energy[qp_, i_, max_] := 
  Sum[If[s == i, 0, (centralangle[qp[[i]], qp[[s]]])^-1], {s, 1, max}];
totalenergy[qp_, max_] := Sum[energy[qp, i, max], {i, 1, max}]/2;

Now test this for 10 particles, and I find the symbolic minimization to be ~6x faster

ClearAll[q, p]
max = 10;
qp = Table[Unique[{q, p}], {max}];
qpFlat = Flatten[qp];
eqn = totalenergy[qp, max];
ceqn = Compile[Evaluate[{#, _Real} & /@ qpFlat], Evaluate[eqn], 
  "RuntimeOptions" -> {"EvaluateSymbolically" -> False}];
Timing[
 NMinimize[eqn, qpFlat];
 ]
Timing[
 NMinimize[ceqn[Sequence @@ qpFlat], qpFlat];
 ]



{0.374402, Null}(symbolic)
{2.246414, Null}(compiled)

Although for 100 particles, I find single evaluations ~100x faster when compiled

ClearAll[q, p]
max = 100;
qp = Table[Unique[{q, p}], {max}];
qpFlat = Flatten[qp];
eqn = totalenergy[qp, max];
ceqn = Compile[Evaluate[{#, _Real} & /@ qpFlat], Evaluate[eqn], 
  "RuntimeOptions" -> {"EvaluateSymbolically" -> False}];

init = RandomReal[1., 2 max];
Timing[
 eqn /. Table[qpFlat[[i]] -> init[[i]], {i, 1, 2 max}];
 ]
Timing[
 ceqn[Sequence @@ init];
 ]



{1.279208, Null}(symbolic)
{0.015600, Null}(compiled)

The minimization for 100 particles is overall 1.5x faster for the compiled function, clearly as a trade-off between whatever is going on here.

|improve this question
$\endgroup$
  • 1
    $\begingroup$ Do you have to use global optimization (NMinimize) rather than local (FindMinimum)? Any information you have such as good initial guesses will help speed up the problem. $\endgroup$ – dr.blochwave Mar 2 '15 at 11:40
  • $\begingroup$ FindMinimum produces rather poor results, and would surely suffer from the same slowing-down that it's buddy does. I'll have another go and see if there's a good initial point I can use though. $\endgroup$ – Crêpo Mar 2 '15 at 12:08
  • 2
    $\begingroup$ NMinimize does not require (and cannot accept) gradients. All its methods are direct search algorithms; there is no more information required than a black box function. If the Nelder-Mead method would be sufficient for you, perhaps you might like to try a compiled implementation of that. $\endgroup$ – Oleksandr R. Mar 2 '15 at 12:29
  • $\begingroup$ Can you explain why NMinimize would be so much slower using a compiled function? $\endgroup$ – Crêpo Mar 2 '15 at 13:00
  • $\begingroup$ Evaluate here simply converts the symbol eqn to the equation it represents. You can see the code has compiled correctly by running the additional code Needs["CompiledFunctionTools`"]; CompilePrint[ceqn]. The second part of my question shows major speed up in single evaluations of the equation, which further evidences correct compilation. $\endgroup$ – Crêpo Mar 2 '15 at 14:17
2
$\begingroup$

Update

In a comment Oleksandr pointed out that the number of iterations may be having an effect. This is evident from running the following code:

max = 10;
qp = Table[Unique[{q, p}], {max}];
qpFlat = Flatten[qp];
eqn = totalenergy[qp, max];
ceqn = Compile[Evaluate[{#, _Real} & /@ qpFlat], Evaluate[eqn], 
   "RuntimeOptions" -> {"EvaluateSymbolically" -> False}];

Block[{i = 0}, 
  r1 = NMinimize[eqn, qpFlat, 
    Method -> {Automatic, RandomSeed -> 1}, 
    EvaluationMonitor :> i++]; 
  i] // AbsoluteTiming

(* Average from 10 runs: 0.399 seconds, 303 evaluations *)

Block[{i = 0}, 
  r2 = NMinimize[ceqn[Sequence @@ qpFlat], qpFlat, 
    Method -> {Automatic, RandomSeed -> 1}, 
    EvaluationMonitor :> i++]; 
  i] // AbsoluteTiming

(* Average from 10 runs: 3.25 seconds, 3056 evaluations *)

For max = 50, the timings are

(* Uncompiled: 15.94 seconds, 521 evaluations *)
(* Compiled: 19.47 seconds, 5000 evaluations *)

Quite why the compiled function needs 10x more evaluations I'm currently not sure...

Original answer

It seems that the value of the WorkingPrecision option in NMinimize has an effect on the performance of the two functions.

Using the OP's definitions for eqn, ceqn and qpFlat, I tried the following. The random seed is set to 1 to ensure the returned parameters are the same.

AbsoluteTiming[
 r1 = NMinimize[eqn, qpFlat, 
    Method -> {"NelderMead", RandomSeed -> 1}, 
    WorkingPrecision -> 10];]
(* 3.119 seconds *)

AbsoluteTiming[
 r2 = NMinimize[ceqn[Sequence @@ qpFlat], qpFlat, 
    Method -> {"NelderMead", RandomSeed -> 1}, 
    WorkingPrecision -> 10];]
(* 0.311 seconds *)

Abs@(First@r1 - First@r2)
(* 6.6*10^-7 : this is the difference in the objective function value *)

Changing back to WorkingPrecision -> MachinePrecision, which is the default, I find that the compiled version takes longer than the uncompiled version, which is what the OP found.

|improve this answer
$\endgroup$
  • $\begingroup$ Are you sure you ran this after eqn had been assigned? The compilation works and prints correctly if you run the above code followed by Needs["CompiledFunctionTools`"]; CompilePrint[ceqn] If it did not, we would see no change in the speed of single evaluations as we do at the bottom of my question. $\endgroup$ – Crêpo Mar 2 '15 at 14:13
  • $\begingroup$ @user5751 have updated the answer having investigated the effect of WorkingPrecision... $\endgroup$ – dr.blochwave Mar 2 '15 at 15:32
  • 1
    $\begingroup$ It has some effect but the difference is not large. Compare e.g. Method -> {Automatic, RandomSeed -> 1}, WorkingPrecision -> MachinePrecision with Method -> {Automatic, RandomSeed -> 1}, WorkingPrecision -> 10, both for the compiled funciton. I think this is consistent with fewer iterations being performed more than anything else. $\endgroup$ – Oleksandr R. Mar 2 '15 at 15:42
  • $\begingroup$ @OleksandrR. I was thinking the same...see my updated answer. $\endgroup$ – dr.blochwave Mar 2 '15 at 18:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.