Description

While investigating an issue with the markov model from GAMS' modlib,
I discovered the model generated differently on different machines.
The cause of this was traced back to an "inaccurate" parameter
calculation (more precisely, a calculation that gives significantly
different results on different machines).

In this model I exaggerate the difference in the results by changing
pn only slightly.  The large relative error or difference arises in
the usual way - by subtracting two quantities of nearly equal value.

Contributor: Steve Dirkse, June 2005

Small Model of Type : GAMS

Category : GAMS Test library

Main file : math1.gms

$title 'Demonstrate an inaccurate calculation' (MATH1,SEQ=228)

$onText
While investigating an issue with the markov model from GAMS' modlib,
I discovered the model generated differently on different machines.
The cause of this was traced back to an "inaccurate" parameter
calculation (more precisely, a calculation that gives significantly
different results on different machines).

In this model I exaggerate the difference in the results by changing
pn only slightly.  The large relative error or difference arises in
the usual way - by subtracting two quantities of nearly equal value.

Contributor: Steve Dirkse, June 2005
$offText

Scalars  b   discount factor             /   .95   /
         e                               /    .1   /
         q   supply                      / 110     /
         d                               /-130     /
         k                               / 342     /
         pn  normal price (us$ per bbl)  /  34.5261889893 /;

scalar p, c, cmin, cmax, mn, mx, relDiff;

p = ( k / (q - d))**10;
c = 1e5 * d * (p-pn);

* cmin=-0.000556809709451045 computed with:
*    22.0 Alfa on LNX, LEX, LNXonLEX, WEX, SOL
$gdxIn math1min
$load cmin
$gdxIn

* cmax=-0.000556624968339747 computed with:
*    22.0 Alfa on WIN, WINonWEX
$gdxIn math1max
$load cmax
$gdxIn

mn=min(cmin,c);
mx=max(cmax,c);
relDiff = (mx-mn)/mx;
option decimals = 8;
display cmin, c, cmax, relDiff;

if {(c < cmin),
   execute_unload 'math1min_new.gdx' c=cmin;
   abort "Computed c < previous minimum, check math1min_new.gdx";
};
if {(c > cmax),
   execute_unload 'math1max_new.gdx' c=cmax;
   abort "Computed c > previous maximum, check math1max_new.gdx";
};