Skip to contents

Optimize an objective function using an adaptive random search algorithm. The function works for both discrete and continuous optimization parameters and allows for box-constraints and sets of allowed values.

Usage

optim_ARS(
  par,
  fn,
  lower = NULL,
  upper = NULL,
  allowed_values = NULL,
  loc_fac = 4,
  no_bounds_sd = par,
  iter = 400,
  iter_adapt = 50,
  adapt_scale = 1,
  max_run = 200,
  trace = TRUE,
  trace_iter = 5,
  new_par_max_it = 200,
  maximize = F,
  parallel = F,
  parallel_type = NULL,
  num_cores = NULL,
  mrgsolve_model = NULL,
  seed = round(runif(1, 0, 1e+07)),
  allow_replicates = TRUE,
  replicates_index = seq(1, length(par)),
  generator = NULL,
  ...
)

Arguments

par

A vector of initial values for the parameters to be optimized over.

fn

A function to be minimized (or maximized), with first argument the vector of parameters over which minimization is to take place. It should return a scalar result.

lower

Lower bounds on the parameters. A vector.

upper

Upper bounds on the parameters. A vector.

allowed_values

A list containing allowed values for each parameter list(par1=c(2,3,4,5,6),par2=c(5,6,7,8)). A vector containing allowed values for all parameters is also allowed c(2,3,4,5,6).

loc_fac

Locality factor for determining the standard deviation of the sampling distribution around the current position of the parameters. The initial standard deviation is normally calculated as (upper - lower)/loc_fac except in cases when there are no upper or lower limits (e.g. when upper=Inf or lower=-Inf).

no_bounds_sd

The standard deviation of the sampling distribution around the current position of the parameters when there are no upper or lower limits (e.g. when upper=Inf or lower=-Inf).

iter

The number of iterations for the algorithm to perform (this is a maximum number, it could be less).

iter_adapt

The number of iterations before adapting (shrinking) the parameter search space.

adapt_scale

The scale for adapting the size of the sampling distribution. The adaptation of the standard deviation of the sampling distribution around the current position of the parameters is done after iter_adapt iteration with no change in the best objective function. When adapting, the standard deviation of the sampling distribution is calculated as (upper - lower)/(loc_fac*ff*adapt_scale) where ff starts at 1 and increases by 1 for each adaptation.

max_run

The maximum number of iterations to run without a change in the best parameter estimates.

trace

Should the algorithm output results intermittently.

trace_iter

How many iterations between each update to the screen about the result of the search.

new_par_max_it

The algorithm randomly chooses samples based on the current best set of parameters. If when drawing these samples the new parameter set has already been tested then a new draw is performed. After new_par_max_it draws, with no new parameter sets, then the algorithm stops.

maximize

Should the function be maximized? Default is to minimize.

parallel

Should we use parallel computations?

parallel_type

Which type of parallelization should be used? Can be "snow" or "multicore". "snow" works on Linux-like systems & Windows. "multicore" works only on Linux-like systems. By default this is chosen for you depending on your operating system. See start_parallel.

num_cores

The number of cores to use in the parallelization. By default is set to the number output from parallel::detectCores(). See start_parallel.

mrgsolve_model

If the computations require a mrgsolve model and you are using the "snow" method then you need to specify the name of the model object created by mread or mcode.

seed

The random seed to use in the algorithm,

allow_replicates

Should the algorithm allow parameters to have the same value?

replicates_index

A vector, the same length as the parameters. If two values are the same in this vector then the parameters may not assume the same value in the optimization.

generator

A user-defined function that generates new parameter sets to try in the algorithm. See examples below.

...

Additional arguments passed to fn and start_parallel.

References

  1. M. Foracchia, A.C. Hooker, P. Vicini and A. Ruggeri, "PopED, a software fir optimal experimental design in population kinetics", Computer Methods and Programs in Biomedicine, 74, 2004.

  2. J. Nyberg, S. Ueckert, E.A. Stroemberg, S. Hennig, M.O. Karlsson and A.C. Hooker, "PopED: An extended, parallelized, nonlinear mixed effects models optimal design tool", Computer Methods and Programs in Biomedicine, 108, 2012.

Examples


## "wild" function , global minimum at about -15.81515
fw <- function(x) 10*sin(0.3*x)*sin(1.3*x^2) + 0.00001*x^4 + 0.2*x+80

# optimization with fewer function evaluations compared to SANN
res1 <- optim_ARS(50, fw,lower = -50, upper=100)
#> Initial OFV = 159.001
#> It.   5 | OFV = 93.7087
#> It.  10 | OFV = 72.1331
#> It.  15 | OFV = 72.1331
#> It.  20 | OFV = 72.1331
#> It.  25 | OFV = 72.1331
#> It.  30 | OFV = 72.1331
#> It.  35 | OFV = 72.1331
#> It.  40 | OFV = 72.1331
#> It.  45 | OFV = 72.1331
#> It.  50 | OFV = 72.1331
#> It.  55 | OFV = 72.1331
#> It.  60 | OFV = 72.1331
#> It.  65 | OFV = 72.1331
#> It.  70 | OFV = 72.1331
#> It.  75 | OFV = 70.2551
#> It.  80 | OFV = 70.2551
#> It.  85 | OFV = 70.2551
#> It.  90 | OFV = 70.2551
#> It.  95 | OFV = 70.2551
#> It. 100 | OFV = 69.6762
#> It. 105 | OFV = 69.6762
#> It. 110 | OFV = 69.6762
#> It. 115 | OFV = 69.6762
#> It. 120 | OFV = 69.6762
#> It. 125 | OFV = 69.4823
#> It. 130 | OFV = 69.4823
#> It. 135 | OFV = 69.4823
#> It. 140 | OFV = 69.4823
#> It. 145 | OFV = 69.4823
#> It. 150 | OFV = 69.4823
#> It. 155 | OFV = 69.4823
#> It. 160 | OFV = 69.4823
#> It. 165 | OFV = 69.4823
#> It. 170 | OFV = 69.4823
#> It. 175 | OFV = 69.4823
#> It. 180 | OFV = 69.4823
#> It. 185 | OFV = 69.4823
#> It. 190 | OFV = 69.4823
#> It. 195 | OFV = 69.4823
#> It. 200 | OFV = 69.4823
#> It. 205 | OFV = 69.4823
#> It. 210 | OFV = 69.4823
#> It. 215 | OFV = 69.4823
#> It. 220 | OFV = 69.4823
#> It. 225 | OFV = 69.4823
#> It. 230 | OFV = 68.5206
#> It. 235 | OFV = 68.5206
#> It. 240 | OFV = 68.5206
#> It. 245 | OFV = 68.244
#> It. 250 | OFV = 68.244
#> It. 255 | OFV = 68.244
#> It. 260 | OFV = 68.244
#> It. 265 | OFV = 68.244
#> It. 270 | OFV = 68.244
#> It. 275 | OFV = 68.244
#> It. 280 | OFV = 68.244
#> It. 285 | OFV = 68.244
#> It. 290 | OFV = 68.244
#> It. 295 | OFV = 68.0162
#> It. 300 | OFV = 68.0162
#> It. 305 | OFV = 68.0162
#> It. 310 | OFV = 68.0162
#> It. 315 | OFV = 68.0162
#> It. 320 | OFV = 68.0162
#> It. 325 | OFV = 68.0162
#> It. 330 | OFV = 68.0162
#> It. 335 | OFV = 68.0162
#> It. 340 | OFV = 68.0162
#> It. 345 | OFV = 68.0162
#> It. 350 | OFV = 68.0162
#> It. 355 | OFV = 68.0162
#> It. 360 | OFV = 68.0162
#> It. 365 | OFV = 67.8963
#> It. 370 | OFV = 67.8963
#> It. 375 | OFV = 67.8963
#> It. 380 | OFV = 67.8963
#> It. 385 | OFV = 67.7936
#> It. 390 | OFV = 67.7936
#> It. 395 | OFV = 67.7936
#> It. 400 | OFV = 67.7936
#> 
#> Total iterations: 400 
#> Elapsed time: 0.159 seconds.
#> 
#> Final OFV =  67.79359 
#> Parameters: -15.18685 
#> 

# often not as good performance when upper and lower bounds are poor
res2 <- optim_ARS(50, fw, lower=-Inf,upper=Inf)
#> Initial OFV = 159.001
#> It.   5 | OFV = 83.4828
#> It.  10 | OFV = 83.4828
#> It.  15 | OFV = 83.4828
#> It.  20 | OFV = 73.7904
#> It.  25 | OFV = 68.0584
#> It.  30 | OFV = 68.0584
#> It.  35 | OFV = 68.0584
#> It.  40 | OFV = 68.0584
#> It.  45 | OFV = 68.0584
#> It.  50 | OFV = 68.0584
#> It.  55 | OFV = 68.0584
#> It.  60 | OFV = 68.0584
#> It.  65 | OFV = 68.0584
#> It.  70 | OFV = 68.0584
#> It.  75 | OFV = 68.0584
#> It.  80 | OFV = 68.0584
#> It.  85 | OFV = 68.0584
#> It.  90 | OFV = 68.0584
#> It.  95 | OFV = 68.0584
#> It. 100 | OFV = 68.0584
#> It. 105 | OFV = 68.0584
#> It. 110 | OFV = 68.0584
#> It. 115 | OFV = 68.0584
#> It. 120 | OFV = 68.0584
#> It. 125 | OFV = 68.0584
#> It. 130 | OFV = 67.9762
#> It. 135 | OFV = 67.9762
#> It. 140 | OFV = 67.9762
#> It. 145 | OFV = 67.9762
#> It. 150 | OFV = 67.9762
#> It. 155 | OFV = 67.9762
#> It. 160 | OFV = 67.9762
#> It. 165 | OFV = 67.9762
#> It. 170 | OFV = 67.9762
#> It. 175 | OFV = 67.9762
#> It. 180 | OFV = 67.9762
#> It. 185 | OFV = 67.9762
#> It. 190 | OFV = 67.9762
#> It. 195 | OFV = 67.9762
#> It. 200 | OFV = 67.9762
#> It. 205 | OFV = 67.9762
#> It. 210 | OFV = 67.9762
#> It. 215 | OFV = 67.9762
#> It. 220 | OFV = 67.9762
#> It. 225 | OFV = 67.9762
#> It. 230 | OFV = 67.9762
#> It. 235 | OFV = 67.9762
#> It. 240 | OFV = 67.9762
#> It. 245 | OFV = 67.9762
#> It. 250 | OFV = 67.9762
#> It. 255 | OFV = 67.9762
#> It. 260 | OFV = 67.9762
#> It. 265 | OFV = 67.9762
#> It. 270 | OFV = 67.9762
#> It. 275 | OFV = 67.9762
#> It. 280 | OFV = 67.9762
#> It. 285 | OFV = 67.9762
#> It. 290 | OFV = 67.6953
#> It. 295 | OFV = 67.6953
#> It. 300 | OFV = 67.6953
#> It. 305 | OFV = 67.6953
#> It. 310 | OFV = 67.6953
#> It. 315 | OFV = 67.6953
#> It. 320 | OFV = 67.6953
#> It. 325 | OFV = 67.6953
#> It. 330 | OFV = 67.6953
#> It. 335 | OFV = 67.6953
#> It. 340 | OFV = 67.6953
#> It. 345 | OFV = 67.6953
#> It. 350 | OFV = 67.6953
#> It. 355 | OFV = 67.6953
#> It. 360 | OFV = 67.6953
#> It. 365 | OFV = 67.6953
#> It. 370 | OFV = 67.6953
#> It. 375 | OFV = 67.6953
#> It. 380 | OFV = 67.6953
#> It. 385 | OFV = 67.6953
#> It. 390 | OFV = 67.6953
#> It. 395 | OFV = 67.6953
#> It. 400 | OFV = 67.6953
#> 
#> Total iterations: 400 
#> Elapsed time: 0.15 seconds.
#> 
#> Final OFV =  67.6953 
#> Parameters: -15.96228 
#> 

# Only integer values allowed
if (FALSE) { 
res_int <- optim_ARS(50, fw, allowed_values = seq(-50,100,by=1))
}

if (FALSE) { 
  #plot of the function and solutions
  require(graphics)
  plot(fw, -50, 50, n = 1000, main = "Minimizing 'wild function'")
  points(-15.81515, fw(-15.81515), pch = 16, col = "red", cex = 1)
  points(res1$par, res1$ofv, pch = 16, col = "green", cex = 1)
  points(res2$par, res2$ofv, pch = 16, col = "blue", cex = 1)
} 

# optim_ARS does not work great for hard to find minima on flat surface:
# Rosenbrock Banana function
# f(x, y) = (a-x)^2 + b(y-x^2)^2
# global minimum at (x, y)=(a, a^2), where f(x, y)=0. 
# Usually a = 1 and b = 100.
if (FALSE) { 
  fr <- function(x,a=1,b=100) {   
    x1 <- x[1]
    x2 <- x[2]
    b*(x2 - x1*x1)^2 + (a - x1)^2
  }
  
  res3 <- optim_ARS(c(-1.2,1), fr,lower = -5, upper = 5)
  
  # plot the surface
  x <- seq(-50, 50, length= 30)
  y <- x
  f <- function(x,y){apply(cbind(x,y),1,fr)}
  z <- outer(x, y, f)
  persp(x, y, z, theta = 30, phi = 30, expand = 0.5, col = "lightblue", ticktype="detailed") -> res
  points(trans3d(1, 1, 0, pmat = res), col = 2, pch = 16,cex=2)
  points(trans3d(res3$par[1], res3$par[1], res3$ofv, pmat = res), col = "green", pch = 16,cex=2)
}

# box constraints
flb <- function(x){
  p <- length(x) 
  sum(c(1, rep(4, p-1)) * (x - c(1, x[-p])^2)^2) 
}
## 25-dimensional box constrained
#optim(rep(3, 25), flb,lower = rep(2, 25), upper = rep(4, 25),method = "L-BFGS-B") 
res_box <- optim_ARS(rep(3, 25), flb,lower = rep(2, 25), upper = rep(4, 25)) 
#> Initial OFV = 3460
#> It.   5 | OFV = 2952.53
#> It.  10 | OFV = 2509.04
#> It.  15 | OFV = 2015.84
#> It.  20 | OFV = 2015.84
#> It.  25 | OFV = 2015.84
#> It.  30 | OFV = 2015.84
#> It.  35 | OFV = 2015.84
#> It.  40 | OFV = 2015.84
#> It.  45 | OFV = 2015.84
#> It.  50 | OFV = 1804.98
#> It.  55 | OFV = 1804.98
#> It.  60 | OFV = 1591.51
#> It.  65 | OFV = 1487.52
#> It.  70 | OFV = 1446.98
#> It.  75 | OFV = 1446.98
#> It.  80 | OFV = 1446.98
#> It.  85 | OFV = 1379.83
#> It.  90 | OFV = 1379.83
#> It.  95 | OFV = 1379.83
#> It. 100 | OFV = 1379.83
#> It. 105 | OFV = 1379.83
#> It. 110 | OFV = 1317.46
#> It. 115 | OFV = 1317.46
#> It. 120 | OFV = 987.338
#> It. 125 | OFV = 987.338
#> It. 130 | OFV = 987.338
#> It. 135 | OFV = 948.102
#> It. 140 | OFV = 948.102
#> It. 145 | OFV = 948.102
#> It. 150 | OFV = 948.102
#> It. 155 | OFV = 948.102
#> It. 160 | OFV = 948.102
#> It. 165 | OFV = 948.102
#> It. 170 | OFV = 948.102
#> It. 175 | OFV = 948.102
#> It. 180 | OFV = 948.102
#> It. 185 | OFV = 927.919
#> It. 190 | OFV = 927.919
#> It. 195 | OFV = 927.919
#> It. 200 | OFV = 927.919
#> It. 205 | OFV = 922.904
#> It. 210 | OFV = 922.904
#> It. 215 | OFV = 922.904
#> It. 220 | OFV = 922.904
#> It. 225 | OFV = 922.904
#> It. 230 | OFV = 842.595
#> It. 235 | OFV = 605.287
#> It. 240 | OFV = 605.287
#> It. 245 | OFV = 605.287
#> It. 250 | OFV = 603.449
#> It. 255 | OFV = 603.449
#> It. 260 | OFV = 603.449
#> It. 265 | OFV = 603.449
#> It. 270 | OFV = 603.449
#> It. 275 | OFV = 603.449
#> It. 280 | OFV = 603.449
#> It. 285 | OFV = 603.449
#> It. 290 | OFV = 538.008
#> It. 295 | OFV = 538.008
#> It. 300 | OFV = 538.008
#> It. 305 | OFV = 538.008
#> It. 310 | OFV = 538.008
#> It. 315 | OFV = 538.008
#> It. 320 | OFV = 538.008
#> It. 325 | OFV = 538.008
#> It. 330 | OFV = 538.008
#> It. 335 | OFV = 538.008
#> It. 340 | OFV = 523.843
#> It. 345 | OFV = 517.833
#> It. 350 | OFV = 517.833
#> It. 355 | OFV = 517.833
#> It. 360 | OFV = 517.833
#> It. 365 | OFV = 517.833
#> It. 370 | OFV = 517.833
#> It. 375 | OFV = 517.833
#> It. 380 | OFV = 517.833
#> It. 385 | OFV = 517.833
#> It. 390 | OFV = 517.833
#> It. 395 | OFV = 517.833
#> It. 400 | OFV = 517.833
#> 
#> Total iterations: 400 
#> Elapsed time: 0.174 seconds.
#> 
#> Final OFV =  517.833 
#> Parameters: 2.052582 2.107788 2.236304 2.054542 2 2.176677 2 2.184549 2 2 2.194875 2.036758 2.098032 2 2 2 2 2 2.166911 2.40921 2 2.28692 2.18428 2.006472 3.260584 
#> 


## Combinatorial optimization: Traveling salesman problem
eurodistmat <- as.matrix(eurodist)

distance <- function(sq) {  # Target function
  sq2 <- embed(sq, 2)
  sum(eurodistmat[cbind(sq2[,2], sq2[,1])])
}

genseq <- function(sq) {  # Generate new candidate sequence
  idx <- seq(2, NROW(eurodistmat)-1)
  changepoints <- sample(idx, size = 2, replace = FALSE)
  tmp <- sq[changepoints[1]]
  sq[changepoints[1]] <- sq[changepoints[2]]
  sq[changepoints[2]] <- tmp
  sq
}

sq <- c(1:nrow(eurodistmat), 1)  # Initial sequence: alphabetic
res3 <- optim_ARS(sq,distance,generator=genseq) # Near optimum distance around 12842
#> Initial OFV = 29625
#> It.   5 | OFV = 29530
#> It.  10 | OFV = 28183
#> It.  15 | OFV = 25418
#> It.  20 | OFV = 25279
#> It.  25 | OFV = 24859
#> It.  30 | OFV = 24859
#> It.  35 | OFV = 23999
#> It.  40 | OFV = 21911
#> It.  45 | OFV = 21285
#> It.  50 | OFV = 20727
#> It.  55 | OFV = 20727
#> It.  60 | OFV = 20673
#> It.  65 | OFV = 20587
#> It.  70 | OFV = 20526
#> It.  75 | OFV = 20526
#> It.  80 | OFV = 20414
#> It.  85 | OFV = 20414
#> It.  90 | OFV = 20414
#> It.  95 | OFV = 19642
#> It. 100 | OFV = 19434
#> It. 105 | OFV = 19398
#> It. 110 | OFV = 19398
#> It. 115 | OFV = 19398
#> It. 120 | OFV = 19398
#> It. 125 | OFV = 19276
#> It. 130 | OFV = 19160
#> It. 135 | OFV = 19160
#> It. 140 | OFV = 18381
#> It. 145 | OFV = 17403
#> It. 150 | OFV = 17187
#> It. 155 | OFV = 17187
#> It. 160 | OFV = 17187
#> It. 165 | OFV = 17187
#> It. 170 | OFV = 16445
#> It. 175 | OFV = 16445
#> It. 180 | OFV = 16445
#> It. 185 | OFV = 16445
#> It. 190 | OFV = 16445
#> It. 195 | OFV = 16445
#> It. 200 | OFV = 16445
#> It. 205 | OFV = 16445
#> It. 210 | OFV = 16445
#> It. 215 | OFV = 16445
#> It. 220 | OFV = 16241
#> It. 225 | OFV = 16241
#> It. 230 | OFV = 16241
#> It. 235 | OFV = 16241
#> It. 240 | OFV = 16241
#> It. 245 | OFV = 16241
#> It. 250 | OFV = 16210
#> It. 255 | OFV = 16210
#> It. 260 | OFV = 16210
#> It. 265 | OFV = 16210
#> It. 270 | OFV = 16210
#> It. 275 | OFV = 16053
#> It. 280 | OFV = 16053
#> It. 285 | OFV = 15494
#> It. 290 | OFV = 15494
#> It. 295 | OFV = 15494
#> It. 300 | OFV = 15494
#> It. 305 | OFV = 15494
#> It. 310 | OFV = 15494
#> It. 315 | OFV = 15494
#> It. 320 | OFV = 15494
#> It. 325 | OFV = 15494
#> It. 330 | OFV = 14903
#> It. 335 | OFV = 14903
#> It. 340 | OFV = 14903
#> It. 345 | OFV = 14903
#> It. 350 | OFV = 14903
#> It. 355 | OFV = 14903
#> It. 360 | OFV = 14798
#> It. 365 | OFV = 14798
#> It. 370 | OFV = 14798
#> It. 375 | OFV = 14798
#> It. 380 | OFV = 14798
#> It. 385 | OFV = 14798
#> It. 390 | OFV = 14798
#> It. 395 | OFV = 14798
#> It. 400 | OFV = 14798
#> 
#> Total iterations: 400 
#> Elapsed time: 0.188 seconds.
#> 
#> Final OFV =  14798 
#> Parameters: 1 19 18 3 11 7 20 10 4 5 12 9 14 2 15 13 8 16 17 6 21 1 
#> 

if (FALSE) { 
  # plot of initial sequence
  # rotate for conventional orientation
  loc <- -cmdscale(eurodist, add = TRUE)$points
  x <- loc[,1]; y <- loc[,2]
  s <- seq_len(nrow(eurodistmat))
  tspinit <- loc[sq,]
  
  plot(x, y, type = "n", asp = 1, xlab = "", ylab = "",
       main = paste("Initial sequence of traveling salesman problem\n",
                    "Distance =",distance(sq)), axes = FALSE)
  arrows(tspinit[s,1], tspinit[s,2], tspinit[s+1,1], tspinit[s+1,2],
         angle = 10, col = "green")
  text(x, y, labels(eurodist), cex = 0.8)
  
  # plot of final sequence from optim_ARS
  tspres <- loc[res3$par,]
  plot(x, y, type = "n", asp = 1, xlab = "", ylab = "",
       main = paste("optim_ARS() 'solving' traveling salesman problem\n",
                    "Distance =",distance(c(1,res3$par,1))),axes = FALSE)
  arrows(tspres[s,1], tspres[s,2], tspres[s+1,1], tspres[s+1,2],
         angle = 10, col = "red")
  text(x, y, labels(eurodist), cex = 0.8)
  
  # using optim
  set.seed(123) # chosen to get a good soln relatively quickly
  (res4 <- optim(sq, distance, genseq, method = "SANN",
                 control = list(maxit = 30000, temp = 2000, trace = TRUE,
                                REPORT = 500))) 
  
  tspres <- loc[res4$par,]
  plot(x, y, type = "n", asp = 1, xlab = "", ylab = "",
       main = paste("optim() 'solving' traveling salesman problem\n",
                    "Distance =",distance(res4$par)),axes = FALSE)
  arrows(tspres[s,1], tspres[s,2], tspres[s+1,1], tspres[s+1,2],
         angle = 10, col = "red")
  text(x, y, labels(eurodist), cex = 0.8)
}  

# one-dimensional function
if (FALSE) { 
  f <- function(x)  abs(x)+cos(x)
  res5 <- optim_ARS(-20,f,lower=-20, upper=20)
  
  curve(f, -20, 20)
  abline(v = res5$par, lty = 4,col="green")
}  

# one-dimensional function
f <- function(x)  (x^2+x)*cos(x) # -10 < x < 10
res_max <- optim_ARS(0,f,lower=-10, upper=10,maximize=TRUE) # sometimes to local maxima
#> Initial OFV = 0
#> It.   5 | OFV = 0.304035
#> It.  10 | OFV = 1.0327
#> It.  15 | OFV = 20.8775
#> It.  20 | OFV = 20.8775
#> It.  25 | OFV = 32.7412
#> It.  30 | OFV = 32.7412
#> It.  35 | OFV = 32.7412
#> It.  40 | OFV = 32.7412
#> It.  45 | OFV = 32.7412
#> It.  50 | OFV = 32.7412
#> It.  55 | OFV = 32.7412
#> It.  60 | OFV = 34.4013
#> It.  65 | OFV = 34.4013
#> It.  70 | OFV = 34.4013
#> It.  75 | OFV = 34.4013
#> It.  80 | OFV = 34.4013
#> It.  85 | OFV = 47.0327
#> It.  90 | OFV = 47.0327
#> It.  95 | OFV = 47.0327
#> It. 100 | OFV = 47.0327
#> It. 105 | OFV = 47.0327
#> It. 110 | OFV = 47.0327
#> It. 115 | OFV = 47.0327
#> It. 120 | OFV = 47.0327
#> It. 125 | OFV = 47.0327
#> It. 130 | OFV = 47.0327
#> It. 135 | OFV = 47.0327
#> It. 140 | OFV = 47.0327
#> It. 145 | OFV = 47.3228
#> It. 150 | OFV = 47.3228
#> It. 155 | OFV = 47.5131
#> It. 160 | OFV = 47.5131
#> It. 165 | OFV = 47.5131
#> It. 170 | OFV = 47.5131
#> It. 175 | OFV = 47.5131
#> It. 180 | OFV = 47.6892
#> It. 185 | OFV = 47.6892
#> It. 190 | OFV = 47.6892
#> It. 195 | OFV = 47.6892
#> It. 200 | OFV = 47.6892
#> It. 205 | OFV = 47.6892
#> It. 210 | OFV = 47.6892
#> It. 215 | OFV = 47.6892
#> It. 220 | OFV = 47.6892
#> It. 225 | OFV = 47.6966
#> It. 230 | OFV = 47.6966
#> It. 235 | OFV = 47.6966
#> It. 240 | OFV = 47.6966
#> It. 245 | OFV = 47.6975
#> It. 250 | OFV = 47.7054
#> It. 255 | OFV = 47.7054
#> It. 260 | OFV = 47.7054
#> It. 265 | OFV = 47.7054
#> It. 270 | OFV = 47.7054
#> It. 275 | OFV = 47.7054
#> It. 280 | OFV = 47.7054
#> It. 285 | OFV = 47.7054
#> It. 290 | OFV = 47.7054
#> It. 295 | OFV = 47.7054
#> It. 300 | OFV = 47.7054
#> It. 305 | OFV = 47.7054
#> It. 310 | OFV = 47.7054
#> It. 315 | OFV = 47.7054
#> It. 320 | OFV = 47.7054
#> It. 325 | OFV = 47.7054
#> It. 330 | OFV = 47.7054
#> It. 335 | OFV = 47.7054
#> It. 340 | OFV = 47.7054
#> It. 345 | OFV = 47.7054
#> It. 350 | OFV = 47.7054
#> It. 355 | OFV = 47.7054
#> It. 360 | OFV = 47.7054
#> It. 365 | OFV = 47.7054
#> It. 370 | OFV = 47.7054
#> It. 375 | OFV = 47.7054
#> It. 380 | OFV = 47.7054
#> It. 385 | OFV = 47.7054
#> It. 390 | OFV = 47.7054
#> It. 395 | OFV = 47.7054
#> It. 400 | OFV = 47.7054
#> 
#> Total iterations: 400 
#> Elapsed time: 0.181 seconds.
#> 
#> Final OFV =  47.7054 
#> Parameters: 6.563397 
#> 

if (FALSE) { 
  res_min <- optim_ARS(0,f,lower=-10, upper=10) # sometimes to local minima
  
  curve(f, -10, 10)
  abline(v = res_min$par, lty = 4,col="green")
  abline(v = res_max$par, lty = 4,col="red")
}


# two-dimensional Rastrigin function
#It has a global minimum at f(x) = f(0) = 0.
if (FALSE) { 
  Rastrigin <- function(x1, x2){
    20 + x1^2 + x2^2 - 10*(cos(2*pi*x1) + cos(2*pi*x2))
  }
  
  
  x1 <- x2 <- seq(-5.12, 5.12, by = 0.1)
  z <- outer(x1, x2, Rastrigin)
  
  res6 <- optim_ARS(c(-4,4),function(x) Rastrigin(x[1], x[2]),lower=-5.12, upper=5.12)
  
  # color scale
  nrz <- nrow(z)
  ncz <- ncol(z)
  jet.colors <-
    colorRampPalette(c("#00007F", "blue", "#007FFF", "cyan",
                       "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000"))
  # Generate the desired number of colors from this palette
  nbcol <- 100
  color <- jet.colors(nbcol)
  # Compute the z-value at the facet centres
  zfacet <- z[-1, -1] + z[-1, -ncz] + z[-nrz, -1] + z[-nrz, -ncz]
  # Recode facet z-values into color indices
  facetcol <- cut(zfacet, nbcol)
  persp(x1, x2, z, col = color[facetcol], phi = 30, theta = 30)
  filled.contour(x1, x2, z, color.palette = jet.colors)
}


## Parallel computation  
## works better when each evaluation takes longer
## here we have added extra time to the computations
## just to show that it works
if (FALSE) { 
  res7 <- optim_ARS(c(-4,4),function(x){Sys.sleep(0.01); Rastrigin(x[1], x[2])},
                    lower=-5.12, upper=5.12)
  res8 <- optim_ARS(c(-4,4),function(x){Sys.sleep(0.01); Rastrigin(x[1], x[2])},
                    lower=-5.12, upper=5.12,parallel = T)
  res9 <- optim_ARS(c(-4,4),function(x){Sys.sleep(0.01); Rastrigin(x[1], x[2])},
                    lower=-5.12, upper=5.12,parallel = T,parallel_type = "snow")
}