15
$\begingroup$

A differential evolution algorithm is given here. I would like to get this kind of animation. I thought I could use NMinimize, given DifferentialEvolution as an option, but it turns out that does not work as I espected.

Is it possible to extract intermediate step in DifferentialEvolution, or do I have to implement algorithm myself?

f[x_, y_] := 
  -20 E^(-0.2 Sqrt[0.5 (x^2 + y^2)]) - E^(0.5 (Cos[2 π x] + Cos[2 π y])) + E + 20

p1 = 
  Plot3D[f[x, y], {x, -5, 5}, {y, -5, 5}, 
  PerformanceGoal -> "Quality",
  ColorFunction -> "WatermelonColors", 
  Mesh -> None, 
  BoxRatios -> {1, 1, 1}];   

p2 = 
   DensityPlot[f[x, y], {x, -5, 5}, {y, -5, 5}, 
     ColorFunction -> "WatermelonColors", 
     PlotPoints -> 200, 
     PerformanceGoal -> "Quality", 
     Frame -> False, 
     PlotRangePadding -> None];

p3 = Plot3D[0, {x, -5, 5}, {y, -5, 5}, PlotStyle -> Texture[p2], Mesh -> None];

Show[p1, p3, PlotRange -> {0, 15}]

enter image description here

When I use StepMonitor to track iterations as follows, it does not work.

{fit, intermediates} = 
  Reap[NMinimize[{f[x, y], -5 <= x <= 5, -5 <= y <= 5}, {x, y}, 
    MaxIterations -> 1000, 
    Method -> {"DifferentialEvolution", "InitialPoints" -> Tuples[Range[-5, 5], 2]}, 
    StepMonitor :> Sow[{x, y}]]];

Table[
  ListPlot[Take[intermediates[[1, i ;; i + 10]]], 
    Frame -> True, ImageSize -> 350, AspectRatio -> 1], 
  {i, 10, 1000, 100}]

EDIT Here is the result when we used @Michael E2 solution. Cool!!

f[x_, y_] := -20 E^(-0.2 Sqrt[0.5 (x^2 + y^2)]) - 
  E^(0.5 (Cos[2 \[Pi] x] + Cos[2 \[Pi] y])) + E + 20

p1 = Plot3D[f[x, y], {x, -5, 5}, {y, -5, 5}, 
   PerformanceGoal -> "Quality", ColorFunction -> "WatermelonColors", 
   Mesh -> None, BoxRatios -> {1, 1, 1}];

p2 = DensityPlot[f[x, y], {x, -5, 5}, {y, -5, 5}, 
   ColorFunction -> "WatermelonColors", PerformanceGoal -> "Quality", 
   Frame -> False, PlotRangePadding -> None];

p3 = Plot3D[-0.5, {x, -5, 5}, {y, -5, 5}, PlotStyle -> Texture[p2], 
   Mesh -> None];

p4 = Show[p1, p3, PlotRange -> {-0.5, 15}]


  {fit, intermediates} = 
   Reap[NMinimize[{f[x, y], -5 <= x <= 5, -5 <= y <= 5}, {x, y}, 
     MaxIterations -> 40, 
     Method -> {"DifferentialEvolution", 
       "InitialPoints" -> Tuples[Range[-5, 5], 2]}, 
     StepMonitor :> 
      Sow[{Optimization`NMinimizeDump`vecs, 
        Optimization`NMinimizeDump`vals}]]]; // Quiet
sim = Table[
   Show[p4, 
    ListPointPlot3D[{Append[#, 0] & /@ intermediates[[1, i, 1]]}, 
     PlotRange -> {{-5, 5}, {-5, 5}, {-5, 5}}, Boxed -> False, 
     PlotStyle -> Directive[AbsolutePointSize[3], Black]]], {i, 40}];

  Export[NotebookDirectory[] <> "DE2DSim.gif", sim, 
  "AnimationRepetitions" -> 50];

enter image description here

p1 = Plot3D[f[x, y], {x, -5, 5}, {y, -5, 5}, 
   PerformanceGoal -> "Quality", ColorFunction -> "WatermelonColors", 
   Mesh -> None, PlotStyle -> Opacity[0.4], BoxRatios -> {1, 1, 1}];

p2 = DensityPlot[f[x, y], {x, -5, 5}, {y, -5, 5}, 
   ColorFunction -> "WatermelonColors", PerformanceGoal -> "Quality", 
   Frame -> False, PlotRangePadding -> None];

p3 = Plot3D[-0.5, {x, -5, 5}, {y, -5, 5}, PlotStyle -> Texture[p2], 
   Mesh -> None];

p4 = Show[p1, p3, PlotRange -> {-0.5, 15}]


{fit, intermediates} = 
   Reap[NMinimize[{f[x, y], -5 <= x <= 5, -5 <= y <= 5}, {x, y}, 
     MaxIterations -> 40, 
     Method -> {"DifferentialEvolution", 
       "InitialPoints" -> Tuples[Range[-5, 5], 2]}, 
     StepMonitor :> 
      Sow[{Optimization`NMinimizeDump`vecs, 
        Optimization`NMinimizeDump`vals}]]]; // Quiet
sim = Table[
   Show[p4, 
    ListPointPlot3D[{Append[#, 0] & /@ intermediates[[1, i, 1]], 
      Partition[
       Flatten@Riffle[intermediates[[1, i, 1]], 
         f @@@ intermediates[[1, i, 1]]], 3]}, Boxed -> False, 
     PlotStyle -> {Directive[AbsolutePointSize[3], Black], 
       Directive[AbsolutePointSize[3], Red]}]], {i, 40}];

Export[NotebookDirectory[] <> "DE2DSim.gif", sim, 
  "AnimationRepetitions" -> 50];

enter image description here

Multicolumn[
 Table[Show[p4, 
   ListPointPlot3D[{Append[#, 0] & /@ intermediates[[1, i, 1]]}, 
    PlotRange -> {{-5, 5}, {-5, 5}, {-5, 5}}, Boxed -> False, 
    PlotStyle -> Directive[AbsolutePointSize[3], Black]]], {i, 1, 30, 
   2}], 5, Appearance -> "Horizontal"]

enter image description here

$\endgroup$
1
  • 2
    $\begingroup$ Note that blocking f (Block[{f}, ...]) isn't necessary. It was just to prevent f from being defined, which is a habit I have with single-lettter symbols on SE, esp. ones I use like f, x, etc. -- thanks for the accept! $\endgroup$
    – Michael E2
    Mar 11, 2019 at 3:01

1 Answer 1

15
$\begingroup$

Here's a way:

Block[{f},
 f[x_, y_] := -20 E^(-0.2 Sqrt[0.5 (x^2 + y^2)]) - 
   E^(0.5 (Cos[2 \[Pi] x] + Cos[2 \[Pi] y])) + E + 20;
 {fit, intermediates} = 
  Reap[NMinimize[{f[x, y], -5 <= x <= 5, -5 <= y <= 5}, {x, y}, 
    MaxIterations -> 30, 
    Method -> {"DifferentialEvolution", 
      "InitialPoints" -> Tuples[Range[-5, 5], 2]}, 
    StepMonitor :> 
     Sow[{Optimization`NMinimizeDump`vecs, 
       Optimization`NMinimizeDump`vals}]]];
 ]

Manipulate[
 Graphics[{
   PointSize[Medium],
   Point[intermediates[[1, n, 1]], 
    VertexColors -> 
     ColorData["Rainbow"] /@ 
      Rescale[intermediates[[1, n, 2]], 
       MinMax[intermediates[[1, All, 2]]]]]
   },
  PlotRange -> 5, Frame -> True],
 {n, 1, Length@intermediates[[1]], 1}
 ]

enter image description here

You can find out about things like Optimization`NMinimizeDump`vecs by inspecting the code for Optimization`NMinimizeDump`CoreDE.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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