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** **M**athematical **P**rogramming **L**anguage (**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?

**N**etwork-**E**nabled **O**ptimization **S**ystem (**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.

You can include your email address in case you want to receive the results. Note that for some solvers you *must* include it. If not, NEOS provides your job with an identifier and a password to have access to the outputs.

## 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.

### Forget about the **dat** file

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
```

## References

- AMPL: https://ampl.com/
- The NEOS Server: https://neos-server.org/neos/
- Fourer, R., Gay, D.M. and Kernighan, B.W. AMPL: A Modeling Language for Mathematical Programming Thompson/Brooks/Cole, 2nd ed., 2003.