In [5]:
### BEGIN SOLUTION
plt.figure()
plt.plot(t,eps)
plt.grid(True)
plt.title('Annual Anthropogenic Emissions')
plt.xlabel('year')
plt.ylabel('GtC/year')
### END SOLUTION
Out[5]:
Text(0, 0.5, 'GtC/year')
In the previous CGI, we constructed annual emissions schedules, $\epsilon(t)$, in units like GtC/year. But there are two good reasons to think about a related quantity, the total amount of carbon that will be released by humans to the atmosphere. One is a reason you already know: too much emission will result in an unlivable planet. Another reason is that there is a finite amount of fossil fuel carbon in the ground. It wasn't always know which of these would ultimately limit anthropogenic emissions (as it turns out, it's the latter). Therefore, one important form of analysis is to evaluate the total anthropogenic emission a given emission scenario represents.
One way to formulate the problem is to use integral calculus,
$$ E(t) = \int_{-\infty}^t \epsilon(t) \ dt \ \ \ \ \ (1) $$If you have an algebraic form for $\epsilon$ (and you're skilled at calculus), you might think about doing this integration analytically, to get an algebraic form for $E(t)$.
We're not going to do that here, though, because there's an easy way to do it in Python. The key is a numpy function called "np.cumsum," which adds up the emissions year after year. But you can see that there's a problem with that: the units would be the same units as $\epsilon(t)$, namely GtC/year. To convert the result of np.cumsum into a total accumulated emission in GtC, you'd need to multiply it by the time interval between steps in the array you're using to represent the emissions. Details are given in the code below.
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import h5io
%matplotlib notebook
To get the file you generated in ScheduledFlows into your current work space, download it to your laptop or desktop, then upload it to the current folder.
Once you've done that, in the cell below, use h5io.read_hdf5 to load the scenario into Python as a dictionary named epsdictionary_fromfile (just as you did in ScheduledFlows).
# Load in your scenario as a dictionary, using h5io.read
### BEGIN SOLUTION
# Assign a name to your file
filename = 'Peaks_in_2030.hdf5'
# This loads the stored dictionary
epsdictionary_fromfile = h5io.read_hdf5(filename)
### END SOLUTION
The cell below extracts time and emissions arrays in the scenario.
# Here we're using "display" to double-check the metadata in your dictionary
display(epsdictionary_fromfile)
# This extracts the dataframe from the dictionary
epsdf = epsdictionary_fromfile['dataframe']
# This extracts the time and emissions from the dataframe
t = np.array(epsdf['time'])
eps = np.array(epsdf['emissions'])
{'dataframe': time emissions
0 1750.000000 1.539637e-02
1 1750.667780 1.565556e-02
2 1751.335559 1.591912e-02
3 1752.003339 1.618711e-02
4 1752.671119 1.645961e-02
.. ... ...
595 2147.328881 1.430219e-07
596 2147.996661 1.272478e-07
597 2148.664441 1.132135e-07
598 2149.332220 1.007270e-07
599 2150.000000 8.961769e-08
[600 rows x 2 columns],
'delta_t_trans': 15,
'emission units': 'GtC/year',
'eps_0': 12.9,
'k': 0.025,
't_0': 2020,
't_peak': 2030,
't_trans': 2039.7295507452766}
In the cell below, plot the emissions we just extracted as a function of time (this will remind us what the scenario looks like, and also verify that the data haven't been corrupted).
### BEGIN SOLUTION
plt.figure()
plt.plot(t,eps)
plt.grid(True)
plt.title('Annual Anthropogenic Emissions')
plt.xlabel('year')
plt.ylabel('GtC/year')
### END SOLUTION
Text(0, 0.5, 'GtC/year')
Eq. (1) says, if you integrate all the emissions up to a certain time, you'll get another function of time that is the accumulated carbon emitted up to that point. Turns out, Numpy has a built-in function, cumsum, that does almost that.
# Find accumulated emissions
E = np.cumsum(eps)
# Graph it
plt.figure()
plt.plot(t,E)
plt.grid()
plt.xlabel('year')
plt.ylabel('GtC')
plt.title('Result of cumsum')
Text(0.5, 1.0, 'Result of cumsum')
What's missing in the above is that cumsum doesn't take into account how long each time step is. How do you find that out? The easy way is to say dt = t[1]-t[0], where t[0] is the time at the beginning of the series, and t[1] is the time after the first time step.
In the cell below, duplicate the integration by cumsum we did in the previous cell, but multiply the result by $dt$. Then graph your result.
# Get the time interval from the array t
### BEGIN SOLUTION
dt = t[1]-t[0]
print('Time step =', dt)
### END SOLUTION
# Find the accumulation by numerical integration of emissions
### BEGIN SOLUTION
E = np.cumsum(eps)*dt
### END SOLUTION
# Graph the cumulative emissions as a function of time
### BEGIN SOLUTION
plt.figure()
plt.plot(t,E)
plt.grid()
plt.xlabel('year')
plt.ylabel('GtC')
plt.title('Cumulative carbon emitted to atmosphere')
### END SOLUTION
Time step = 0.6677796327212491
Text(0.5, 1.0, 'Cumulative carbon emitted to atmosphere')
It is thought that known reserves of fossil carbon (mostly in the form of coal) tally up to around 4000 GtC (see https://www.nature.com/articles/nature14016). Hopefully, your cumulative total is less than that amount (otherwise, we have an unrealistic scenario!). In the cell below, calculate the percentage of known reserves of fossil carbon does your schedule leaves in the ground.
Some hints on how to do this:
### BEGIN SOLUTION
remaining = (1-E[-1]/4000)*100
print(remaining)
### END SOLUTION
77.91807282729626
Almost done! To double-check everything is OK, repeat the "Three steps for refreshing and saving your code," and press the "Validate" button (as usual).