Merge pull request #281 from khaeru/issue/278

Correct logic of Scenario.years_active()

RELEASE_NOTES.md

@@ -10,6 +10,7 @@
 
 ## All changes
 
+- [#281](https://github.com/iiasa/message_ix/pull/281): Test and improve logic of `years_active` and `vintage_and_active_years`.
 - [#269](https://github.com/iiasa/message_ix/pull/269): Enforce 'year'-indexed columns as integers.
 - [#256](https://github.com/iiasa/message_ix/pull/256): Update to use :obj:`ixmp.config` and improve CLI.
 - [#255](https://github.com/iiasa/message_ix/pull/249): Add :mod:`message_ix.testing.nightly` and `nightly` CLI command group for slow-running tests.

message_ix/core.py

@@ -1,6 +1,6 @@
 import collections
 from functools import lru_cache
-import itertools
+from itertools import product
 
 import ixmp
 from ixmp.utils import as_str_list, pd_read, pd_write, isscalar, logger
@@ -56,6 +56,9 @@ def _year_as_int(self, name, df):
 
         if len(year_idx):
             return df.astype({col_name: 'int' for _, col_name in year_idx})
+        elif name == 'year':
+            # The 'year' set itself
+            return df.astype(int)
         else:
             return df
 
@@ -289,50 +292,56 @@ def vintage_and_active_years(self, ya_args=None, in_horizon=True):
 
         Parameters
         ----------
-        ya_args : tuple of (node, technology, year_vtg), optional
-            Arguments to :meth:`ixmp.Scenario.years_active`.
+        ya_args : tuple of (node, tec, yr_vtg), optional
+            Arguments to :meth:`years_active`.
         in_horizon : bool, optional
-            Restrict years returned to be within the current model horizon.
+            Only return years within the model horizon
+            (:obj:`firstmodelyear` or later).
 
         Returns
         -------
         pandas.DataFrame
-            with columns, "year_vtg" and "year_act", in which each row is a
+            with columns 'year_vtg' and 'year_act', in which each row is a
             valid pair.
         """
-        horizon = self.set('year')
         first = self.firstmodelyear
 
+        # Prepare lists of vintage (yv) and active (ya) years
         if ya_args:
             if len(ya_args) != 3:
                 raise ValueError('3 arguments are required if using `ya_args`')
-            years_active = self.years_active(*ya_args)
-            combos = itertools.product([ya_args[2]], years_active)
+            ya = self.years_active(*ya_args)
+            yv = ya[0:1]  # Just the first element, as a list
         else:
-            combos = itertools.product(horizon, horizon)
+            # Product of all years
+            yv = ya = self.set('year')
 
-        combos = [(int(y1), int(y2)) for y1, y2 in combos]
+        # Predicate for filtering years
+        def _valid(elem):
+            yv, ya = elem
+            return (yv <= ya) and (not in_horizon or (first <= ya))
 
-        def valid(y_v, y_a):
-            ret = y_v <= y_a
-            if in_horizon:
-                ret &= y_a >= first
-            return ret
-
-        year_pairs = [(y_v, y_a) for y_v, y_a in combos if valid(y_v, y_a)]
-        v_years, a_years = zip(*year_pairs)
-        return pd.DataFrame({'year_vtg': v_years, 'year_act': a_years})
+        # - Cartesian product of all yv and ya.
+        # - Filter only valid years.
+        # - Convert to data frame.
+        return pd.DataFrame(
+            filter(_valid, product(yv, ya)),
+            columns=['year_vtg', 'year_act'])
 
     def years_active(self, node, tec, yr_vtg):
         """Return years in which *tec* of *yr_vtg* can be active in *node*.
 
+        The :ref:`parameters <params-tech>` ``duration_period`` and
+        ``technical_lifetime`` are used to determine which periods are partly
+        or fully within the lifetime of the technology.
+
         Parameters
         ----------
         node : str
             Node name.
         tec : str
             Technology name.
-        yr_vtg : str
+        yr_vtg : int or str
             Vintage year.
 
         Returns
@@ -348,11 +357,16 @@ def years_active(self, node, tec, yr_vtg):
 
         # Duration of periods
         data = self.par('duration_period')
-        # Cumulative sum of period duration for periods after the vintage year
-        data['age'] = data.where(data.year >= yv)['value'].cumsum()
-
-        # Return years where the age is less than or equal to the lifetime
-        return data.where(data.age <= lt)['year'].dropna().astype(int).tolist()
+        # Cumulative sum for periods including the vintage period
+        data['age'] = data.where(data.year >= yv, 0)['value'].cumsum()
+
+        # Return periods:
+        # - the tec's age at the end of the *prior* period is less than or
+        #   equal to its lifetime, and
+        # - at or later than the vintage year.
+        return data.where(data.age.shift(1, fill_value=0) < lt) \
+                   .where(data.year >= yv)['year'] \
+                   .dropna().astype(int).tolist()
 
     @property
     def firstmodelyear(self):

message_ix/model/MESSAGE/parameter_def.gms

@@ -141,6 +141,8 @@ Parameter
 ;
 
 ***
+* .. _params-tech:
+*
 * Parameters of the `Technology` section
 * --------------------------------------
 *

tests/test_core.py

@@ -5,6 +5,7 @@
 from numpy import testing as npt
 import pandas as pd
 import pandas.util.testing as pdt
+import pytest
 
 from message_ix import Scenario
 from message_ix.testing import make_dantzig
@@ -85,19 +86,15 @@ def test_add_spatial_hierarchy(test_mp):
 
 def test_vintage_and_active_years(test_mp):
     scen = Scenario(test_mp, *msg_args, version='new')
-    scen.add_horizon({'year': ['2000', '2010', '2020'],
-                      'firstmodelyear': '2010'})
+
+    years = [2000, 2010, 2020]
+    scen.add_horizon({'year': years, 'firstmodelyear': 2010})
     obs = scen.vintage_and_active_years()
     exp = pd.DataFrame({'year_vtg': (2000, 2000, 2010, 2010, 2020),
                         'year_act': (2010, 2020, 2010, 2020, 2020)})
     pdt.assert_frame_equal(exp, obs, check_like=True)  # ignore col order
 
-
-def test_vintage_and_active_years_with_lifetime(test_mp):
-    scen = Scenario(test_mp, *msg_args, version='new')
-    years = ['2000', '2010', '2020']
-    scen.add_horizon({'year': years,
-                      'firstmodelyear': '2010'})
+    # Add a technology, its lifetime, and period durations
     scen.add_set('node', 'foo')
     scen.add_set('technology', 'bar')
     scen.add_par('duration_period', pd.DataFrame({
@@ -137,6 +134,22 @@ def test_vintage_and_active_years_with_lifetime(test_mp):
                         'year_act': (2020,)})
     pdt.assert_frame_equal(exp, obs, check_like=True)  # ignore col order
 
+    # Advance the first model year
+    scen.add_cat('year', 'firstmodelyear', years[-1], is_unique=True)
+
+    # Empty data frame: only 2000 and 2010 valid year_act for this node/tec;
+    # but both are before the first model year
+    obs = scen.vintage_and_active_years(ya_args=('foo', 'bar', years[0]),
+                                        in_horizon=True)
+    pdt.assert_frame_equal(
+        pd.DataFrame(columns=['year_vtg', 'year_act']),
+        obs)
+
+    # Exception is raised for incorrect arguments
+    with pytest.raises(ValueError,
+                       match='3 arguments are required if using `ya_args`'):
+        scen.vintage_and_active_years(ya_args=('foo', 'bar'))
+
 
 def test_cat_all(test_mp):
     scen = Scenario(test_mp, *msg_args)
@@ -169,23 +182,63 @@ def test_add_cat_unique(test_mp):
 
 
 def test_years_active(test_mp):
-    scen = Scenario(test_mp, *msg_multiyear_args)
-    years = scen.years_active('seattle', 'canning_plant', '2020')
+    test_mp.add_unit('year')
+    scen = Scenario(test_mp, *msg_args, version='new')
+    scen.add_set('node', 'foo')
+    scen.add_set('technology', 'bar')
+
+    # Periods of uneven length
+    years = [1990, 1995, 2000, 2005, 2010, 2020, 2030]
+
+    # First period length is immaterial
+    duration = [1900, 5, 5, 5, 5, 10, 10]
+    scen.add_horizon({'year': years, 'firstmodelyear': years[-1]})
+    scen.add_par('duration_period',
+                 pd.DataFrame(zip(years, duration), columns=['year', 'value']))
+
+    # 'bar' built in period '1995' with 25-year lifetime:
+    # - is constructed in 1991-01-01.
+    # - by 1995-12-31, has operated 5 years.
+    # - operates until 2015-12-31. This is within the period '2020'.
+    scen.add_par('technical_lifetime', pd.DataFrame(dict(
+        node_loc='foo',
+        technology='bar',
+        unit='year',
+        value=25,
+        year_vtg=years[1]), index=[0]))
+
+    result = scen.years_active('foo', 'bar', years[1])
+
+    # Correct return type
     assert isinstance(years, list)
     assert isinstance(years[0], int)
-    npt.assert_array_equal(years, [2020, 2030])
+
+    # Years 1995 through 2020
+    npt.assert_array_equal(result, years[1:-1])
 
 
 def test_years_active_extend(test_mp):
     scen = Scenario(test_mp, *msg_multiyear_args)
-    scen = scen.clone(keep_solution=False)
+
+    # Existing time horizon
+    years = [2010, 2020, 2030]
+    result = scen.years_active('seattle', 'canning_plant', years[1])
+    npt.assert_array_equal(result, years[1:])
+
+    # Add years to the scenario
+    years.extend([2040, 2050])
     scen.check_out()
-    scen.add_set('year', ['2040', '2050'])
+    scen.add_set('year', years[-2:])
     scen.add_par('duration_period', '2040', 10, 'y')
     scen.add_par('duration_period', '2050', 10, 'y')
-    years = scen.years_active('seattle', 'canning_plant', '2020')
-    npt.assert_array_equal(years, [2020, 2030, 2040])
-    scen.discard_changes()
+
+    # technical_lifetime of seattle/canning_plant/2020 is 30 years.
+    # - constructed in 2011-01-01.
+    # - by 2020-12-31, has operated 10 years.
+    # - operates until 2040-12-31.
+    # - is NOT active within the period '2050' (2041-01-01 to 2050-12-31)
+    result = scen.years_active('seattle', 'canning_plant', '2020')
+    npt.assert_array_equal(result, years[1:-1])
 
 
 def test_new_timeseries_long_name64(test_mp):