option display_precision 9; # model unit_commit_lin.mod # data DATA/toy.dat # data DATA/u10t24.dat # data DATA/u20t24.dat # data DATA/u30t24.dat # data DATA/u40t24.dat # data DATA/u50t24.dat # data DATA/u60t24.dat # data DATA/u70t24.dat # data DATA/u80t24.dat # data DATA/u90t24.dat # data DATA/u100t24.dat model run.mod data run.dat option solver cplexamp; option show_stats 1; option cplex_options 'timing 1 mipdisplay 1 threads 1'; # option cplex_options 'threads=1 mipgap=1.e-12 integrality=1.e-9 absmipgap=1.e-12 qcpconvergetol=1.e-12'; # rel_boundtol=1.e-6'; # option solver gurobi_ampl; # option gurobi_options 'threads 1'; param epsilon := 1.e-6; param nCUTS; set TMP ordered; let {u in UNITS} nBRK[u] := 6; solve; for {u in UNITS} { for {i in 1..nBRK[u]} { printf "%g\t", p_med[u,i]; } printf "\n"; } let nCUTS := 0; for {u in UNITS} { printf "\n\nUNIT %s\n", u; let TMP := setof {i in 1..nBRK[u]} p_med[u,i]; for {t in 1..T : p[u,t] > 0 and p[u,t] not in TMP and abs(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t] - fuel[u,t])/(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t]) > epsilon} { printf "quadratic costs: %g\n", a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t]; printf "linearized costs: %g\n", fuel[u,t]; printf "error: %g\n", abs(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t] - fuel[u,t])/(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t]); let TMP := TMP union {p[u,t]}; let nCUTS := nCUTS + 1; } # display TMP; let nBRK[u] := card(TMP); let {i in 1..nBRK[u]} p_med[u,i] := member(i,TMP); # display p_med; } for {u in UNITS} { printf "\n\nUNIT %s, CURRENT ABS ERROR %g, CURRENT FUEL COST %g\n", u, sum {t in 1..T} abs(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t] - fuel[u,t]), sum {t in 1..T} (a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t]); } printf "\n\nCURRENT TOTAL ABS ERROR %g, CURRENT FUEL COST %g\n", sum {u in UNITS, t in 1..T} abs(a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t] - fuel[u,t]), sum {u in UNITS, t in 1..T} (a[u] * y[u,t] + b[u] * p[u,t] + c[u] * p[u,t] * p[u,t]); display sum{u in UNITS} nBRK[u], nCUTS, cost; display p; display on; display off; display y; display cost;