Modelling the renewable transition: Scenarios and pathways for a decarbonized future using pymedeas, a new open-source energy systems model

This paper reviews different approaches to modelling the energy transition towards a zero carbon economy. It identifies a number of limitations in current approaches such as a lack of consideration of out-of-equilibrium situations (like an energy transition) and non-linear feedbacks. To tackle those issues, the new open source integrated assessment model pymedeas is introduced, which allows the exploration of the design and planning of appropriate strategies and policies for decarbonizing the energy sector at World and EU level. The main novelty of the new open-source model is that it addresses the energy transition by considering biophysical limits, availability of raw materials, and climate change impacts. This paper showcases the model capabilities through several simulation experiments to explore alternative pathways for the renewable transition. In the selected scenarios of this work, future shortage of fossil fuels is found to be the most influential factor of the simulations system evolution. Changes in efficiency and climate change damages are also important determinants influencing model outcomes.

This paper reviews different approaches to modelling the energy transition towards a zero carbon economy. It identifies a number of limitations in current approaches such as a lack of consideration of out-of-equilibrium situations (like an energy transition) and non-linear feedbacks. To tackle those issues, the new open source integrated assessment model pymedeas is introduced, which allows the exploration of the design and planning of appropriate strategies and policies for decarbonizing the energy sector at World and EU level. The main novelty of the new open-source model is that it addresses the energy transition by considering biophysical limits, availability of raw materials, and climate change impacts. This paper showcases the model capabilities through several simulation experiments to explore alternative pathways for the renewable transition. In the selected scenarios of this work, future shortage of fossil fuels is found to be the most influential factor of the simulations system evolution. Changes in efficiency and climate change damages are also important determinants influencing model outcomes.

Introduction
Today's societal challenges require new tools and models that consider, in an integrative way, aspects such as climate change impacts [1], impacts and vulnerabilities due to resource limitations [2], biophysical resources management [3] and human impacts on ecosystems [4].
Models can be classified following the fields/areas they analyse or the features they can represent. For instance, for economic analysis the most used models are Computational General Equilibrium (CGE) models and the Dynamic Stochastic General Equilibrium (DSGE) models. In the field of energy market modelling, bottom-up dynamic partial equilibrium models like the MARket ALlocation models and TIMES, are widely used, and combine technical engineering and economic approaches [5]. If the features they represent are analysed, eight aspects could be considered [6]: (i) complexity, non-linearity, non-ergodicity and deep uncertainty, (ii) the importance of time, (iii) agents' heterogeneity and behavioural elements, (iv) interdisciplinary aspects (v) role of institutions and social context, (vi) ethical and philosophical aspects, (vii) finance and (vii) multiple equilibria/disequilibrium. In the case of CGE or the DSGE for instance, they do not have the possibility to represent out-of-equilibrium situations (like an economic crisis or energy transition), non-linear feedbacks or other system characteristics related to complexity. In this line [6], classified eleven models following these eight aspects and four general types: econometric, system dynamics, agent-based and Stock-Flow Consistent (SFC) models. In the econometric type they distinguish those that are Keynesian and post-Keynesian [6]. identified eleven models that met such criteria and in the case of pymedeas, it is a system dynamics, econometric (using an input-output approach) model framed in a post-Keynesian approach.
Although there are different models combining System Dynamics and input-output analysis in the literature [7], most of them put the emphasis on only one side: System Dynamics [8,9] or input-output analysis [10,11]. More balanced presentations are rare and specific for applications, such as ecosystem services analysis [12,13] or other more focused applications. From the point of view of climate impacts and socio-economy aspects, Integrated Assessment Models (IAM) have been used to explore future scenarios under Climate Change IPCC projections [14]. They use different models for each of 5 Shared Socio-economic Pathways (SSP), some of them are IMAGE [15], MESSAGE [16] or POLES [17] (a review of IAMs for IPCC scenarios could be found at the IAMC website https://www.iamconsortium.org, and at [18]). However, SSPs do not consider resource limitation scenarios as a main hypothesis. The resources role (minerals, fossil fuels availability/cost) in the projections is considered, for instance, in the World 6 model [19]. World models have their origin in the limits to growth work [20], all World models are global with no geographical disaggregation. Thus, the current different simulation tools here mentioned (economy, energy, IAMs) have a constellation of different aspects, scattered through different models, but they do not consider the possibility to jointly analyse the resource limitations' role, the impact of climate and the feedbacks with the environment, the economy evolution (in terms of efficiency) and also how these different aspects are affecting lower geographical levels (region and country).
Thus, this approach of linking the two methodologies for energyeconomic analyses at different geographical levels represents a step forward to the current models' state of the art because pymedeas is a model that mixes economy and energy-resources analysis using a novel coupled Input-Output energy-resources methodology. The novelty of pymedeas within the existing modelling approaches lies also in the new scenarios defined to frame the model runs. Different scenarios, which define initial hypotheses on key parameters and system evolution hypotheses, can be fed to the model as input data. The scenarios act as a framework that helps understand the links between exogenous variables and endogenous variables of the model like GHG emissions and the energy-economy evolution and thus, climate change impacts on the economy. The key exogenous variables that make up each scenario include those related to the expected rate of deployment/implementation of RES, transport electrification, the efficiency improvement, population growth or economic growth, amongst others.
Here we introduce the detail of pymedeas, a new open-source piece of software (written in Python 3), for exploring the transition to a low carbon socio-economy [21,22]. This process relies mainly on the commitment of the European Union to reach at least 32% renewables by 2030 and reduce its GHG by 40% until 2030 and by 80-90% until 2050, compared to 1990 [23,24]. This requires intense emission reduction [25], particularly in the energy sector, changing it from fossil fuel-based to RES based, in a process called the Renewable Transition.
The pymedeas models have three main objectives. The first is to identify the key physical variables for the energy transition and their relationships with economic indicators, socio-economic variables and a range of different environmental impacts. The second objective is to quantify the implications of the emerging challenges within the implementation of a transition to a low-carbon economy. These may include the impact of technological parameters (e.g. RES capacity factors, Energy Return On Energy Investment or energy intensities evolution), and how to overcome possible drawbacks and provide solutions. Finally, the third objective is to suggest policies and strategies to face such challenges when drafting the roadmap to a European future socio-economic transition to a sustainable energy system.
The approach to energy systems modelling builds on four pillars: open source, transparency, user friendly and community -based. The model was originally developed with the Vensim® software [26]  guaranteed by the MIT license; all code can be downloaded for free from the project's website (https://medeas.eu/model/medeas-model) or from a public Gitlab repository (https://gitlab.com/MEDEAS/pymede as_models). The transparency aspect is guaranteed by a detailed and extensive documentation of the code, also available in the website. The MIT license, the Git repository (contributions from any third party are welcome) and a users' forum in the projects' website are expected to encourage the development of a sense of a community around the tool. Finally, user friendliness is achieved by providing a GUI for plotting simulation results and allowing the visualization and comparison of outputs from independent simulations.
The pymedeas model can be framed in the set of energy system models within the Open Energy Modelling (OpenMod) initiative (https ://www.openmod-initiative.org). Currently, there is a growing interest in opening up energy models for increasing transparency [27,28] so initiatives as OpenMod are gaining users and the number of models there, has been enlarged. In this sense, pymedeas supports the open modelling initiatives, not only by distributing the software under an open-source license but also the MEDEAS database management system (MEDEAS database: https://medeas.cmima.csic.es).
This paper is structured in four sections. In materials and methods, the structure and features of the pymedeas software, and the scenarios and hypotheses are described. In addition, the benchmark procedure used to compare simulation results from different scenarios is presented. In the Results section the simulation results under the different hypotheses are interpreted, focusing on characteristic output variables. Finally, the last section is dedicated to discuss the model and the main results in the framework of the state of the art of the current modelling tools for energy systems analysis.

Methods
The pymedeas models can simulate at three geographical levels: World, the European Union and Austria (pymedeas_w, pymedeas_eu and pymedeas_at respectively). Although in this work we will only use World and EU geographical levels, the user can also take the region-based (pymedeas_eu) and country level (pymedeas_at) models and adapt them to different regions and countries. The model design and structure accounts for aggregated variables in an hybrid bottom-up top-down approach, and the default simulation period spans from 1995 to 2050. This simulation period can be shortened or extended according to the needs of the analysis. The pymedeas models at the different geographical levels are nested in a one-way approach (the model in a larger geographical level, called parent model, includes the one with a smaller geographical level, also referred to as child model). Therefore, a coherent simulation of the parent model must be run in order to generate the boundary conditions to be used by the child model. Both models have the same structure and characteristics, and they both include the following aspects: a) resources limitation, b) climate change impacts, c) energy efficiency and d) dynamic evolution of the EROI [29,30]. The EROI of a system corresponds to the ratio between the final usable energy (exergy) and the amount of energy used to obtain such energy. The model structure, software, scenarios and main hypotheses to compare the behaviour of the model are introduced bellow.

The model's general structure
The model takes into account the effects of biophysical limits, resource limitations and climate change impacts on the economy [26,31,32]. Fig. 1 shows the interrelations between the 7 modules included in the model, represented by boxes, the main characteristics of which are [26]: � Economy and population: the economy is modelled [33] assuming non-clearing markets (i.e. not forcing general equilibrium), demand-led growth and supply constraints. This assumes that demand matters in both the long as well as the short run, so that a competitive market economy has no natural or automatic tendency towards full employment [34,35]. The economic structure is captured by the integration of IOT (35 industrial sectors and households) following the WIOD (www.wiod.org), for each region (World and EU). MEDEAS builds on the open-access database WIOD [36], to remain consistent with the open-source aspect. The link between the monetary and the energy one is performed through the modelling of sectoral final energy intensities [37]. � Energy: this module includes renewable and non-renewable energy (extraction reserves and resources) and availability, taking into account biophysical and temporal constraints. In total, 5 final energy carriers (electricity, heat, solids, gases and liquids) and a diversity of energy technologies are modelled. A net energy approach is applied [38] endogenously and dynamically accounting for the EROI of both individual technologies and that of the system. The demand of energy is affected by the variation of the EROI of the system [38]. � Energy infrastructures: represents the infrastructures of power plants to generate electricity and heat. � Climate: this module projects the climate change levels (rising temperature) due to the anthropogenic GHG emissions, which will produce economic costs [39,40]. This module is an adapted version of the C-Roads climatic model [41]. � Materials: the module tracks the material requirements for the construction, operation and maintenance of the energy infrastructures. The extraction demand is subsequently compared with the levels of available reserves and resources. � Land-use: this part of the model accounts for the land requirements of the RES. � Social and environmental impacts: this module contextualizes the implications for human societies in terms of well-being for each simulation. It translates the "biophysical" results of the simulations into metrics related with social and environmental impacts.
The model evolution and integration across modules requires, for each time step, the energy demand and supply to be balanced dynamically. The energy demand is driven endogenously by the economic activity simulated by the evolution of IOT structure [42] and the energy intensities (which translate the monetary flows into final energy by carrier) and also by the expected GDP evolution introduced exogenously by the user (in the input file). The IOT provide the basis to account for all direct and indirect demands from each sector, including energy demand-based on the sector-level energy intensities. WIOD includes a detailed, sector-level disaggregation of energy demand in its environmental accounts [43]. When final energy demand is lower than final energy availability, the expected demand (external, desired GDP) is fulfilled and economic growth is possible. On the other hand, if demand exceeds supply, then final demand adapts to available final energy. Thus, the demand of the next time-step is estimated taking as reference the consumption of the last time-step, constrained by an energy scarcity function to account for fuel supply limitation.
The final energy consumption by fuel is covered by a mix of technologies (energy infrastructures module), in the consumption of primary energy. In the current version of the model, material availability does not directly constrain the deployment of technologies given the uncertainty in available metrics of reserves and resources.
The extent of the primary energy consumption by fuel translates into a certain amount of GHG emissions, which affect climate change which, in turn, feeds back into human societies through a certain level of impacts on e.g. adaptation costs, increase in energy needs. Climate effects are translated to the economy by a cost damage function [41], which will increase the energy consumption depending on the atmospheric GHG concentrations in each time step.

The pymedeas code
The pymedeas code is the result of translating to Python the Vensim® versions of the models [26], using the PySD library [45]. The code can be downloaded either from the project website (https://medeas.eu/mod el/medeas-model) or from the public git repository (https://medeas.eu/ model/medeas-model). The version of the code used in this work corresponds to v0.3.1.
When running simulations, either the pymedeas_w.py or the pyme-deas_eu.py (depending on the users' choice) is loaded as a PySD Model object [45]. By design, the Model class maintains only a single state of the system in memory, and therefore the future state of the model is calculated strictly based upon their current state. The integration is carried out using the Euler method.
The user interacts with the models through a command line interface that wraps the main functionalities of the PySD. This way, users can modify model parameters, the time-step, and the final simulation date (among others) without any previous knowledge of Python. The outputs of the simulations are stored in the form of a Pandas DataFrame. The three model files (pymedeas_w.py, pymedeas_eu.py and pymedeas_at. py) share elements in common listed in Appendix A. pymedeas comes with a plot tool (plot_tool.py) which is a GUI aimed at displaying simulation results obtained with the pymedeas models (for more description please refer to Appendix A).
The model input data is provided with the downloaded package, but can also be visualized/downloaded from the MEDEAS database management system linked in the project website or at www.medeas.cmima. csic.es.

Scenarios
The scenarios used here are Business as Usual (BAU), which projects into the future the current trends of all input variables, and Optimal Level Transition (OLT), which increases significantly the deployment of RES. OLT, represents a scenario of economic growth and reduction of energy use through improvements in energy efficiency. In addition, OLT increases recycling rates of minerals and boosts the transition to low carbon energy production [23,46].
The differences between the two scenarios start to manifest from 2020, which is the date in which most policies start to take action in the OLT scenario. The two scenarios were previously explored from the GHG emissions perspective [47], where the carbon budget (the cumulative carbon emissions) [48] for the world and the EU were used to estimate decarbonization pathways until 2050. Thus, taking the policy of the European Union to reduce emissions by 80-90% with respect to their 1990 value by 2050 [23,46], the emissions from 2020 to 2050 can be estimated with a further constraint on cumulative greenhouse gases released until 2050. The scenarios used here take a decarbonization hypotheses [37] as starting point, but, the model trajectories are not constrained by the cumulative GHG emissions till 2050. OLT sets a target on the annual GHG emissions, rather than cumulative emissions, following EU policy. In this way, the BAU scenario does not follow Paris Agreement while OLT does. The interest of BAU is to warn about the consequences of not complying with the current international agreements on climate change [49].
The BAU and OLT scenarios take the SSP [51] as a reference. SSPs are scenario frameworks for facilitating the integrative analysis of future climate change impacts, vulnerabilities, adaptation and mitigation. They are based in five narratives, describing socio-economic and  political evolution, including aspects such as fossil-fuel development, regional rivalry or inequality. A summary of the characteristics for the two scenarios is shown in Table 1. Both BAU and OLT follow SSP2 for population and GDP per capita projections, while the other variables and hypotheses differ from SSP2. In BAU, RES deployment growth follows current trends, while in OLT it assumes significant, very rapid, deployment, to reduce GHG emissions. OLT RES rapid deployment is supposed by a medium technological development (compared to slow development of BAU).
In each scenario the general framework is setup; however, particular conditions can be considered exogenously in each scenario, such as the assumptions for the evolution of economic indicators, technological efficiencies or biophysical constraints. Some of the multiple choices of the model are explained in the next subsection.

Hypotheses used
The behaviour of the model is evaluated by comparing a set of simulation results obtained by turning on one at a time, the following model hypotheses in both scenarios: � Evolution of Energy Efficiencies (EFFI): defines the way the energy intensities of the different industrial sectors evolve over time. Historical values of the energy intensities are used between 1995 and 2009. The evolution of the energy efficiency using the energy intensities is already explored in previous literature [34,37]. By default, in the BAU scenario, energy intensities from 2009 evolve following historical trends while in OLT the trends from 2009 are imposed. Modifying the default state of this parameter corresponds to keeping the energy intensity constant at its value in 2009 until 2050, for BAU and OLT. It is assumed fixed trends in the evolution of Energy Efficiencies (using the energy intensity) are the simplest assumption to consider. This behaviour is also consistent with the IO tables, which in their historical yearly evolution do not show high variability of sectoral flows. � EROI feedback on energy demand (EROI): the EROI accounts for net energy and corresponds to the energy delivered from a process divided by the energy required to get it over its lifetime [52,53]. The pymedeas model endogenously calculates the EROI of the whole system over time [38] Using the default setting (for both BAU and OLT), the variations of the EROI of the system over time do not affect the energy demand. Modifying the default state corresponds to making the energy demand sensitive to the variation of the EROI at every time step. � Environmental and Climate damage (CLIMA): The modelling of the damages caused by environmental/climate changes is performed through the integration of an energy loss function that reduces the overall net energy delivered to the society, based on growing atmospheric CO2 concentrations and land-use change. The severity of the damage is regulated by selecting among the four RCPs [54], which provide different scenarios of the temperature change due to GHG emissions. With the default setting, climate change does reduce the net energy delivered to society, RCP 8.5 is used for BAU and RCP 2.6 for OLT. This is also consistent with the damage functions used in the literature [55]. � Fossil fuels scarcity (SCAR): limitations in the availability of fossil fuels caused by resource constraints e.g. Peak-Oil [56]. The default setting imposes maximum depletion curves [57], for the different fuels obtained from recent literature, therefore putting limits to production. A Hubbert curve is a bell-shaped curve that shows an approximation for the production rate of a limited resource over time. Modifying the default setting removes any limitations for the production of non-renewable fuels, making them infinitely available.
The four hypotheses introduced here help to analyse the role of the combined impact of climate/environment, economic evolution, technological change and resources in the future evolution of the socioeconomy. Such analysis are aligned with recent research focused in the trade-off between sustainability and over-consumption [58]. Overall there are 4 hypotheses that can take two different values (default, non-default) each. This results in a total of 16 simulations per scenario (see Table B1), hence 32 simulations per geographical level and 64 model simulations considering the two nested models. Consistency is kept between the scenarios of the parent and the nested models (World and EU), meaning that, the EU simulations for BAU correspond to the World BAU, and the same for OLT. Each simulation is given an identifier following the convention shown in Fig. 2.

Comparison methodology
The results of each of the simulations are compared to those of the simulation with default settings (S1111). Although it is theoretically possible that Europe might opt to undertake the renewable transition regardless of the situation at global scale, this situation will not be considered and therefore the simulations World-BAU and EU-OLT will not be compared. This choice is taken to reduce the number of scenarios' combinations to study, as the aim of this work is to analyse the model hypotheses that have more influence on the outputs.
The Root Mean Square Difference (RMSD) is used to quantify the difference between the simulation with default settings and the other simulations for twelve output variables only, to facilitate the interpretation of the results. The twelve output variables considered are shown in Table 2 (note that temperature change is only endogenously calculated in the World model):

Results
In this section, we show some of the simulation results for the scenarios and hypotheses presented above. For each scenario (BAU and OLT) a set of simulations is compared with the simulation with default settings (S1111) (see Fig. 2 and Table B1 (Appendix B) for more detail on the keys of the simulation codes). Thus, based on the default simulation (all hypothesis switches activated) we switch off the four hypotheses one at time (section 3.1) or in combination (section 3.2). The aim is to identify the most influential hypothesis, compared to the reference (section 3.1) and the impacts and induced changes by the combinations of switches (section 3.2). Fig. 3 shows the values of the default simulation with the World and EU models (left column: World, right column: EU) until 2050 and in both scenarios for variables: GDP, share of blue water, temperature change (World) and CO2 emissions for EU. For each model, the differences between the two scenarios are particularly important for the environmental variables water and temperature. The World model (left column) shows systemic changes when the system achieves the maximum of fossil fuel production (Peak Oil or 'Hubbert peak') [59] in the BAU scenario (green line). These systemic changes have important effects on the economy (shown through the decline of GDP). In comparison, for the OLT scenario in the World model (red line),

Fig. 2.
Graphical representation of the simulation identifiers. "S" for simulation, each of the four digits corresponds to one hypothesis, number 1 corresponds to the default value, and number 0 implies that the default state has been altered.
GDP keeps growing while temperature change has slowed. One can see similar effects in the EU simulations, with a higher impact on the economy due to the lower EU GDP historical growth. Table C1, C2, C3 and C4 in Appendix C show the RMSD over the whole simulated period (1995-2050) for simulations performed for BAU and OLT at the two geographical levels (64 simulations in total) compared with the corresponding default simulation (S1111).
To analyse the RMSD shown in Appendix C, the analysis is split in: 1. Altering the default state of one hypothesis at a time: those that compare the results of altering the default scenario of one hypotheses at a time in each scenario and each geographical level. This corresponds to simulations S1110, S1101, S1011, and S0111. 2. Altering the default state of multiple hypotheses simultaneously: these are simulations with different combinations of hypotheses, including those in which the default state of more than one hypothesis are altered. This will give the rest of simulations in Tables C1-C4 and is shown by means of the maximum combination variable differences in Table 3.

Altering the default state of one hypothesis at a time
To begin with, the following model-scenario combinations are considered: World-BAU, World-OLT, EU-BAU and EU-OLT. The combinations shown below for all four model-scenario duets include the simulation with the default hypothesis (S1111) vs. S0111, S1011, S1101 and S1110. To interpret the results, we focus only on the variable 'GDP' ( Table 2).

World-BAU
The values of the RMSD calculated for the selected output variables obtained from each simulation with respect to the values of the same variables from the simulation with default hypothesis (S1111) are presented in Table C1, Appendix C. With regards to GDP, it can be noted that climate impacts are modulated by the effect of resource limitation (Fig. 4). If there is energy scarcity, then no big differences in the system behaviour due to climate change impacts are observed. This can be observed by comparing S1111 and S1101 (simulation neglecting the impacts of climate change on net energy availability) and small differences are detected. Among all hypotheses, the fossil fuel scarcity is the one that gives the largest differences, followed by the energy intensities and the EROI of the system.

World-OLT
This corresponds to the scenario in which the transition is undertaken with maximum efforts and at the highest pace under a plausible future. The results for each of the twelve output variables in terms of the RMSD calculated with respect to simulation S1111, is shown in Table C2 in Appendix C. For the sake of example, and as it was done in the previous section (World-BAU), in Fig. 5 the results of the different simulations are compared using GDP. The energy intensity of 2009 in the future (S0111) presents worse results in terms of GDP, while activating the link between EROI and available net energy (S1011) results in a lower GDP growth due to the to the net energy decrease of RES over time.

EU-BAU
The results at EU level are coherent with those at global scale. Indeed, the effects of limiting the availability of fossil fuels (default state for SCAR hypothesis) at World scale translates into energy scarcity at EU level. Indeed, only when the default hypothesis is altered (S1110: unrestricted availability of fossil fuels) can GDP continue to grow after 2020 (Fig. 6d).

EU-OLT
Results at EU level are in accordance with the results at World level, also for the OLT scenario. Fig. 7 shows GDP stagnation of the World translates into a decrease of European GDP, although less dramatically, mainly due to the lower historical GDP growth for the region. From such decline follows stagnation and a slight increase at the end of the period.

Altering the default state of multiple hypotheses simultaneously
To evaluate the effects of altering the default state of more than one of the four hypotheses at a time, the simulation with the maximum RMSD (with respect to S1111) is highlighted for each variable from the selected set of outputs (i.e. 12 for the World and 10 for the EU, since variables temperature_change and total_extraction_nre_ej are not endogenously calculated in the EU model). As in the previous section, the results are interpreted individually for each model-scenario couple (i.e. World-BAU, World-OLT, EU-BAU and EU-OLT).

World-BAU
As in the previous subsection, in this case the most significant differences between the simulations with default hypothesis and those altering the default hypotheses are observed when assuming limitless availability of fossil fuels (the right-most figure in the simulation identifier set to 0) ( Table 3). The only exception is the variable representing the energy intensity of the whole system (tfes_intensity_ej_tdollar), which is the ratio between the energy consumption (in EJ) and the GDP (in T$). This exception is due to the fact that both in S1111 and S0011 resources are limited, but in the later the energy demand is higher (constant energy intensities of 2009 over time). Therefore, in S0011 the resource limitation hits before they do in S1111, sending the GDP down, which results in higher values of tfes_intensity_ej_tdollar. On the contrary, when removing the non-renewable energy limitation (as in S0000), even when the increase of energy demand caused by keeping energy intensities constant is considered, it is compensated with infinite availability of non-renewable energy sources, and the ratio between energy consumption and GDP remains relatively constant and smaller than in S0011, hence more similar to that of S1111. As expected, CO2 emissions and temperature changes give maximum differences with respect to S1111 in simulation S0000, where all hypotheses have been Fig. 4. GDP in selected World BAU simulations: a) S1111 vs. S0111, b) S1111 vs. S1011, c) S1111vs. S1101 and d) S1111 vs. S1110. In the name of the simulation 1 means default and 0 non-default. The first number after the S corresponds to Energy efficiency evolution (EFFI), the second corresponds to the effects of EROI on net energy availability, the third is the climate change impacts (CLIMA) and the fourth is the availability of fossil fuels (SCAR).
altered. In such simulation, although there is no limit for fossil fuels availability the GDP is disconnected from the climate damage produced by GHG emissions.

World-OLT
Very similar results to those of World-BAU are obtained in this case.
The only notable difference occurs for variables GDP and per capita GDP, which have the highest RMSD when including the impacts of EROI on the net energy available. This result can be understood by the fact that, in OLT, more RES is implemented compared to BAU and thus, the resulting lower EROI of the system has a bigger impact even in the case of unlimited supply or non-renewable energy.

Fig. 5. GDP in selected World
OLT simulations: a) S1111 vs. S0111, b) S1111 vs. S1011, c) S1111vs. S1101 and d) S1111 vs. S1110. In the name of the simulation 1 means default and 0 non-default. The first number after the S corresponds to Energy efficiency evolution (EFFI), the second corresponds to the effects of EROI on net energy availability, the third is the climate change impacts (CLIMA) and the fourth is the availability of fossil fuels (SCAR).  6. GDP in selected EU BAU simulations: a) S1111 vs. S0111, b) S1111 vs. S1011, c) S1111vs. S1101 and d) S1111 vs. S1110. In the name of the simulation 1 means default and 0 non-default. The first number after the S corresponds to Energy efficiency evolution (EFFI), the second corresponds to the effects of EROI on net energy availability, the third is the climate change impacts (CLIMA) and the fourth is the availability of fossil fuels (SCAR).

EU-BAU
At EU scale results are again similar than those of the couple World-BAU. In this case, the main differences are observed for variables total land requirements and GDP and per capita GDP. This is mainly due to the role played by the energy intensities and the imports in the child EU model.

EU-OLT
Again, in this regional simulation (Table 3) the differences from the corresponding parent model (World-OLT) are related to the energy intensities, which as in the EU-BAU are related to the import function and the scarcity function. Both will affect the economic activity. As a general result, it can be seen that the most restrictive hypothesis is the limited availability of fossil fuels, followed, in some cases, by the choice of the evolution of energy intensities. The climate change impacts on the economy are important when they act together with resource limitations. The role of EROI is important when a massive deployment of RES takes place and it also contributes significantly in the OLT scenarios for the EU, due to constraints in imports and the scarcity function on the evolution of the system.
Overall, results show the importance of non-renewable resources in the future evolution of the economy. In particular, simulation results show the influence of fossil fuel supply in the economy and their impact (through GHG emissions) on global temperature. Following the Business as Usual path (non-compliant with Paris Agreement) will conduct to not only environmental but also to potentially catastrophic economic consequences, which is an scenario that must be avoided by all means [49,60].

Discussion
This paper introduces a new open-source, python programmed, modelling tool for energy-economy-environment analysis. Following the two categorizations of [61], pymedeas is a Mathematical (mechanistic) and Energy Economics model with a flexible bottom-up top-down approach depending on the sector. The most used (and known) bottom-up energy models are MARKAL (MARket ALlocation) family of models [5] and their evolution [5] TIMES, developed by the International Energy Agency's Energy Technology Systems Analysis Program (IEA-ETSAP), and the MESSAGE (Model for Energy Supply Strategy Alternatives and their General Environmental Impacts) family of models [62]. The open source energy model systems, OSeMOSYS framework [63] is a widely used energy modelling system, which is a bottom-up model system that has a global domain and can be applied to different geographical levels and with different levels of complexity. pymedeas complements such kind of bottom-up approaches by providing the top-down approach, using aggregated variables.
There are a number of novelties [26]: First, this tool combines system dynamics and input-output analysis Fig. 7. GDP in selected EU OLT simulations: a) S1111 vs. S0111, b) S1111 vs. S1011, c) S1111vs. S1101 and d) S1111 vs. S1110. In the name of the simulation 1 means default and 0 non-default. The first number after the S corresponds to Energy efficiency evolution (EFFI), the second corresponds to the effects of EROI on net energy availability, the third is the climate change impacts (CLIMA) and the fourth is the availability of fossil fuels (SCAR).

Table 3
Simulations with the highest RMSD per variable, compared to the reference simulation S1111 for World and EU under BAU and OLT scenarios. The code follows the reference key shown in Appendix B, Table B1. In the name of the simulation 1 means default and 0 non-default. The first number after the S corresponds to Energy efficiency evolution (EFFI), the second corresponds to the effects of EROI on net energy availability, the third is the climate change impacts (CLIMA) and the fourth is the availability of fossil fuels (SCAR). S0000 S0000 S0000 gdp S1100 S1000 S0100 S0000 gdppc S1100 S1000 S0100 S0000 percent_res_vs_tpes S1000 S0000 S0000 S0000 real_tfec S0000 S0000 S0000 S0000 share_blue_water_use_vs_ar S1000 S1000 S0000 S0000 temperature_change S0000 S0000 --total_land_requirements_renew_mha S1000 S1000 S0000 S0000 tpe_from_res_ej S0000 S0000 S0000 S0000 [7,64] to study the evolution of the energy system at two geographical levels under environmental and biophysical constraints. The system dynamics approach allows the analysis of non-linear feedback in the system whereas the input-output approach provides very detailed representation of socio-economic variables albeit using fixed coefficients (e. g. Refs. [42]). The input-output matrix at the core of pymedeas enables the analyst to assess environmental, economic and social impacts resulting from each scenario. The pymedeas model takes the data from the open-access database World Input-Output tables and Data (WIOD, www.wiod.org), which provides the social and environmental accounts to conduct such analysis [36]. Secondly, a set of scenarios and hypotheses were developed to enable projections on the future evolution of the system. The MEDEAS scenarios used are Business as Usual, based on current trends, and Optimal Level Transition, supposing economic growth while achieving a reduction of energy, improvements in energy efficiency, increase in recycling rates of minerals, and low carbon energy production. Both are based on the shared socio-economic pathways also used by the IPCC. The set of hypotheses and the differences show how the fossil fuels availability constraint is the most important factor modifying the models' output. This is particularly significant depending on the scenario chosen. For BAU, the consequences of fossil resources scarcity [2,43] are considerable compared to OLT [53,65,66]. GHG emissions act at the end of the simulation period where the climate change damage on the economy is playing a visible role. The other hypotheses have minor effects compared to resources' limitation. In general, the propagation of 2009 energy intensities is the second most impactful factor, followed by EROI activation in the case of the OLT scenarios. The massive requirement for energy during RES deployment reduces the net energy availability for the functioning of the system. It should be noted here that pymedeas also gives estimations of raw materials needs for each scenario but the analysis of the models' output in this aspect is out of the scope of this paper.
Thirdly: The modularity of the model in both structure and functions, gives the user multiple options to make detailed investigations on a range of subjects. This allows the analyst to include for example job effects of energy policies at the low carbon transition. Current studies usually show positive job effects associated with higher renewable energy deployment compared to BAU among those using input-output approaches (e.g. Refs. [67][68][69][70]). The effects of renewable energy expansion on employment (see arrow to Social Module in Fig. 1) are especially critical to foster public acceptance of the transition [71,72]. Thus, future expansion of model capabilities could take advantage of such modularity, for instance, by considering the implementation of an agent-based model in the existing pymedeas social module to analyse the job effects.
Fourthly, pymedeas also adds a new original approach focusing not only on the demand side. The main energy system model approaches are oriented to analyse the demand side of the energy sector and thus, the behaviour of the model is driven by the energy needed to keep the economic activity. These equilibrium models solve the trade-off between the demand side (which try to maximize utility) and the supply side (which maximizes profit) [61] by using energy prices. The pymedeas model takes a different and original approach, firstly by using energy (instead of prices) in a dynamical equilibrium between supply and demand. The pymedeas model considers not only the systems' demand driven evolution and a (optionally unlimited) supply, but also situations where scarcity of a fuel (or the constraints imposed by climate change, allowing the user the option to keep fossil resources underground) have an important role in the behaviour of the system. Moreover, the pymedeas design uses a (one way) nested approach by a system dynamics import function. This also implies a nested approach for the IOT, linked at all geographical levels. To illustrate the model's capabilities, we have constructed a comparative set of simulations to show how the model responds to some key features.
These set of novelties aim to give a complementary approach to the existent literature on the Shared Socioeconomic Pathways (SSP) [51]. These scenarios provide different socio-economic projections up to 2100, based on a set of five narratives: sustainability (SSP1) [22], middle of the road (SSP2) [73], regional rivalry (SSP3) [74], inequality (SSP4) [75] and fossil-fuelled development (SSP5) [76]. In this work it is considered SSP2 scenario as a reference, with additional hypotheses.
The new aspects and hypothesis implemented in pymedeas allow to analyse the impact of the economy (energy intensities evolution) RES deployment (EROI impact) climate change and energy scarcity on the evolution of the system. The implications of using pymedeas, besides their openness and transparency, are the possibility to analyse how the different climate and sustainability policies are affected by resources limitations (energy/materials) and the climate impacts itself as a sink of energy/resources. On the other side, it allows to analyse the impacts of RES massive deployment and the role played by the socio-economy (need of growth). As this study shows, all these factors can be analysed together or (by switching on/off the functions) separately.

Limitations
The current versions of the MEDEAS models have some limitations that are expected to be overcome in future releases [26]. The main limitation in the Economy module resides in the fact that the IOT are static (historical values for year 2009); hence even in scenarios of severe regime changes such as those simulated in the current work, the economic structure remains unchanged. Another limitation, in this case in the Energy and infrastructures module, lies in the value of the EROI for the different energy sources, which is only dynamically calculated for the RES technologies that generate electricity. The EROIs of the rest of energy sources remain constant over time. Hence it is not possible to evaluate the effects of the progressive decline of the EROI of fossil fuels that has historically occurred and is expected to continue in the future [79].
The Energy and infrastructures module does not explicitly model electricity grids, given that these infrastructures are regional/national by definition. However, an estimate of the additional grids per MW of variable RES (overgrids and inter-regional grids) to be constructed to integrate the renewable variable electricity generation is made. Thereafter the additional material requirements associated to such grid developments are computed, which affects the EROI of the RES variables for electricity generation.
Regarding the climate module, only the emissions of CO 2 are simulated, while those of the other GHG (N 2 O, PFC HFC, SF 6 and CH 4 ) are obtained from projections (RCPs). Similarly, the models do not include any carbon sequestration technology, although they allow to indirectly store carbon through an afforestation program covering 345 MHa. The potential CO 2 sequestration obtained by enforcing such program (1.5 GtC/year over 50 years) is discussed in Ref. [78]. The inclusion of carbon sequestration technologies in the code requires a detailed analysis on the impacts of these technologies in the EROI of RES [65] and the other energy sources, which could be a subject for a future work.
Despite it being possible to simulate environmental impacts such as land and water use and CO 2 emissions, others such as land degradation, biodiversity loss, noise and visual impacts, pollution of water and air cannot be evaluated.
Finally, it should be noted that the pymedeas code has been tested in November 2018 in a workshop with stakeholders (MEDEAS lab). Some of the stakeholders' suggestions helped to improve the GUI and model capabilities, other are part of the limitations listed above. Based on the mentioned above limitations, future research could take several directions. One possible development, besides the already mentioned before, could be implementing a social behaviour module (i.e. agentbased model) and also including an open source model for the electricity production sector in pymedeas regional (EU) or country models. As all community-based codes, pymedeas aims to evolve as the users of the model enlarge. Thus, pymedeas aims to contribute to the energy research area in their urgent need for transparency and open-source codes and databases [59].

Conclusions
The pymedeas models are useful tools that can easily be adapted, not only by modifying the input scenarios, but also changing the model equations and relationships based on the particular needs of the user. This flexibility is attributable to the fact that the models are written in Python, which is a popular high-level programming language with English-like syntax and semantics and with a gentle learning curve. The strategy followed by pymedeas is totally aligned with the open-source software initiative [60], which include the rights to study, change, and distribute the software to anyone and for any purpose. Any suggestions or contributions to the code can be shared through the MEDEAS User Forum (https://medeas.eu/forum). This allows for greater transparency and usability, and opens this tool to a wide audience, from specialists to policy makers and from consultancy to teaching at undergraduate level.