analysis.py

import collections
import numpy as np
import matplotlib.pyplot as plt
import sqlite3 as lite
import sys
from itertools import cycle
from matplotlib import cm
from pyne import nucname
from collections import Counter

if len(sys.argv) < 2:
    print('Usage: python analysis.py [cylus_output_file]')


def cursor(file_name):
    """Connects and returns a cursor to an sqlite output file

    Parameters
    ----------
    file_name: str
        name of the sqlite file

    Returns
    -------
    sqlite cursor3
    """
    con = lite.connect(file_name)
    con.row_factory = lite.Row
    return con.cursor()


def agent_ids(cur, archetype):
    """Gets all agentids from Agententry table for wanted archetype

        agententry table has the following format:
            SimId / AgentId / Kind / Spec /
            Prototype / ParentID / Lifetime / EnterTime

    Parameters
    ----------
    cur: cursor
        sqlite cursor3
    archetype: str
        agent's archetype specification

    Returns
    -------
    agentids: list
        list of all agentId strings
    """
    agents = cur.execute("SELECT agentid FROM agententry WHERE spec "
                         "LIKE '%" + archetype + "%' COLLATE NOCASE"
                         ).fetchall()

    return list(str(agent['agentid']) for agent in agents)


def prototype_id(cur, prototype):
    """Returns agentid of a prototype

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    prototype: str
        name of prototype

    Returns
    -------
    agent_id: list
        list of prototype agentids as strings
    """
    ids = cur.execute('SELECT agentid FROM agententry '
                      'WHERE prototype = "' +
                      str(prototype) + '" COLLATE NOCASE').fetchall()

    return list(str(agent['agentid']) for agent in ids)


def institutions(cur):
    """Returns prototype and agentids of institutions

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    sqlite query result (list of tuples)
    """
    return cur.execute('SELECT prototype, agentid FROM agententry '
                       'WHERE kind = "Inst"').fetchall()


def timestep_to_years(init_year, timestep):
    """Returns list of years in simulation

    Parameters
    ----------
    init_year: int
        initial year in simulation
    timestep: np.array
        timestep of simulation (months)

    Returns
    -------
    array of years
    """

    return init_year + (timestep / 12)


def exec_string(specific_search, search, request_colmn):
    """Generates sqlite query command to select things and
        inner join resources and transactions.

    Parameters
    ----------
    specific_search: list
        list of items to specify search
        This variable will be inserted as sqlite
        query arugment following the search keyword
    search: str
        criteria for specific_search search
        This variable will be inserted as sqlite
        query arugment following the WHERE keyword
    request_colmn: str
        column (set of values) that the sqlite query should return
        This variable will be inserted as sqlite
        query arugment following the SELECT keyword

    Returns
    -------
    str
        sqlite query command.
    """
    if len(specific_search) == 0:
        raise Exception('Cannot create an exec_string with an empty list')
    if isinstance(specific_search[0], str):
        specific_search = ['"' + x + '"' for x in specific_search]

    query = ("SELECT " + request_colmn +
             " FROM resources INNER JOIN transactions"
             " ON transactions.resourceid = resources.resourceid"
             " WHERE (" + str(search) + ' = ' + str(specific_search[0])
             )
    for item in specific_search[1:]:
        query += ' OR ' + str(search) + ' = ' + str(item)
    query += ')'

    return query


def simulation_timesteps(cur):
    """Returns simulation start year, month,
    duration and timesteps (in numpy linspace).

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    init_year: int
        start year of simulation
    init_month: int
        start month of simulation
    duration: int
        duration of simulation
    timestep: list
        linspace up to duration: ar array with
        [start of the sequence, end of the sequence]
    """
    info = cur.execute('SELECT initialyear, initialmonth, '
                       'duration FROM info').fetchone()
    init_year = info['initialyear']
    init_month = info['initialmonth']
    duration = info['duration']
    timestep = np.linspace(0, duration - 1, num=duration)

    return init_year, init_month, duration, timestep


def timeseries(specific_search, duration, kg_to_tons):
    """returns a timeseries list from specific_search data.

    Parameters
    ----------
    specific_search: list
        list of data to be created into timeseries
        list[0] = time
        list[1] = value, quantity
    duration: int
        duration of the simulation
    kg_to_tons: bool
        if True, list returned has units of tons
        if False, list returned as units of kilograms

    Returns
    -------
    timeseries list of commodities stored in specific_search
    """
    value = 0
    value_timeseries = []
    array = np.array(specific_search)
    if len(specific_search) > 0:
        for i in range(0, duration):
            value = sum(array[array[:, 0] == i][:, 1])
            if kg_to_tons:
                value_timeseries.append(value * 0.001)
            else:
                value_timeseries.append(value)
    return value_timeseries


def timeseries_cum(specific_search, duration, kg_to_tons):
    """returns a timeseries list from specific_search data.

    Parameters
    ----------
    specific_search: list
        list of data to be created into timeseries
        list[0] = time
        list[1] = value, quantity
    multiplyby: int
        integer to multiply the value in the list by for
        unit conversion from kilograms
    kg_to_tons: bool
        if True, list returned has units of tons
        if False, list returned as units of kilograms

    Returns
    -------
    timeseries of commodities in kg or tons
    """
    value = 0
    value_timeseries = []
    array = np.array(specific_search)
    if len(specific_search) > 0:
        for i in range(0, duration):
            value += sum(array[array[:, 0] == i][:, 1])
            if kg_to_tons:
                value_timeseries.append(value * 0.001)
            else:
                value_timeseries.append(value)
    return value_timeseries


def isotope_transactions(resources, compositions):
    """Creates a dictionary with isotope name, mass, and time

    Parameters
    ----------
    resources: list of tuples
        resource data from the resources table
        (times, sum(quantity), qualid)
    compositions: list of tuples
        composition data from the compositions table
        (qualid, nucid, massfrac)

    Returns
    -------
    transactions: dictionary
        dictionary with "key=isotope, and
        value=list of tuples (time, mass_moved)"
    """
    transactions = collections.defaultdict(list)
    for res in resources:
        for comp in compositions:
            if res['qualid'] == comp['qualid']:
                transactions[comp['nucid']].append((res['time'],
                                                    res['sum(quantity)'] *
                                                    comp['massfrac']))

    return transactions


def facility_commodity_flux(cur, agentids,
                            facility_commodities, is_outflux,
                            is_cum=True):
    """Returns dictionary of commodity in/outflux from agents

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    agentids: list
        list of agentids
    facility_commodities: list
        list of commodities
    is_outflux: bool
        gets outflux if True, influx if False
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    commodity_region: dictionary
        dictionary with "key=commodity, and
        value=timeseries list of masses in kg"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    commodity_region = collections.OrderedDict()
    for comm in facility_commodities:
        query = (exec_string(agentids, 'receiverid',
                             'time, sum(quantity), qualid') +
                 ' and (commodity = "' + str(comm) +
                 '") GROUP BY time')
        # outflux changes receiverid to senderid
        if is_outflux:
            query = query.replace('receiverid', 'senderid')

        res = cur.execute(query).fetchall()
        if is_cum:
            commodity_region[comm] = timeseries_cum(res, duration, True)
        else:
            commodity_region[comm] = timeseries(res, duration, True)

    return commodity_region


def commodity_flux_region(cur, agentids, commodities,
                          is_outflux, is_cum=True):
    """Returns dictionary of timeseries of all the commodity outflux,
        that is either coming in/out of the agent
        separated by region

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    agentids: list
        list of agentids
    commodities: list
        list of commodities to include
    is_outflux: bool
        gets outflux from agent if True
        gets influx to agent if False
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    commodity_region: dictionary
        dictionary with "key=region, and
        value= timeseries list of masses in kg"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    commodity_region = collections.OrderedDict()
    commodities = ['"' + x + '"' for x in commodities]
    query = ('SELECT time, sum(quantity), parentid '
             'FROM transactions '
             'INNER JOIN resources '
             'ON resources.resourceid = '
             'transactions.resourceid '
             'INNER JOIN agententry '
             'ON agententry.agentid = transactions.SENDERID '
             'WHERE (commodity = ' +
             ' OR commodity = '.join(commodities) + ') AND ('
             'receiverid = ' +
             ' OR receiverid = '.join(agentids) + ') GROUP BY '
             'time, parentid')
    if is_outflux:
        query = query.replace('receiverid', 'senderid')
        query = query.replace('SENDERID', 'RECEIVERID')
    resources = cur.execute(query).fetchall()
    govs = cur.execute('SELECT agentid, prototype FROM agententry '
                       'WHERE kind = "Inst"').fetchall()
    for gov in govs:
        from_gov = [(x['time'], x['sum(quantity)'])
                    for x in resources if x['parentid'] == gov['agentid']]
        if is_cum:
            commodity_region[gov['prototype']] = timeseries_cum(
                from_gov, duration, True)
        else:
            commodity_region[gov['prototype']] = timeseries(
                from_gov, duration, True)
    return commodity_region


def facility_commodity_flux_isotopics(
        cur,
        agentids,
        facility_commodities,
        is_outflux,
        is_cum=True):
    """Returns timeseries isotoptics of commodity in/outflux
    from agents

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    agentids: list
        list of agentids
    facility_commodities: list
        list of commodities
    is_outflux: bool
        gets outflux if True, influx if False
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    isotope_timeseries: dictionary
        dictionary with "key=isotope, and
        value=timeseries list of masses in kg"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    isotope_timeseries = collections.defaultdict(list)
    for comm in facility_commodities:
        query = ('SELECT time, sum(quantity)*massfrac, nucid '
                 'FROM transactions INNER JOIN resources '
                 'ON resources.resourceid = transactions.resourceid '
                 'LEFT OUTER JOIN compositions '
                 'ON compositions.qualid = resources.qualid '
                 'WHERE (receiverid = ' +
                 ' OR receiverid = '.join(agentids) +
                 ') AND (commodity = "' + str(comm) +
                 '") GROUP BY time, nucid')
        # outflux changes receiverid to senderid
        if is_outflux:
            query = query.replace('receiverid', 'senderid')

        res = cur.execute(query).fetchall()
        for time, amount, nucid in res:
            isotope_timeseries[nucname.name(nucid)].append((time, amount))
    for key in isotope_timeseries:
        if is_cum:
            isotope_timeseries[key] = timeseries_cum(
                isotope_timeseries[key], duration, True)
        else:
            isotope_timeseries[key] = timeseries(
                isotope_timeseries[key], duration, True)
    return isotope_timeseries


def stockpiles(cur, facility, is_cum=True):
    """gets inventory timeseries in a fuel facility

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        name of facility
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    pile: dictionary
        dictionary with "key=agent type, and
        value=timeseries list of stockpile"
    """
    pile = collections.OrderedDict()
    agentid = agent_ids(cur, facility)
    query = exec_string(agentid, 'agentid', 'timecreated, quantity, qualid')
    query = query.replace('transactions', 'agentstateinventories')
    stockpile = cur.execute(query).fetchall()
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    if is_cum:
        stock_timeseries = timeseries_cum(stockpile, duration, True)
    else:
        stock_timeseries = timeseries(stockpile, duration, True)
    pile[facility] = stock_timeseries

    return pile


def swu_timeseries(cur, is_cum=True):
    """returns dictionary of swu timeseries for each enrichment plant

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    swu: dictionary
        dictionary with "key=Enrichment (facility number), and
        value=swu timeseries list"
    """
    swu = collections.OrderedDict()
    agentid = agent_ids(cur, 'Enrichment')
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    for num in agentid:
        swu_data = cur.execute('SELECT time, value '
                               'FROM timeseriesenrichmentswu '
                               'WHERE agentid = ' + str(num)).fetchall()
        if is_cum:
            swu_timeseries = timeseries_cum(swu_data, duration, False)
        else:
            swu_timeseries = timeseries(swu_data, duration, False)

        swu['Enrichment_' + str(num)] = swu_timeseries

    return swu


def power_capacity(cur):
    """Gets dictionary of power capacity by calling capacity_calc

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    ------
    power: dictionary
        "dictionary with key=government, and
        value=timeseries list of installed capacity"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    insts = institutions(cur)

    # get power cap values
    entry_exit = cur.execute('SELECT max(value), timeseriespower.agentid, '
                             'parentid, entertime, entertime + lifetime'
                             ' FROM agententry '
                             'INNER JOIN timeseriespower '
                             'ON agententry.agentid = timeseriespower.agentid '
                             'GROUP BY timeseriespower.agentid').fetchall()

    return capacity_calc(insts, timestep, entry_exit)


def power_capacity_of_region(cur, region_name):
    """Gets dictionary of power capacity of a specific region

    Parameters
    ----------
    cur: sqlite cursor
    region_name: str
        name of region to serach for

    Returns
    -------
    power: dictionary
        "dictionary with key=government and
        value=timeseries list of installed capacity"
    """
    parentid = cur.exectue('SELECT agentid FROM agententry WHERE '
                           'Prototype LIKE "%' + region_name + '%" '
                           'AND Kind = "Inst"').fetchone()

    entry_exit = cur.execute('SELECT max(value), timeseriespower.agentid, '
                             'parentid, entrytime, entertime + lifetime'
                             ' FROM agententry '
                             'INNER JOIN timeseriespower '
                             'ON agententry.agentid = timeseriespower.agentid '
                             'GROUP BY timeseriespower.agentid '
                             'WHERE parentid = %i' % parentid[0]).fetchall()


def deployments(cur):
    """Gets dictionary of reactors deployed over time
    by calling reactor_deployments

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    ------
    deployment_government: dictionary
        "dictionary with key=government, and
        value=timeseries list of number of reactors"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    insts = institutions(cur)

    # get power cap values
    entry = cur.execute('SELECT max(value), timeseriespower.agentid, '
                        'parentid, entertime FROM agententry '
                        'INNER JOIN timeseriespower '
                        'ON agententry.agentid = timeseriespower.agentid '
                        'GROUP BY timeseriespower.agentid').fetchall()

    exit_step = cur.execute('SELECT max(value), timeseriespower.agentid, '
                            'parentid, exittime FROM agentexit '
                            'INNER JOIN timeseriespower '
                            'ON agentexit.agentid = timeseriespower.agentid'
                            ' INNER JOIN agententry '
                            'ON agentexit.agentid = agententry.agentid '
                            'GROUP BY timeseriespower.agentid').fetchall()
    return reactor_deployments(insts, timestep, entry, exit_step)


def fuel_usage_timeseries(cur, fuels, is_cum=True):
    """Calculates total fuel usage over time

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    fuels: list
        list of fuel commodity names (eg. uox, mox) as string
        to consider in fuel usage.
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    fuel_usage: dictionary
        dictionary with "key=fuel (from fuels),
        value=timeseries list of fuel amount [kg]"
    """
    fuel_usage = collections.OrderedDict()
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    for fuel in fuels:
        temporary_fuels = [fuel]
        fuel_quantity = cur.execute(exec_string(temporary_fuels, 'commodity',
                                                'time, sum(quantity)') +
                                    ' GROUP BY time').fetchall()
        quantity_timeseries = []
        try:
            if is_cum:
                quantity_timeseries = timeseries_cum(
                    fuel_quantity, duration, True)
            else:
                quantity_timeseries = timeseries(
                    fuel_quantity, duration, True)
            fuel_usage[fuel] = quantity_timeseries
        except KeyError:
            print(str(fuel) + ' has not been used.')

    return fuel_usage


def nat_u_timeseries(cur, is_cum=True):
    """Finds natural uranium supply from source
        Since currently the source supplies all its capacity,
        the timeseriesenrichmentfeed is used.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    timeseries: function
        calls a function that returns timeseries list of natural U
        demand from enrichment [MTHM]
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    # Get Nat U feed to enrichment from timeseriesenrichmentfeed
    feed = cur.execute('SELECT time, sum(value) '
                       'FROM timeseriesenrichmentfeed '
                       'GROUP BY time').fetchall()
    if is_cum:
        return timeseries_cum(feed, duration, True)
    else:
        return timeseries(feed, duration, True)


def trade_timeseries(cur, sender, receiver,
                     is_prototype, do_isotopic,
                     is_cum=True):
    """Returns trade timeseries between two prototypes' or facilities
    with or without isotopics

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    sender: str
        name of sender as facility type or prototype name
    receiver: str
        name of receiver as facility type or prototype name
    is_prototype: bool
        if True, search sender and receiver as prototype,
        if False, as facility type from spec.
    do_isotopic: bool
        if True, perform isotopics (takes significantly longer)
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    trades: dictionary
        if do_isotopic:
            dictionary with "key=isotope, and
                        value=timeseries list
                        of mass traded between
                        two prototypes"
        else:
            dictionary with "key=string, sender to receiver,
                        value=timeseries list of mass traded
                        between two prototypes"

    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    isotope_timeseries = collections.defaultdict(list)
    trades = collections.defaultdict()

    if is_prototype:
        sender_id = prototype_id(cur, sender)
        receiver_id = prototype_id(cur, receiver)
    else:
        sender_id = agent_ids(cur, sender)
        receiver_id = agent_ids(cur, receiver)

    if do_isotopic:
        trade = cur.execute('SELECT time, sum(quantity)*massfrac, nucid '
                            'FROM transactions INNER JOIN resources ON '
                            'resources.resourceid = transactions.resourceid '
                            'LEFT OUTER JOIN compositions '
                            'ON compositions.qualid = resources.qualid '
                            'WHERE (senderid = ' +
                            ' OR senderid = '.join(sender_id) +
                            ') AND (receiverid = ' +
                            ' OR receiverid = '.join(receiver_id) +
                            ') GROUP BY time, nucid').fetchall()
    else:
        trade = cur.execute('SELECT time, sum(quantity), qualid '
                            'FROM transactions INNER JOIN resources ON '
                            'resources.resourceid = transactions.resourceid'
                            ' WHERE (senderid = ' +
                            ' OR senderid = '.join(sender_id) +
                            ') AND (receiverid = ' +
                            ' OR receiverid = '.join(receiver_id) +
                            ') GROUP BY time').fetchall(
        )
    if do_isotopic:
        for time, amount, nucid in trade:
            isotope_timeseries[nucname.name(nucid)].append((time, amount))
        for key in isotope_timeseries:
            if is_cum:
                isotope_timeseries[key] = timeseries_cum(
                    isotope_timeseries[key], duration, True)
            else:
                isotope_timeseries[key] = timeseries(
                    isotope_timeseries[key], duration, True)
        return isotope_timeseries
    else:
        key_name = str(sender)[:5] + ' to ' + str(receiver)[:5]
        if is_cum:
            trades[key_name] = timeseries_cum(trade, duration, True)
        else:
            trades[key_name] = timeseries(trade, duration, True)
        return trades


def final_stockpile(cur, facility):
    """get final stockpile in a fuel facility

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        name of facility

    Returns
    -------
    mthm_stockpile: str
        MTHM value of stockpile
    """
    agentid = agent_ids(cur, facility)
    mthm_stockpile = ''
    for agent in agentid:
        count = 1
        name = cur.execute('SELECT prototype FROM agententry'
                           'WHERE agentid = ' + str(agent)).fetchone()

        mthm_stockpile += 'The Stockpile in ' + str(name[0]) + ' : \n \n'
        stkpile = cur.execute('SELECT sum(quantity), inventoryname, qualid'
                              ' FROM agentstateinventories'
                              ' INNER JOIN resources'
                              ' ON resources.resourceid'
                              ' = agentstateinventories.resourceid'
                              ' WHERE agentstateinventories.agentid'
                              ' = """ + str(agent) + """ GROUP BY'
                              ' inventoryname').fetchall()
        for stream in stkpile:
            masses = cur.execute('SELECT qualid, nucid, massfrac '
                                 'FROM compositions '
                                 'WHERE qualid = ' +
                                 str(stream['qualid'])).fetchall()

            mthm_stockpile += ('Stream ' + str(count) +
                               ' Total = ' + str(stream['sum(quantity)']) +
                               ' kg \n')
            for isotope in masses:
                mthm_stockpile += (str(isotope['nucid']) + ' = ' +
                                   str(isotope['massfrac'] *
                                       stream['sum(quantity)']) +
                                   ' kg \n')
            mthm_stockpile += '\n'
            count += 1
        mthm_stockpile += '\n'
    mthm_stockpile += '\n'

    return mthm_stockpile


def fuel_into_reactors(cur, is_cum=True):
    """Finds timeseries of mass of fuel received by reactors

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    is_cum: bool
        gets cumulative timeseris if True, monthly value if False

    Returns
    -------
    timeseries list of fuel into reactors [tons]
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    fuel = cur.execute('SELECT time, sum(quantity) FROM transactions '
                       'INNER JOIN resources ON '
                       'resources.resourceid = transactions.resourceid '
                       'INNER JOIN agententry ON '
                       'transactions.receiverid = agententry.agentid '
                       'WHERE spec LIKE "%Reactor%" '
                       'GROUP BY time').fetchall()

    if is_cum:
        return timeseries_cum(fuel, duration, True)
    else:
        return timeseries(fuel, duration, True)


def u_util_calc(cur):
    """Returns fuel utilization factor of fuel cycle

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    u_util_timeseries: numpy array
        Timeseries of Uranium utilization factor
    Prints simulation average Uranium Utilization
    """
    # timeseries of natural uranium
    u_supply_timeseries = np.array(nat_u_timeseries(cur))

    # timeseries of fuel into reactors
    fuel_usage = np.array(fuel_into_reactors(cur))

    # timeseries of Uranium utilization
    u_util_timeseries = np.nan_to_num(fuel_usage / u_supply_timeseries)
    print('The Average Fuel Utilization Factor is: ')
    print(sum(u_util_timeseries) / len(u_util_timeseries))

    return u_util_timeseries


def plot_uranium_utilization(cur):
    """Plots uranium utilization factor of fuel cycle

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    none
    """

    u_util_timeseries = u_util_calc(cur)
    plt.plot(u_util_timeseries, label='Uranium utilization')
    plt.xlabel('time [months]')
    plt.ylabel('Uranium Utilization')
    plt.legend()
    plt.show()


def commodity_origin(cur, commodity, prototypes, is_cum=True):
    """Returns dict of where a commodity is from

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    commodity: str
        name of commodity
    prototypes: list
        list of prototypes that provide the commodity

    Returns
    -------
    prototype_trades: dictionary
        "dictionary with key=prototype name, and
        value=timeseries list of commodity sent from prototypes"
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    query = ('SELECT time, sum(quantity) FROM transactions '
             'INNER JOIN resources ON resources.resourceid = '
             'transactions.resourceid WHERE commodity = "' +
             str(commodity) + '" AND senderid '
             '= 9999 GROUP BY time')
    prototype_trades = collections.OrderedDict()
    for agent in prototypes:
        agent_id = prototype_id(cur, agent)
        from_agent = cur.execute(query.replace(
            '9999', ' OR senderid = '.join(agent_id))).fetchall()
        if is_cum:
            prototype_trades[agent] = timeseries_cum(
                from_agent, duration, True)
        else:
            prototype_trades[agent] = timeseries(from_agent, duration, True)
    return prototype_trades


def commodity_per_institution(cur, commodity, timestep=10000):
    """Outputs outflux of commodity per institution
        before timestep

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    commodity: str
        commodity to search for

    Returns
    -------
    institution_output: dictionary
        key = institution
        value = timeseries list of outflux of commodity
    """

    insts = institutions(cur)
    institution_output = collections.OrderedDict()
    for inst in insts:
        inst_id = inst[1]
        inst_name = inst[0]
        facilities = cur.execute('SELECT agentid FROM agententry '
                                 'WHERE parentid = ' + str(inst_id)).fetchall()
        facilities_collected = []
        for fac in facilities:
            facilities_collected.append(fac[0])
        query = exec_string(facilities_collected, 'senderid', 'sum(quantity)')
        query += ' AND commodity = "' + commodity + \
            '" and time < ' + str(timestep)
        institution_output[inst_name] = cur.execute(query).fetchone()[0]

    return institution_output


def waste_mass_series(isotopes, mass_timeseries, duration):
    """Given an isotope, mass and time list, creates a dictionary
       With key as isotope and time series of the isotope mass.

    Parameters
    ----------
    isotopes: list
        list with all the isotopes from resources table
    mass_timeseries: list
        a list of lists.  each outer list corresponds to a different isotope
        and contains tuples in the form (time,mass)
        for the isotope transaction.
    duration: integer
        simulation duration

    Returns
    -------
    waste_mass: dictionary
        dictionary with "key=isotope, and
        value=mass timeseries of each unique isotope"   """
    waste_mass = {}
    for isotope in isotopes:
        postion = [i for i, x in enumerate(isotopes) if x == isotope][0]
        mass = [item[1] for item in mass_timeseries[postion]]
        waste_mass[isotope] = mass
    return waste_mass


def waste_timeseries(isotopes, mass_timeseries, duration):
    """Given an isotope, mass and time list, creates a dictionary
       With key as isotope and time series of the isotope mass.

    Parameters
    ----------
    isotopes: list
        list with all the isotopes from resources table
    mass_timeseries: list
        a list of lists.  each outer list corresponds to a different isotope
        and contains tuples in the form
        (time,mass) for the isotope transaction.
    duration: integer
        simulation duration

    Returns
    -------
    waste_time: dictionary
        dictionary with "key=isotope, and
        value=mass timeseries of each unique isotope"   """
    waste_time = {}
    for isotope in isotopes:
        postion = [i for i, x in enumerate(isotopes) if x == isotope][0]
        time = [item[0] for item in mass_timeseries[postion]]
        waste_time[isotope] = time
    return waste_time


def capacity_calc(insts, timestep, entry_exit):
    """Adds and subtracts capacity over time for plotting

    Parameters
    ----------
    insts: list
        list of insts (countries)
    timestep: np.linspace
        list of timestep from 0 to simulation time
    entry_exit: list
        power_cap, agentid, parentid, entertime, exittime
        of all entered reactors

    Returns
    -------
    power: dictionary
        "dictionary with key=government, and
        value=timeseries list capacity"
    """
    power = collections.OrderedDict()
    for inst in insts:
        capacity = []
        cap = 0
        for t in timestep:
            for agent in entry_exit:
                if (agent['entertime'] == t and
                        agent['parentid'] == inst['agentid']):
                    cap += agent['max(value)'] * 0.001
                if (agent['entertime + lifetime'] == t and
                        agent['parentid'] == inst['agentid']):
                    cap -= agent['max(value)'] * 0.001
            capacity.append(cap)
        power[inst['prototype']] = np.asarray(capacity)

    return power


def reactor_deployments(insts, timestep, entry, exit_step):
    """Adds and subtracts number of reactors deployed over time
    for plotting

    Parameters
    ----------
    insts: list
        list of insts (countries)
    timestep: np.linspace
        list of timestep from 0 to simulation time
    entry: list
        power_cap, agentid, parentid, entertime
        of all entered reactors

    exit_step: list
        power_cap, agentid, parenitd, exittime
        of all decommissioned reactors

    Returns
    -------
    deployment: dictionary
        "dictionary with key=government, and
        value=timeseries number of reactors"
    """
    deployment = collections.OrderedDict()
    for inst in insts:
        num_reactors = []
        count = 0
        for t in timestep:
            for enter in entry:
                if (enter['entertime'] == t and
                        enter['parentid'] == inst['agentid']):
                    count += 1
            for dec in exit_step:
                if (dec['exittime'] == t and
                        dec['parentid'] == inst['agentid']):
                    count -= 1
            num_reactors.append(count)
        deployment[inst['prototype']] = np.asarray(num_reactors)

    return deployment


def multiple_line_plots(dictionary, timestep,
                        xlabel, ylabel, title,
                        outputname, init_year):
    """Creates multiple line plots of timestep vs dictionary

    Parameters
    ----------
    dictionary: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    timestep: numpy linspace
        timestep of simulation
    xlabel: str
        xlabel of plot
    ylabel: str
        ylabel of plot
    title: str
        title of plot
    init_year: int
        initial year of simulation

    Returns
    -------
    plot: plot
        multiple line plots
    """
    # set different colors for each bar
    color_index = 0
    # for every country, create bar chart with different color
    for key in dictionary:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)

        plt.plot(timestep_to_years(init_year, timestep),
                 dictionary[key],
                 label=label)
        color_index += 1
        if sum(sum(dictionary[k]) for k in dictionary) > 1000:
            ax = plt.gca()
            ax.get_yaxis().set_major_formatter(
                plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))
        plt.ylabel(ylabel)
        plt.title(title)
        plt.xlabel(xlabel)
        plt.legend(loc=(1.0, 0), prop={'size': 10})
        plt.grid(True)
        plt.savefig(label + '_' + outputname + '.png',
                    format='png',
                    bbox_inches='tight')
        plt.close()


def combined_line_plot(dictionary, timestep,
                       xlabel, ylabel, title,
                       outputname, init_year,
                       colormap=cm.viridis):
    """Creates a combined line plot of timestep vs dictionary

    Parameters
    ----------
    dictionary: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    timestep: numpy linspace
        timestep of simulation
    xlabel: str
        xlabel of plot
    ylabel: str
        ylabel of plot
    title: str
        title of plot
    init_year: int
        initial year of simulation

    Returns
    -------
    plot: plot
        produces a combined line plot
    """
    # set different colors for each bar
    color_index = 0
    plt.figure()
    # for every country, create bar chart with different color
    for key in dictionary:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)

        plt.plot(timestep_to_years(init_year, timestep),
                 dictionary[key],
                 label=label,
                 color=colormap(float(color_index) / len(dictionary)))
        color_index += 1

    if sum(sum(dictionary[k]) for k in dictionary) > 1000:
        ax = plt.gca()
        ax.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))
    plt.ylabel(ylabel)
    plt.title(title)
    plt.xlabel(xlabel)
    plt.legend(loc=(1.0, 0), prop={'size': 10})
    plt.grid(True)
    plt.savefig(label + '_' + outputname + '.png',
                format='png',
                bbox_inches='tight')
    plt.close()


def double_axis_bar_line_plot(dictionary1, dictionary2, timestep,
                              xlabel, ylabel1, ylabel2,
                              title, outputname, init_year):
    """Creates a double-axis plot of timestep vs dictionary

    It is recommended that a non-cumulative timeseries is on dictionary1.

    Parameters
    ----------
    dictionary1: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    dictionary2: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    timestep: numpy linspace
        timestep of simulation
    xlabel: str
        xlabel of plot
    ylabel: str
        ylabel of plot
    title: str
        title of plot
    init_year: int
        initial year of simulation

    Returns
    -------
    plot:plot
        double-axis bar-line plot
    """
    # set different colors for each bar

    fig, ax1 = plt.subplots()
    # for every country, create bar chart with different color
    color1 = 'r'
    color2 = 'b'
    for key in dictionary1:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)
        if sum(dictionary1[key]) == 0:
            print(label + ' has no values')
        else:
            ax1.bar(timestep_to_years(init_year, timestep),
                    dictionary1[key],
                    label=label,
                    color=color1)
    ax1.set_xlabel(xlabel)
    ax1.set_ylabel(ylabel1, color=color1)
    ax1.tick_params('y', colors=color1)
    if sum(sum(dictionary1[k]) for k in dictionary1) > 1000:
        ax1 = plt.gca()
        ax1.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))
    ax2 = ax1.twinx()

    lines = ['-', '--', '-.', ':']
    linecycler = cycle(lines)
    for key in dictionary2:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)

        if sum(dictionary2[key]) == 0:
            print(label + ' has no values')
        else:
            ax2.plot(timestep_to_years(init_year, timestep),
                     dictionary2[key],
                     label=label,
                     color=color2,
                     linestyle=next(linecycler))
    ax2.set_ylabel(ylabel2, color=color2)
    ax2.tick_params('y', colors=color2)

    if sum(sum(dictionary2[k]) for k in dictionary2) > 1000:
        ax2 = plt.gca()
        ax2.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))

    plt.title(title)
    plt.grid(True)
    plt.savefig(label + '_' + outputname + '.png',
                format='png',
                bbox_inches='tight')
    plt.close()


def double_axis_line_line_plot(dictionary1, dictionary2, timestep,
                               xlabel, ylabel1, ylabel2,
                               title, outputname, init_year):
    """Creates a double-axis plot of timestep vs dictionary

    Parameters
    ----------
    dictionary1: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    dictionary2: dictionary
        dictionary with "key=description of timestep, and
        value=list of timestep progressions"
    timestep: numpy linspace
        timestep of simulation
    xlabel: str
        xlabel of plot
    ylabel: str
        ylabel of plot
    title: str
        title of plot
    init_year: int
        initial year of simulation

    Returns
    -------
    plot: plot
        double-axis plot
    """
    # set different colors for each bar
    lines = ['-', '--', '-.', ':']
    linecycler = cycle(lines)
    fig, ax1 = plt.subplots()
    top = True
    color1 = 'r'
    color2 = 'b'
    # for every country, create bar chart with different color
    for key in dictionary1:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)
        if top:
            lns = ax1.plot(timestep_to_years(init_year, timestep),
                           dictionary1[key],
                           label=label,
                           color=color1,
                           linestyle=next(linecycler))
            top = False
        else:
            lns += ax1.plot(timestep_to_years(init_year, timestep),
                            dictionary1[key],
                            label=label,
                            color=color1,
                            linestyle=next(linecycler))
    ax1.set_xlabel(xlabel)
    ax1.set_ylabel(ylabel1, color=color1)
    ax1.tick_params('y', colors=color1)
    if sum(sum(dictionary1[k]) for k in dictionary1) > 1000:
        ax1 = plt.gca()
        ax1.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))
    ax2 = ax1.twinx()

    linecycler = cycle(lines)

    for key in dictionary2:
        # label is the name of the nuclide (converted from ZZAAA0000 format)
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)

        lns += ax2.plot(timestep_to_years(init_year, timestep),
                        dictionary2[key],
                        label=label,
                        color=color2,
                        linestyle=next(linecycler))
    ax2.set_ylabel(ylabel2, color=color2)
    ax2.tick_params('y', colors=color2)

    if sum(sum(dictionary2[k]) for k in dictionary2) > 1000:
        ax2 = plt.gca()
        ax2.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))

    plt.title(title)
    labs = [l.get_label() for l in lns]
    plt.legend(lns, labs, loc=0, prop={'size': 10})
    plt.grid(True)
    plt.savefig(label + '_' + outputname + '.png',
                format='png',
                bbox_inches='tight')
    plt.close()


def stacked_bar_chart(dictionary, timestep,
                      xlabel, ylabel, title,
                      outputname, init_year,
                      colormap=cm.viridis):
    """Creates stacked bar chart of timstep vs dictionary

    Parameters
    ----------
    dictionary: dictionary
        dictionary with value: timeseries data
    timestep: numpy linspace
        list of timestep (x axis)
    xlabel: str
        xlabel of plot
    ylabel: str
        ylabel of plot
    title: str
        title of plot
    init_year: int
        simulation start year

    Returns
    -------
    plot : plot
        plot of stacked bar chart of timstep vs dictionary
    """
    # set different colors for each bar
    color_index = 0
    top_index = True
    prev = np.zeros(1)
    plots = []
    # for every country, create bar chart with different color
    for key in dictionary:
        if isinstance(key, str) is True:
            label = key.replace('_government', '')
        else:
            label = str(key)
        # very first country does not have a 'bottom' argument
        if sum(dictionary[key]) == 0:
            print(label + ' has no values')
        elif top_index is True:
            plot = plt.bar(x=timestep_to_years(init_year, timestep),
                           height=dictionary[key],
                           width=0.5,
                           color=colormap(
                float(color_index) / len(dictionary)),
                edgecolor='none',
                label=label)
            prev = dictionary[key]
            top_index = False
            plots.append(plot)

        # All curves except the first have a 'bottom'
        # defined by the previous curve
        else:
            plot = plt.bar(x=timestep_to_years(init_year, timestep),
                           height=dictionary[key],
                           width=0.5,
                           color=colormap(
                float(color_index) / len(dictionary)),
                edgecolor='none',
                bottom=prev,
                label=label)
            prev = np.add(prev, dictionary[key])
            plots.append(plot)

        color_index += 1

    # plot
    if sum(sum(dictionary[k]) for k in dictionary) > 1000:
        ax = plt.gca()
        ax.get_yaxis().set_major_formatter(
            plt.FuncFormatter(lambda x, loc: "{:,}".format(int(x))))
    plt.ylabel(ylabel)
    plt.title(title)
    plt.xlabel(xlabel)
    axes = plt.gca()
    handles, labels = ax.get_legend_handles_labels()
    if len(dictionary) > 1:
        plt.legend(handles[::-1], labels[::-1], loc=(1.0, 0), )
    plt.grid(True)
    plt.savefig(outputname + '.png', format='png', bbox_inches='tight')
    plt.close()


def plot_power(cur):
    """Gets capacity vs time for every country
        in stacked bar chart.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    """
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    power = power_capacity(cur)
    stacked_bar_chart(power, timestep,
                      'Years', 'Net_Capacity [GWe]',
                      'Net Capacity vs Time',
                      'power_plot', init_year)

    deploys = deployments(cur)
    stacked_bar_chart(deploys, timestep,
                      'Years', 'Number of Reactors',
                      'Number of Reactors vs Time',
                      'num_plot', init_year)


def plot_in_out_flux(cur, facility, influx_bool,
                     title, is_cum=False, is_tot=False):
    """Plots timeseries influx/ outflux from facility name in kg.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        facility name
    influx_bool: bool
        if true, calculates influx,
        if false, calculates outflux
    title: str
        title of the multi line plot
    outputname: str
        filename of the multi line plot file
    is_cum: Boolean:
        true: add isotope masses over time
        false: do not add isotope masses at each timestep

    Returns
    -------
    plot : plot
        plot of in flux of isotopes from facility
    """

    agentids = prototype_id(cur, facility)

    if influx_bool is True:
        resources = cur.execute(exec_string(agentids,
                                            'transactions.receiverId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()
    else:
        resources = cur.execute(exec_string(agentids,
                                            'transactions.senderId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()

    compositions = cur.execute('SELECT qualid, nucid, massfrac '
                               'FROM compositions').fetchall()

    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    transactions = isotope_transactions(resources, compositions)

    time_mass = []
    time_waste = {}
    for key in transactions.keys():

        time_mass.append(transactions[key])
        time_waste[key] = transactions[key]

    waste_mass = waste_mass_series(transactions.keys(),
                                   time_mass,
                                   duration)

    if not is_cum and not is_tot:
        keys = []
        for key in waste_mass.keys():
            keys.append(key)

        for element in range(len(keys)):
            time_and_mass = np.array(time_waste[keys[element]])
            time = [item[0] for item in time_and_mass]
            mass = [item[1] for item in time_and_mass]
            plt.plot(time, mass, linestyle=' ', marker='.',
                     markersize=1, label=nucname.name(keys[0]))

        plt.legend(loc='upper left')
        plt.title(title)
        plt.xlabel('time [months]')
        plt.ylabel('mass [kg]')
        plt.xlim(left=0.0)
        plt.ylim(bottom=0.0)
        plt.show()

    elif is_cum and not is_tot:
        value = 0
        keys = []
        for key in waste_mass.keys():
            keys.append(key)

        for element in range(len(waste_mass.keys())):
            placeholder = []
            value = 0
            key = keys[element]

            for index in range(len(waste_mass[key])):
                value += waste_mass[key][index]
                placeholder.append(value)
            waste_mass[key] = placeholder

        times = []
        nuclides = []
        masstime = {}
        for element in range(len(keys)):
            time_and_mass = np.array(time_waste[keys[element]])
            time = [item[0] for item in time_and_mass]
            mass = [item[1] for item in time_and_mass]
            nuclide = nucname.name(keys[element])
            mass_cum = np.cumsum(mass)
            times.append(time)
            nuclides.append(str(nuclide))
            masstime[nucname.name(keys[element])] = mass_cum
        mass_sort = sorted(
            masstime.items(), key=lambda e: e[1][-1], reverse=True)
        nuclides = [item[0] for item in mass_sort]
        masses = [item[1] for item in mass_sort]
        plt.stackplot(times[0], masses, labels=nuclides)
        plt.legend(loc='upper left')
        plt.title(title)
        plt.xlabel('time [months]')
        plt.ylabel('mass [kg]')
        plt.xlim(left=0.0)
        plt.ylim(bottom=0.0)
        plt.show()

    elif not is_cum and is_tot:
        keys = []
        for key in waste_mass.keys():
            keys.append(key)

        total_mass = np.zeros(len(waste_mass[keys[0]]))
        for element in range(len(keys)):
            for index in range(len(waste_mass[keys[0]])):
                total_mass[index] += waste_mass[keys[element]][index]

        total_mass[total_mass == 0] = np.nan
        plt.plot(total_mass, linestyle=' ', marker='.', markersize=1)
        plt.title(title)
        plt.xlabel('time [months]')
        plt.ylabel('mass [kg]')
        plt.xlim(left=0.0)
        plt.ylim(bottom=0.0)
        plt.show()

    elif is_cum and is_tot:
        value = 0
        keys = []
        for key in waste_mass.keys():
            keys.append(key)

        times = []
        nuclides = []
        masstime = {}
        for element in range(len(keys)):
            time_and_mass = np.array(time_waste[keys[element]])
            time = [item[0] for item in time_and_mass]
            mass = [item[1] for item in time_and_mass]
            nuclide = nucname.name(keys[element])
            mass_cum = np.cumsum(mass)
            times.append(time)
            nuclides.append(str(nuclide))
            masstime[nucname.name(keys[element])] = mass_cum
        mass_sort = sorted(
            masstime.items(), key=lambda e: e[1][-1], reverse=True)
        nuclides = [item[0] for item in mass_sort]
        masses = [item[1] for item in mass_sort]
        plt.stackplot(times[0], masses, labels=nuclides)
        plt.legend(loc='upper left')
        plt.title(title)
        plt.xlabel('time [months]')
        plt.ylabel('mass [kg]')
        plt.xlim(left=0.0)
        plt.ylim(bottom=0.0)
        plt.show()


def entered_power(cur):
    """Returns dictionary of power entered into simulation.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor

    Returns
    -------
    power: dictionary
        key: 'power'
        value: timeseries of power entered (non-cumulative)
    """
    power = {}
    entered = cur.execute('SELECT entertime, max(value) FROM '
                          'agententry INNER JOIN timeseriespower '
                          'ON agententry.agentid = timeseriespower.agentid '
                          'WHERE spec LIKE "%reactor%" '
                          'GROUP BY agententry.agentid').fetchall()
    init_year, init_month, duration, timestep = simulation_timesteps(cur)
    power['power'] = timeseries(entered, duration, False)
    return power


def source_throughput(cur, duration, frac_prod, frac_tail):
    """Calculates throughput required for nat_u source before enrichment
    by calculating the average mass of fuel gone into reactors over
    simulation. Assuming natural uranium is put as feed

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    duration: int
        duration of simulation
    frac_prod: float
        mass fraction of U235 in fuel after enrichment in decimals
    frac_tail: float
        mass fraction of U235 in tailings after enrichment in decimals

    Returns
    -------
    throughput: float
        appropriate nat_u throughput for source
    """
    avg_fuel_used = fuel_into_reactors(cur)[-1] * 1000 / duration
    feed_factor = (frac_prod - frac_tail) / (0.00711 - frac_tail)
    print('Throughput should be at least: ' +
          str(feed_factor * avg_fuel_used) + ' [kg]')
    return feed_factor * avg_fuel_used


def plot_in_flux_cumulative(
        cur,
        facility,
        title):
    """Plots timeseries influx/ outflux from facility name in kg.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        facility name
    influx_bool: bool
        if true, calculates influx,
        if false, calculates outflux
    title: str
        title of the multi line plot
    outputname: str
        filename of the multi line plot file
    is_cum: Boolean:
        true: add isotope masses over time
        false: do not add isotope masses at each timestep

    Returns
    -------
    """

    masstime = cumulative_mass_timeseries(cur, facility, flux='in')[0]
    times = cumulative_mass_timeseries(cur, facility, flux='in')[1]
    mass_sort = sorted(masstime.items(), key=lambda e: e[
                       1][-1], reverse=True)
    nuclides = [item[0] for item in mass_sort]
    masses = [item[1] for item in mass_sort]
    plt.stackplot(times[0], masses, labels=nuclides)
    plt.legend(loc='upper left')
    plt.title(title)
    plt.xlabel('time [months]')
    plt.ylabel('mass [kg]')
    plt.xlim(left=0.0)
    plt.ylim(bottom=0.0)
    plt.show()


def plot_out_flux_cumulative(
        cur,
        facility,
        title):
    """Plots timeseries influx/ outflux from facility name in kg.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        facility name
    influx_bool: bool
        if true, calculates influx,
        if false, calculates outflux
    title: str
        title of the multi line plot
    outputname: str
        filename of the multi line plot file
    is_cum: Boolean:
        true: add isotope masses over time
        false: do not add isotope masses at each timestep

    Returns
    -------
    """

    masses = cumulative_mass_timeseries(cur, facility, flux='out')
    masstime = masses[0]
    times = masses[1]

    mass_sort = sorted(masstime.items(), key=lambda e: e[
                       1][-1], reverse=True)
    nuclides = [item[0] for item in mass_sort]
    masses = [item[1] for item in mass_sort]
    plt.stackplot(times[0], masses, labels=nuclides)
    plt.legend(loc='upper left')
    plt.title(title)
    plt.xlabel('time [months]')
    plt.ylabel('mass [kg]')
    plt.xlim(left=0.0)
    plt.ylim(bottom=0.0)
    plt.show()


def plot_in_flux_basic(
        cur,
        facility,
        title):
    """Plots timeseries influx/ outflux from facility name in kg.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        facility name
    influx_bool: bool
        if true, calculates influx,
        if false, calculates outflux
    title: str
        title of the multi line plot
    outputname: str
        filename of the multi line plot file
    is_cum: Boolean:
        true: add isotope masses over time
        false: do not add isotope masses at each timestep

    Returns
    -------
    """
    masstime = mass_timeseries(cur, facility, flux='in')[0]
    times = mass_timeseries(cur, facility, flux='in')[1]
    nuclides = [item[0] for item in masstime]
    masses = [item[1] for item in masstime]
    mass_sort = sorted(masstime.items(), key=lambda e: e[
        1][-1], reverse=True)
    nuclides = [item[0] for item in mass_sort]
    masses = [item[1] for item in mass_sort]
    for i in range(len(times)):
        plt.plot(times[i], masses[i], label=nuclides[i])
    plt.legend(loc='upper left')
    plt.title(title)
    plt.xlabel('time [months]')
    plt.ylabel('mass [kg]')
    plt.xlim(left=0.0)
    plt.ylim(bottom=0.0)
    plt.show()


def plot_out_flux_basic(
        cur,
        facility,
        title):
    """Plots timeseries influx/ outflux from facility name in kg.

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facility: str
        facility name
    influx_bool: bool
        if true, calculates influx,
        if false, calculates outflux
    title: str
        title of the multi line plot
    outputname: str
        filename of the multi line plot file
    is_cum: Boolean:
        true: add isotope masses over time
        false: do not add isotope masses at each timestep

    Returns
    -------
    plot : plot
        plot of out flux of isotopes
    """
    masstime = mass_timeseries(cur, facility, flux='out')[0]
    times = mass_timeseries(cur, facility, flux='out')[1]
    nuclides = [item[0] for item in masstime]
    masses = [item[1] for item in masstime]
    mass_sort = sorted(masstime.items(), key=lambda e: e[
        1][-1], reverse=True)
    nuclides = [item[0] for item in mass_sort]
    masses = [item[1] for item in mass_sort]
    for i in range(len(times)):
        plt.plot(times[i], masses[i], label=nuclides[i])
    plt.legend(loc='upper left')
    plt.title(title)
    plt.xlabel('time [months]')
    plt.ylabel('mass [kg]')
    plt.xlim(left=0.0)
    plt.ylim(bottom=0.0)
    plt.show()


def plot_net_flux(
        cur,
        facility,
        title):
    """Plots net flux of all isotopes over the duration of the simulation.

    Parameters
    ----------
    cur : sqlite cursor
        sqlite cursor
    facility : str
        name of facility
    title : str
        title of plot

    Returns
    -------
    plot : plot
        plot of net flux of isotopes
    """
    masstime_in = mass_timeseries(cur, facility, flux='in')[0]
    times_in = mass_timeseries(cur, facility, flux='in')[1]
    masstime_out = mass_timeseries(cur, facility, flux='out')[0]
    times_out = mass_timeseries(cur, facility, flux='out')[1]
    mass_sort_in = sorted(masstime_in.items(), key=lambda e: e[
        1][-1], reverse=True)
    mass_sort_out = sorted(masstime_out.items(), key=lambda e: e[
        1][-1], reverse=True)
    nuclides_in = [item[0] for item in mass_sort_in]
    masses_in = [item[1] for item in mass_sort_in]
    nuclides_out = [item[0] for item in mass_sort_out]
    masses_out = np.negative([item[1] for item in mass_sort_out])
    plt.stackplot(times_in[0], masses_in, labels=nuclides_in)
    plt.stackplot(times_out[0], masses_out, labels=nuclides_out)
    plt.legend(loc='upper left')
    plt.title(title)
    plt.xlabel('time [months]')
    plt.ylabel('mass [kg]')
    plt.xlim(left=0.0)
    plt.show()


def mass_timeseries(cur, facility, flux):
    """Returns dictionary of mass timeseries of each isotope at a facility.

    Parameters
    ----------
    cur : sqlite cursor
        sqlite cursor
    facility : str
        name of facility
    flux : str
        direction of flux

    Returns
    -------
    masstime : dict
        dictionary of isotopes and their mass series
    times : list
        list of times in the simulation
    """
    agentids = prototype_id(cur, facility)

    if flux == 'in':
        resources = cur.execute(exec_string(agentids,
                                            'transactions.receiverId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()
    else:
        resources = cur.execute(exec_string(agentids,
                                            'transactions.senderId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()

    compositions = cur.execute('SELECT qualid, nucid, massfrac '
                               'FROM compositions').fetchall()

    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    transactions = isotope_transactions(resources, compositions)

    time_mass = []
    time_waste = {}
    for key in transactions.keys():

        time_mass.append(transactions[key])
        time_waste[key] = transactions[key]

    waste_mass = waste_mass_series(transactions.keys(),
                                   time_mass,
                                   duration)
    keys = []
    for key in waste_mass.keys():
        keys.append(key)

    times = []
    nuclides = []
    masstime = {}
    for element in range(len(keys)):
        time_and_mass = np.array(time_waste[keys[element]])
        time = [item[0] for item in time_and_mass]
        mass = [item[1] for item in time_and_mass]
        mass_per_time = Counter()
        for time, mass in time_and_mass:
            mass_per_time.update({time: mass})
        time = list(mass_per_time.keys())
        mass_val = list(mass_per_time.values())
        nuclide = nucname.name(keys[element])
        for j in np.arange(0, int(time[-1]+1)):
            if j not in time:
                time.insert(j, j)
                mass_val.insert(j, 0)
        times.append(time)
        nuclides.append(str(nuclide))
        masstime[nucname.name(keys[element])] = mass_val
    return masstime, times


def cumulative_mass_timeseries(cur, facility, flux):
    """Returns dictionary of the cumulative mass
       timeseries of each isotope at a facility.

    Parameters
    ----------
    cur : sqlite cursor
        sqlite cursor
    facility : str
        name of facility
    flux : str
        direction of flux

    Returns
    -------
    masstime : dict
        dictionary of isotopes and their mass series
    times : list
        list of times in the simulation
    """
    agentids = prototype_id(cur, facility)

    if flux == 'in':
        resources = cur.execute(exec_string(agentids,
                                            'transactions.receiverId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()
    else:
        resources = cur.execute(exec_string(agentids,
                                            'transactions.senderId',
                                            'time, sum(quantity), '
                                            'qualid') +
                                ' GROUP BY time, qualid').fetchall()

    compositions = cur.execute('SELECT qualid, nucid, massfrac '
                               'FROM compositions').fetchall()

    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    transactions = isotope_transactions(resources, compositions)

    time_mass = []
    time_waste = {}
    for key in transactions.keys():

        time_mass.append(transactions[key])
        time_waste[key] = transactions[key]

    waste_mass = waste_mass_series(transactions.keys(),
                                   time_mass,
                                   duration)
    keys = []
    for key in waste_mass.keys():
        keys.append(key)

    times = []
    nuclides = []
    masstime = {}
    for element in range(len(keys)):
        time_and_mass = np.array(time_waste[keys[element]])
        time = [item[0] for item in time_and_mass]
        mass = [item[1] for item in time_and_mass]
        mass_per_time = Counter()
        for time, mass in time_and_mass:
            mass_per_time.update({time: mass})
        time = list(mass_per_time.keys())
        mass_val = list(mass_per_time.values())
        nuclide = nucname.name(keys[element])
        for j in np.arange(0, int(time[-1]+1)):
            if j not in time:
                time.insert(j, j)
                mass_val.insert(j, 0)
        mass_cum = list(np.cumsum(mass_val))
        times.append(time)
        nuclides.append(str(nuclide))
        masstime[nucname.name(keys[element])] = mass_cum
    return masstime, times


def plot_cumulative_swu(cur, facilities=[]):
    """Returns dictionary of swu timeseries for each enrichment plant

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facilities : list
        list of facilities to plot

    Returns
    -------
    swu_dict: dictionary
        dictionary with "key=Enrichment (facility number), and
        value=swu timeseries list"
    """
    swu_dict = {}
    agentid = agent_ids(cur, 'Enrichment')
    if len(facilities) != 0:
        agentid = facilities
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    for num in agentid:
        swu_data = cur.execute('SELECT time, value '
                               'FROM timeseriesenrichmentswu '
                               'WHERE agentid = ' + str(num)).fetchall()
        swu_timeseries = timeseries_cum(swu_data, duration, False)
        swu_dict['Enrichment_' + str(num)] = swu_timeseries

    keys = []
    for key in swu_dict.keys():
        keys.append(key)
    swu = []
    facilities = []
    swu_time = {}
    for element in range(len(keys)):
        swu_cum = np.array(swu_dict[keys[element]])
        facility = keys[element]
        swu.append(swu_cum)
        facilities.append(facility)
        swu_time[keys[element]] = swu_cum
    swu_sort = sorted(swu_time.items(), key=lambda e: e[
        1][-1], reverse=True)
    facilities = [item[0] for item in swu_sort]
    swus = [item[1] for item in swu_sort]
    times = np.arange(0, duration, 1)
    plt.stackplot(times, swus, labels=facilities)
    plt.legend(loc='upper left')
    plt.xlabel('Time [months]')
    plt.ylabel('SWU')
    plt.title('Cumulative SWU by Facility')
    plt.show()


def plot_swu(cur, facilities=[]):
    """Returns dictionary of swu timeseries for each enrichment plant

    Parameters
    ----------
    cur: sqlite cursor
        sqlite cursor
    facilities : list
        list of facilities to plot

    Returns
    -------
    swu_dict: dictionary
        dictionary with "key=Enrichment (facility number), and
        value=swu timeseries list"
    """
    swu_dict = {}
    if len(facilities) != 0:
        agentid = facilities
    agentid = agent_ids(cur, 'Enrichment')
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    for num in agentid:
        swu_data = cur.execute('SELECT time, value '
                               'FROM timeseriesenrichmentswu '
                               'WHERE agentid = ' + str(num)).fetchall()
        swu_timeseries = timeseries(swu_data, duration, False)
        swu_dict['Enrichment_' + str(num)] = swu_timeseries

    keys = []
    for key in swu_dict.keys():
        keys.append(key)
    swu = []
    facilities = []
    swu_time = {}
    for element in range(len(keys)):
        swu_list = np.array(swu_dict[keys[element]])
        facility = keys[element]
        swu.append(swu_list)
        facilities.append(facility)
        swu_time[keys[element]] = swu_list
    swu_sort = sorted(swu_time.items(), key=lambda e: e[
        1][-1], reverse=True)
    facilities = [item[0] for item in swu_sort]
    swus = [item[1] for item in swu_sort]
    times = np.arange(0, duration, 1)
    plt.stackplot(times, swus, labels=facilities)
    plt.legend(loc='upper left')
    plt.xlabel('Time [months]')
    plt.ylabel('SWU')
    plt.title('SWU by Facility')
    plt.show()


def plot_cumulative_power(cur, reactors):
    """Plots cumulative power of reactor fleet over the simulation duration.

    Parameters
    ----------
    cur : sqlite cursor
        sqlite cursor
    reactors : list
        list of reactors to plot

    Returns
    -------
    plot : plot
        plot of cumulative powers
    """
    power_dict = {}
    agentid = agent_ids(cur, 'Reactor')
    if len(reactors) != 0:
        agentid = reactors
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    for num in agentid:
        power_data = cur.execute('SELECT time, value '
                                 'FROM timeseriespower '
                                 'WHERE agentid = ' + str(num)).fetchall()
        power_timeseries = timeseries_cum(power_data, duration, False)
        power_dict['Reactor_' + str(num)] = power_timeseries

    keys = []
    for key in power_dict.keys():
        keys.append(key)
    power = []
    reactors = []
    power_time = {}
    for element in range(len(keys)):
        power_level = np.array(power_dict[keys[element]])
        reactor = keys[element]
        power.append(power_level)
        reactors.append(reactor)
        power_time[keys[element]] = power_level
    power_sort = sorted(power_time.items(), key=lambda e: e[
        1][-1], reverse=True)
    reactors = [item[0] for item in power_sort]
    powers = [item[1] for item in power_sort]
    times = np.arange(0, duration, 1)
    plt.stackplot(times, powers, labels=reactors)
    plt.legend(loc='upper left')
    plt.xlabel('Time [months]')
    plt.ylabel('Power [MWe]')
    plt.title('Power: cumulative')
    plt.show()


def plot_power_reactor(cur, reactors):
    """Plots power of reactor fleet over the simulation duration.

    Parameters
    ----------
    cur :  mlite cursor
        sqlite cursor
    reactors : list
        list of reactors to plot

    Returns
    -------
    plot : plot
        plot of reactor powers
    """
    power_dict = {}
    agentid = agent_ids(cur, 'Reactor')
    if len(reactors) != 0:
        agentid = reactors
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    for num in agentid:
        power_data = cur.execute('SELECT time, value '
                                 'FROM timeseriespower '
                                 'WHERE agentid = ' + str(num)).fetchall()
        power_timeseries = timeseries(power_data, duration, False)
        power_dict['Reactor_' + str(num)] = power_timeseries

    keys = []
    for key in power_dict.keys():
        keys.append(key)
    power = []
    reactors = []
    power_time = {}
    for element in range(len(keys)):
        power_list = np.array(power_dict[keys[element]])
        reactor = keys[element]
        power.append(power_list)
        reactors.append(reactor)
        power_time[keys[element]] = power_list
    power_sort = sorted(power_time.items(), key=lambda e: e[
        1][-1], reverse=True)
    reactors = [item[0] for item in power_sort]
    powers = [item[1] for item in power_sort]
    times = np.arange(0, duration, 1)
    plt.stackplot(times, powers, labels=reactors)
    plt.legend(loc='upper left')
    plt.xlabel('Time [months]')
    plt.ylabel('Power [MWe]')
    plt.title('Reactor Power')
    plt.show()


def powerseries_reactor(cur, reactors):
    """Returns power of reactor fleet over the simulation duration.

    Parameters
    ----------
    cur :  mlite cursor
        sqlite cursor
    reactors : list
        list of reactors to plot

    Returns
    -------
    power_dict : dict
        dictionary of reactor fleet powers
    """
    power_dict = {}
    agentid = agent_ids(cur, 'Reactor')
    if len(reactors) != 0:
        agentid = reactors
    init_year, init_month, duration, timestep = simulation_timesteps(cur)

    for num in agentid:
        power_data = cur.execute('SELECT time, value '
                                 'FROM timeseriespower '
                                 'WHERE agentid = ' + str(num)).fetchall()
        power_timeseries = timeseries(power_data, duration, False)
        power_dict['Reactor_' + str(num)] = power_timeseries
    return power_dict


def total_isotope_used(cur, facility):
    """Returns dictionary of total masses of isotopes mined

    Parameters
    ----------
    cur :  cursor
        sqlite cursor
    facility : str
        str of mine facility

    Returns
    -------
    total_isotopes_mined : dict
        dictionary of isotopes mined and the total mass mined
    """
    flux = 'out'
    isotope_masses_used = cumulative_mass_timeseries(cur, facility, flux)[0]

    total_isotope = [item[1][-1] for item in isotope_masses_used.items()]
    nuclides = [item[0] for item in isotope_masses_used.items()]
    total_mass_used = {}
    for i in range(len(total_isotope)):
        nuclide = nuclides[i]
        mass = total_isotope[i]
        total_mass_used[nuclide] = mass
    return total_mass_used