**a**which has not been declared. For the first error the macro

**%x%**causes the cursor to point to

**parameterb**.

In case of a syntax error, the GAMS/IDE will usually place the cursor on the correct place. This makes it easy to fix errors quickly without even consulting the listing file. Here is a case where the GAMS/IDE gets confused due to substitution of a macro:

In this case the cursor should go to symbol **a** which has not been declared. For the first error the macro **%x%** causes the cursor to point to **parameterb**.

Follow up on http://yetanothermathprogrammingconsultant.blogspot.com/2008/05/cutting-stock-problem-example.html.

*>I decided to have a go at the challenge set in this paper - hope you don't mind.>*

> [length, quantity] pairs:[100,100],[200,2],[300,2],[493,250],[590,2],

> [630,2],[780,2],[930,2]

> [100,100,100,100,100,630]2 * [100,100,930]2 * [100,300,780]2 *

>[100,493,590]118 * [200,493,493]

>

>

>(100*100+200*2+300*2+493*250+590*2+630*2+780*2+930*2)/1200=116.76

>

>2 * 1 + 2 * 1 = 102200: 118 * 1 = 118300: 2 * 1 = 2493: 118 * 2 + 2 * 1 + 12 * 1 =

>250590: 2 * 1 = 2630: 2 * 1 = 2780: 2 * 1 = 2930: 2 * 1 = 2

See also: http://people.sabanciuniv.edu/ertekg/papers/2006/alp_et_al_ims2006.pdf.

I can find a better solution using the column generation model. The following input was used:

`*-----------------------------------------------------`

* Data

*-----------------------------------------------------

scalar r 'raw width' /1200/;

table demanddata(i,*)

width demand

width1 100 100

width2 200 2

width3 300 2

width4 493 250

width5 590 2

width6 630 2

width7 780 2

width8 930 2

;

The result is a solution (not proven optimal) of 131:

`---- 174 PARAMETER pat pattern usage`

p4 p5 p9 p10 p11 p15 Total

width1 1 2 100

width2 1 2

width3 1 2

width4 2 1 2 250

width5 2 2

width6 1 2

width7 1 2

width8 1 2

Count 75 1 2 2 2 49 131

Mosek does not return the best integer solution found sofar when interrupted. Below we interrupted a MIP run, while integer solutions were found by the solver (the model was no longer intermediate integer infeasible).

I don't know if this is caused by a misbehavior of the link or by Mosek itself. Work around: when stopped on a RESLIM limit, GAMS/Mosek will report the best integer solution. Don't press the Interrupt button if you don't want to loose your solution.

Note also that the iteration count is reported as zero, which is another bug.

`**** SOLVER STATUS 8 USER INTERRUPT `

**** MODEL STATUS 14 NO SOLUTION RETURNED

**** OBJECTIVE VALUE 3885.0000

RESOURCE USAGE, LIMIT 6.618 3600.000

ITERATION COUNT, LIMIT 0 10000000

MOSEK Link May 1, 2008 22.7.1 WIN 3927.4700 VIS x86/MS Windows

M O S E K version 5.0.0.79 (Build date: Jan 29 2008 13:41:45)

Copyright (C) MOSEK ApS, Fruebjergvej 3, Box 16

DK-2100 Copenhagen, Denmark

http://www.mosek.com

No solution returned

I don't know if this is caused by a misbehavior of the link or by Mosek itself. Work around: when stopped on a RESLIM limit, GAMS/Mosek will report the best integer solution. Don't press the Interrupt button if you don't want to loose your solution.

Note also that the iteration count is reported as zero, which is another bug.

The expression smin(s, xxx) over an empty set s may or may not return -INF. The following fragment from a real world scheduling model with a problem in its input data illustrates inconsistencies in how GAMS evaluates this.

This should give the same result, however when we run it we see:

This is of course highly undesirable. In general I would try to prevent such situations by not using SMIN or SMAX over an empty set. As GAMS does not flag these occurences you will need to add explicit tests for this.

####
Update

This problem has been fixed in GAMS release 24,3,1 (july 2014).

```
set si /i1,i2/;
parameter duedateadj(si) /i1 10, i2 10/;
set k /k1,k2/;
parameter duration(si,k) / i1.k1 1 /;
parameters
minduration(si)
mintime(si) "calculate in one step"
mintime2(si) "calculate in two steps"
;
mintime(si) = duedateadj(si)-smin(k$duration(si,k),duration(si,k));
minduration(si) = smin(k$duration(si,k),duration(si,k));
mintime2(si) = duedateadj(si)-minduration(si);
display mintime, mintime2;
```

This should give the same result, however when we run it we see:

```
---- 17 PARAMETER mintime calculate in one step
i1 9.000, i2 -INF
---- 17 PARAMETER mintime2 calculate in two steps
i1 9.000, i2 10.000
```

This is of course highly undesirable. In general I would try to prevent such situations by not using SMIN or SMAX over an empty set. As GAMS does not flag these occurences you will need to add explicit tests for this.

The date format is coded as the number of days since jan 1, 1900. However, it is one off compared to Excel, Access and compared to the date format in Delphi. It is also one off compared to the IDE Charting facility:

The result display is:

The chart is:

The dates on the x-axis are wrong (click to enlarge). The work around is to add one to each date before exporting (and subtract one after importing).

set

j /x,y/

i /i1,i2,i3/

;

parameter data(i,j);

* may 30, may 31, june 1

data(i,'x') = jdate(2008,5,30)+ord(i)-1;

data(i,'y') = ord(i);

* check dates

parameter check(i,*);

check(i,'year') = gyear(data(i,'x'));

check(i,'month') = gmonth(data(i,'x'));

check(i,'day') = gday(data(i,'x'));

option check:0:1:1;

display data,check;

execute_unload 'data',data;

The result display is:

`---- 20 PARAMETER data `

x y

i1 39597.000 1.000

i2 39598.000 2.000

i3 39599.000 3.000

---- 20 PARAMETER check

year month day

i1 2008 5 30

i2 2008 5 31

i3 2008 6 1

The chart is:

The dates on the x-axis are wrong (click to enlarge). The work around is to add one to each date before exporting (and subtract one after importing).

I received a very large spreadsheet to be converted to GAMS. To quickly get a feeling how many formulas there are and where they are located, this code colors all cells with formulas. The color depends on the type of the cell (number, string, logical or error).

An earlier version used HasFormula on each individual cell. As the spreadsheet I am working on is very large, that approach was very slow. This version returns complete ranges and is much faster.

`Sub color()`

ActiveSheet.UsedRange.Style = "Normal"

Call colorsub(xlNumbers, "Accent1")

Call colorsub(xlTextValues, "Accent2")

Call colorsub(xlErrors, "Accent3")

Call colorsub(xlLogical, "Accent4")

End Sub

Sub colorsub(SpecialCells As Integer, styletype As String)

Dim r As Range

On Error GoTo err

Set r = ActiveSheet.UsedRange.SpecialCells(xlCellTypeFormulas, SpecialCells)

r.Style = styletype

err:

End Sub

An earlier version used HasFormula on each individual cell. As the spreadsheet I am working on is very large, that approach was very slow. This version returns complete ranges and is much faster.

I would highly prefer that the listing file is self-contained. In many cases solvers stop without issuing proper messages to the listing file. Here is an example with PATHNLP:

There is no indication why PATHNLP stopped.

` S O L V E S U M M A R Y`

MODEL TIME OBJECTIVE LL

TYPE NLP DIRECTION MAXIMIZE

SOLVER PATHNLP FROM LINE 160

**** SOLVER STATUS 4 TERMINATED BY SOLVER

**** MODEL STATUS 7 INTERMEDIATE NONOPTIMAL

**** OBJECTIVE VALUE 47.2365

RESOURCE USAGE, LIMIT 35.697 1000.000

ITERATION COUNT, LIMIT 0 100000

EVALUATION ERRORS 0 0

PATH-NLP May 1, 2008 22.7.2 WEX 3906.4799 WEI x86_64/MS Windows

NLP size: 829 rows, 881 cols, 10145 non-zeros, 1.39% dense.

MCP size: 1708 rows/cols, 26626 non-zeros, 0.91% dense.

There is no indication why PATHNLP stopped.

When asked for an example of a small difficult MIP, the following model may be of interest:

http://www.amsterdamoptimization.com/models/etc/ran14x18.gms. Some notes on this model from a posting on sci.op-research:

*The story of the instance ran14x18 *

We tried CPLEX 5.0 on this instance, obtaining a solution with objective value 3712 but

without proving optimality. The same solution was also found in some runs of the genetic

algorithm from reference 1.

In response to my posting on sci.op-research, some people tried exact solvers on this

instance (fctp, mps), but without success. I also contacted CPLEX, but they also failed

to prove optimality.

On Friday, 13 November 1998, I announced this web page on sci.op-research. Some

hours later I received good news: Jeff proved optimality for ran14x18. Some time later,

a parallel CPLEX version also solved this problem to optimality. I thank the people cited

below for their efforts and comments concerning this problem!

To get an impression of the difficulty of this instance, have a look at the discussion

about the problem, which has remained unsolved for more than 10 months.

* *

Jeff Linderoth (*..@mcs.anl.gov**)*

*[Tue, 18 Aug 1998] Probably (repeat probably) a MIP solver that uses flow cover *

inequalities would be able to solve the instances. (At least it should do better than

CPLEX for problems of this class). If you don’t have too many instances, then you can

send me the MPS files, and I would be happy to try the solver MINTO out on them for you.

[Fri, 28 Aug 1998] For the problem ran14x18.mps, I ran MINTO until it ran out of memory,

and was only able to improve the lower bound to 3617.994472. When I get the opportunity,

I will perform some longer runs and also try out PARINO (which is like a parallel

version of MINTO) on the problem. I also tried XPRESS and CPLEX and MINTO was much

better than either of these packages (due to flow cover inequalities). How did bc-opt

do? I’ll email more results later. (If you don’t hear from me, send me an email to

remind me). Anyway, sorry I wasn’t able to prove optimality for you, but hopefully the

results are somewhat useful.

[Fri, 13 Nov 1998] I proved the optimality of the solution 3712 for the ran14x18 FCTP

instance. In order to do this I ran my PARINO parallel MIP code for about 52 hours on

32 Pentium II (300MHz) processors. 7624677 nodes of the branch and bound tree were

evaluated in this time. PARINO adds "flow cover inequalities" which are usually quite

useful in helping to improve the LP bounds for these types of problems. Also, I used a

parallel pseudocost-type branching rule which is often better than the CPLEX default or

CPLEX’s strong branching.

*Mark Wiley (..@lindo.com) *

[Tue, 18 Aug 1998] I saw your recent posting regarding solving randomly generated fixed

charge transportation models. As Jeff Linderoth mentioned in a subsequent posting, flow

cover cuts should help significantly on this type of problem. Our LINGO and What’sBest

packages include this enhancement. If you have some examples in MPS format we could try

them out.

[Fri, 21 Aug 1998] Thanks for sending the model. Unfortunately, flow cover cuts (or at

least our implementation) don’t do as much for your models as I would have liked. Almost

immediately LINGO finds a feasible integer solution of 4151 and a bound of 3359.81.

However, after that it makes virtually no progress. I’d be curious to know how MINTO

performs. I wish I could have given you better news. We will keep your model in our IP

test set. Perhaps a future version will be able to solve it quickly. It is an

interesting problem, small but tough.

Good luck.

[Fri, 23 Oct 1998] We tinkered with some cut options and let it run a bit longer and

found a better solution but failed to solve to optimality. Our best solution was 3764

with a bound of 3522.

*Ed Klotz (..@cplex.com) *

[Mon, 2 Feb 1998] I tried a few runs with your problem, without success. I tried this on

a very fast machine, and it did not help. I think the fundamental difficulty with this

problem is that is very easy for the optimizer to move fractional solutions in the LP

subproblems from one variable to the next with little or no degradation in the objective

value. Adding cuts, or removing the symmetry of the problem, can help.

[Tue, 3 Feb 1998] Thinking about it [the FCTP model], it’s not very surprising that this

is difficult for the branch and bound to solve. When you relax the integrality

requirement on the y_ij, the problem essentially reduces to a regular transportion

problem (except that the cost coefficients are no longer cij). So the convex hull of the

LP relaxation is quite different from the convex hull of the integer feasible solutions.

This characteristic tends to lead to difficult MIPs. And, I am not the only one who

thinks this is the case. Here’s what Nemhauser and Wolsey say about this problem on page

498 of "Integer and Combinatorial Optimization." "Unfortunately, the bounds obtained

from these relaxations are frequently very poor primarily because they do not accurately

represent the fixed costs..." I do not think we will be able to get CPLEX 5.0 to solve

this problem quickly just by changing parameters. It may be possible to come up with a

priority order file, but I think the best strategy will be to tighten the formulation

by adding cuts.

[Mon, 9 Nov 1998] While CPLEX 6.0 did not perform significantly better on this problem,

we have a development code that established a bound of 3578.6418 in two hours on a fast

machine. The associated integer solution of 3774 is therefore within 5.46 percent of

optimal. We obtained these results with default settings except for setting the

heuristicfreq parameter to 5 (i.e. we applied a rounding heuristic on every 5th node).

*We also did a run with variableselect 3 (strong branching), heuristicfreq 5 that we let *

go for several days. It eventually found a solution of 3714 and a bound of 3655. So,

CPLEX found a solution within 2 percent of optimal, although it took a long time to do

so. I suspect that with additional experimentation we can get the development version

of CPLEX to find a 2 percent solution in a more reasonable amount of time.

*Irv Lustig (..@cplex.com) *

[Wed, 25 Nov 1998] We are in the process of developing a new version of CPLEX, so we

decided to try the ran14x18 problem with our new version on a large number of parallel

processors. We used an SGI Origin system with 64 processors. Our parallel development

code solved it with the following output:

Integer optimal solution (0.0001/0): Objective = 3.7120000000e+03

Current MIP best bound = 3.7116289184e+03 (gap = 0.371082, 0.010%)

Solution time = 11517.98 sec.

Iterations = 500288158 Nodes = 8982253 (101132)

This corresponds to 3 hours, 12 minutes wall clock time.

By default, we solve the problem to a .01% gap, but we could decrease this and the time

would only go up a little bit.

Just ran the model on a (cheap) quad core PC with Cplex 11 using default settings (we used*threads 4*):

So we need 630 seconds on this machine to prove optimality (optcr=0).

http://www.amsterdamoptimization.com/models/etc/ran14x18.gms. Some notes on this model from a posting on sci.op-research:

We tried CPLEX 5.0 on this instance, obtaining a solution with objective value 3712 but

without proving optimality. The same solution was also found in some runs of the genetic

algorithm from reference 1.

In response to my posting on sci.op-research, some people tried exact solvers on this

instance (fctp, mps), but without success. I also contacted CPLEX, but they also failed

to prove optimality.

On Friday, 13 November 1998, I announced this web page on sci.op-research. Some

hours later I received good news: Jeff proved optimality for ran14x18. Some time later,

a parallel CPLEX version also solved this problem to optimality. I thank the people cited

below for their efforts and comments concerning this problem!

To get an impression of the difficulty of this instance, have a look at the discussion

about the problem, which has remained unsolved for more than 10 months.

Jeff Linderoth (

inequalities would be able to solve the instances. (At least it should do better than

CPLEX for problems of this class). If you don’t have too many instances, then you can

send me the MPS files, and I would be happy to try the solver MINTO out on them for you.

[Fri, 28 Aug 1998] For the problem ran14x18.mps, I ran MINTO until it ran out of memory,

and was only able to improve the lower bound to 3617.994472. When I get the opportunity,

I will perform some longer runs and also try out PARINO (which is like a parallel

version of MINTO) on the problem. I also tried XPRESS and CPLEX and MINTO was much

better than either of these packages (due to flow cover inequalities). How did bc-opt

do? I’ll email more results later. (If you don’t hear from me, send me an email to

remind me). Anyway, sorry I wasn’t able to prove optimality for you, but hopefully the

results are somewhat useful.

[Fri, 13 Nov 1998] I proved the optimality of the solution 3712 for the ran14x18 FCTP

instance. In order to do this I ran my PARINO parallel MIP code for about 52 hours on

32 Pentium II (300MHz) processors. 7624677 nodes of the branch and bound tree were

evaluated in this time. PARINO adds "flow cover inequalities" which are usually quite

useful in helping to improve the LP bounds for these types of problems. Also, I used a

parallel pseudocost-type branching rule which is often better than the CPLEX default or

CPLEX’s strong branching.

[Tue, 18 Aug 1998] I saw your recent posting regarding solving randomly generated fixed

charge transportation models. As Jeff Linderoth mentioned in a subsequent posting, flow

cover cuts should help significantly on this type of problem. Our LINGO and What’sBest

packages include this enhancement. If you have some examples in MPS format we could try

them out.

[Fri, 21 Aug 1998] Thanks for sending the model. Unfortunately, flow cover cuts (or at

least our implementation) don’t do as much for your models as I would have liked. Almost

immediately LINGO finds a feasible integer solution of 4151 and a bound of 3359.81.

However, after that it makes virtually no progress. I’d be curious to know how MINTO

performs. I wish I could have given you better news. We will keep your model in our IP

test set. Perhaps a future version will be able to solve it quickly. It is an

interesting problem, small but tough.

Good luck.

[Fri, 23 Oct 1998] We tinkered with some cut options and let it run a bit longer and

found a better solution but failed to solve to optimality. Our best solution was 3764

with a bound of 3522.

[Mon, 2 Feb 1998] I tried a few runs with your problem, without success. I tried this on

a very fast machine, and it did not help. I think the fundamental difficulty with this

problem is that is very easy for the optimizer to move fractional solutions in the LP

subproblems from one variable to the next with little or no degradation in the objective

value. Adding cuts, or removing the symmetry of the problem, can help.

[Tue, 3 Feb 1998] Thinking about it [the FCTP model], it’s not very surprising that this

is difficult for the branch and bound to solve. When you relax the integrality

requirement on the y_ij, the problem essentially reduces to a regular transportion

problem (except that the cost coefficients are no longer cij). So the convex hull of the

LP relaxation is quite different from the convex hull of the integer feasible solutions.

This characteristic tends to lead to difficult MIPs. And, I am not the only one who

thinks this is the case. Here’s what Nemhauser and Wolsey say about this problem on page

498 of "Integer and Combinatorial Optimization." "Unfortunately, the bounds obtained

from these relaxations are frequently very poor primarily because they do not accurately

represent the fixed costs..." I do not think we will be able to get CPLEX 5.0 to solve

this problem quickly just by changing parameters. It may be possible to come up with a

priority order file, but I think the best strategy will be to tighten the formulation

by adding cuts.

[Mon, 9 Nov 1998] While CPLEX 6.0 did not perform significantly better on this problem,

we have a development code that established a bound of 3578.6418 in two hours on a fast

machine. The associated integer solution of 3774 is therefore within 5.46 percent of

optimal. We obtained these results with default settings except for setting the

heuristicfreq parameter to 5 (i.e. we applied a rounding heuristic on every 5th node).

go for several days. It eventually found a solution of 3714 and a bound of 3655. So,

CPLEX found a solution within 2 percent of optimal, although it took a long time to do

so. I suspect that with additional experimentation we can get the development version

of CPLEX to find a 2 percent solution in a more reasonable amount of time.

[Wed, 25 Nov 1998] We are in the process of developing a new version of CPLEX, so we

decided to try the ran14x18 problem with our new version on a large number of parallel

processors. We used an SGI Origin system with 64 processors. Our parallel development

code solved it with the following output:

Integer optimal solution (0.0001/0): Objective = 3.7120000000e+03

Current MIP best bound = 3.7116289184e+03 (gap = 0.371082, 0.010%)

Solution time = 11517.98 sec.

Iterations = 500288158 Nodes = 8982253 (101132)

This corresponds to 3 hours, 12 minutes wall clock time.

By default, we solve the problem to a .01% gap, but we could decrease this and the time

would only go up a little bit.

Just ran the model on a (cheap) quad core PC with Cplex 11 using default settings (we used

` S O L V E S U M M A R Y`

MODEL m OBJECTIVE cost

TYPE MIP DIRECTION MINIMIZE

SOLVER CPLEX FROM LINE 168

**** SOLVER STATUS 1 NORMAL COMPLETION

**** MODEL STATUS 1 OPTIMAL

**** OBJECTIVE VALUE 3712.0000

RESOURCE USAGE, LIMIT 630.895 3600.000

ITERATION COUNT, LIMIT 48023779 10000000

So we need 630 seconds on this machine to prove optimality (optcr=0).

Usually I do Excel programming in VBA, or Delphi and VBScript via Excel Com Automation. The Delphi tools from http://www.axolot.com/components/xlsrwii20.htm are reputed to be very fast, bypassing Excel altogether. Unfortunately it cannot read the spreadsheets I am interested in. We'll try to get the spreadsheet to them pending clients permission.

The following options are useful to make pictures of MIP performance:

- miptrace
*filename*: write trace file to the indicated file - miptracenode
*n*: write a trace record each*n*nodes (default:*n*=100) - miptracetime
*r*: write a trace record each*r*seconds (default:*r*=1)

Example output file:

* miptrace file trace.csv: ID = Cplex 11.0.1

* fields are lineNum, seriesID, node, seconds, bestFound, bestBound

1, S, 0, 0.075, na, 4053.06

2, N, 0, 0.084, na, 4053.06

3, N, 0, 0.09, 1740.37, 3852.17

4, N, 0, 0.095, 1740.37, 3804.84

5, N, 0, 0.101, 1740.37, 3798.37

6, N, 0, 0.107, 1740.37, 3781.77

7, N, 0, 0.113, 1740.37, 3781.53

8, N, 0, 0.149, 1826.58, 3779.84

9, N, 100, 0.229, 2769.62, 3327.77

10, N, 200, 0.304, 2891.26, 3135.02

11, N, 300, 0.391, 2891.26, 2910.57

12, E, 320, 0.397, 2891.26, 2891.26

* miptrace file trace.csv closed

It is easy to make a plot of this in Excel:

These options are actually available in all or most MIP solvers under GAMS.

Can this be done? This question comes up regularly.

There are a few interesting papers by Richard Howitt of UC Davis in this respect:

There are a few interesting papers by Richard Howitt of UC Davis in this respect:

Here is my list of features I would like to see added or improved to GAMS/CPLEX:

- Support for user cuts (redundant constraints that can help the solver)
- Support for lazy constraints (constraints that are only included when violated)
- Better language support for piecewise linear expressions
- Proper language support for indicator constraints
- Sensitivity analysis results to GDX file
- Sensitivity analysis for bounds
- Performance improvement for QP/QCP related models
- Provide unbounded ray
- Fixing a number of reported bugs (some of them reported on in this blog)

See also http://www.ampl.com/BOOKLETS/amplcplex110userguide.pdf (the AMPL link has a number of these features implemented) and http://www.gams.com/dd/docs/solvers/cplex.pdf.

Cplex needed a long time to solve the root node of an MIQCP model. I tried Mosek, but that gave

This is with version:

This is with version:

I sent it to GAMS for further investigation. This problem seems only to happen with the 64 bit version. The work around would be to use the 32 bit version.

MOSEK Link May 1, 2008 22.7.2 WEX 3927.4799 WEI x86_64/MS Windows

M O S E K version 5.0.0.79 (Build date: Jan 29 2008 13:47:58)

The tool SQL2GMS can be used to import large CSV files into GAMS. This is especially useful if the CSV file does not fit in Excel. Excel before Office 2007 could only handle up to 65,536 rows (in Office 2007 this is increased to 1,048,576). Examples:

- example 8.4 and 19.4 in http://www.amsterdamoptimization.com/pdf/sql2gms.pdf
- a practical example: http://www.amsterdamoptimization.com/models/blog/rawcsvdata3-2.gms
- technical info on the ODBC Text Driver: http://support.microsoft.com/kb/178717, http://msdn.microsoft.com/en-us/library/ms714091%28VS.85%29.aspx

To export GAMS models with semi-continuous variables to MPS files the tool GAMS/Convert can be used. However note that semi-continuous variables with lower bound = upper bound are exported incorrectly. Consider the small example:

*x=0*.

Convert generates the following MPS file:

*sc1* corresponds to GAMS variable *x*. However, this representation excludes *sc1=0* as *sc1* is fixed to one using FX. The correct MPS representation should look like:

*writemps* to create an MPS file. This will generate correct MPS files with fixed semi-continuous variables.

`semicont variable x;`

variable y;

equation e;

y.lo=0;

x.fx=1;

e.. y =e= x;

model m /e/;

solve m minimizing y using mip;

The solution is Convert generates the following MPS file:

`ROWS`

N obj

E c1

COLUMNS

sc1 c1 -1

x2 obj 1

x2 c1 1

RHS

BOUNDS

FX bnd sc1 1

ENDATA

In this MPS file column `ROWS`

N obj

E c1

COLUMNS

sc1 c1 -1

x2 obj 1

x2 c1 1

RHS

BOUNDS

LO bnd sc1 1

SC bnd sc1 1

ENDATA

Work around: use the Cplex option
In some cases GAMS/Cplex reports erroneously that zero iterations have been performed. Example: run THAI model from the model library with iterlim=5:

S O L V E S U M M A R Y

MODEL thainavy OBJECTIVE obj

TYPE MIP DIRECTION MINIMIZE

SOLVER CPLEX FROM LINE 84

**** SOLVER STATUS 4 TERMINATED BY SOLVER

**** MODEL STATUS 14 NO SOLUTION RETURNED

**** OBJECTIVE VALUE 0.0000

RESOURCE USAGE, LIMIT 0.037 1000.000

ITERATION COUNT, LIMIT 0 5

The SOLVER STATUS is also wrong, this should be 2 ITERATION INTERRUPT. Turns out I have reported this already in january 2006, but it does not seem to get fixed.

Correct values for these settings are often important when running Cplex from a more complex algorithm formulated in GAMS.

> Suppose I want to maximise, say:

> a1 + a2 + a3 + b1 + b2 + b3 + c1 + c2 + c3 < 30

>

> ...but with the extra condition that I must have two items from each

> of the following sets (at least two members of each set must be

> greater than zero): {a1, a2, a3}, {b1, b2, b3}, {c1, c2, c3}

>

> Is there an easy way to ask lp_solve for this condition, please?

>

One way of achieving this is to add the constraints:

a(i) ≥ 0.001·p(i)

b(i) ≥ 0.001·q(i)

c(i) ≥ 0.001·r(i)

sum(i,p(i)) ≥ 2

sum(i,q(i)) ≥ 2

sum(i,r(i)) ≥ 2

p(i),q(i),r(i) binary variables

Note that we do not need to model the additional implications

p(i)=0 => a(i)=0

because we exploit the "at least" part of the condition. If we want exactly two nonzeroes for each group, this implication would need to be added.

>

> Can anybody help me to formulate this ? just not able to figure this out

>

> x is a continuous variable. Can take any positive value including 0.

> w is a 0/1 binary variable

>

> if x > 0, I want w = 0

> if x = 0, I want w = 1

This should work:

x ≤ (1-w)·M where M=x.up

x ≥ 0.001·(1-w)

The update button (File|Options|Execute|Update) does not always work very well with Vista (access denied or some other failure). This has to do with missing administrator rights. Workaround: start **cmd** from start menu using Shift-Ctrl-Enter, cd to gams system directory and run GAMSINST.

It is quite amazing how fast these three alternatives for Cholesky are. Here is a model to time the performance and accuracy. Some results are listed below.

Method | Time | Accuracy | Time | Accuracy | Time | Accuracy |
---|---|---|---|---|---|---|

n=50 | n=100 | n=200 | ||||

Nonlinear Equations | 0.148 | 9.6e-8 | 0.917 | 1.2e-4 | 8.991 | 0.621 |

External Solver | 0.024 | 1.6e-13 | 0.025 | 1.9e-10 | 0.056 | 1.4e-5 |

GAMS Algorithm | 0.034 | 1.9e-13 | 0.482 | 6.3e-11 | 7.571 | 7.7e-6 |

I am surprised that the GAMS algorithm is competitive. Notice the accuracy issues when using nonlinear equations due to ill-conditioning of the system of non-linear equations (small changes in a(i,j) lead to large changes in L(i,j)).

Here is some info about the external solver:

`c:\>cholesky`

CHOLESKY: matrix decomposition A=LL^T

Usage

> cholesky gdxin i a gdxout L

where

gdxin : name of gdxfile with matrix

i : name of set used in matrix

a : name of 2 dimensional parameter inside gdxin

gdxout : name of gdxfile for results (factor L)

L : name of 2 dimensional parameter inside gdxout

Calculates the Cholesky decomposition A=LL^T of a symmetric

positive definite matrix A=a(i,j) where i and j are aliased sets.

L will contain the Cholesky factor L(i,j)

C:\>

Three ways to do a Cholesky decomposition in GAMS:

- Write A=LL
^{T}as nonlinear constraints. This has the advantage it can be part of an NLP model and maintains the factorization during the solve. Of course, this makes the NLP (much) more difficult. Example gams model. - Program a Cholesky algorithm directly in GAMS. Example gams model. This only works on parameters
- Call an external Cholesky solver. This is probably better for large matrices. Example gams model.

This seems to work reasonably well: http://www.hj-gym.dk/~hj/writer2latex/index.html (called from OpenOffice by File|Export|LaTeX 2e).

I just provided this example:

which may be surprising to non-expert GAMS users.

Here is one that may even stump hardcore GAMS modelers:

Question: can you predict the result of the display statement?

if(1,

$set name "hello"

else

$set name "world"

);

display "%name%";

which may be surprising to non-expert GAMS users.

Here is one that may even stump hardcore GAMS modelers:

scalar s /NA/;

s$(s>1) = 3.14;

display s;

Question: can you predict the result of the display statement?

>When I run REAP I get a execution error division by EPS message and I

>can't figure out why. There must be something I'm missing. There are $

>control statements to prevent this from happening and its never been a

>problem before. I've tried placing the $ control conditions on either

>side of the assignment symbol and still get the error.

Yes, EPS is numerically the same as zero. As soon as you start calculating

with it, it is considered the same as zero. However when used as a $

condition it is considered to be nonzero. So

1 scalar a /EPS/;

2 scalar b;

3 b$a=1/a;

4 display a,b;

gives:

**** Exec Error at line 3: division by eps

Indeed this fix is correct:

scalar a /EPS/;

a$(a=EPS) = 0;

scalar b;

b$a=1/a;

display a,b;

this will make a equal to zero if it was EPS and the expression b$a works

again as expected.

>can't figure out why. There must be something I'm missing. There are $

>control statements to prevent this from happening and its never been a

>problem before. I've tried placing the $ control conditions on either

>side of the assignment symbol and still get the error.

Yes, EPS is numerically the same as zero. As soon as you start calculating

with it, it is considered the same as zero. However when used as a $

condition it is considered to be nonzero. So

1 scalar a /EPS/;

2 scalar b;

3 b$a=1/a;

4 display a,b;

gives:

**** Exec Error at line 3: division by eps

Indeed this fix is correct:

scalar a /EPS/;

a$(a=EPS) = 0;

scalar b;

b$a=1/a;

display a,b;

this will make a equal to zero if it was EPS and the expression b$a works

again as expected.

I always thought that this would work under Windows:

> touch f1.txt f2.txt f3.txtHowever it fails miserably:

> ren f*.txt gg*.txt

A ("somewhat" complicated) solution is presented here.C:\projects\test\test1>touch f1.txt f2.txt f3.txt

C:\projects\test\test1>ren f*.txt gg*.txt

A duplicate file name exists, or the file cannot be found.

A duplicate file name exists, or the file cannot be found.

C:\projects\test\test1>dir

Volume in drive C has no label.

Volume Serial Number is 7563-3993

Directory of C:\projects\test\test1

05/16/2008 08:42 PM <DIR> .

05/16/2008 08:42 PM <DIR> ..

05/16/2008 08:42 PM 0 f2.txt

05/16/2008 08:42 PM 0 f3.txt

05/16/2008 08:42 PM 0 gg.txt

3 File(s) 0 bytes

2 Dir(s) 159,707,815,936 bytes free

C:\projects\test\test1>

A question on the GAMS list was:

>I would like to write the results of a GAMS program into a text file

>whose name should depend on the date and time of the run. How could

>that be done?

>

>Ex: output-05.15.08-09.39.txt

>The above is just a sample and it does not have to be exactly the same.

The suggestion to use `myfile /filename_%system.date%_%system.time%/` does not work for several reasons: quotes are needed to handle embedded blanks and the ‘/’ characters in the system date cause problems when used inside a file name. Here is a different approach that actually works (I ran it before suggesting it):

The generated include file will look like:$onecho > createfile.awk

BEGIN { print "file f / \"f_" strftime("%Y-%m-%d_%H.%M.%S") ".txt\"/" }

$offecho

$call "awk -f createfile.awk > ftime.inc"

$include ftime.inc

put f "hello"/;

For more info on AWK see gawk.htmlfile f / "f_2008-05-16_14.28.59.txt"/

Have not seen this error before. This does not bode well for this model.

Trouble polishing MIP start.

For large MIQCP/QCP models there is an issue with GAMS and with GAMS/CPLEX. It is simply very slow in handling the quadratics. There are two places where there is likely a complexity problem (e.g. a quadratic term in run time vs linear in the linear nonzeroes):

- In GAMS itself. Besides generating the model itself as usual, GAMS calls a program called QMAKER that exports the quadratic structure in the form of a GDX file. This QMAKER program is not suited for large models. E.g. for a very large MIQCP, GAMS can generate the model in 6.7 seconds, and then QMAKER needs 906.6 seconds to handle the quadratic terms.
- In the GAMS/CPLEX link there is a similar problem. Setting up the quadratic constraints is very slow for large models. I don't have timings because they are not provided.

I was hit again by this long standing bug. In some cases GAMS/Cplex will terminate with

**** SOLVER STATUS 4 TERMINATED BY SOLVER

**** MODEL STATUS 14 NO SOLUTION RETURNED

**** OBJECTIVE VALUE 0.0000

but the log file nor the listing file will contain any hints about why Cplex stopped. It may be related to multiple threads being used and running out of memory, but this is just a guess. I suspect the link does not capture correctly the return code from the cplex library call. This can happen with many GAMS versions including 22.1 (march 2006), 22.2, 22.3, 22.4, 22.5, 22.6, and 22.7.

**** SOLVER STATUS 4 TERMINATED BY SOLVER

**** MODEL STATUS 14 NO SOLUTION RETURNED

**** OBJECTIVE VALUE 0.0000

but the log file nor the listing file will contain any hints about why Cplex stopped. It may be related to multiple threads being used and running out of memory, but this is just a guess. I suspect the link does not capture correctly the return code from the cplex library call. This can happen with many GAMS versions including 22.1 (march 2006), 22.2, 22.3, 22.4, 22.5, 22.6, and 22.7.

On sci.op-research the following message was posted:

*> There is a worked example of how to solve a basic Cutting Stock Problem at> http://docs.google.com/View?docid=dfkkh8qj_64rrpz86db> It uses simple methods, and the open source program lp_solve. Briefly, > the solution shown is:>> * create a list of cut combinations (not permutations - there would be > too many of those to work with in most cases)>> * from this list, generate a set of linear programming (simplex) expressions, > and output them to a file>> * this file becomes the input to lp_solve, which calculates an optimal solution>> * most of the cut combinations in the solution will have a score of zero. > For those that have a score greater than zero, look them up on the list > created in step 1*

A clever way to solve these problems is using a Delayed Column Generation approach. Such a model is not too difficult to code in GAMS. See http://www.amsterdamoptimization.com/pdf/colgen.pdf. A nice write up on the cutting stock problem is available from http://en.wikipedia.org/wiki/Cutting_stock_problem.

The expression \(z = x \cdot \delta \) where \(x\) is a **continuous** variable and \(\delta\) is a **binary** variable can also be linearized fairly easily. The exact form depends on the bounds of \(x\). Assume \(x\) has lower bound \(x^{lo}\) and upper bound \(x^{up}\), then we can form the inequalities:

For the special case \(x^{lo}=0\) this reduces to:

The last restriction is essentially a big-M constraint:

\[\boxed{\begin{split}\min\{0,x^{lo}\} &\le z \le \max\{0,x^{up}\} \\

x^{lo} \cdot \delta &\le z \le x^{up} \cdot \delta \\

x − x^{up} \cdot (1 − \delta) &\le z \le x − x^{lo} \cdot (1 − \delta)\end{split}}\]

where \(M_1\), \(M_2\) are chosen as tightly as possible while not excluding \(z=0\) when \(\delta=0\).

\[\boxed{x − M_1 \cdot (1 − \delta) \le z \le x + M_2 \cdot (1 − \delta)}\]

For the special case \(x^{lo}=0\) this reduces to:

The construct \(z=x· \delta\) can be used to model an OR condition: "\(z = 0\) OR \(z = x\)".

\[\boxed{\begin{split}& z \in [0,x^{up}]\\ & z ≤ x^{up} · \delta\\ & z ≤ x\\ & z ≥ x − x^{up} · (1 − \delta)\end{split}}\]

Added a few regression models from the SciStat.org nonlinear regression datasets. They are solved by the linear and non-linear regression solvers GAMS/LS and GAMS/NLS.

- Dugongs implemented as dugongs.gms This is actual a linear regression problem, so we solve with LS. The log/log fit is better, see also dugongs.png and dugongsplot.gms. What is a Dugong? See: http://en.wikipedia.org/wiki/Dugong.
- Dolphin in gams formulated as dolphin.gms. If we add a constant term the problem converges to a singular solution. The problem may be related to an outlier visible in the plot dolphin.png (created by dolphinplot.gms).
- brunhild is implemented as brunhilda.gms. The first non-linear fit. The seconds one is better. The third regression is a linearized version of model 2.
- troutpcb is implemented as troutpcb.gms. The linear regression problem provides initial values for the non-linear regression problem.
- chloride is implemented as chloride.gms. See also chloride.png, chlorideplot.gms
- rabbit is implemented as rabbit.gms. See also rabbit.png, rabbitplot.gms
- chlorine: chlorine.gms

I used GNUPLOT to plot results, e.g.:

I see often MINLP models where the only nonlinear expression is the multiplication of two binary variables. Such models can be linearized quite easily and are in most cases better solved as a linear MIP model. A good formulation for \(z=x_1 \cdot x_2\) where \(x_1\), \(x_2\) are **binary variables** is as follows:

Another weaker formulation that is sometimes used is:

Notice that we relaxed \(z\) to be continuous between 0 and 1: \(z\) is automatically integer; a property we can use to reduce the number of discrete variable. Note that modern high performance solvers may reintroduce these variables as binary as they can exploit this in their algorithms.

\[\boxed{\begin{split}0 &≤ z ≤ 1\\ z &≤ x_1\\ z &≤ x_2\\ z &≥ x_1 + x_2 − 1\end{split}}\]

Another weaker formulation that is sometimes used is:

This formulation has fewer equations but has two major disadvantages:

\[\boxed{\begin{split}z &\in \{0,1\}\\ z & ≤ \frac{x_1+x_2}{2}\\z &≥ x_1 + x_2 − 1\end{split}}\]

- \(z\) needs to be integer
- it is not as tight as the first formulation. E.g. \((x_1,x_2,z)=(0,1,\frac{1}{2})\) is infeasible in the first formulation, while it is allowed in relaxations of the second formulation.

An application of this reformulation is a QP (Quadratic Program) with binary variables:

\[\boxed{\begin{split} \min\> & x'Qx\\ & Ax=b\\ & x \in \{0,1\} \end{split}}\]

This can be modeled as:

`variable y(i,j);`

set ij(i,j);

ij(i,j)$( ord(i)<ord(j) ) = yes;

y.lo(ij)=0;

y.up(ij)=1;

obj.. z =e= sum(i, q(i,i)*x(i)) + sum(ij(i,j), (q(i,j)+q(j,i))*y(i,j) );

ymul1(ij(i,j)).. y(i,j) =L= x(i);

ymul2(ij(i,j)).. y(i,j) =L= x(j);

ymul3(ij(i,j)).. y(i,j) =G= x(i)+x(j)-1;

The more general product \(z = \prod_i x_i\) with \(x_i\in \{0,1\}\) can be linearized in the same way:

where \(n\) is the cardinality of the set \(I\). Again if we want we can relax \(z\) to \(z \in [0,1]\).

\[\boxed{\begin{split} &z ≤ x_i\\ & z ≥ \sum_i x_i – (n-1)\\& z \in \{0,1\}\end{split}}\]

The GAMS construct $ifthen has more issues. This bug in $ifthen only happens in release 22.6. The following code does not work in GAMS 22.6:

This has been fixed in GAMS 22.7. If you need to work with 22.6, a work around is to add $else everywhere:

In my case I have added a $ statement that makes sure the user can only run this code with GAMS 22.7:

This has been fixed in GAMS 22.7. If you need to work with 22.6, a work around is to add $else everywhere:

In my case I have added a $ statement that makes sure the user can only run this code with GAMS 22.7:

My first blog starts with a down note. I was working last night on this very large GAMS model with many model settings implemented as $set/$setglobal macros. In this context I want to nest some $ifthen constructs:

$set stochastic 0

$set erslivestockmodel 1

$ifthen %stochastic%==1

$ ifthen %erslivestockmodel%==1

b(i,"Exports",k) = tempbx(k);

b(i,"Season Average Price",k) = tempbp(k);

b(i,"Feed Use",k) = tempbdf(k);

b(i,"Food",k) = tempbdi(k);

cfixed(i,"Exports",k) = 0;

cfixed(i,"Season Average Price",k) = 0;

cfixed(i,"Feed Use",k) = 0;

cfixed(i,"Food",k) = 0;

$ else

$ include grains.inc

$ endif

$endif

This does not work as nesting is not allowed. You don't get an error message either indicating nesting is not allowed. In some cases you may get a possibly misleading error message, in the worst case it just behaves incorrectly. The documentation in the McCarl User Guide does not mention this limitation. Looks like $ifthen will be less useful in practice than I had hoped.

The workaround is to use a run-time if statement:

$set stochastic 0

$set erslivestockmodel 1

scalars

stochastic /%stochastic%/

erslivestockmodel /%erslivestockmodel%/

;

if(stochastic=1,

if(erslivestockmodel=1,

b(i,"Exports",k) = tempbx(k);

b(i,"Season Average Price",k) = tempbp(k);

b(i,"Feed Use",k) = tempbdf(k);

b(i,"Food",k) = tempbdi(k);

cfixed(i,"Exports",k) = 0;

cfixed(i,"Season Average Price",k) = 0;

cfixed(i,"Feed Use",k) = 0;

cfixed(i,"Food",k) = 0;

else

$ include grains.inc

);

);

This has additional advantages: all branches are syntax-checked even if they are not executed, and we check the type of the parameters using the scalar initializations. Of course we lose some of the compile time features using this approach.

Another workaround would be to place the inner $ifthen in an include file:

$set stochastic 0Conclusion: GAMS users be aware when using nested $ifthen constructs. They are not implemented (correctly) and no error message is produced to warn you about this.

$set erslivestockmodel 1

$ifthen %stochastic%==1

$include gamsbugworkaround.inc

$endif

Subscribe to:
Posts (Atom)