Endogenous if-then-else

> I am trying to define a nonlinear function where the variable in the
> objective is dependent on a parameter (x) which varies across
> i=1,2,...,N cases and two variables, a cut-off value (cutoff) and a
> coefficient (coeff) which transforms x above the value of cutoff. The
> logic, in pseudocode is:
>
> if x[i] > cutoff {
>     x2[i] = x[i] * coeff
>     }
> else {
>     x2[i] = x[i]
>     }
>
> The objective function would be something like:
>
> minimize Objective: sum {i in N} [complicated function goes here] * x2[i] ;
>
> How, if at all, can this be implemented in AMPL?
>

You could try:
e{i in I}: x2[i] = if (x[i]>cutoff) then x[i]*coeff else x[i];
However this is a discontinuous function so beware that many NLP solvers assume smooth functions and may get into trouble near the jump. A smooth approximation would be better in that case (there are several possible functional forms for such a smooth approximation; one would need to invest some time to try a few of them; an example is shown here). If the model was otherwise linear, one could look into a piecewise linear formulation (AMPL has very good support for this; this approach is less interesting if the model is otherwise nonlinear as it would become an MINLP).

Indexing by a variable

> I've:
>
> param N_MAX = 1000;
> var N >= 1, <= N_MAX;
>
> var x {0 .. N_MAX} >= 0;
>
> I want to do *something like*:
>
> subject to {i in  {(N/2)-1, (N/2)+1} }: #calculation of indices
> simplified for sake of discussion
>     x[i] = 1000; #RHS simplified for sake of discussion
>
> But of course I can't use variables for indexing. Neither can I do it
> explicitly as there are too many possible values that N can take. So,
> I'm looking for suggestions.
>
> Sincerely,

param N_MAX := 10;

set I := {1..N_MAX};
set J := {0..N_MAX};
set IJ{i in I} := setof{j in J: ((j >= i/2-1) and (j <= i/2+1))}(j);
display IJ;

var N >= 1, <= N_MAX;
param xup := 10000;
var x {J} >= 0, <= xup;
var b {I} binary;

minimize obj: sum{j in J} x[j];
s.t. e1: sum{i in I} i*b[i] = N;
s.t. e2: sum{i in I} b[i] = 1;
s.t. e3{i in I, j in IJ[i]}:  x[j] <= 1000*b[i] + xup*(1-b[i]);
s.t. e4{i in I, j in IJ[i]}:  x[j] >= 1000*b[i];

solve;
display N,b,x;
end;

Usually it is best to introduce a binary variable b[i] with b[i]=1 <=> i=N. Then you can introduce constraints based on b[i]. xup is here an upper bound on x[j]. Choose this number as tight as possible: it is a big-M constant. As usual it is a good idea to move as much as possible complexity away from the constraints. In this case we use an intermediate set IJ to handle a more difficult concept. Such a set can be tested and debugged before a complete model is ready to run.

AMPL/Cplex: No solver needed

 

Hello,

Anyone can help?
I have run the following model using AMPL code with CPLEX 11.2 solver
in window xp. The program run well and take a couple seconds for the
small size of problem. But once the problem size is more than 200 (of
j and k), the solve time increases dramatically. When the problem size
of set j and set k is 900 ,  it took so long ( longer than 10 hrs) to
find the solution. I wonder, since the model looks simple with only 2
variables per constraint. I'm not sure did I write the AMPL code wrong
or there is something wrong about this. (However, I checked the
solution from CPLEX after it done. The solution is correct). Thank you
very much.

Best regards,
Alex

The Mathematical Model is:

Objective: Min Delta
S.t. (1) Xj *(Bj + Sum(i) M(ICij + TC1ij + TC2ijk)) <= Delta,     For all j, k

(2) Sum(i) Xj = 1

(3) Xj = {0, 1}

The AMPL code is:

set stag ;
set sup;
set w;

param tc1{i in sup,j in stag}>=0;
param b{j in stag}>=0;
param ic {i in sup,j in stag}>=0;
param M integer;
param tc2{i in sup,j in stag,k in w}>=0;

var x{j in stag} binary;
var Delta;

minimize cost: Delta;
subject to delt{k in w, j in stag}:  x[j]*(b[j]+ sum{i in sup} (M*(ic[i,j]+ tc1[i,j]+ tc2[i,j,k]) )) <= Delta;
subject to c: sum{j in stag} x[j] = 1;

Ten hours is a long time for 900 binary variables in this rather simple model. The problem has many rows as a result of the way the max is formulated. This construct where we minimize the max is sometimes difficult for MIP solvers. Note that some MIP models are just very difficult to solve (example). However in this case, I believe this problem can be solved extremely fast without even using Cplex. First we calculate a[j,k] such that we can simplify the problem to:

minimize cost: Delta;
subject to delt{k in w, j in stag}:  x[j]*a[j,k] <= Delta;
subject to c: sum{j in stag} x[j] = 1;

The first constraint can be interpreted as: if (x[j]=1) then Delta = max_k a[j,k]. The complete model can now be restated as:

param a{j in stag, k in w} := b[j] + sum{i in sup} (M*(ic[i,j]+ tc1[i,j]+ tc2[i,j,k]));
param d{j in stag} := max{k in w} a[j,k];
param dmin := min{j in stag} d[j];
var x{j in stag} binary := if (d[j]=dmin) then 1 else 0;
var Delta := dmin;
display x,Delta;

(I left as an exercise to make sure there is only one x[j]=1 in case multiple d[j]’s assume the value dmin). If you really want to call Cplex (especially after shelling out big bucks for it), the following should be very fast:

var x{j in stag} binary;
var Delta;

minimize cost: Delta;
subject to c: sum{j in stag} x[j] = 1;
subject to c2: sum{j in stag} x[j]*d[j] = dmin;
subject to c3: Delta = dmin;

option solver cplex;
solve;

This formulation automatically deals with the problem of multiple d[j]’s being equal to the minimum.

AMPL vs. GAMS on radiotherapy treatment planning

 

hi all,
i would just like to solicit your advice regarding whether to get an
academic license for AMPL or GAMS. i'm currently doing a master's
thesis on cancer radiotherapy treatment planning (RTP), and i've yet
to decide between AMPL or GAMS for prototyping the models. some of my
decision criteria include language syntax, user interface, software
cost and available literature on RTP.
i think that the syntax for AMPL is better compared to GAMS because
the former more closely resembles the algebraic model formulations.
however, there seem to be more literature on RTP in GAMS and there's
also an integrated development environment for it. thanks for your
help.

 

Both systems are very good in what they are doing and widely used, so you cannot really go wrong with either choice. I would probably suggest to add extra points for the modeling system that is used by your colleagues and collaborators. That makes exchanging models and data easier and also is easier when discussing problems, tricks, issues etc.

Bob Fourer (AMPL) answered:

It's hard to find someone who can give equally expert advice on two competing
systems, as once you become familiar with one of them you don't usually have
much incentive to keep learning about the other.  But here are a few comments
from my hardly unbiased view.

AMPL was designed with the idea of being much closer to mathematical notation
and generally much more natural to use than GAMS, and it's superior on that
score.  A GAMS model typically relies on more special conventions and
reformulations than its AMPL counterpart; a case in point is the often extensive
use of the GAMS $ operator to impose various conditions.  Also, the IDE
notwithstanding, GAMS is fundamentally more of a batch system whereas AMPL
offers a more flexible option of interactively exploring models and results.
Finally while in certain areas GAMS is established through long use, still I see
modelers in these areas choosing AMPL, particularly when they are undertaking
new projects that do not depend on existing GAMS models.

In my opinion AMPL and GAMS are closer in practice than suggested here (e.g. where you use $ in GAMS, one would use : in AMPL). I actually slightly prefer the GAMS syntax when doing real work, as it is a little bit more compact and it is obvious where a summation ends (in AMPL this is based on operator priority, in GAMS a sum is visually bracketed by parentheses).

AMPL NLP question

> I want to formulate some equations like ln(z(k)) <= D, 
> z(k) is a variable. Is there some function or keyword used as logarithm? 
> Or is there some method to calculate logarithm? Thanks 

You can use the log() function. Note that if D is a parameter the constraint can better be formulated as z[k] <= exp(D). This makes it a simple linear constraint that can be converted into a bound.

The difference can be illustrated with a small example:

ampl: var x := 1;
ampl: param D := 3;
ampl: maximize obj: x;
ampl: e: log(x) <= D;
ampl: solve;
MINOS 5.5: optimal solution found.
11 iterations, objective 20.08553692
Nonlin evals: constrs = 42, Jac = 41.
-----------------------------------------
ampl: var x := 1;
ampl: param D := 3;
ampl: maximize obj: x;
ampl: e: x <= exp(D);
ampl: solve;
MINOS 5.5: optimal solution found.
1 iterations, objective 20.08553692

AMPL has a very good presolver, but not so good this reformulation is applied automatically.

Converting MathProg (and AMPL) models to GAMS

The tool glpk2gams can convert models written in MathProg (subset of AMPL) so we can quickly assess how GAMS solvers perform on such a model.

If the model only has a .mod file, you can just use a single argument:

D:\glpk\glpk2gams\Debug>glpk2gams.exe diet.mod
Glpk Model filename: diet.mod
Gams filename: diet.gms
  
Reading model section from diet.mod...
Reading data section from diet.mod...
99 lines were read
Generating nb...
Generating cost...
Model has been successfully generated
  
Rows : 10  Columns : 20  Nonzero elements : 179
Writing diet.gms
  
D:\glpk\glpk2gams\Debug>

If the model also has a .dat file, you need to give all three command line arguments:

L:\glpk>glpk2gams egypt2.mod egypt2.gms egypt2.dat
Glpk Model filename: egypt2.mod
Gams filename: egypt2.gms
Data filename: egypt2.dat
 
Reading model section from egypt2.mod...
egypt2.mod:264: warning: unexpected end of file; missing end statement inserted
264 lines were read
Reading data section from egypt2.dat...
egypt2.dat:276: warning: unexpected end of file; missing end statement inserted
276 lines were read
Generating Psi...
Generating mbd...
Generating mbdb...
Generating mb...
Generating cc...
Generating ap...
Generating al...
Generating ai...
Model has been successfully generated
  
Rows : 285  Columns : 351  Nonzero elements : 1336
Writing egypt2.gms
 
L:\glpk>

It is generating scalar GAMS code, so this is not suited to convert models that have to be maintained in GAMS.

AMPL defined variables

AMPL has a construct called defined variables that can sometimes reduce the size of an NLP by a large amount. Here are some statistics for a client model (confidential so I cannot share it):

GAMS:

BLOCKS OF EQUATIONS           6     SINGLE EQUATIONS       26,192
BLOCKS OF VARIABLES           5     SINGLE VARIABLES       26,201
NON ZERO ELEMENTS       119,031     NON LINEAR N-Z         45,220
DERIVATIVE POOL               8     CONSTANT POOL           1,541
CODE LENGTH             359,381

AMPL:

Substitution eliminates 26180 variables.
Adjusted problem:
10 variables, all nonlinear
1 constraint, all linear; 10 nonzeros 
1 nonlinear objective; 10 nonzeros.

GAMS has some similar facilities built in the LGO solver, but these are not available as a general feature.