In [1]:
%pylab inline
Populating the interactive namespace from numpy and matplotlib
In [2]:
from brightway2 import *
from bw2temporalis import *
import arrow
import numpy as np

The basic approach

After banging my head against the unspoken assumptions and terse methodology of the ESPA paper, I finally decided to go with a simpler and easier to program alternative, what I calling (at least for now) temporal graph traversal. The idea is that we use the existing machinery in Brightway2 for graph traversal to create a supply chain graph of all the possibly important nodes and vertices, and then traverse this graph, accounting for the temporal differences between inputs and outputs.

There is a more complete explanation in the brightway2-temporalis documentation.

Simple example

Let's look at a simple example of a wooden house. The numbers here are wildly incorrect, and should be politely ignored. Here is a graph of our supply chain:

In [3]:
from IPython.core.display import Image 
Image(filename='images/example-system.png') 
Out[3]:

And here is how we define this system as an inventory dataset. Note the use of the new keys temporal distribution and absolute date. You don't, however, need to read through this in detail.

In [4]:
EXAMPLE_DATABASE = {
    ("example", "replanting"): {
        "type": "process",
        "name": "replanting",
        "unit": "hectare",
        "exchanges": [
            {
                "type": "biosphere",
                "amount": -100,
                "input": ("biosphere", "co2"),
                "temporal distribution": [
                    (0, -2.5),  # Use default of years
                    (1, -5),
                    (2, -10),
                    (3, -10),
                    (4, -5)
                ] + [(x + 5, -2.5) for x in range(675 / 25)]
            }
        ]
    },
    ("example", "harvest"): {
        "type": "process",
        "name": "forest harvest",
        "unit": "kilogram",
        "exchanges": [
            {
                "type": "biosphere",
                "amount": 1.5,
                "input": ("biosphere", "co2"),
                "comment": "flux from soil carbon",
                "temporal distribution": [(x, 0.25) for x in range(0, 30, 5)]
            }, {
                "type": "technosphere",
                "amount": 0.002,
                "input": ("example", "replanting"),
            }
        ]
    },
    ("example", "sawing"): {
        "type": "process",
        "name": "sawing",
        "unit": "kilogram",
        "exchanges": [
            {
                "type": "biosphere",
                "amount": 0.05,
                "input": ("biosphere", "co2"),
            },
            {
                "type": "technosphere",
                "amount": 1,
                "input": ("example", "harvest"),
                "temporal distribution": [(-2, 1)]  # Two years to dry
            },
            {
                "type": "technosphere",
                "amount": 1e-6,
                "input": ("example", "sawmill"),
            }
        ]
    },
    ("example", "sawmill"): {
        "type": "process",
        "name": "sawmill",
        "unit": "plant",
        "absolute date": "1970-1-1",
        "exchanges": [
            {
                "type": "biosphere",
                "amount": 1e5,
                "input": ("biosphere", "co2"),
            }
        ]
    },
    ("example", "house"): {
        "type": "process",
        "name": "house",
        "unit": "house",
        "exchanges": [
            {
                "type": "technosphere",
                "amount": 1e4,
                "input": ("example", "sawing"),
                "temporal distribution": [
                    (0, 9e3),
                    (10, 250),
                    (20, 250),
                    (30, 250),
                    (40, 250)
                ]
            }
        ]
    }
}

Register our example database, and write it to disk:

In [5]:
example_db = Database("example")
if "example" not in databases:
    example_db.register()
example_db.write(EXAMPLE_DATABASE)
example_db.process()

Check to make sure our temporal distribution values are correct.

In [6]:
check_temporal_distribution_totals("example")
Out[6]:
True

Extended dynamic characterization function: CO2

Extended dynamic characterization factors have two qualities: their total amounts are a function of time, and their impact is spread throughout time, as opposed to happening at only one time.

We can understand impacts defined as a function of time with the concept of climate change urgency - emissions now are worse than those in the future. Instead of a single global warming potential, we spread the radiative forcing throughout time, with a time horizon of 100 years.

The numbers in this example are completely artificial.

In [7]:
def cf_function(tp):
    """Importance of CO2 halves every twenty years from 2010"""
    from bw2temporalis import data_point
    from datetime import timedelta
    import arrow

    year_in_seconds = timedelta(seconds=365.24 * 24 * 60 * 60)
    base_cf = min(1., 0.5 ** ((arrow.get(tp.dt) - arrow.get(2010, 1, 1)).days / 365.24 / 20))
    # Assume median atmospheric residence time of 50 years
    dropoff = lambda x: 0.5 ** (x / 50.)
    TIME_HORIZON = 100
    return [data_point(
        tp.dt + x * year_in_seconds, 
        tp.flow, 
        tp.ds, 
        dropoff(x) * base_cf * tp.amount
        ) for x in range(0, TIME_HORIZON, 10)]
In [8]:
starting_points = [
    data_point(arrow.get(2000, 1, 1), 'foo', 'bar', 1),
    data_point(arrow.get(2020, 1, 1), 'foo', 'bar', 1),
    data_point(arrow.get(2040, 1, 1), 'foo', 'bar', 1),
]

results = [cf_function(x) for x in starting_points]

figure(figsize=(8, 6))

for ds, label, style in zip(results, ("Emission at 2000", "Emission at 2020", "Emission at 2040"), ('b-', 'g-', 'r-')):
    plot_date([x.dt.datetime for x in ds], [x.amount for x in ds], style, lw=3, alpha=0.5, label=label)
xlabel("Year")
ylabel(r"CF")
legend()
xlim(results[0][0].dt.datetime, results[0][-1].dt.datetime)
print "Graph of example extended dynamic CFs"
Graph of example extended dynamic CFs

Defining a dynamic LCIA method

Now that we have our function, we need to define, register, and write our dynamic LCIA method:

In [9]:
NAME = "Dynamic GWP (CO2)"
dynamic_method = DynamicIAMethod(NAME)
if NAME not in dynamic_methods:
    dynamic_method.register()
dynamic_method.write({("biosphere", "co2"): """def %s(dt):
    from datetime import timedelta
    import arrow
    import collections

    nt = collections.namedtuple('nt', ['dt', 'amount'])
    dt = arrow.get(dt)
    year_in_seconds = timedelta(seconds=365.24 * 24 * 60 * 60)
    base_cf = min(1., 0.5 ** ((arrow.get(dt) - arrow.get(2010, 1, 1)).days / 365.24 / 20))
    # Assume median atmospheric residence time of 50 years
    dropoff = lambda x: 0.5 ** (x / 50.)
    TIME_HORIZON = 100
    return [nt(
        dt + x * year_in_seconds, 
        dropoff(x) * base_cf
        ) for x in range(TIME_HORIZON)]    
"""})

We also need to create our "worst case" static LCIA method.

In [10]:
if ("ipcc gwp", "worst case") in methods:
        del methods[("ipcc gwp", "worst case")]
dynamic_method.to_worst_case_method(("ipcc gwp", "worst case"))
Out[10]:
<bw2data.method.Method at 0x10b911dd0>
In [11]:
Method(("ipcc gwp", "worst case")).load()
Out[11]:
[[('biosphere', 'co2'), 54.47693046453665]]

Temporal graph traversal

We are now ready to do dynamic LCA calculations.

In [12]:
dynamic_lca = DynamicLCA(
    demand={("example", "house"): 1}, 
    dynamic_method='Dynamic GWP (CO2)', 
    worst_case_method=(u'ipcc gwp', u'worst case'), 
    now=arrow.now()
)
dynamic_lca.calculate()

Plot the results

In [13]:
figure(figsize=(12, 8))
xs, ys = dynamic_lca.timeline.characterize_dynamic('Dynamic GWP (CO2)', stepped=True)
plot_date([x.datetime for x in xs], ys, 'r-', lw=2, alpha=0.5, label="Dynamic")

# Format graph
xlabel("year")
ylabel(r"GWP")

print "Plot of dynamic GWP"
Plot of dynamic GWP
In [13]: