Example MSFCli.exe

Example of using the Microsoft solver with an mps file

C:\projects\test>"c:\Program Files (x86)\Microsoft Solver Foundation\1.0.6.4960\Bin\MSFCli.exe" +type mps +verbose +log indus89.mps
===== Processing .\indus89.mps =====
Algorithm: Primal
CostingExact: SteepestEdge
CostingDouble: SteepestEdge
Variables: 6570 -> 6366 + 2363
Rows: 2727 -> 2564
NonZeros: 39490
Eliminated Slack Variables: 201
Phase 1 double Pivots: 4639
Phase 2 double Pivots: 1864
double Factorings: 138
Rational Factorings: 1
Degenerate pivots: 1367 (21.02 %)
{
x1 -> 114873.65553186,
x2 -> 1522.85610775999,
x3 -> 3460.00034188404,
x4 -> 479.152767198711,
x5 -> 22532.3511432455,
x6 -> 2490.55460455224,
....
x6530 -> 1,
x6550 -> 1,
x6570 -> 1
}
Zero variables: 3949
Positive variables: 2621
Negative variables: 0
Objective value: 'obj' Min Optimal -114873.65553186

===== Finished .\indus89.mps: 00:00:05.8776890 Optimal =====

The solver seems to feature some kind of rational arithmetic.

Not sure if this is related, but the Microsoft solver finds a correct solution to the Neumaier model http://www.gams.com/modlib/libhtml/badmip.htm.

Microsoft Solver Foundation

.NET based solvers:

• Simplex solver (LP,MIP)
• Interior Point solver (LP,QP)
• Constraint solver (CSP)
• Unconstrained non-linear solver (NLP)

Can be used as library via any .NET language, or input can be specified as MPS and a simple equation oriented format called OML. In addition, Excel bindings are available.

http://code.msdn.microsoft.com/solverfoundation

After downloading, documentation and samples are available in My Documents\Microsoft Solver Foundation.

GBD feasibility cut

>Dear Mr. Kalvelagen,
>
>I am trying to do a comparison between the performance
>of GBD and OA in the case of convex and nonconvex minlp
>problems.
>
>I am using the GBD-code from the script for MINLP methods
>that you have written and posted on your page, and comparing
>it with DICOPT-OA. My application is on transportation
>network design problem. The theory as highlighted in Floudas
>1995 about the difference between the two methods seems to
>match for the most part in my problem. However, I am facing
>infeasibilities in the sub-problem.
>
>Your GBD code on the portfolio optimization generates support
>functions only from the feasible primal. Can you give me an
>idea (if that is not much of a hassle for you) how a
>feasibility primal subproblem can be implemented in your
>GBD code to generate additional supports (from the infeasible
>subproblems), or if you have provided an example of such a
>case somewhere else I’d be grateful if you can direct me to that?
>
>I do find your webpage and the relevant files very useful when
>testing theoretical concepts through implementation in GAMS.
>
>Many thanks in advance for your time!

You can add feasibility cuts by forming a feasibility subproblem.

I do something like that for the dual subproblem in http://amsterdamoptimization.com/pdf/benders.pdf

A straightforward discussion can be found in Ullmann's Modeling and Simulation By Fritz Ullmann page 113, and 116.

Of course it may be possibly to bypass the problem by formulating an elastic version such that the subproblems are always feasible.

Conditional generation

> OK thanx for this information.
> Did you have time to look at my "corrected" constraints? I'm really sorry for
> having written an error ...
> The correct conditions are now written:
>------------------------------------------
> if s=1 then t > b
> if s=-1 then t < -b
> if s=0 then -b <= t <= b
> where:
> variables: t real in [-1,1] b real in [0,1]
> parameter: s integer in {-1,0,1}
> ------------------------------------------
>I just can't see how to write these conditions as linear constraints in AMPL ...
>Thanx again for your help or any advice which could guide me towards the solution ...

param s = 1;
var t >= -1, <= 1;
var b >= 0, <= 1;

minimize obj: b; subject to
e1{1..1:s=1}: t >= b;
e2{1..1:s=-1}: t <= -b;
e3{1..1:s=0}: -b <= t;
e4{1..1:s=0}: t <= b;


Conditional generation of equations is an important feature in formulating real-world models. You may want to consider hiring a consultant with experience in AMPL modeling. That would speed up the implementation considerably.

GAMS Command Prompt

Many Windows programs that have a command line interface install in the Start menu a Console Prompt that have path and environment automatically setup correctly. Examples are Lahey and the Microsoft SDK's. This would be useful for GAMS also for users that use the command line.

$load x=x.l not allowed (GAMS bug)

When a variable x is stored in a GDX file we can read it as a parameter as follows:

parameter x(k,i,t);
execute_load 'merged.gdx',x=x.l;

Unfortunately at compile time we don't have the same facility:

parameter x(k,i,t);
$gdxin merged.gdx
$load x=x.l

yields

5923  parameter x(k,i,t);
GDXIN   C:\projects\china\input\merged.gdx
5925  $load x=x.l
****            $490
**** LINE     19 INPUT       C:\projects\china\input\report.gms
**** 490  Cannot load/unload this suffix

On solving a large system of nonlinear equation

When trying to solve a very large system of nonlinear F(x)=0 I encountered a number of problems. Large meant here about 224k variables and equations. Solving as CNS gave "pivot too small". An NLP formulation:

min dummy

F(x)=0

caused lack of progress and stopped with an infeasible solution. In the end I could solve the model reliably using a min absolute deviations formulation where explicit slacks are omitted for linear equations, i.e.:

min ∑ s_i

Ax−b=0

G(x)+s = 0

A standard variable splitting technique is used to model the absolute value:

min ∑ v_i+w_i

Ax−b=0

G(x)+v−w=0

v,w≥0

sub models

For large models it is often convenient to build models in separate modules. The module can exchange data via save/restart files or gdx files. A driver module will call the individual sum modules.

However this does not work well with the IDE. For instance the log of a submodel if called by

$call =gams submodel.gms

will be invisible and lost. There are few thing we can do to improve the situation. The settings I use are as follows:

* flags to run submodels better under the IDE
$setglobal ide "ide=%gams.ide% lo=%gams.lo% errorlog=%gams.errorlog% errmsg=1"
 
$call =gams submodel.gms %ide%
$if errorlevel 1 $abort submodel.gms has issues

This will solve most problems of calling gams via $call when using the IDE, and also works good when not running under the IDE.

The one remaining and somewhat annoying problem (bug) is that clicking on black lines of the submodel log, will cause the IDE to go to listing file of the main model. I.e. the wrong listing file is loaded. Clicking on the red lines work fine: the IDE will go the correct source file.

I think it would be an improvement if GAMS has better built-in support for executing submodels (compile time and execution time) where we don't have to worry about all these parameters and with better interaction with the IDE.

Restarts

> I am dealing with a class of large scale LPs. The class has nice
> structure though, let's say it has 1000 constraints. The objective
> function and the first 999 constraints are identical among all
> instances of the problems in this class. So everytime that user wants
> to solve a new problem it is only different from other problems in the
> 1000th constraint.
>
> I was wondering if this feature can be exploited to improve the speed
> of solving problems in this class. Is it possible that something be
> pre-computed and stored about the problem from the fixed part
> (objective and 999 constraints) and then everytime that user inserts
> the 1000'th constraint, that knowledge be used efficiently?

Absolutely. 1000 constraints is very small, but in general simplex based LP solvers have excellent restart facilities using the concept of a basis. Using the GAMS modeling language this even completely automated. I.e. the fragment:

set scenario /s1*s5/;
 
parameter growthq_s(scenario)   /
 s1 2
 s2 2.2
 s3 2.3
 s4 2.5
 s5 2.7
/;
 
parameter report(scenario,*);
 
loop(scenario,
   growthq = growthq_s(scenario);
   solve wsisn maximizing cps using lp;
   report(scenario,'time') = wsisn.resusd;
   report(scenario,'iterations') = wsisn.iterusd;
);
 
display report;

solves a series of related LP's: one coefficient in the model
changes here. Iteration count and cpu time are recorded and
they are displayed as:

 ----   3642 PARAMETER report
 
          time  iterations
 
 s1       0.630    3704.000
 s2       0.108      76.000
 s3       0.070      76.000
 s4       0.066      76.000
 s5       0.069      76.000
 

so indeed the first solve is the most expensive one.

Bug in gams

The fragment below demonstrates a bug in GAMS. From early analysis it looks like this is fairly difficult to encounter in a real model, so it may be less alarming than it looks like at first sight.

   1   SET NUM /1*2/;
   2   VARIABLE X(NUM);
   3   X.L(NUM) = ORD(NUM);
   4
   5   X.FX('2') $(X.LO('2') EQ X.UP('2')) = 234;
   6   DISPLAY X.LO,X.L,X.UP;
 
COMPILATION TIME     =        0.000 SECONDS      3 Mb  WEX228-228 Jul 26, 2008
 
----      6 VARIABLE X.Lo
 2 234.000
----      6 VARIABLE X.L
 1   1.000,    2 234.000
----      6 VARIABLE X.Up
 2 234.000

GAMS/Xpress documentation bug

loadmipsol in xpress.pdf is documented as having a default of 1. That is hard to believe. Indeed the def file indicates this is wrong.

Smooth approximations

Some examples of using the logistic function in smooth approximations for ifthen and max/min:

• ifthen(x1>x2, y, z) = z + (y-z)*0.5*(1+tanh(k*(x1-x2)))
• ifthen(x1<x2, y, z) = y + (z-y)*0.5*(1+tanh(k*(x1-x2)))
• max(x,y) = ifthen(x>y,x,y) = y + (x-y)*0.5*(1+tanh(k*(x-y)))
• min(x,y) = ifthen(x<y,x,y) = x + (y-x)*0.5*(1+tanh(k*(x-y)))

Some graphs with k=100:

See also: http://yetanothermathprogrammingconsultant.blogspot.com/2008/08/smooth-approximation-for-ifthen.html.

Bug in ifthen (GAMS)

Consider the models:

   1  variable x;
  2  equation e1,e2;
  3
  4  scalar s /1/;
  5
  6  e1.. x$(s>0) =e= 0;
  7  e2.. ifthen(s>0,x,0) =e= 0;
  8
  9  model m1 /e1/;
 10  solve m1 using cns;
 11
 12  model m2 /e2/;
 13  solve m2 using cns;

The models are the same. The ifthen is exogenous but GAMS does not recognize this. As a result the fist model solves fine and the second model leads to numeric problems:

** Error in Square System: Pivot too small.

Turns out the gradients are not correct for ifthen:

---- e2  =E= 
e2..  (0)*x =E= 0 ; (LHS = 0)

Advice: don't use ifthen in equations as the derivatives of this function are wrong. If possible use $ conditions (exogenous if) or a smooth approximation (endogenous if). Smooth approximations for ifthen are discussed in this blog.

Arithmetic on indices

>Is it possible to write an objective function as
>"obj.. PROFIT =e= sum((l,m,dig) $(ord(l) le p and ord(m) le p and ord(dig) le digmax ), x('(1+p*(ord(l)-1))','(1+p*(ord(m)-1))',dig)*ord(dig)); "
>if I only want to add particular x's to my objective function (indexed as 1+p*(ord(l)-1) and 1+p*(ord(m)-1) where p is a constant )?
>In this case, gams gives me domain violation error for first and second index and i couldn't find a way to solve it
>(the domain is not violated mathematically for example for p=3 and digmax= 9).
>Is there a way to do arithmetic operations in indexes?

Indices are always strings in GAMS. But you can do what you want as follows:

alias (l,i), (m,j); 
obj.. PROFIT =e=
  sum((l,m,dig,i,j) $(ord(l) le p and ord(m) le p and ord(dig) le digmax and ord(i)=(1+p*(ord(l)-1)) and ord(j)=(1+p*(ord(m)-1))),
       x(i,j,dig)*ord(dig));

This is somewhat messy, but we can do this better. You can simplify this by introducing extra sets:

alias (l,i), (m,j);  
set l2(l);     l2(l)$(ord(l)<=p) = yes;
set m2(m);     m2(m)$(ord(m)<=p) = yes;
set dig2(dig); dig2(dig)$(ord(dig)<=digmax) = yes;
set mapl(l,i); mapl(l2(l),i)$(ord(i)=(1+p*(ord(l)-1))) = yes;
set mapm(m,j); mapm(m2(m),j)$(ord(j)=(1+p*(ord(m)-1))) = yes;
display mapl,mapm;
...
equation obj;
obj.. PROFIT =e= sum((mapl(l,i),mapm(m,j),dig2(dig)) , x(i,j,dig)*ord(dig));

This way you can display the intermediate sets and debug the model more easily. The equation becomes much simpler. 

For large, complex models this strategy is very helpful: move complexities from the equations to set operations. Sets can be debugged easily using display or using the gdx viewer. In addition they can be debugged before the complete model is ready to be generated.

Privacy.pdf not found

> You told me to look at http://www.gams.com/docs/privacy.pdf but I get
> 403 Forbidden. Please help!!!!

Sorry. GAMS moved their web pages and this caused problems here and there. I'll email a copy.

GAMS bug in redeclaration

GAMS allows redeclarations of symbols as long as they don't cause a contradiction. I.e. the following is valid GAMS:

set i;
set i /1,2/;
set i(*);

The GAMS compiler gets confused if i is used as domain:

   1  set i;
   2  parameter p(i);
   3  set i /1,2/;
   4  set i(*);
****         $184
**** 184  Domain list redefined

Work-around: make sure line 1 and 4 use both i or i(*).