Update costs for pre-y₀ vintages in periods y₀ and later

[khaeru]

@OFR-IIASA noted today in Slack that:

For historic vintages the the costs are not updated for the activity years in the model horizon.

The expectation is that the code should provide these values, by default.

He also provided this code as a work-around/starting point for a PR that would resolve the issue:

fix_cost = fix_cost[fix_cost["technology"].isin(model_tec_set)]

# Ensure that fix_costs are also added for activity years of vintage years
# prior to the firstmodelyear and remove any entries for activity years for
# which the technology is not active.
tmp_fc = pd.DataFrame()
remove = []
for node in fix_cost.node_loc.unique().tolist():
    for tec in fix_cost.loc[fix_cost.node_loc == node].technology.unique().tolist():
        yva_found = False
        tmp = fix_cost.loc[(fix_cost.node_loc == node) & (fix_cost.technology == tec)].copy()
        try:
            yva = scen.vintage_and_active_years([node, tec])
            yva_found = True
        except Exception:
            yva = scen.vintage_and_active_years()
        if yva_found == False or max(yva.year_vtg) < max(scen.set("year")):
            inp = scen.par("input", filters={"node_loc": node, "technology": tec})[["year_vtg", "year_act"]]
            inp = inp.loc[inp.year_act >= scen.firstmodelyear]
            out = scen.par("output", filters={"node_loc": node, "technology": tec})[["year_vtg", "year_act"]]
            out = out.loc[out.year_act >= scen.firstmodelyear]
            if not inp.empty:
                yva = yva.merge(inp, on=["year_vtg", "year_act"], how="right")
            elif not out.empty:
                yva = yva.merge(out, on=["year_vtg", "year_act"], how="right")
            else:
                yva = scen.par("relation_activity", filters={"node_loc": node, "technology": tec})[["year_rel", "year_act"]].drop_duplicates().rename(columns={"year_rel": "year_vtg"})
                if yva.empty:
                    yva = scen.par("relation_total_capacity", filters={"node_rel": node, "technology": tec})[["year_rel"]].drop_duplicates().rename(columns={"year_rel": "year_vtg"})
                    yva["year_act"] = yva["year_vtg"]
        if yva.empty:
            print(f"No data found for technology {tec} in node {node}. Continuing")
            remove.append([node, tec])
            continue
        
        tmp = tmp.merge(yva, on=["year_vtg", "year_act"], how="right").sort_values(by=["year_act", "year_vtg"]).bfill()
        tmp_fc = pd.concat([tmp_fc, tmp])

fix_cost = tmp_fc.copy()

[OFR-IIASA]

I have written some code (based on what was provided above) that addresses some of the short comings of the cost tool in the case that is being used to update cost to a message_ix scenario.

• Fix cost are adjusted to be applied only to technology specific activity years.
• Fix_cost are applied to also adjust the activity years of historic vintages.
• Fix and inv_cost are only added for regions in which the technologies exist.
• Before adding parameters, the old entries are removed.
  The code can be copied to an ipynb or a py file, which can be used as a post-processing step, after the cost tool has been run and can later be integrated into the cost tool itself.

import message_ix
import ixmp
mp = ixmp.Platform()

from message_ix_models.model.material.data_util import (
    gen_te_projections,
)


base = message_ix.Scenario(mp, "SSP_SSP2_v3.0", "baseline_DEFAULT_step_14")
scen = base.clone("SSP_SSP2_v3.0", "baseline_DEFAULT_step_14_fix_cost_test", keep_solution=False)


inv_cost, fix_cost = gen_te_projections(
            scen, "SSP2", method="gdp", module="energy", reduction_year=2150
        )


import pandas as pd
# Ensure that fix_costs are also added for activity years of vintage years
# prior to the firstmodelyear and remove any entries for activity years for
# which the technology is not active.

par_input = scen.par("input", filters={"technology": fix_cost.technology.unique()})
par_output = scen.par("output", filters={"technology": fix_cost.technology.unique()})
par_relation_activty = scen.par("relation_activity", filters={"technology": fix_cost.technology.unique()})
par_relation_total_capacity = scen.par("relation_total_capacity", filters={"technology": fix_cost.technology.unique()})

tmp_fc = pd.DataFrame()
remove = []
for node in fix_cost.node_loc.unique().tolist():
    for tec in fix_cost.loc[fix_cost.node_loc == node].technology.unique().tolist():
        print(f"Adapting fix_cost for {tec} in node {node}.")
        
        # Filter out node/technology specific fix-cost
        tmp = fix_cost.loc[(fix_cost.node_loc == node) & (fix_cost.technology == tec)].copy()

        # Combine input and output vintage and activity years
        # Filter out node and technology specific input parameter
        inp = par_input.loc[(par_input.technology == tec) & (par_input.node_loc == node)].copy()[["year_vtg", "year_act"]]
        # Filter all entires for year_act which fall into the model-horizon
        inp = inp.loc[inp.year_act >= scen.firstmodelyear]
        
        # Filter out node and technology specific output parameter
        out = par_output.loc[(par_output.technology == tec) & (par_output.node_loc == node)].copy()[["year_vtg", "year_act"]]
        # Filter all entires for year_act which fall into the model-horizon
        out = out.loc[out.year_act >= scen.firstmodelyear]

        if not inp.empty and not out.empty:
            print("    Using input and output parameter to derive vintage and active years")
            yva = inp.merge(out, on=["year_vtg", "year_act"], how="inner")
        elif not inp.empty and out.empty:
            print("    Using only input parameter to derive vintage and active years")
            yva = inp.copy()
        elif inp.empty and not out.empty:
            print("    Using only output parameter to derive vintage and active years")
            yva = out.copy()
        else:
            print("    Using relation_activity parameter to derive vintage and active years")
            yva = par_relation_activty.loc[(par_relation_activty.technology == tec) & (par_relation_activty.node_rel == node)].copy()[["year_rel", "year_act"]].drop_duplicates().rename(columns={"year_rel": "year_vtg"})
        if yva.empty:
            print("    Using relation_total_capacity parameter to derive vintage and active years")
            yva = par_relation_total_capacity.loc[(par_relation_total_capacity.technology == tec) & (par_relation_total_capacity.node_rel == node)].copy()[["year_rel"]].drop_duplicates().rename(columns={"year_rel": "year_vtg"})
            yva["year_act"] = yva["year_vtg"]
        if yva.empty:
            print(f"No data found for technology {tec} in node {node}. Continuing")
            remove.append([node, tec])
            continue

        # For some technologies, like nuc_lc, the costs are projected only for future vintages,
        # but the last vintage in the model is prior to the furstmodelyear. In order to backfill
        # values and outer merge is done, the values are backfilled, and then an "right" merge
        # is performed to only include historic vintages
        if max(yva.year_vtg) < min(tmp.year_vtg):
            tmp = tmp.merge(yva, on=["year_vtg", "year_act"], how="outer").sort_values(by=["year_act", "year_vtg"]).bfill()
            tmp = tmp.merge(yva, on=["year_vtg", "year_act"], how="right").sort_values(by=["year_act", "year_vtg"])

        else:
            # Backfill activity years within the model horizon
            tmp = tmp.merge(yva, on=["year_vtg", "year_act"], how="right").sort_values(by=["year_act", "year_vtg"]).bfill()
            # Forwardfill any years for historic vintages outside the model-year vintage
            # This can come from the lifetime extensions algorithm applied as part of the capacity calibration workflow
            tmp.loc[tmp.year_vtg < scen.firstmodelyear] = tmp.loc[tmp.year_vtg < scen.firstmodelyear].ffill()

        assert tmp.loc[tmp.value.isnull()].empty

        tmp_fc = pd.concat([tmp_fc, tmp])
fix_cost = tmp_fc.copy()

for items in remove:
    print(f"Removing technology {items[1]} from node {items[0]}")
    inv_cost = inv_cost.loc[~((inv_cost.node_loc == items[0]) & (inv_cost.technology == items[1]))]

with scen.transact():
    tec = inv_cost.technology.unique().tolist()
    inv_cost_remove = scen.par("inv_cost", filters={"technology": tec})
    inv_cost_remove = inv_cost_remove.loc[inv_cost_remove.year_vtg >= scen.firstmodelyear]
    scen.remove_par("inv_cost", inv_cost_remove)

    tec = fix_cost.technology.unique().tolist()
    fix_cost_remove = scen.par("fix_cost", filters={"technology": tec})
    fix_cost_remove = fix_cost_remove.loc[fix_cost_remove.year_act >= scen.firstmodelyear]
    scen.remove_par("fix_cost", fix_cost_remove)

    scen.add_par("inv_cost", inv_cost)
    scen.add_par("fix_cost", fix_cost)

scen.solve()