I want to position several balls (8 in the example below) on a plane in such a way that distance between all balls is approximately 0.5. Precisely my objective is to minimize the highest deviance from 0.5 of distance between any pair of balls.

The interesting feature of this problem is that visualization of 8 balls on a plane is simple and intuitive - it is important because it allows me to visually assess the solution quality - not having to rely on objective function value only.

On the other hand the optimization involves 16 decision variables for 8 balls and the objective function is not smooth and has many local minima - a perfect case to compare different optimization algorithms.

The code that generates the plots for BFGS optimization algorithm is given below:

```
library(animation)
b.init <- function(balls) {
best.fit <<- +Inf
pos <- runif(2 * balls)
names(pos) <- paste(c("x", "y"), rep(1:balls, each = 2),
sep = "")
return(pos)
}
b.plot <- function(pos) {
par(mar = c(2, 2, 1.5, 0.5))
# extract x and y coordinates of balls
x <- pos[c(T, F)]
y <- pos[c(F, T)]
# plot balls
plot(x, y, cex = 5, pch = 20, main = b.fit(pos),
xlim = c(min(x) - 0.1, max(x) + 0.1),
ylim = c(min(y) - 0.1, max(y) + 0.1))
# and lines connecting them
for (i in 1:(length(x) - 1)) {
for (j in (i + 1):length(x)) {
lines(x[c(i, j)], y[c(i, j)])
}
}
}
b.fit <- function(pos) {
# our goal is to make distance between all balls equal to 0.5
max(abs(dist(t(matrix(pos, 2))) - 0.5))
}
f <- function(pos) {
fit <- b.fit(pos)
# update plot if improvement is found
if (fit < best.fit) {
best.fit <<- fit
b.plot(pos)
}
return(fit)
}
run <- function() {
set.seed(1)
pos <- b.init(8)
for (i in 1:20) {
b.plot(pos)
}
pos <- optim(pos, f, method = "BFGS", control = list(maxit = 1000))$par
for (i in 1:20) {
b.plot(pos)
}
}
ani.options(interval = 0.05)
saveGIF(run())
```

It generates the following animation:It is interesting to run this simulation for different initial positions of balls as BFGS might converge to different local minima or to try different optimization algorithms either from optim function (like SANN) or any other optiization routine (function f keeps track of best solution found independently from optimization procedure).

Very interesting post! Just a tiny reminder, maybe the code in the last line should be "saveGIF(run())", rather than "saveGIF(runa())"?

ReplyDeleteOf course. Corrected.

Delete