R interface to NOMAD

Description

snomadr is an R interface to NOMAD (Nonsmooth Optimization by Mesh Adaptive Direct Search, Abramson, Audet, Couture and Le Digabel (2011)), an open source software C++ implementation of the Mesh Adaptive Direct Search (MADS, Le Digabel (2011)) algorithm designed for constrained optimization of blackbox functions.

NOMAD is designed to find (local) solutions of mathematical optimization problems of the form

1
2
3
4
5
   min     f(x)
x in R^n

s.t.       g(x) <= 0 
           x_L <=  x   <= x_U

where f(x): R^n --> R^k is the objective function, and g(x): R^n --> R^m are the constraint functions. The vectors x_L and x_U are the bounds on the variables x. The functions f(x) and g(x) can be nonlinear and nonconvex. The variables can be integer, continuous real number, binary, and categorical.

Kindly see http://www.gerad.ca/nomad/Project/Home.html and the references below for details.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
snomadr(eval.f, 
        n, 
        bbin = NULL,  
        bbout = NULL, 
        x0 = NULL, 
        lb = NULL, 
        ub = NULL, 
        nmulti = 0,
        random.seed = 0,
        opts = list(),
        print.output = TRUE, 
        information = list(), 
        snomadr.environment = new.env(), 
        ... ) 

Arguments

eval.f

function that returns the value of the objective function

n

the number of variables

bbin

types of variables. Variable types are 0 (CONTINUOUS), 1 (INTEGER), 2 (CATEGORICAL), 3 (BINARY)

bbout

types of output of eval.f. See the NOMAD User Guide http://www.gerad.ca/NOMAD/Downloads/user_guide.pdf

x0

vector with starting values for the optimization. If it is provided and nmulti is bigger than 1, x0 will be the first initial point for multiple initial points

lb

vector with lower bounds of the controls (use -1.0e19 for controls without lower bound)

ub

vector with upper bounds of the controls (use 1.0e19 for controls without upper bound)

nmulti

when it is missing, or it is equal to 0 and x0 is provided, snomadRSolve will be called to solve the problem. Otherwise, smultinomadRSolve will be called

random.seed

when it is not missing and not equal to 0, the initial points will be generated using this seed when nmulti > 0

opts

list of options for NOMAD, see the NOMAD user guide http://www.gerad.ca/NOMAD/Downloads/user_guide.pdf. Options can also be set by nomad.opt which should be in the folder where R (snomadr) is executed. Options that affect the solution and their defaults and some potential values are

"MAX_BB_EVAL"=10000

"INITIAL_MESH_SIZE"=1

"MIN_MESH_SIZE"="r1.0e-10"

"MIN_POLL_SIZE"="r1.0e-10"

Note that the "r..." denotes relative measurement (relative to lb and ub)

Note that decreasing the maximum number of black box evaluations ("MAX_BB_EVAL") will terminate search sooner and may result in a less accurate solution. For complicated problems you may want to increase this value. When experimenting it is desirable to set "DISPLAY_DEGREE"=1 (or a larger integer) to get some sense for how the algorithm is progressing

print.output

when FALSE, no output from snomadr is displayed on the screen. If the NOMAD option "DISPLAY_DEGREE"=0, is set, there will also be no output from NOMAD. Higher integer values for "DISPLAY_DEGREE"= provide successively more detail

information

is a list. snomadr will call snomadRInfo to return the information about NOMAD according to the values of "info", "version" and "help".

"info"="-i": display the usage and copyright of NOMAD

"version"="-v": display the version of NOMAD you are using

"help"="-h": display all options

You also can display a specific option, for example, "help"="-h x0", this will tell you how to set x0

snomadr.environment

environment that is used to evaluate the functions. Use this to pass additional data or parameters to a function

...

arguments that will be passed to the user-defined objective and constraints functions. See the examples below

Details

snomadr is used in the crs package to numerically minimize an objective function with respect to the spline degree, number of knots, and optionally the kernel bandwidths when using crs with the option cv="nomad" (default). This is a constrained mixed integer combinatoric problem and is known to be computationally ‘hard’. See frscvNOMAD and krscvNOMAD for the functions called when cv="nomad" while using crs.

However, the user should note that for simple problems involving one predictor exhaustive search may be faster and potentially more accurate, so please bear in mind that cv="exhaustive" can be useful when using crs.

Naturally, exhaustive search is also useful for verifying solutions returned by snomadr. See frscv and krscv for the functions called when cv="exhaustive" while using crs.

Value

The return value contains a list with the inputs, and additional elements

call

the call that was made to solve

status

integer value with the status of the optimization

message

more informative message with the status of the optimization

iterations

number of iterations that were executed, if multiple initial points are set, this number will be the sum for each initial point.

objective

value if the objective function in the solution

solution

optimal value of the controls

Author(s)

Zhenghua Nie <niez@mcmaster.ca>

References

Abramson, M.A. and C. Audet and G. Couture and J.E. Dennis Jr. and S. Le Digabel (2011), “The NOMAD project”. Software available at http://www.gerad.ca/nomad.

Le Digabel, S. (2011), “Algorithm 909: NOMAD: Nonlinear Optimization With The MADS Algorithm”. ACM Transactions on Mathematical Software, 37(4):44:1-44:15.

See Also

optim, nlm, nlminb

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
## Not run: 
## List all options
snomadr(information=list("help"="-h"))

## Print given option,  for example,  MESH_SIZE
snomadr(information=list("help"="-h MESH_SIZE"))

## Print the version of NOMAD
snomadr(information=list("version"="-v"))

## Print usage and copyright
snomadr(information=list("info"="-i"))

## This is the example found in
## NOMAD/examples/basic/library/single_obj/basic_lib.cpp

eval.f <- function ( x ) {

    f <- c(Inf, Inf, Inf);
    n <- length (x);

    if ( n == 5 && ( is.double(x) || is.integer(x) ) ) {
        f[1] <- x[5];
        f[2] <- sum ( (x-1)^2 ) - 25;
        f[3] <- 25 - sum ( (x+1)^2 );
    }  

    return ( as.double(f) );
}

## Initial values
x0 <- rep(0.0, 5 )

bbin <-c(1, 1, 1, 1, 1)
## Bounds
lb <- rep(-6.0,5 )
ub <- c(5.0, 6.0, 7.0, 1000000, 100000)

bbout <-c(0, 2, 1)
## Options
opts <-list("MAX_BB_EVAL"=500,
            "MIN_MESH_SIZE"=0.001,
            "INITIAL_MESH_SIZE"=0.1,
            "MIN_POLL_SIZE"=0.0001)

snomadr(eval.f=eval.f,n=5,  x0=x0, bbin=bbin, bbout=bbout, lb=lb, ub=ub, opts=opts)


## How to transfer other parameters into eval.f
##
## First example: supply additional arguments in user-defined functions
##

## objective function and gradient in terms of parameters
eval.f.ex1 <- function(x, params) { 
    return( params[1]*x^2 + params[2]*x + params[3] ) 
}

## Define parameters that we want to use
params <- c(1,2,3)

## Define initial value of the optimization problem
x0 <- 0

## solve using snomadr 
snomadr( n          =1, 
        x0          = x0, 
        eval.f      = eval.f.ex1, 
        params      = params )


##
## Second example: define an environment that contains extra parameters
##

## Objective function and gradient in terms of parameters
## without supplying params as an argument
eval.f.ex2 <- function(x) { 
    return( params[1]*x^2 + params[2]*x + params[3] ) 
}

## Define initial value of the optimization problem
x0 <- 0

## Define a new environment that contains params
auxdata        <- new.env()
auxdata$params <- c(1,2,3)

## pass The environment that should be used to evaluate functions to snomadr 
snomadr(n                  =1, 
        x0                 = x0, 
        eval.f             = eval.f.ex2, 
        snomadr.environment = auxdata )

## Solve using algebra
cat( paste( "Minimizing f(x) = ax^2 + bx + c\n" ) )
cat( paste( "Optimal value of control is -b/(2a) = ", -params[2]/(2*params[1]), "\n" ) )
cat( paste( "With value of the objective function f(-b/(2a)) = ",
           eval.f.ex1( -params[2]/(2*params[1]), params ), "\n" ) )

## The following example is NOMAD/examples/advanced/multi_start/multi.cpp
## This will call smultinomadRSolve to resolve the problem.  
eval.f.ex1 <- function(x, params) { 
    M<-as.numeric(params$M)
    NBC<-as.numeric(params$NBC)

    f<-rep(0, M+1)
    x<-as.numeric(x)

    x1 <- rep(0.0, NBC)
    y1 <- rep(0.0, NBC)

    x1[1]<-x[1]
    x1[2]<-x[2]
    y1[3]<-x[3]
    x1[4]<-x[4]
    y1[4]<-x[5]

    epi <- 6
    for(i in 5:NBC){
        x1[i]<-x[epi]
        epi <- epi+1
        y1[i]<-x[epi]
        epi<-epi+1
    }
    constraint <- 0.0
    ic <- 1
    f[ic]<-constraint
    ic <- ic+1

    constraint <- as.numeric(1.0)
    distmax <- as.numeric(0.0)
    avg_dist <- as.numeric(0.0)
    dist1<-as.numeric(0.0)

    for(i in 1:(NBC-1)){
        for (j in (i+1):NBC){
            dist1 <- as.numeric((x1[i]-x1[j])*(x1[i]-x1[j])+(y1[i]-y1[j])*(y1[i]-y1[j]))
            
            if((dist1 > distmax)) {distmax <- as.numeric(dist1)}
            if((dist1[1]) < 1) {constraint <- constraint*sqrt(dist1)}
            else if((dist1) > 14) {avg_dist <- avg_dist+sqrt(dist1)}
        }
    }

    if(constraint < 0.9999) constraint <- 1001.0-constraint
    else constraint = sqrt(distmax)+avg_dist/(10.0*NBC)

    f[2] <- 0.0
    f[M+1] <- constraint 


    return(as.numeric(f) ) 
}

## Define parameters that we want to use
params<-list()
NBC <- 5
M <- 2
n<-2*NBC-3

params$NBC<-NBC
params$M<-M
x0<-rep(0.1, n)
lb<-rep(0, n)
ub<-rep(4.5, n)

eval.f.ex1(x0, params)

bbout<-c(2, 2, 0)
nmulti=5
bbin<-rep(0, n)
## Define initial value of the optimization problem

## Solve using snomadRSolve
snomadr(n            = as.integer(n), 
        x0           = x0, 
        eval.f       = eval.f.ex1, 
        bbin         = bbin, 
        bbout        = bbout, 
        lb           = lb, 
        ub           = ub, 
        params       = params )

## Solve using smultinomadRSolve, if x0 is provided,  x0 will
## be the first initial point,  otherwise,  the program will
## check best_x.txt,  if it exists,  it will be read in as
## the first initial point. Other initial points will be
## generated by uniform distribution.
## nmulti represents the number of mads will run.
##
snomadr(n            = as.integer(n), 
        eval.f       = eval.f.ex1, 
        bbin         = bbin, 
        bbout        = bbout, 
        lb           = lb, 
        ub           = ub, 
        nmulti = as.integer(nmulti), 
        print.output = TRUE, 
        params       = params )

## End(Not run)