import os
import gams

def get_model_text():
    return """
$title Electricity Generation

*###############################################################################
*                                DATA
*###############################################################################
set
    g       set of technologies
    /coal   technology running on hard coal
     nuc    nuclear
     ccgt   combined cycle gas turbine/
    f       set of fuels
    /hcoa   hard coal
     uran
     ng      natural gas/
    m(g,f)  mapping from technologies to
    /coal.hcoa
     nuc.uran
     ccgt.ng/
;

parameter
    pf(f)   fuel price  in $
    /hcoa   2
     uran   1
     ng     3/
    c(g)    other unit cost in $
    /coal   1
     nuc    1
     ccgt   2/
    eta(g)  heat efficiency
    /coal   0.35
     nuc    0.4
     ccgt   0.5/
    theta(g) availability in hours
    /coal   20
     nuc    25
     ccgt   15/
    cap(g)  installed capacity in MW
    /coal   1
     nuc    1
     ccgt   1/
;

scalar
    d       electricity demand in MWH       /50/
    pe0     electricity price               /8/
    elas    price elasticity of demand      /-0.5/
    dem_a   intercept demand function
    dem_b   slope demand function
;

dem_b = -elas*d/pe0;
dem_a = (1-elas)*d;

*###############################################################################
*                            MODEL
*###############################################################################
Variable
    WELFARE welfare (QCP objective)
    COST    total generation cost in $ (LP objective)
;

Positive Variable
    X(g)    generation in MWh
    DEM     demand in MWh
;

Equation
    obj_qp          objective QP = WELFARE definition in $
    obj_lp          objective LP = COST definition in $
    mkt_E_LP        market clearing electricity in MWh - LP version
    mkt_E_QP        market clearing electricity in MWh - QP version
    mkt_CAP(g)      market clearing capacity in MWh
;


obj_qp..
    WELFARE         =E=     dem_a/dem_b*DEM - 1/(2*dem_b)*(DEM*DEM)
                            - COST
;

obj_lp..
    COST            =E=     sum(g, X(g)*c(g)) + sum(m(g,f), pf(f)*X(g)/eta(g))
;

mkt_E_lp..
    sum(g, X(g))    =E=     d
;

mkt_E_qp..
    sum(g, X(g))    =G=     DEM
;


mkt_CAP(g)..
    cap(g)*theta(g) =G=     X(g)
;

model elec_LP /obj_lp, mkt_E_lp, mkt_CAP/;
model elec_QP /obj_qp, obj_lp, mkt_E_qp, mkt_CAP/;

* Set some initial values
X.L(g) = 1;"""

def print_model_stat(mi, objective):
    """ Prints solution statistics
    :param mi: <gams.GamsModelInstance>
    :param objective: <string> name of objevtive variable
    """
    print("Modelstatus: %d " % mi.model_status)
    print("Solvestatus: %d" % mi.solver_status)
    if objective is not None:
        print("Objective: %f" % mi.sync_db.get_variable(objective).find_record().level)


# create model space and prepare model
ws = gams.GamsWorkspace(os.getcwd())
cp = ws.add_checkpoint()
job_ = ws.add_job_from_string(get_model_text())
job_.run(checkpoint=cp)
gdxdata = job_.out_db

# Instantiate LP model with fuel prices as modifier
print("----- LP model with fuel price modifier")
mi_lp = cp.add_modelinstance()
pf = mi_lp.sync_db.add_parameter("pf", 1, "fuel prices")
mi_lp.instantiate("elec_LP use LP min COST", gams.GamsModifier(pf))
gdxdata.get_symbol("pf").copy_symbol(pf)   # copy orginal data from checkpoint
mi_lp.solve()
print_model_stat(mi_lp, "COST")

# Instantiate QP model with fuel prices as modifier
print("----- QCP model with fuel price modifier")
mi_qp1 = cp.add_modelinstance()
pf = mi_qp1.sync_db.add_parameter("pf", 1, "fuel prices")
mi_qp1.instantiate("elec_QP use QCP max WELFARE", gams.GamsModifier(pf))
gdxdata.get_symbol("pf").copy_symbol(pf)   # copy orginal data from checkpoint
mi_qp1.solve()
print_model_stat(mi_qp1, "WELFARE")

# Instantiate QP model with fuel prices AND demand intercept as modifier
print("----- QCP model with fuel price and demand intercept modifier")
mi_qp2 = cp.add_modelinstance()
pf = mi_qp2.sync_db.add_parameter("pf", 1, "fuel prices")
dem_a = mi_qp2.sync_db.add_parameter("dem_a", 0, "demand_intercept")
mi_qp2.instantiate("elec_QP use QCP max WELFARE",
                   [gams.GamsModifier(pf), gams.GamsModifier(dem_a)])
gdxdata.get_symbol("pf").copy_symbol(pf)   # copy orginal data from checkpoint
gdxdata.get_symbol("dem_a").copy_symbol(dem_a)
mi_qp2.solve()
print_model_stat(mi_qp2, "WELFARE")

# Instantiate QP model with fuel prices AND demand intercept as modifier
print("----- QCP model with fuel price, demand intercept, and slope modifier")
mi_qp3 = cp.add_modelinstance()
pf = mi_qp3.sync_db.add_parameter("pf", 1, "fuel prices")
dem_a = mi_qp3.sync_db.add_parameter("dem_a", 0, "demand_intercept")
dem_b = mi_qp3.sync_db.add_parameter("dem_b", 0, "demand_intercept")
mi_qp3.instantiate("elec_QP use QCP max WELFARE",
                   [gams.GamsModifier(pf),
                   gams.GamsModifier(dem_a),gams.GamsModifier(dem_b)])
gdxdata.get_symbol("pf").copy_symbol(pf)   # copy orginal data from checkpoint
gdxdata.get_symbol("dem_a").copy_symbol(dem_a)
gdxdata.get_symbol("dem_b").copy_symbol(dem_b)
mi_qp3.solve()
print_model_stat(mi_qp3, "WELFARE")