Optimizing using AMPL & NEOS Server

In this session we will learn how to model and solve optimization problems of different types (linear, nonlinear, continuous and/or mixed integer) using AMPL as modeling language and some off-the-shelf solvers through the NEOS Server.

What is AMPL?

A Mathematical Programming Language (AMPL) is an algebraic modeling language broadly used to manipulate data and implement mathematical optimization models and strategies. It was developed by Robert Fourer, David Gay, and Brian Kernighan at Bell Laboratories in the 80s, and their book, AMPL: A Modeling Language for Mathematical Programming, is a complete guide for modelers at all levels of experience. Both the book and all the examples it contains are publicly available in The AMPL Book.

AMPL sintax is similar to usual notation of mathematical optimization models in which you describe data, variables, objectives and constraints. It also provides commands for debugging models and analyzing results. In order to solve these models, more than 40 open source and commercial solvers are supported. One of its advantages is that AMPL jobs (of reasonable sizes) can be implemented and solved “in the cloud” in collaboration with the NEOS Server project.

What is the NEOS Server?

Network-Enabled Optimization System (NEOS) Server is a free internet-based service for solving numerical optimization problems, which is hosted by the Wisconsin Institute for Discovery at the University of Wisconsin in Madison.

NEOS provides access to more than 60 state-of-the-art solvers in more than a dozen optimization categories. Jobs can be submitted directly to NEOS, via its web page or email, or indirectly, using for instance Pyomo or R, using the rneos package. An appropriate remote computer on which to run it is sought (holding the job in a queue if necessary). When a remote computer is available, it reads and executes the submitted files, runs the requested solver, and gathers the resulting output. The Server reports the output on the submission website and optionally emails a copy to the user. Most of the solvers in NEOS admit and AMPL input.

Build a simple NEOS job using AMPL

First, we need to create three files using the AMPL sintax. These files contain, respectively, the optimization model (mod file), the data to be used by the model (dat file) and the commands to be executed (run file). Second, we select a suitable solver from the NEOS pool to solve the optimization problem at hand. Finally, the AMPL files are uploaded and executed in a server.

To illustrate the process we will use the well-known “diet problem”, which finds a mix of foods that satisfies some requirements on the amounts of various vitamins minimizing the total cost. This problem is fully described in Chapter 2 in the AMPL book.

Let $J$ be the set of foods and $I$ the set of vitamins, $x_j$, the amount of food of type $j$, $c_j$ its cost and $v_{ij}$ the amount of vitamin $i$ in food $j$, for $i\in I$ and $j\in J$. The diet problem is stated as:

$$\begin{array}{ll} \min \,\, \displaystyle\sum_{j \in J} c_jx_j & \\ \text{s.t.} & \\ l_i \leq \displaystyle\sum_{j \in J} v_{ij}x_j \leq u_i & \,\, i\in I, \\ a_j \leq x_j \leq b_j & \,\, j\in J, \end{array}$$

where $l_i$ and $u_i$ denote the lower and upper bounds for the amount of vitamin $i$, for $i\in I$.

The dat file

The parameters to be used in the optimization model are included in this file. We will learn how to declare numbers, vectors and matrices.

Numbers

To declare a numerical parameter we use the keyword param, followed by a name, the “:=” symbols, its value and “;”.

param n:=8; #number of foods
param m:=6; #number of vitamins


Vectors

To declare a numerical vector we use the keyword param, followed by a name, the “:=” symbols, a column of indeces and values and “;”. One can even declare more than one vector at once, using different names for each column.

param c:=
1 3.19
2 2.59
3 2.29
4 2.89
5 1.89
6 1.99
7 1.99
8 2.49
;

param: a b:=
1 2	10
2 2	10
3 2	10
4 2	10
5 2	10
6 2	10
7 2	10
8 2	10
;

param: l u:=
1 700   20000
2 700   20000
3 700   20000
4 700   20000
5 0     50000
6 16000 24000
;



Matrices

To declare a numerical matrix we use the keyword param, followed by “:“, a list of indices for the columns of the matrix (either characters or numbers), “:=” symbols, a column of indeces for the rows, the matrix itself and “;”.

param v:   1  2 3 4 5 6 7 8:=
1   60	8     8     40	15	    70	25	    60
2   20	0     10	40	35	    30	50	    20
3   10	20    15	35	15	    15	25	    15
4   15	20    10	10	15	    15	15	    10
5   938	2180  945	278	1182    896	1329    1397
6   295	770   440	430	315	    400	370	    450
;


The mod file

Parameters and Sets

The list of the parameters’ names included in the dat file, which are used by the optimization model, must be included in this file, together with some requirements if necessary (positive, integer, etc.) and the indexing sets. To declare a one-dimensional set we use the keyword set, followed by a name, the “:=” symbols, the elements of the set and “;”.

param n > 0, integer;
param m > 0, integer;

set J := 1..n;
set I := 1..m;

param c{J} >= 0;
param l{I} >= 0;
param u{i in I} >= l[i];
param a{J} >= 0;
param b{j in J} >= a[j];
param v{I,J} >=0 ;



Variables

Variables are declared with the keyword var, a name, an indexing set if necessary, the variable (binary, integer) type and/or bounds and “;”.

var x {j in J} >= a[j], <=b[j];


Objective function

First, we have to specify the optimization sense with the keywords minimize or maximize followed by a name for the objective and “:“. Second, the expression of the function to optimize is typed using the appropriate commands of AMPL.

minimize Cost:  sum {j in J} c[j] * x[j];


Constraints

To introduce constraints in the model we use the keyword subject to, a name, an indexing set if necessary, “:“, and the expression of the constraint using “<=“, “>=” or “==” to set the sense. Finish the line with “;”.

subject to Diet_ub {i in I}:
sum {j in J} v[i,j] * x[j] <= u[i];

subject to Diet_lb {i in I}:
sum {j in J} v[i,j] * x[j] >= l[i];



The run file

To solve the problem the keyword solve must be included in this file, together with commands that, for instance, print the solution of the problem.

solve;
display Cost;
display x;


Submit the job to NEOS

To solve our instance of the diet problem, we first need to select a linear programming solver from the list, then choose the AMPL input option, upload your three files (data, mod and run) and submit.

Build a complex NEOS job using AMPL

In order to illustrate some additional tools of AMPL and NEOS, let us consider the Lasso Regression optimization problem defined as

$$\displaystyle\min_\beta \,\, |Y - X\beta|_2 + \lambda|\beta|_1,$$

where $\lambda\geq 0$ is the penalty parameter. In this case, this is an unconstrained nonlinear convex optimization problem, so we can choose any of the solvers from the corresponding list in NEOS.

The dat file is not really needed. The parameters used by the optimization model can be declared directly in the mod file as we did in the dat file. However, if we proceed this way, those values cannot be modified to, for instance, solve the optimization problem for different values of $\lambda$. Declaring the parameters as default is a possible way to avoid this drawback. Using these tools, the mod file for the Lasso Regression problem is:

# Problem size
param m default 0, integer; # number of rows
param n default 0, integer; # number of columns

# Sets
set M := 1..m;
set N := 1..n;

# Problem data
param X {M, N} default 0;
param Y {M} default 0;

# Penalty parameter
param Lambda >= 0, default 1;

# Variables (beta = xp - xn)
var xp{N} >= 0;
var xn{N} >= 0;

# Objective function
minimize LassoReg: sum{i in M}( sum{j in N}( X[i,j]*(xp[j]-xn[j]) - Y[i]) )^2
+ Lambda*sum{j in N}( xp[j] + xn[j] );


Then, in the run file we can declare the parameters values using the keyword let. If we were working in local, we could read the values from a file and assign them to the parameters also in the run. For our example, we randomly generate $Y$ and $X$, using the Normal and Uniform distributions and a for loop, and choose $m=50$ and $n=128$. The run file is:

let m:=50;
let n:=128;

for {i in 1..m}{

let Y[i]:=Uniform(-3,3);

for{j in 1..n}{
let X[i,j]:=Normal(0,5);

}
}

solve;

display {j in M}(xp[j]-xn[j]);


Solving many problems in one submission

If we are interested in analyzing the behaviour of the solution for different values of the penalty parameter $\lambda$, we can solve as many cases as we want in the same NEOS submission. The set of values for $\lambda$ are defined in the run file, as well as the corresponding assignment for every case within a for loop. The run file in this case is:

let m:=50;
let n:=128;

for {i in 1..m}{

let Y[i]:=Uniform(-3,3);

for{j in 1..n}{
let X[i,j]:=Normal(0,5);

}
}

set LambdaSet := {1, 10, 100, 1000, 10000};

for {l in LambdaSet}{

let Lambda := l;
solve;
display {j in N}(xp[j]-xn[j]);

}


Analyzing the output

We could be interested in analyzing the number of nonzero elements in $\beta$ for each value of $\lambda$. Declaring a new parameter to keep that number and using an if statement, those numbers can be easily saved. In this case, the run file is:

let m:=50;
let n:=128;

for {i in 1..m}{

let Y[i]:=Uniform(-3,3);

for{j in 1..n}{
let X[i,j]:=Normal(0,5);

}
}

set LambdaSet := {1, 10, 100, 1000, 10000};

param nnzs{LambdaSet} default 0;

for {l in LambdaSet}{

let Lambda := l;
solve;

let nnzs[l] := 0;
for {i in 1..n}{
if (abs(xp[i]-xn[i])>0) then let nnzs[l] := nnzs[l]+1;
}

}

display nnzs;


Initial guesses (to help the solvers) and variable fixing

It is possible to define a starting solution using also the let command with the variables defined in your models. We can even fix the value for a certain variable using fix command instead, and unfix it using unfix.

Let us illustrate the fix/unfix commands by means of the Uncapacitated Facility Location Problem (UFLP):

$$\begin{array}{ll} \min \,\, \displaystyle\sum_{i=1}^n\sum_{j=1}^m c_{ij}y_{ij} + \sum_{i =1}^n f_ix_i & \\ \text{s.t.} & \\ \displaystyle\sum_{i=1}^n y_{ij} = 1& \,\, j=1,\ldots,m, \\ y_{ij} \leq x_i & \,\, i=1,\ldots,n,\,j=1,\ldots,m, \\ x_i,y_{ij}\in\lbrace0,1\rbrace & \,\, i=1,\ldots,n,\,j=1,\ldots,m. \\ \end{array}$$

Let us solve $n$ UFLPs, imposing for each of them that facility $i$ must be open, i.e. $x_i=1$ ($i=1,\ldots,n$).

The mod file states as:

# Problem size
param m default 0, integer; # number of demand points
param n default 0, integer; # number of facilities

# Sets
set M := 1..m;
set N := 1..n;

# Problem data
param c {N, M} default 0;
param f {N} default 0;

# Variables
var y{N,M}, binary;
var x{N}, binary;

# Objective function
minimize costs: sum{i in N, j in M} c[i,j]*y[i,j] + sum{i in N} f[i]*x[i];

s.t.

c1 {j in M}: sum{i in N} y[i,j] == 1;
c2 {i in N, j in M}: y[i,j] <= x[i];


Whereas the run file states as:

let m:=50;
let n:=100;

for {j in 1..m}{
let f[j]:=Uniform(1,10);
for{i in 1..n}{
let c[i,j]:=Uniform(1,10);
}
}

for {i in N}{

fix x[i] := 1;
solve;
display costs;
unfix x[i];

} # end for