! Waste water pooling problem.
 Decide the routing of wastewater sources through
 treatment facilities, so as to minimize the
 costs of the linkages used.
 Three levels:
    source -> treatment -> sink.
 There are constraints on the level of
several pollutant concentrations;
!Given parameters:
   fs(s) = flow out of source s,
   qcs(c,s) = concentration of pollutant c in flow from source s,
   rct(c,t) = fraction of pollutant c removed from 
              the stream passing through treatment plant t,
 Decision variables:
   fst(s,t) = flow from source s directly to treatment t,
   fse(s,e) = flow from source s directly to sink e,
   fte(t,e) = flow from treatment t directly to sink e,
   ftt(t,t2) = flow from treatment t directly to treatment t2,
   ft(t) = total flow out of treatment t,
   qct(c,t) = concentration of pollutant c in effluent of treatment t,
   yt(t) = 1 if we use treatment t, else 0,
   yst(s,t) = 1 if we ship from source s to treatment t, else 0;
! Keywords: Water, Waste Treatment, Network Design, Pooling;

!Best performance is obtained by settings:
   LINGO  | Options...  |  Global  |  off
   LINGO  | Options...  |  Global  |  Multistart Solver Attempts | 3
 
Reference: Misener, R. and C. Floudas(2010), 
"Global Optimization of Large-Scale Generalized Pooling Problems:
 Quadratically Constrained MINLP Models", Industrial & Engineering Chemistry Research,
 ;
! Define the various primitive entitities;
SETS:
 source: fs ;
 treatment: costt, costyt, yt, ft, ftmax ;
 sink: femax ;
 concentration;
! and combinations of entitities, e.g., links, etc.;
 st(source, treatment): cst, cyst, fst, dst, yst ;
 se(source, sink): cse, cyse, fse, dse, yse, fsemax, fsemax1;
 te(treatment, sink): cte, cyte, fte, dte, yte;
 tt(treatment, treatment): ctt, cytt, ftt, dtt, ytt;
 ct(concentration, treatment): rct, qct, qctsub;
 cs(concentration, source): qcs ;
 ce(concentration, sink): qce, qcemax ;
ENDSETS
DATA: source = 1..7; treatment = 1..4; sink = 1..1; concentration = 1..3; fmin = 0.2 ; cfix = 124.6 ; cvar = 3603.4 ; v = 3600 ; fs = 20 50 47.5 28 100 30 25 ; ! Amount at each source; qcs= 100 800 400 1200 500 50 1000 ! Three qualities; 500 1750 80 1000 700 100 50 ! at the sources; 500 2000 100 400 250 50 150 ; qceMAX = 5 5 10 ; costt = 1102.9 2895.2 1102.9 1102.9 ; ! Variable cost, each treatment; costyt = 13972 36676 13972 13972 ; ! Fixed cost, each treatment; dse = 150 135 100 90 40 70 45 ; dte = 75 80 160 190 ; dtt = 0 80 210 190 80 0 130 100 210 130 0 110 190 100 110 0 ; dst = 75 150 280 245 55 125 260 215 30 115 240 220 55 140 245 245 55 40 150 150 40 120 230 230 30 60 175 165 ; rct = 0.99 0 .1 .7 .9 .87 .99 .2 .95 .90 0 .3 ; ! Simple upper bounds on qct(i,j):; qctsub = 6 1000 500 300 10 230 10 700 20 200 300 400 ; ! 12 1200 1080 360 175 227.5 17.5 1400 100 200 2000 1400; ! Optimal qct; ! qct = 4.598752 800.0000 417.9658 222.5000 5.000000 227.5000 6.154791 633.3333 9.956739 200.0000 259.5370 309.1667; ! Optimal solution has Obj=1086187.0; ! Variable counts for this data set: fse(s,e): 7, fte(t,e): 4, ftt(t,t2): 16, ft(t): 4 qct(c,t) : 12 ; ENDDATA ! Compute some bounds on flow; ! We (arbitrarily?) set the max flow into a sink = total flow from all sources; fssum = @SUM(source: fs) ; @FOR(sink(s): femax(s) = fssum) ; @FOR(source(s): @FOR(sink(e): fsemax1(s, e) = @MIN(concentration(c): qcemax(c,e)*fssum/qcs(c, s)) ; fsemax(s, e) = @if(fsemax1(s, e) #GT# fs(s), fs(s), fsemax1(s, e)) ; )) ; ! Minimize sum of variable + fixed costs; min = z; z = @sum(se: cse * fse + cyse * yse) + ! Source to sink costs; @sum(st: cst * fst + cyst * yst) + ! Source to treatment costs; @sum(te: cte * fte + cyte * yte) + ! Treatment to sink costs; @sum(tt: ctt * ftt + cytt * ytt) + ! Treatment to treatment costs; @sum(treatment(t): costt(t) * ft(t) + costyt(t) * yt(t) ) ; ! Treatment costs; ! 1.0 Flow balance constraints; ! 1.1 Flow balance at source s. fs(s) is given, = flows to end sinks + flows to treatment facilities; @FOR(source(s): [f1] fs(s) = @SUM(sink(e): fse(s, e)) + @SUM(treatment(t): fst(s, t)) ) ; ! 1.2 Flow balance at treatment t; @FOR(treatment(t): ! Set flow throught treatment t, ft(t) = flow into t from sources and other treatments; [f2] @SUM(source(s): fst(s, t)) + @SUM(treatment(t2)| t #ne# t2: ftt(t2, t)) = ft(t) ; ! = Total flow out of treatment t; [f3] ft(t) = @SUM(treatment(t2)|t #ne# t2: ftt(t, t2)) + @SUM(sink(e): fte(t, e)) ; ) ; [f4] @sum(source(s): fs(s)) = @sum(se(s,e): fse(s,e)) + @sum(te(t,e): fte(t,e)) ; ! Redundant? implied by f1+f2+f3?; ! 1.3 Flow balance at sink e; @FOR(sink(e): [f5] @sum(source(s): fse(s, e)) + @sum(treatment(t): fte(t, e)) <= femax(e)) ; ! Redundant? implied by f1 + f3 - f2? ; ! 2.0 Force the y variables, based on the flow, f, variables; ! 2.1 Source to treatments; @FOR(st(s,t): cyst = dst * cfix ; cst = dst * cvar/v ; @BIN(yst) ; ! Lower and upper limits on flow from source s to treatment t; [c4] fmin * yst(s,t) <= fst(s,t) ; [c5] fssum *yst(s,t) >= fst(s,t) ; ) ; ! 2.2 Source to sinks; @FOR(se(s,e) : cyse(s,e) = dse(s,e) * cfix ; cse(s,e) = dse(s,e) * cvar/v ; @BIN(yse(s,e)) ; ! Lower and upper limits on flow, fse(sk,e), from source s to sink e; [cfn] fmin * yse(s,e) <= fse(s,e) ; [cfx] fsemax * yse(s,e) >= fse(s,e) ; ) ; ! 2.3 Treatments; @FOR( treatment(t): @BIN(yt(t)) ; [c6] fmin *yt(t) <= ft(t) ; [c7] fssum*yt(t) >= ft(t) ; ); ! 2.4 Treatments to treatments; @FOR(tt(t1,t2): cytt(t1,t2) = dtt(t1,t2) * cfix ; ctt(t1,t2) = dtt(t1,t2) * cvar/v ; @BIN(ytt(t1,t2)) ; ! Lower and upper limits on flow from treatment t1 to treatment t1; [c2] fmin * ytt(t1,t2) <= ftt(t1,t2) ; [c3] fssum * ytt(t1,t2) >= ftt(t1,t2) ; ) ; @FOR(treatment(t): ! Can ship from t to t2, or t2 to t, but not both; @FOR(treatment(t2)| t #ne# t2: [L1] ytt(t, t2) + ytt(t2, t) <= 1 )) ; @FOR(treatment(t): ytt(t, t) = 0 ; ftt(t, t) = 0 ) ; ! 2.5 Treatments to sinks; @FOR(te: cyte = dte * cfix ; cte = dte * cvar/v ; @BIN(yte) ; [L2] fmin * yte <= fte ; [L3] fssum * yte >= fte ; ) ; ! 3.0 Quality constraints; ! Compute qct(c,t) flow out of treatment t; @FOR(concentration(c): @FOR(treatment(t): !(1-rct(c,t))*(pollutant c flowing into t = (pollutant c flowing out of t); [cmpq1] (1 - rct(c, t))*(@SUM(source(s):qcs(c, s) *fst(s, t)) + @sum(treatment(t2)|t #ne# t2: qct(c, t2)*ftt(t2, t))) = ! Nonlinear; qct(c,t)*@sum(treatment(t2)|t #ne# t2: ftt(t, t2)) ! Nonlinear; + qct(c,t)*@sum(sink(e): fte(t, e)) ; ! Nonlinear; [cmpq2] (1 - rct(c, t))*(@SUM(source(s):qcs(c, s)*fst(s, t)) + @sum(treatment(t2)|t #ne# t2: qct(c, t2)*ftt(t2, t))) = ! Nonlinear (Redundant? cmq1 & cmq2 equivalent | definition of ft?); qct(c,t)*ft(t) ; ! Nonlinear; [qub] qct(c, t) <= (1 - rct(c, t))*@MAX(source(s): qcs(c, s)); ! Some enforced bounds on the qct; @BND(0, qct(c,t), qctsub(c,t)); )) ; ! Here is the driving quality constraint. For each quality c, the fraction concentration going into sink e cannot exceed qcemax(c,e); @FOR(concentration(c): @FOR(sink(e): [qsink] qcemax(c, e)* ( @SUM(source(s): fse(s,e)) + @SUM(treatment(t): fte(t,e))) >= @SUM(source(s): qcs(c, s)*fse(s,e)) + @sum(treatment(t): qct(c, t)*fte(t,e)); ! Nonlinear; )) ;