gastrans.gms : Gas Transmission Problem - Belgium

Description

The problem of distributing gas through a network of pipelines is formulated as
a cost minimization subject to nonlinear flow-pressure relations, material
balances, and pressure bounds. The Belgian gas network is used as an example.

First, we model a straight-forward NLP formulation that can be solved fine
by todays NLP solvers.
Afterwards, the 3-stage approach by Wolf & Smeers is implemented (line 160ff).

de Wolf, D, and Smeers, Y, The Gas Transmission Problem Solved by
and Extension of the Simplex Algorithm. Management Science 46, 11
(2000), 1454-1465.

Keywords: nonlinear programming, discontinuous derivatives, engineering, distribution
problem, gas transmission, network problem

Small Model of Type : NLP

Category : GAMS Model library

Main file : gastrans.gms

\$title Gas Transmission Problem - Belgium (GASTRANS,SEQ=217)

\$onText
The problem of distributing gas through a network of pipelines is formulated as
a cost minimization subject to nonlinear flow-pressure relations, material
balances, and pressure bounds. The Belgian gas network is used as an example.

First, we model a straight-forward NLP formulation that can be solved fine
by todays NLP solvers.
Afterwards, the 3-stage approach by Wolf & Smeers is implemented (line 160ff).

de Wolf, D, and Smeers, Y, The Gas Transmission Problem Solved by
and Extension of the Simplex Algorithm. Management Science 46, 11
(2000), 1454-1465.

Keywords: nonlinear programming, discontinuous derivatives, engineering, distribution
problem, gas transmission, network problem
\$offText

\$eolCom //

Set
i 'towns'
/ Anderlues, Antwerpen, Arlon,   Berneau,   Blaregnies
Brugge,    Dudzele,   Gent,    Liege,     Loenhout
Mons,      Namur,     Petange, Peronnes,  Sinsin
Voeren,    Wanze,     Warnand, Zeebrugge, Zomergem   /
a          'arcs' / 1*24 /
as(a)      'active arcs'
ap(a)      'passive arcs'
aij(a,i,i) 'arc description';

Alias (i,j);

Set nc 'node data column headers'
/ slo 'supply lower bound (mill M3 per day)'
sup 'supply upper bound (mill M3 per day)'
plo 'pressure lower bound           (bar)'
pup 'pressure upper bound           (bar)'
c   'cost                    (\$ per MBTU)' /;

Table Ndata(i,nc) 'node data'
slo       sup   plo    pup   c
Anderlues     0     1.2       0   66.2   0
Antwerpen  -inf    -4.034    30   80.0   0
Arlon      -inf    -0.222     0   66.2   0
Berneau       0     0         0   66.2   0
Blaregnies -inf   -15.616    50   66.2   0
Brugge     -inf    -3.918    30   80.0   0
Dudzele       0     8.4       0   77.0   2.28
Gent       -inf    -5.256    30   80     0
Liege      -inf    -6.385    30   66.2   0
Loenhout      0     4.8       0   77.0   2.28
Mons       -inf    -6.848     0   66.2   0
Namur      -inf    -2.120     0   66.2   0
Petange    -inf    -1.919    25   66.2   0
Peronnes      0     0.96      0   66.2   1.68
Sinsin        0     0         0   63.0   0
Voeren     20.344  22.012    50   66.2   1.68
Wanze         0     0         0   66.2   0
Warnand       0     0         0   66.2   0
Zeebrugge   8.870  11.594     0   77.0   2.28
Zomergem      0     0         0   80.0   0   ;

Set ac 'arc data column headers' / D   'diameter (mm)'
L   'length   (km)'
act 'type indicator (1 = active)' /;

D     L  act
1.Zeebrugge. Dudzele      890.0   4.0
2.Zeebrugge. Dudzele      890.0   4.0
3.Dudzele  . Brugge       890.0   6.0
4.Dudzele  . Brugge       890.0   6.0
5.Brugge   . Zomergem     890.0  26.0
6.Loenhout . Antwerpen    590.1  43.0
7.Antwerpen. Gent         590.1  29.0
8.Gent     . Zomergem     590.1  19.0
9.Zomergem . Peronnes     890.0  55.0
10.Voeren   . Berneau      890.0   5.0    1
11.Voeren   . Berneau      395.0   5.0    1
12.Berneau  . Liege        890.0  20.0
13.Berneau  . Liege        395.0  20.0
14.Liege    . Warnand      890.0  25.0
15.Liege    . Warnand      395.0  25.0
16.Warnand  . Namur        890.0  42.0
17.Namur    . Anderlues    890.0  40.0
18.Anderlues. Peronnes     890.0   5.0
19.Peronnes . Mons         890.0  10.0
20.Mons     . Blaregnies   890.0  25.0
21.Warnand  . Wanze        395.5  10.5
22.Wanze    . Sinsin       315.5  26.0    1
23.Sinsin   . Arlon        315.5  98.0
24.Arlon    . Petange      315.5   6.0     ;

Scalar
T   'gas temperature                (K)' / 281.15  /
e   'absolute rugosity             (mm)' /   0.05  /
den 'density of gas relative to air (-)' /   0.616 /
z   'compressibility factor         (-)' /   0.8   /;

Parameter
lam(a,i,j)    'lambda constant (page 1464)'
c2(a,i,j)     'pipe constant   (page 1463)'
arep(a,i,j,*) 'arc report';

ap(a) = not as(a);

arep(aij,'lam') = lam(aij);
arep(aij,'c2')  = c2(aij);

option  arep:6:3:1;
display arep, as, ap;

Variable
f(a,i,j) 'arc flow        (1e6 SCM)'
s(i)     'supply - demand (1e6 SCM)'
pi(i)    'squared pressure'
sc       'supply cost';

Equation
flowbalance(i)   'flow conservation'
weymouthp(a,i,j) 'flow pressure relationship - passive'
weymoutha(a,i,j) 'flow pressure relationship - active'
defsc            'definition of supply cost';

flowbalance(i).. sum(aij(a,i,j), f(aij)) =e= sum(aij(a,j,i), f(aij)) + s(i);

weymouthp(aij(ap,i,j)).. signpower(f(aij),2) =e= c2(aij)*(pi(i) - pi(j));

weymoutha(aij(as,i,j))..       - sqr(f(aij)) =g= c2(aij)*(pi(i) - pi(j));

defsc.. sc =e= sum(i, ndata(i,'c')*s(i));

* supply, demand, pressure, and flow bounds
s.lo(i)  = ndata(i,'slo');
s.up(i)  = ndata(i,'sup');
pi.lo(i) = sqr(ndata(i,'plo'));
pi.up(i) = sqr(ndata(i,'pup'));
f.lo(aij(as,i,j)) = 0;

* initialize flow to avoid getting trapped at signpower(0,2)
f.l(aij) = uniform(-1,1);

Model gastrans / flowbalance, weymouthp, weymoutha, defsc /;

solve gastrans using nlp min sc;

display f.l;

* to run the Wolf & Smeers approach, remove the following \$exit
\$exit

Parameter
frep 'flow report'
sdp  'supply demand and pressure';

frep('NLP',aij,'Flow')  =  f.l(aij);
sdp('NLP',i,'Supply')   =  s.l(i)\$(s.l(i) > 0);
sdp('NLP',i,'Demand')   = -s.l(i)\$(s.l(i) < 0);
sdp('NLP',i,'Pressure') =  sqrt(pi.l(i));
sdp('NLP','','Obj')     =  sc.l;

Parameter
flow(a,i,j,*)    'flow bounds'
pirange(a,i,j,*) 'squared pressure bounds';

flow(aij(ap,i,j),'min') = -sqrt(c2(aij)*(pi.up(j) - pi.lo(i)));
flow(aij(ap,i,j),'max') =  sqrt(c2(aij)*(pi.up(i) - pi.lo(j)));

pirange(aij(ap,i,j),'min') = pi.lo(i) - pi.up(j);
pirange(aij(ap,i,j),'max') = pi.up(i) - pi.lo(j);

option  flow:3:3:1, pirange:3:3:1;
display flow, pirange;

Equation
defh              'definition of Smeers obj'
defsig(a,i,j)     'definition of flow direction'
weymouthp2(a,i,j) 'flow pressure relationship - passive'
flo
fup
pilo
piup;

Variable
sig(a,i,j) 'flow direction (-1 or +1 for passive elements)'
b(a,i,j)   'flow direction ( 1 = i to j)'
h          'Smeers obj';

Binary Variable b;

weymouthp2(aij(ap,i,j)).. sig(aij)*sqr(f(aij)) =e= c2(aij)*(pi(i) - pi(j));

defh.. h =e= sum(aij(a,i,j), abs(f(aij))*sqr(f(aij))/3/c2(aij));

defsig(aij(ap,i,j)).. sig(aij) =e= 2*b(aij) - 1;

flo(aij(ap,i,j))..  f(aij) =g= flow(aij,'min')*(1 - b(aij));

fup(aij(ap,i,j))..  f(aij) =l= flow(aij,'max')* b(aij);

pilo(aij(ap,i,j)).. pi(i) - pi(j) =g= pirange(aij,'min')*(1 - b(aij));

piup(aij(ap,i,j)).. pi(i) - pi(j) =l= pirange(aij,'max')*b(aij);

* drop the previous solution
pi.l(i)  = 0;
s.l(i)   = 0;
f.l(aij) = uniform(-1,1);

Model
one   / defh,  flowbalance /
two   / defsc, flowbalance, weymouthp2, weymoutha /
three / defsc, flowbalance, weymouthp2, weymoutha, defsig, pilo, piup, flo, fup /;

option limRow = 0, limCol = 0;

solve one using dnlp min h; // provides good initial point

solve two using  nlp min sc;

* assignmenst below fix known solution
* b.fx(*aij)         = 1;
* b.fx(aij('7',i,j)) = 0;
* b.fx(aij('8',i,j)) = 0;

solve three using minlp min sc;

frep('Smeers',aij,'Flow')  =  f.l(aij);
sdp('Smeers',i,'Supply')   =  s.l(i)\$(s.l(i) > 0);
sdp('Smeers',i,'Demand')   = -s.l(i)\$(s.l(i) < 0);
sdp('Smeers',i,'Pressure') =  sqrt(pi.l(i));
sdp('Smeers','','Obj')     =  sc.l;

option  frep:6:4:1;
display frep, sdp;
GAMS Development Corp.
GAMS Software GmbH

General Information and Sales
U.S. (+1) 202 342-0180
Europe: (+49) 221 949-9170