diff --git a/.gitignore b/.gitignore index ac252a84..10243def 100644 --- a/.gitignore +++ b/.gitignore @@ -28,7 +28,8 @@ gurobi.log /data/.nfs* /data/Industrial_Database.csv /data/retro/tabula-calculator-calcsetbuilding.csv -/data +/data/nuts* + *.org *.nc diff --git a/LICENSE.txt b/LICENSE.txt index 9cecc1d4..dc10fd32 100644 --- a/LICENSE.txt +++ b/LICENSE.txt @@ -1,674 +1,20 @@ - GNU GENERAL PUBLIC LICENSE - Version 3, 29 June 2007 - - Copyright (C) 2007 Free Software Foundation, Inc. - Everyone is permitted to copy and distribute verbatim copies - of this license document, but changing it is not allowed. - - Preamble - - The GNU General Public License is a free, copyleft license for -software and other kinds of works. - - The licenses for most software and other practical works are designed -to take away your freedom to share and change the works. By contrast, -the GNU General Public License is intended to guarantee your freedom to -share and change all versions of a program--to make sure it remains free -software for all its users. We, the Free Software Foundation, use the -GNU General Public License for most of our software; it applies also to -any other work released this way by its authors. You can apply it to -your programs, too. - - When we speak of free software, we are referring to freedom, not -price. Our General Public Licenses are designed to make sure that you -have the freedom to distribute copies of free software (and charge for -them if you wish), that you receive source code or can get it if you -want it, that you can change the software or use pieces of it in new -free programs, and that you know you can do these things. - - To protect your rights, we need to prevent others from denying you -these rights or asking you to surrender the rights. Therefore, you have -certain responsibilities if you distribute copies of the software, or if -you modify it: responsibilities to respect the freedom of others. - - For example, if you distribute copies of such a program, whether -gratis or for a fee, you must pass on to the recipients the same -freedoms that you received. You must make sure that they, too, receive -or can get the source code. And you must show them these terms so they -know their rights. - - Developers that use the GNU GPL protect your rights with two steps: -(1) assert copyright on the software, and (2) offer you this License -giving you legal permission to copy, distribute and/or modify it. - - For the developers' and authors' protection, the GPL clearly explains -that there is no warranty for this free software. For both users' and -authors' sake, the GPL requires that modified versions be marked as -changed, so that their problems will not be attributed erroneously to -authors of previous versions. - - Some devices are designed to deny users access to install or run -modified versions of the software inside them, although the manufacturer -can do so. This is fundamentally incompatible with the aim of -protecting users' freedom to change the software. The systematic -pattern of such abuse occurs in the area of products for individuals to -use, which is precisely where it is most unacceptable. Therefore, we -have designed this version of the GPL to prohibit the practice for those -products. If such problems arise substantially in other domains, we -stand ready to extend this provision to those domains in future versions -of the GPL, as needed to protect the freedom of users. - - Finally, every program is threatened constantly by software patents. -States should not allow patents to restrict development and use of -software on general-purpose computers, but in those that do, we wish to -avoid the special danger that patents applied to a free program could -make it effectively proprietary. To prevent this, the GPL assures that -patents cannot be used to render the program non-free. - - The precise terms and conditions for copying, distribution and -modification follow. - - TERMS AND CONDITIONS - - 0. Definitions. - - "This License" refers to version 3 of the GNU General Public License. - - "Copyright" also means copyright-like laws that apply to other kinds of -works, such as semiconductor masks. - - "The Program" refers to any copyrightable work licensed under this -License. Each licensee is addressed as "you". "Licensees" and -"recipients" may be individuals or organizations. - - To "modify" a work means to copy from or adapt all or part of the work -in a fashion requiring copyright permission, other than the making of an -exact copy. The resulting work is called a "modified version" of the -earlier work or a work "based on" the earlier work. - - A "covered work" means either the unmodified Program or a work based -on the Program. - - To "propagate" a work means to do anything with it that, without -permission, would make you directly or secondarily liable for -infringement under applicable copyright law, except executing it on a -computer or modifying a private copy. Propagation includes copying, -distribution (with or without modification), making available to the -public, and in some countries other activities as well. - - To "convey" a work means any kind of propagation that enables other -parties to make or receive copies. Mere interaction with a user through -a computer network, with no transfer of a copy, is not conveying. - - An interactive user interface displays "Appropriate Legal Notices" -to the extent that it includes a convenient and prominently visible -feature that (1) displays an appropriate copyright notice, and (2) -tells the user that there is no warranty for the work (except to the -extent that warranties are provided), that licensees may convey the -work under this License, and how to view a copy of this License. If -the interface presents a list of user commands or options, such as a -menu, a prominent item in the list meets this criterion. - - 1. Source Code. - - The "source code" for a work means the preferred form of the work -for making modifications to it. "Object code" means any non-source -form of a work. - - A "Standard Interface" means an interface that either is an official -standard defined by a recognized standards body, or, in the case of -interfaces specified for a particular programming language, one that -is widely used among developers working in that language. - - The "System Libraries" of an executable work include anything, other -than the work as a whole, that (a) is included in the normal form of -packaging a Major Component, but which is not part of that Major -Component, and (b) serves only to enable use of the work with that -Major Component, or to implement a Standard Interface for which an -implementation is available to the public in source code form. A -"Major Component", in this context, means a major essential component -(kernel, window system, and so on) of the specific operating system -(if any) on which the executable work runs, or a compiler used to -produce the work, or an object code interpreter used to run it. - - The "Corresponding Source" for a work in object code form means all -the source code needed to generate, install, and (for an executable -work) run the object code and to modify the work, including scripts to -control those activities. However, it does not include the work's -System Libraries, or general-purpose tools or generally available free -programs which are used unmodified in performing those activities but -which are not part of the work. For example, Corresponding Source -includes interface definition files associated with source files for -the work, and the source code for shared libraries and dynamically -linked subprograms that the work is specifically designed to require, -such as by intimate data communication or control flow between those -subprograms and other parts of the work. - - The Corresponding Source need not include anything that users -can regenerate automatically from other parts of the Corresponding -Source. - - The Corresponding Source for a work in source code form is that -same work. - - 2. Basic Permissions. - - All rights granted under this License are granted for the term of -copyright on the Program, and are irrevocable provided the stated -conditions are met. This License explicitly affirms your unlimited -permission to run the unmodified Program. The output from running a -covered work is covered by this License only if the output, given its -content, constitutes a covered work. This License acknowledges your -rights of fair use or other equivalent, as provided by copyright law. - - You may make, run and propagate covered works that you do not -convey, without conditions so long as your license otherwise remains -in force. You may convey covered works to others for the sole purpose -of having them make modifications exclusively for you, or provide you -with facilities for running those works, provided that you comply with -the terms of this License in conveying all material for which you do -not control copyright. Those thus making or running the covered works -for you must do so exclusively on your behalf, under your direction -and control, on terms that prohibit them from making any copies of -your copyrighted material outside their relationship with you. - - Conveying under any other circumstances is permitted solely under -the conditions stated below. Sublicensing is not allowed; section 10 -makes it unnecessary. - - 3. Protecting Users' Legal Rights From Anti-Circumvention Law. - - No covered work shall be deemed part of an effective technological -measure under any applicable law fulfilling obligations under article -11 of the WIPO copyright treaty adopted on 20 December 1996, or -similar laws prohibiting or restricting circumvention of such -measures. - - When you convey a covered work, you waive any legal power to forbid -circumvention of technological measures to the extent such circumvention -is effected by exercising rights under this License with respect to -the covered work, and you disclaim any intention to limit operation or -modification of the work as a means of enforcing, against the work's -users, your or third parties' legal rights to forbid circumvention of -technological measures. - - 4. Conveying Verbatim Copies. - - You may convey verbatim copies of the Program's source code as you -receive it, in any medium, provided that you conspicuously and -appropriately publish on each copy an appropriate copyright notice; -keep intact all notices stating that this License and any -non-permissive terms added in accord with section 7 apply to the code; -keep intact all notices of the absence of any warranty; and give all -recipients a copy of this License along with the Program. - - You may charge any price or no price for each copy that you convey, -and you may offer support or warranty protection for a fee. - - 5. Conveying Modified Source Versions. - - You may convey a work based on the Program, or the modifications to -produce it from the Program, in the form of source code under the -terms of section 4, provided that you also meet all of these conditions: - - a) The work must carry prominent notices stating that you modified - it, and giving a relevant date. - - b) The work must carry prominent notices stating that it is - released under this License and any conditions added under section - 7. This requirement modifies the requirement in section 4 to - "keep intact all notices". - - c) You must license the entire work, as a whole, under this - License to anyone who comes into possession of a copy. This - License will therefore apply, along with any applicable section 7 - additional terms, to the whole of the work, and all its parts, - regardless of how they are packaged. This License gives no - permission to license the work in any other way, but it does not - invalidate such permission if you have separately received it. - - d) If the work has interactive user interfaces, each must display - Appropriate Legal Notices; however, if the Program has interactive - interfaces that do not display Appropriate Legal Notices, your - work need not make them do so. - - A compilation of a covered work with other separate and independent -works, which are not by their nature extensions of the covered work, -and which are not combined with it such as to form a larger program, -in or on a volume of a storage or distribution medium, is called an -"aggregate" if the compilation and its resulting copyright are not -used to limit the access or legal rights of the compilation's users -beyond what the individual works permit. Inclusion of a covered work -in an aggregate does not cause this License to apply to the other -parts of the aggregate. - - 6. Conveying Non-Source Forms. - - You may convey a covered work in object code form under the terms -of sections 4 and 5, provided that you also convey the -machine-readable Corresponding Source under the terms of this License, -in one of these ways: - - a) Convey the object code in, or embodied in, a physical product - (including a physical distribution medium), accompanied by the - Corresponding Source fixed on a durable physical medium - customarily used for software interchange. - - b) Convey the object code in, or embodied in, a physical product - (including a physical distribution medium), accompanied by a - written offer, valid for at least three years and valid for as - long as you offer spare parts or customer support for that product - model, to give anyone who possesses the object code either (1) a - copy of the Corresponding Source for all the software in the - product that is covered by this License, on a durable physical - medium customarily used for software interchange, for a price no - more than your reasonable cost of physically performing this - conveying of source, or (2) access to copy the - Corresponding Source from a network server at no charge. - - c) Convey individual copies of the object code with a copy of the - written offer to provide the Corresponding Source. This - alternative is allowed only occasionally and noncommercially, and - only if you received the object code with such an offer, in accord - with subsection 6b. - - d) Convey the object code by offering access from a designated - place (gratis or for a charge), and offer equivalent access to the - Corresponding Source in the same way through the same place at no - further charge. You need not require recipients to copy the - Corresponding Source along with the object code. If the place to - copy the object code is a network server, the Corresponding Source - may be on a different server (operated by you or a third party) - that supports equivalent copying facilities, provided you maintain - clear directions next to the object code saying where to find the - Corresponding Source. Regardless of what server hosts the - Corresponding Source, you remain obligated to ensure that it is - available for as long as needed to satisfy these requirements. - - e) Convey the object code using peer-to-peer transmission, provided - you inform other peers where the object code and Corresponding - Source of the work are being offered to the general public at no - charge under subsection 6d. - - A separable portion of the object code, whose source code is excluded -from the Corresponding Source as a System Library, need not be -included in conveying the object code work. - - A "User Product" is either (1) a "consumer product", which means any -tangible personal property which is normally used for personal, family, -or household purposes, or (2) anything designed or sold for incorporation -into a dwelling. In determining whether a product is a consumer product, -doubtful cases shall be resolved in favor of coverage. For a particular -product received by a particular user, "normally used" refers to a -typical or common use of that class of product, regardless of the status -of the particular user or of the way in which the particular user -actually uses, or expects or is expected to use, the product. A product -is a consumer product regardless of whether the product has substantial -commercial, industrial or non-consumer uses, unless such uses represent -the only significant mode of use of the product. - - "Installation Information" for a User Product means any methods, -procedures, authorization keys, or other information required to install -and execute modified versions of a covered work in that User Product from -a modified version of its Corresponding Source. The information must -suffice to ensure that the continued functioning of the modified object -code is in no case prevented or interfered with solely because -modification has been made. - - If you convey an object code work under this section in, or with, or -specifically for use in, a User Product, and the conveying occurs as -part of a transaction in which the right of possession and use of the -User Product is transferred to the recipient in perpetuity or for a -fixed term (regardless of how the transaction is characterized), the -Corresponding Source conveyed under this section must be accompanied -by the Installation Information. But this requirement does not apply -if neither you nor any third party retains the ability to install -modified object code on the User Product (for example, the work has -been installed in ROM). - - The requirement to provide Installation Information does not include a -requirement to continue to provide support service, warranty, or updates -for a work that has been modified or installed by the recipient, or for -the User Product in which it has been modified or installed. Access to a -network may be denied when the modification itself materially and -adversely affects the operation of the network or violates the rules and -protocols for communication across the network. - - Corresponding Source conveyed, and Installation Information provided, -in accord with this section must be in a format that is publicly -documented (and with an implementation available to the public in -source code form), and must require no special password or key for -unpacking, reading or copying. - - 7. Additional Terms. - - "Additional permissions" are terms that supplement the terms of this -License by making exceptions from one or more of its conditions. -Additional permissions that are applicable to the entire Program shall -be treated as though they were included in this License, to the extent -that they are valid under applicable law. If additional permissions -apply only to part of the Program, that part may be used separately -under those permissions, but the entire Program remains governed by -this License without regard to the additional permissions. - - When you convey a copy of a covered work, you may at your option -remove any additional permissions from that copy, or from any part of -it. (Additional permissions may be written to require their own -removal in certain cases when you modify the work.) You may place -additional permissions on material, added by you to a covered work, -for which you have or can give appropriate copyright permission. - - Notwithstanding any other provision of this License, for material you -add to a covered work, you may (if authorized by the copyright holders of -that material) supplement the terms of this License with terms: - - a) Disclaiming warranty or limiting liability differently from the - terms of sections 15 and 16 of this License; or - - b) Requiring preservation of specified reasonable legal notices or - author attributions in that material or in the Appropriate Legal - Notices displayed by works containing it; or - - c) Prohibiting misrepresentation of the origin of that material, or - requiring that modified versions of such material be marked in - reasonable ways as different from the original version; or - - d) Limiting the use for publicity purposes of names of licensors or - authors of the material; or - - e) Declining to grant rights under trademark law for use of some - trade names, trademarks, or service marks; or - - f) Requiring indemnification of licensors and authors of that - material by anyone who conveys the material (or modified versions of - it) with contractual assumptions of liability to the recipient, for - any liability that these contractual assumptions directly impose on - those licensors and authors. - - All other non-permissive additional terms are considered "further -restrictions" within the meaning of section 10. If the Program as you -received it, or any part of it, contains a notice stating that it is -governed by this License along with a term that is a further -restriction, you may remove that term. If a license document contains -a further restriction but permits relicensing or conveying under this -License, you may add to a covered work material governed by the terms -of that license document, provided that the further restriction does -not survive such relicensing or conveying. - - If you add terms to a covered work in accord with this section, you -must place, in the relevant source files, a statement of the -additional terms that apply to those files, or a notice indicating -where to find the applicable terms. - - Additional terms, permissive or non-permissive, may be stated in the -form of a separately written license, or stated as exceptions; -the above requirements apply either way. - - 8. Termination. - - You may not propagate or modify a covered work except as expressly -provided under this License. Any attempt otherwise to propagate or -modify it is void, and will automatically terminate your rights under -this License (including any patent licenses granted under the third -paragraph of section 11). - - However, if you cease all violation of this License, then your -license from a particular copyright holder is reinstated (a) -provisionally, unless and until the copyright holder explicitly and -finally terminates your license, and (b) permanently, if the copyright -holder fails to notify you of the violation by some reasonable means -prior to 60 days after the cessation. - - Moreover, your license from a particular copyright holder is -reinstated permanently if the copyright holder notifies you of the -violation by some reasonable means, this is the first time you have -received notice of violation of this License (for any work) from that -copyright holder, and you cure the violation prior to 30 days after -your receipt of the notice. - - Termination of your rights under this section does not terminate the -licenses of parties who have received copies or rights from you under -this License. If your rights have been terminated and not permanently -reinstated, you do not qualify to receive new licenses for the same -material under section 10. - - 9. Acceptance Not Required for Having Copies. - - You are not required to accept this License in order to receive or -run a copy of the Program. Ancillary propagation of a covered work -occurring solely as a consequence of using peer-to-peer transmission -to receive a copy likewise does not require acceptance. However, -nothing other than this License grants you permission to propagate or -modify any covered work. These actions infringe copyright if you do -not accept this License. Therefore, by modifying or propagating a -covered work, you indicate your acceptance of this License to do so. - - 10. Automatic Licensing of Downstream Recipients. - - Each time you convey a covered work, the recipient automatically -receives a license from the original licensors, to run, modify and -propagate that work, subject to this License. You are not responsible -for enforcing compliance by third parties with this License. - - An "entity transaction" is a transaction transferring control of an -organization, or substantially all assets of one, or subdividing an -organization, or merging organizations. If propagation of a covered -work results from an entity transaction, each party to that -transaction who receives a copy of the work also receives whatever -licenses to the work the party's predecessor in interest had or could -give under the previous paragraph, plus a right to possession of the -Corresponding Source of the work from the predecessor in interest, if -the predecessor has it or can get it with reasonable efforts. - - You may not impose any further restrictions on the exercise of the -rights granted or affirmed under this License. For example, you may -not impose a license fee, royalty, or other charge for exercise of -rights granted under this License, and you may not initiate litigation -(including a cross-claim or counterclaim in a lawsuit) alleging that -any patent claim is infringed by making, using, selling, offering for -sale, or importing the Program or any portion of it. - - 11. Patents. - - A "contributor" is a copyright holder who authorizes use under this -License of the Program or a work on which the Program is based. The -work thus licensed is called the contributor's "contributor version". - - A contributor's "essential patent claims" are all patent claims -owned or controlled by the contributor, whether already acquired or -hereafter acquired, that would be infringed by some manner, permitted -by this License, of making, using, or selling its contributor version, -but do not include claims that would be infringed only as a -consequence of further modification of the contributor version. For -purposes of this definition, "control" includes the right to grant -patent sublicenses in a manner consistent with the requirements of -this License. - - Each contributor grants you a non-exclusive, worldwide, royalty-free -patent license under the contributor's essential patent claims, to -make, use, sell, offer for sale, import and otherwise run, modify and -propagate the contents of its contributor version. - - In the following three paragraphs, a "patent license" is any express -agreement or commitment, however denominated, not to enforce a patent -(such as an express permission to practice a patent or covenant not to -sue for patent infringement). To "grant" such a patent license to a -party means to make such an agreement or commitment not to enforce a -patent against the party. - - If you convey a covered work, knowingly relying on a patent license, -and the Corresponding Source of the work is not available for anyone -to copy, free of charge and under the terms of this License, through a -publicly available network server or other readily accessible means, -then you must either (1) cause the Corresponding Source to be so -available, or (2) arrange to deprive yourself of the benefit of the -patent license for this particular work, or (3) arrange, in a manner -consistent with the requirements of this License, to extend the patent -license to downstream recipients. "Knowingly relying" means you have -actual knowledge that, but for the patent license, your conveying the -covered work in a country, or your recipient's use of the covered work -in a country, would infringe one or more identifiable patents in that -country that you have reason to believe are valid. - - If, pursuant to or in connection with a single transaction or -arrangement, you convey, or propagate by procuring conveyance of, a -covered work, and grant a patent license to some of the parties -receiving the covered work authorizing them to use, propagate, modify -or convey a specific copy of the covered work, then the patent license -you grant is automatically extended to all recipients of the covered -work and works based on it. - - A patent license is "discriminatory" if it does not include within -the scope of its coverage, prohibits the exercise of, or is -conditioned on the non-exercise of one or more of the rights that are -specifically granted under this License. You may not convey a covered -work if you are a party to an arrangement with a third party that is -in the business of distributing software, under which you make payment -to the third party based on the extent of your activity of conveying -the work, and under which the third party grants, to any of the -parties who would receive the covered work from you, a discriminatory -patent license (a) in connection with copies of the covered work -conveyed by you (or copies made from those copies), or (b) primarily -for and in connection with specific products or compilations that -contain the covered work, unless you entered into that arrangement, -or that patent license was granted, prior to 28 March 2007. - - Nothing in this License shall be construed as excluding or limiting -any implied license or other defenses to infringement that may -otherwise be available to you under applicable patent law. - - 12. No Surrender of Others' Freedom. - - If conditions are imposed on you (whether by court order, agreement or -otherwise) that contradict the conditions of this License, they do not -excuse you from the conditions of this License. If you cannot convey a -covered work so as to satisfy simultaneously your obligations under this -License and any other pertinent obligations, then as a consequence you may -not convey it at all. For example, if you agree to terms that obligate you -to collect a royalty for further conveying from those to whom you convey -the Program, the only way you could satisfy both those terms and this -License would be to refrain entirely from conveying the Program. - - 13. Use with the GNU Affero General Public License. - - Notwithstanding any other provision of this License, you have -permission to link or combine any covered work with a work licensed -under version 3 of the GNU Affero General Public License into a single -combined work, and to convey the resulting work. The terms of this -License will continue to apply to the part which is the covered work, -but the special requirements of the GNU Affero General Public License, -section 13, concerning interaction through a network will apply to the -combination as such. - - 14. Revised Versions of this License. - - The Free Software Foundation may publish revised and/or new versions of -the GNU General Public License from time to time. Such new versions will -be similar in spirit to the present version, but may differ in detail to -address new problems or concerns. - - Each version is given a distinguishing version number. If the -Program specifies that a certain numbered version of the GNU General -Public License "or any later version" applies to it, you have the -option of following the terms and conditions either of that numbered -version or of any later version published by the Free Software -Foundation. If the Program does not specify a version number of the -GNU General Public License, you may choose any version ever published -by the Free Software Foundation. - - If the Program specifies that a proxy can decide which future -versions of the GNU General Public License can be used, that proxy's -public statement of acceptance of a version permanently authorizes you -to choose that version for the Program. - - Later license versions may give you additional or different -permissions. However, no additional obligations are imposed on any -author or copyright holder as a result of your choosing to follow a -later version. - - 15. Disclaimer of Warranty. - - THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY -APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT -HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY -OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, -THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR -PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM -IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF -ALL NECESSARY SERVICING, REPAIR OR CORRECTION. - - 16. Limitation of Liability. - - IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING -WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS -THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY -GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE -USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF -DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD -PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS), -EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF -SUCH DAMAGES. - - 17. Interpretation of Sections 15 and 16. - - If the disclaimer of warranty and limitation of liability provided -above cannot be given local legal effect according to their terms, -reviewing courts shall apply local law that most closely approximates -an absolute waiver of all civil liability in connection with the -Program, unless a warranty or assumption of liability accompanies a -copy of the Program in return for a fee. - - END OF TERMS AND CONDITIONS - - How to Apply These Terms to Your New Programs - - If you develop a new program, and you want it to be of the greatest -possible use to the public, the best way to achieve this is to make it -free software which everyone can redistribute and change under these terms. - - To do so, attach the following notices to the program. It is safest -to attach them to the start of each source file to most effectively -state the exclusion of warranty; and each file should have at least -the "copyright" line and a pointer to where the full notice is found. - - {one line to give the program's name and a brief idea of what it does.} - Copyright (C) {year} {name of author} - - This program is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - This program is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . - -Also add information on how to contact you by electronic and paper mail. - - If the program does terminal interaction, make it output a short -notice like this when it starts in an interactive mode: - - {project} Copyright (C) {year} {fullname} - This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'. - This is free software, and you are welcome to redistribute it - under certain conditions; type `show c' for details. - -The hypothetical commands `show w' and `show c' should show the appropriate -parts of the General Public License. Of course, your program's commands -might be different; for a GUI interface, you would use an "about box". - - You should also get your employer (if you work as a programmer) or school, -if any, to sign a "copyright disclaimer" for the program, if necessary. -For more information on this, and how to apply and follow the GNU GPL, see -. - - The GNU General Public License does not permit incorporating your program -into proprietary programs. If your program is a subroutine library, you -may consider it more useful to permit linking proprietary applications with -the library. If this is what you want to do, use the GNU Lesser General -Public License instead of this License. But first, please read -. +MIT License + +Copyright 2017-2021 The PyPSA-Eur Authors + +Permission is hereby granted, free of charge, to any person obtaining a copy of +this software and associated documentation files (the "Software"), to deal in +the Software without restriction, including without limitation the rights to +use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of +the Software, and to permit persons to whom the Software is furnished to do so, +subject to the following conditions: + +The above copyright notice and this permission notice shall be included in all +copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS +FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR +COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER +IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN +CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. \ No newline at end of file diff --git a/README.md b/README.md index 31aaa0ef..f4cc34cc 100644 --- a/README.md +++ b/README.md @@ -65,6 +65,6 @@ the additional sectors. # Licence The code in PyPSA-Eur-Sec is released as free software under the -[GPLv3](http://www.gnu.org/licenses/gpl-3.0.en.html), see LICENSE.txt. +[MIT License](https://opensource.org/licenses/MIT), see `LICENSE.txt`. However, different licenses and terms of use may apply to the various input data. diff --git a/Snakefile b/Snakefile index db53f7c8..0b72ae4b 100644 --- a/Snakefile +++ b/Snakefile @@ -1,4 +1,7 @@ +from snakemake.remote.HTTP import RemoteProvider as HTTPRemoteProvider +HTTP = HTTPRemoteProvider() + configfile: "config.yaml" @@ -20,7 +23,6 @@ subworkflow pypsaeur: snakefile: "../pypsa-eur/Snakefile" configfile: "../pypsa-eur/config.yaml" - rule all: input: SDIR + '/graphs/costs.pdf' @@ -156,6 +158,7 @@ rule build_energy_totals: co2="data/eea/UNFCCC_v23.csv", swiss="data/switzerland-sfoe/switzerland-new_format.csv", idees="data/jrc-idees-2015", + district_heat_share='data/district_heat_share.csv', eurostat=input_eurostat output: energy_name='resources/energy_totals.csv', @@ -169,16 +172,37 @@ rule build_energy_totals: rule build_biomass_potentials: input: - jrc_potentials="data/biomass/JRC Biomass Potentials.xlsx" + enspreso_biomass=HTTP.remote("https://cidportal.jrc.ec.europa.eu/ftp/jrc-opendata/ENSPRESO/ENSPRESO_BIOMASS.xlsx", keep_local=True), + nuts2="data/nuts/NUTS_RG_10M_2013_4326_LEVL_2.geojson", # https://gisco-services.ec.europa.eu/distribution/v2/nuts/download/#nuts21 + regions_onshore=pypsaeur("resources/regions_onshore_elec_s{simpl}_{clusters}.geojson"), + nuts3_population="../pypsa-eur/data/bundle/nama_10r_3popgdp.tsv.gz", + swiss_cantons="../pypsa-eur/data/bundle/ch_cantons.csv", + swiss_population="../pypsa-eur/data/bundle/je-e-21.03.02.xls", + country_shapes=pypsaeur('resources/country_shapes.geojson') output: - biomass_potentials_all='resources/biomass_potentials_all.csv', - biomass_potentials='resources/biomass_potentials.csv' + biomass_potentials_all='resources/biomass_potentials_all_s{simpl}_{clusters}.csv', + biomass_potentials='resources/biomass_potentials_s{simpl}_{clusters}.csv' threads: 1 resources: mem_mb=1000 - benchmark: "benchmarks/build_biomass_potentials" + benchmark: "benchmarks/build_biomass_potentials_s{simpl}_{clusters}" script: 'scripts/build_biomass_potentials.py' +if config["sector"]["biomass_transport"]: + rule build_biomass_transport_costs: + input: + transport_cost_data=HTTP.remote("publications.jrc.ec.europa.eu/repository/bitstream/JRC98626/biomass potentials in europe_web rev.pdf", keep_local=True) + output: + biomass_transport_costs="resources/biomass_transport_costs.csv", + threads: 1 + resources: mem_mb=1000 + benchmark: "benchmarks/build_biomass_transport_costs" + script: 'scripts/build_biomass_transport_costs.py' + build_biomass_transport_costs_output = rules.build_biomass_transport_costs.output +else: + build_biomass_transport_costs_output = {} + + rule build_ammonia_production: input: usgs="data/myb1-2017-nitro.xls" @@ -322,7 +346,7 @@ rule prepare_sector_network: transport_name='resources/transport_data.csv', traffic_data_KFZ = "data/emobility/KFZ__count", traffic_data_Pkw = "data/emobility/Pkw__count", - biomass_potentials='resources/biomass_potentials.csv', + biomass_potentials='resources/biomass_potentials_s{simpl}_{clusters}.csv', heat_profile="data/heat_load_profile_BDEW.csv", costs=CDIR + "costs_{planning_horizons}.csv", profile_offwind_ac=pypsaeur("resources/profile_offwind-ac.nc"), @@ -351,7 +375,8 @@ rule prepare_sector_network: solar_thermal_total="resources/solar_thermal_total_elec_s{simpl}_{clusters}.nc", solar_thermal_urban="resources/solar_thermal_urban_elec_s{simpl}_{clusters}.nc", solar_thermal_rural="resources/solar_thermal_rural_elec_s{simpl}_{clusters}.nc", - **build_retro_cost_output + **build_retro_cost_output, + **build_biomass_transport_costs_output output: RDIR + '/prenetworks/elec_s{simpl}_{clusters}_lv{lv}_{opts}_{sector_opts}_{planning_horizons}.nc' threads: 1 resources: mem_mb=2000 diff --git a/config.default.yaml b/config.default.yaml index 05a4a751..fc53e4a0 100644 --- a/config.default.yaml +++ b/config.default.yaml @@ -73,7 +73,7 @@ electricity: # regulate what components with which carriers are kept from PyPSA-Eur; # some technologies are removed because they are implemented differently -# (e.g. battery or H2 storage) or have different year-dependent costs +# (e.g. battery or H2 storage) or have different year-dependent costs # in PyPSA-Eur-Sec pypsa_eur: Bus: @@ -100,28 +100,28 @@ energy: biomass: year: 2030 - scenario: Med + scenario: ENS_Med classes: solid biomass: - - Primary agricultural residues - - Forestry energy residue - - Secondary forestry residues - - Secondary Forestry residues sawdust - - Forestry residues from landscape care biomass + - Agricultural waste + - Fuelwood residues + - Secondary Forestry residues - woodchips + - Sawdust + - Residues from landscape care - Municipal waste not included: - - Bioethanol sugar beet biomass - - Rapeseeds for biodiesel - - sunflower and soya for Biodiesel - - Starchy crops biomass - - Grassy crops biomass - - Willow biomass - - Poplar biomass potential - - Roundwood fuelwood - - Roundwood Chips & Pellets + - Sugar from sugar beet + - Rape seed + - "Sunflower, soya seed " + - Bioethanol barley, wheat, grain maize, oats, other cereals and rye + - Miscanthus, switchgrass, RCG + - Willow + - Poplar + - FuelwoodRW + - C&P_RW biogas: - - Manure biomass potential - - Sludge biomass + - Manure solid, liquid + - Sludge solar_thermal: @@ -142,8 +142,16 @@ existing_capacities: sector: - central: true - central_fraction: 0.6 + district_heating: + potential: 0.6 # maximum fraction of urban demand which can be supplied by district heating + # increase of today's district heating demand to potential maximum district heating share + # progress = 0 means today's district heating share, progress = 1 means maximum fraction of urban demand is supplied by district heating + progress: + 2020: 0.0 + 2030: 0.3 + 2040: 0.6 + 2050: 1.0 + district_heating_loss: 0.15 bev_dsm_restriction_value: 0.75 #Set to 0 for no restriction on BEV DSM bev_dsm_restriction_time: 7 #Time at which SOC of BEV has to be dsm_restriction_value transport_heating_deadband_upper: 20. @@ -152,7 +160,6 @@ sector: ICE_upper_degree_factor: 1.6 EV_lower_degree_factor: 0.98 EV_upper_degree_factor: 0.63 - district_heating_loss: 0.15 bev_dsm: true #turns on EV battery bev_availability: 0.5 #How many cars do smart charging bev_energy: 0.05 #average battery size in MWh @@ -179,7 +186,7 @@ sector: agriculture_machinery_fuel_efficiency: 0.7 # fuel oil per use agriculture_machinery_electric_efficiency: 0.3 # electricity per use shipping_average_efficiency: 0.4 #For conversion of fuel oil to propulsion in 2011 - shipping_hydrogen_liquefaction: true # whether to consider liquefaction costs for shipping H2 demands + shipping_hydrogen_liquefaction: false # whether to consider liquefaction costs for shipping H2 demands shipping_hydrogen_share: # 1 means all hydrogen FC 2020: 0 2025: 0 @@ -227,7 +234,8 @@ sector: co2_vent: true SMR: true co2_sequestration_potential: 200 #MtCO2/a sequestration potential for Europe - co2_sequestration_cost: 20 #EUR/tCO2 for transport and sequestration of CO2 + co2_sequestration_cost: 10 #EUR/tCO2 for sequestration of CO2 + co2_network: false cc_fraction: 0.9 # default fraction of CO2 captured with post-combustion capture hydrogen_underground_storage: true use_fischer_tropsch_waste_heat: true @@ -237,6 +245,7 @@ sector: electricity_grid_connection: true # only applies to onshore wind and utility PV gas_distribution_grid: true gas_distribution_grid_cost_factor: 1.0 #multiplies cost in data/costs.csv + biomass_transport: false # biomass transport between nodes conventional_generation: # generator : carrier OCGT: gas @@ -274,10 +283,23 @@ industry: MWh_elec_per_tNH3_electrolysis: 1.17 # from https://doi.org/10.1016/j.joule.2018.04.017 Table 13 (air separation and HB) NH3_process_emissions: 24.5 # in MtCO2/a from SMR for H2 production for NH3 from UNFCCC for 2015 for EU28 petrochemical_process_emissions: 25.5 # in MtCO2/a for petrochemical and other from UNFCCC for 2015 for EU28 - HVC_primary_fraction: 1.0 #fraction of current non-ammonia basic chemicals produced via primary route + HVC_primary_fraction: 1. # fraction of today's HVC produced via primary route + HVC_mechanical_recycling_fraction: 0. # fraction of today's HVC produced via mechanical recycling + HVC_chemical_recycling_fraction: 0. # fraction of today's HVC produced via chemical recycling + HVC_production_today: 52. # MtHVC/a from DECHEMA (2017), Figure 16, page 107; includes ethylene, propylene and BTX + MWh_elec_per_tHVC_mechanical_recycling: 0.547 # from SI of https://doi.org/10.1016/j.resconrec.2020.105010, Table S5, for HDPE, PP, PS, PET. LDPE would be 0.756. + MWh_elec_per_tHVC_chemical_recycling: 6.9 # Material Economics (2019), page 125; based on pyrolysis and electric steam cracking + chlorine_production_today: 9.58 # MtCl/a from DECHEMA (2017), Table 7, page 43 + MWh_elec_per_tCl: 3.6 # DECHEMA (2017), Table 6, page 43 + MWh_H2_per_tCl: -0.9372 # DECHEMA (2017), page 43; negative since hydrogen produced in chloralkali process + methanol_production_today: 1.5 # MtMeOH/a from DECHEMA (2017), page 62 + MWh_elec_per_tMeOH: 0.167 # DECHEMA (2017), Table 14, page 65 + MWh_CH4_per_tMeOH: 10.25 # DECHEMA (2017), Table 14, page 65 hotmaps_locate_missing: false reference_year: 2015 - + # references: + # DECHEMA (2017): https://dechema.de/dechema_media/Downloads/Positionspapiere/Technology_study_Low_carbon_energy_and_feedstock_for_the_European_chemical_industry-p-20002750.pdf + # Material Economics (2019): https://materialeconomics.com/latest-updates/industrial-transformation-2050 costs: lifetime: 25 #default lifetime @@ -339,7 +361,7 @@ solving: plotting: map: - boundaries: [-11, 30, 34, 71] + boundaries: [-11, 30, 34, 71] color_geomap: ocean: white land: whitesmoke @@ -424,6 +446,7 @@ plotting: lines: k transmission lines: k H2: m + H2 liquefaction: m hydrogen storage: m battery: slategray battery storage: slategray @@ -470,6 +493,7 @@ plotting: hot water storage: '#BBBBBB' hot water charging: '#BBBBBB' hot water discharging: '#999999' + CO2 pipeline: '#999999' CHP: r CHP heat: r CHP electric: r @@ -510,5 +534,6 @@ plotting: shipping oil: "#6495ED" shipping oil emissions: "#6495ED" electricity distribution grid: '#333333' + solid biomass transport: green H2 for industry: "#222222" H2 for shipping: "#6495ED" diff --git a/data/district_heat_share.csv b/data/district_heat_share.csv new file mode 100644 index 00000000..5afd65c8 --- /dev/null +++ b/data/district_heat_share.csv @@ -0,0 +1,34 @@ +country,share to satisfy heat demand (residential) in percent,capacity[MWth] +AT,14,11200 +BG,16,6162 +BA,8, +HR,6.3,2221 +CZ,40, +DK,65, +FI,38,23390 +FR,5, +DE,13.8, +HU,7.92875588637399,8549 +IS,90,8079000 +IE,0.8, +IT,3,8727 +LV,73,2254 +LT,56, +MK,23.7745607009008,636 +NO,4,3400 +PL,42,54912 +PT,0.070754716981132,34 +RS,25,5821 +SI,8.86,1739 +ES,0.251589260787732,1273 +SE,50.4, +UK,2, +BY,70, +EE,52,5406 +KO,3,207 +RO,23,9962 +SK,54,15000 +NL,4,9800 +CH,4,2792 +AL,0, +ME,0, diff --git a/doc/data.csv b/doc/data.csv index 8e316281..01fa04b3 100644 --- a/doc/data.csv +++ b/doc/data.csv @@ -2,6 +2,7 @@ description,file/folder,licence,source JRC IDEES database,jrc-idees-2015/,CC BY 4.0,https://ec.europa.eu/jrc/en/potencia/jrc-idees urban/rural fraction,urban_percent.csv,unknown,unknown JRC biomass potentials,biomass/,unknown,https://doi.org/10.2790/39014 +JRC ENSPRESO biomass potentials,remote,CC BY 4.0,https://data.jrc.ec.europa.eu/dataset/74ed5a04-7d74-4807-9eab-b94774309d9f EEA emission statistics,eea/UNFCCC_v23.csv,EEA standard re-use policy,https://www.eea.europa.eu/data-and-maps/data/national-emissions-reported-to-the-unfccc-and-to-the-eu-greenhouse-gas-monitoring-mechanism-16 Eurostat Energy Balances,eurostat-energy_balances-*/,Eurostat,https://ec.europa.eu/eurostat/web/energy/data/energy-balances Swiss energy statistics from Swiss Federal Office of Energy,switzerland-sfoe/,unknown,http://www.bfe.admin.ch/themen/00526/00541/00542/02167/index.html?dossier_id=02169 @@ -24,3 +25,6 @@ Comparative level investment,comparative_level_investment.csv,Eurostat,https://e Electricity taxes,electricity_taxes_eu.csv,Eurostat,https://appsso.eurostat.ec.europa.eu/nui/show.do?dataset=nrg_pc_204&lang=en Building topologies and corresponding standard values,tabula-calculator-calcsetbuilding.csv,unknown,https://episcope.eu/fileadmin/tabula/public/calc/tabula-calculator.xlsx Retrofitting thermal envelope costs for Germany,retro_cost_germany.csv,unkown,https://www.iwu.de/forschung/handlungslogiken/kosten-energierelevanter-bau-und-anlagenteile-bei-modernisierung/ +District heating most countries,jrc-idees-2015/,CC BY 4.0,https://ec.europa.eu/jrc/en/potencia/jrc-idees,, +District heating missing countries,district_heat_share.csv,unkown,https://www.euroheat.org/knowledge-hub/country-profiles,, + diff --git a/doc/index.rst b/doc/index.rst index 1bf307f5..0c6f929d 100644 --- a/doc/index.rst +++ b/doc/index.rst @@ -134,7 +134,7 @@ it. Licence ======= -The code in PyPSA-Eur-Sec is released as free software under the `GPLv3 -`_, see +The code in PyPSA-Eur-Sec is released as free software under the +`MIT license `_, see `LICENSE `_. However, different licenses and terms of use may apply to the various input data. diff --git a/doc/release_notes.rst b/doc/release_notes.rst index 0c65c8c2..8a1ed562 100644 --- a/doc/release_notes.rst +++ b/doc/release_notes.rst @@ -8,6 +8,8 @@ Future release .. note:: This unreleased version currently requires the master branches of PyPSA, PyPSA-Eur, and the technology-data repository. +* With this release, we change the license from copyleft GPLv3 to the more + liberal MIT license with the consent of all contributors. * Extended use of ``multiprocessing`` for much better performance (from up to 20 minutes to less than one minute). * Compatibility with ``atlite>=0.2``. Older versions of ``atlite`` will no longer work. @@ -60,17 +62,38 @@ Future release These are included in the environment specifications of PyPSA-Eur. * Consistent use of ``__main__`` block and further unspecific code cleaning. * Distinguish costs for home battery storage and inverter from utility-scale battery costs. +* Add option to regionally resolve CO2 storage and add CO2 pipeline transport because geological storage potential, + CO2 utilisation sites and CO2 capture sites may be separated. + The CO2 network is built from zero based on the topology of the electricity grid (greenfield). + Pipelines are assumed to be bidirectional and lossless. + Furthermore, neither retrofitting of natural gas pipelines (required pressures are too high, 80-160 bar vs <80 bar) + nor other modes of CO2 transport (by ship, road or rail) are considered. + The regional representation of CO2 is activated with the config setting ``sector: co2_network: true`` but is deactivated by default. + The global limit for CO2 sequestration now applies to the sum of all CO2 stores via an ``extra_functionality`` constraint. * Added option for hydrogen liquefaction costs for hydrogen demand in shipping. This introduces a new ``H2 liquid`` bus at each location. It is activated via ``sector: shipping_hydrogen_liquefaction: true``. * The share of shipping transformed into hydrogen fuel cell can be now defined for different years in the ``config.yaml`` file. The carbon emission from the remaining share is treated as a negative load on the atmospheric carbon dioxide bus, just like aviation and land transport emissions. * The transformation of the Steel and Aluminium production can be now defined for different years in the ``config.yaml`` file. * Include the option to alter the maximum energy capacity of a store via the ``carrier+factor`` in the ``{sector_opts}`` wildcard. This can be useful for sensitivity analyses. Example: ``co2 stored+e2`` multiplies the ``e_nom_max`` by factor 2. In this example, ``e_nom_max`` represents the CO2 sequestration potential in Europe. +* Add option to regionally disaggregate biomass potential to individual nodes + (currently given per country, then distributed by population density within) + and allow the transport of solid biomass. + The transport costs are determined based on the `JRC-EU-Times Bioenergy report `_ + in the new optional rule ``build_biomass_transport_costs``. + Biomass transport can be activated with the setting ``sector: biomass_transport: true``. +* Use `JRC ENSPRESO database `_ to + spatially disaggregate biomass potentials to PyPSA-Eur regions based on overlaps with NUTS2 regions from ENSPRESO + (proportional to area) (`#151 `_). * Compatibility with ``xarray`` version 0.19. * Added option to include emissions and energy demands of agriculture, forestry and fishing sector via the letter ``A`` in the ``{sector_opts}`` wildcard. Demands are separated into electricity, heat and oil for machinery. Fuel-switching for machinery from oil to electricity can be set exogenously in the ``config.yaml`` `#147 `_. +* Separate basic chemicals into HVC, chlorine, methanol and ammonia [`#166 `_]. +* Add option to specify reuse, primary production, and mechanical and chemical recycling fraction of platics [`#166 `_]. +* Include today's district heating shares in myopic optimisation and add option to specify exogenous path for district heating share increase under ``sector: district_heating:`` [`#149 `_]. +* The myopic option can now be used together with different clustering for the generators and the network. The existing renewable capacities are split evenly among the regions in every country [`#144 `_]. PyPSA-Eur-Sec 0.5.0 (21st May 2021) =================================== diff --git a/doc/spatial_resolution.rst b/doc/spatial_resolution.rst index 1be9f3ad..83a33f73 100644 --- a/doc/spatial_resolution.rst +++ b/doc/spatial_resolution.rst @@ -44,11 +44,13 @@ Hydrogen network: nodal. Methane network: single node for Europe, since future demand is so low and no bottlenecks are expected. -Solid biomass: single node for Europe, until transport costs can be -incorporated. +Solid biomass: choice between single node for Europe and nodal where biomass +potential is regionally disaggregated (currently given per country, +then distributed by population density within) +and transport of solid biomass is possible. CO2: single node for Europe, but a transport and storage cost is added for -sequestered CO2. +sequestered CO2. Optionally: nodal, with CO2 transport via pipelines. Liquid hydrocarbons: single node for Europe, since transport costs for liquids are low. diff --git a/doc/supply_demand.rst b/doc/supply_demand.rst index 31563c43..3ab00d8e 100644 --- a/doc/supply_demand.rst +++ b/doc/supply_demand.rst @@ -183,7 +183,7 @@ Solid biomass provides process heat up to 500 Celsius in industry, as well as fe Solid biomass supply ===================== -Only wastes and residues from the JRC biomass dataset. +Only wastes and residues from the JRC ENSPRESO biomass dataset. Oil product demand diff --git a/scripts/add_existing_baseyear.py b/scripts/add_existing_baseyear.py index 47c31c0e..bb35e378 100644 --- a/scripts/add_existing_baseyear.py +++ b/scripts/add_existing_baseyear.py @@ -28,7 +28,7 @@ def add_build_year_to_new_assets(n, baseyear): # Give assets with lifetimes and no build year the build year baseyear for c in n.iterate_components(["Link", "Generator", "Store"]): - assets = c.df.index[~c.df.lifetime.isna() & c.df.build_year.isna()] + assets = c.df.index[~c.df.lifetime.isna() & c.df.build_year==0] c.df.loc[assets, "build_year"] = baseyear # add -baseyear to name @@ -60,7 +60,7 @@ def add_existing_renewables(df_agg): } for tech in ['solar', 'onwind', 'offwind']: - + carrier = carriers[tech] df = pd.read_csv(snakemake.input[f"existing_{tech}"], index_col=0).fillna(0.) @@ -112,9 +112,9 @@ def add_power_capacities_installed_before_baseyear(n, grouping_years, costs, bas Parameters ---------- n : pypsa.Network - grouping_years : + grouping_years : intervals to group existing capacities - costs : + costs : to read lifetime to estimate YearDecomissioning baseyear : int """ @@ -155,6 +155,11 @@ def add_power_capacities_installed_before_baseyear(n, grouping_years, costs, bas # assign clustered bus busmap_s = pd.read_csv(snakemake.input.busmap_s, index_col=0, squeeze=True) busmap = pd.read_csv(snakemake.input.busmap, index_col=0, squeeze=True) + + inv_busmap = {} + for k, v in busmap.iteritems(): + inv_busmap[v] = inv_busmap.get(v, []) + [k] + clustermaps = busmap_s.map(busmap) clustermaps.index = clustermaps.index.astype(int) @@ -192,24 +197,54 @@ def add_power_capacities_installed_before_baseyear(n, grouping_years, costs, bas capacity = capacity[capacity > snakemake.config['existing_capacities']['threshold_capacity']] if generator in ['solar', 'onwind', 'offwind']: - - rename = {"offwind": "offwind-ac"} - p_max_pu=n.generators_t.p_max_pu[capacity.index + ' ' + rename.get(generator, generator) + '-' + str(baseyear)] - - n.madd("Generator", - capacity.index, - suffix=' ' + generator +"-"+ str(grouping_year), - bus=capacity.index, - carrier=generator, - p_nom=capacity, - marginal_cost=costs.at[generator, 'VOM'], - capital_cost=costs.at[generator, 'fixed'], - efficiency=costs.at[generator, 'efficiency'], - p_max_pu=p_max_pu.rename(columns=n.generators.bus), - build_year=grouping_year, - lifetime=costs.at[generator, 'lifetime'] - ) + suffix = '-ac' if generator == 'offwind' else '' + name_suffix = f' {generator}{suffix}-{baseyear}' + + if 'm' in snakemake.wildcards.clusters: + + for ind in capacity.index: + + # existing capacities are split evenly among regions in every country + inv_ind = [i for i in inv_busmap[ind]] + + # for offshore the spliting only inludes coastal regions + inv_ind = [i for i in inv_ind if (i + name_suffix) in n.generators.index] + + p_max_pu = n.generators_t.p_max_pu[[i + name_suffix for i in inv_ind]] + p_max_pu.columns=[i + name_suffix for i in inv_ind ] + + n.madd("Generator", + [i + name_suffix for i in inv_ind], + bus=ind, + carrier=generator, + p_nom=capacity[ind] / len(inv_ind), # split among regions in a country + marginal_cost=costs.at[generator,'VOM'], + capital_cost=costs.at[generator,'fixed'], + efficiency=costs.at[generator, 'efficiency'], + p_max_pu=p_max_pu, + build_year=grouping_year, + lifetime=costs.at[generator,'lifetime'] + ) + + else: + + p_max_pu = n.generators_t.p_max_pu[capacity.index + name_suffix] + + n.madd("Generator", + capacity.index, + suffix=' ' + generator +"-"+ str(grouping_year), + bus=capacity.index, + carrier=generator, + p_nom=capacity, + marginal_cost=costs.at[generator, 'VOM'], + capital_cost=costs.at[generator, 'fixed'], + efficiency=costs.at[generator, 'efficiency'], + p_max_pu=p_max_pu.rename(columns=n.generators.bus), + build_year=grouping_year, + lifetime=costs.at[generator, 'lifetime'] + ) + else: n.madd("Link", @@ -268,7 +303,7 @@ def add_heating_capacities_installed_before_baseyear(n, baseyear, grouping_years df.fillna(0., inplace=True) # convert GW to MW - df *= 1e3 + df *= 1e3 cc = pd.read_csv(snakemake.input.country_codes, index_col=0) @@ -327,7 +362,7 @@ def add_heating_capacities_installed_before_baseyear(n, baseyear, grouping_years efficiency = cop[heat_pump_type][nodes[name]] else: efficiency = costs.at[costs_name, 'efficiency'] - + for i, grouping_year in enumerate(grouping_years): if int(grouping_year) + default_lifetime <= int(baseyear): @@ -378,7 +413,7 @@ def add_heating_capacities_installed_before_baseyear(n, baseyear, grouping_years build_year=int(grouping_year), lifetime=costs.at[name_type + ' gas boiler', 'lifetime'] ) - + n.madd("Link", nodes[name], suffix=f" {name} oil boiler-{grouping_year}", @@ -410,7 +445,8 @@ if __name__ == "__main__": simpl='', clusters=45, lv=1.0, - sector_opts='Co2L0-168H-T-H-B-I-solar3-dist1', + opts='', + sector_opts='Co2L0-168H-T-H-B-I-solar+p3-dist1', planning_horizons=2020, ) diff --git a/scripts/build_biomass_potentials.py b/scripts/build_biomass_potentials.py index f02c9093..68d87808 100644 --- a/scripts/build_biomass_potentials.py +++ b/scripts/build_biomass_potentials.py @@ -1,55 +1,194 @@ import pandas as pd - -rename = {"UK" : "GB", "BH" : "BA"} +import geopandas as gpd -def build_biomass_potentials(): +def build_nuts_population_data(year=2013): - config = snakemake.config['biomass'] - year = config["year"] - scenario = config["scenario"] + pop = pd.read_csv( + snakemake.input.nuts3_population, + sep=r'\,| \t|\t', + engine='python', + na_values=[":"], + index_col=1 + )[str(year)] + + # only countries + pop.drop("EU28", inplace=True) - df = pd.read_excel(snakemake.input.jrc_potentials, - "Potentials (PJ)", - index_col=[0,1]) + # mapping from Cantons to NUTS3 + cantons = pd.read_csv(snakemake.input.swiss_cantons) + cantons = cantons.set_index(cantons.HASC.str[3:]).NUTS + cantons = cantons.str.pad(5, side='right', fillchar='0') - df.rename(columns={"Unnamed: 18": "Municipal waste"}, inplace=True) - df.drop(columns="Total", inplace=True) - df.replace("-", 0., inplace=True) + # get population by NUTS3 + swiss = pd.read_excel(snakemake.input.swiss_population, skiprows=3, index_col=0).loc["Residents in 1000"] + swiss = swiss.rename(cantons).filter(like="CH") - column = df.iloc[:,0] - countries = column.where(column.str.isalpha()).pad() - countries = [rename.get(ct, ct) for ct in countries] - countries_i = pd.Index(countries, name='country') - df.set_index(countries_i, append=True, inplace=True) + # aggregate also to higher order NUTS levels + swiss = [swiss.groupby(swiss.index.str[:i]).sum() for i in range(2, 6)] - df.drop(index='MS', level=0, inplace=True) + # merge Europe + Switzerland + pop = pd.DataFrame(pop.append(swiss), columns=["total"]) + + # add missing manually + pop["AL"] = 2893 + pop["BA"] = 3871 + pop["RS"] = 7210 + + pop["ct"] = pop.index.str[:2] + + return pop - # convert from PJ to MWh - df = df / 3.6 * 1e6 - df.to_csv(snakemake.output.biomass_potentials_all) +def enspreso_biomass_potentials(year=2020, scenario="ENS_Low"): + """ + Loads the JRC ENSPRESO biomass potentials. + + Parameters + ---------- + year : int + The year for which potentials are to be taken. + Can be {2010, 2020, 2030, 2040, 2050}. + scenario : str + The scenario. Can be {"ENS_Low", "ENS_Med", "ENS_High"}. + + Returns + ------- + pd.DataFrame + Biomass potentials for given year and scenario + in TWh/a by commodity and NUTS2 region. + """ - # solid biomass includes: - # Primary agricultural residues (MINBIOAGRW1), - # Forestry energy residue (MINBIOFRSF1), - # Secondary forestry residues (MINBIOWOOW1), - # Secondary Forestry residues – sawdust (MINBIOWOO1a)', - # Forestry residues from landscape care biomass (MINBIOFRSF1a), - # Municipal waste (MINBIOMUN1)', + glossary = pd.read_excel( + str(snakemake.input.enspreso_biomass), + sheet_name="Glossary", + usecols="B:D", + skiprows=1, + index_col=0 + ) + + df = pd.read_excel( + str(snakemake.input.enspreso_biomass), + sheet_name="ENER - NUTS2 BioCom E", + usecols="A:H" + ) - # biogas includes: - # Manure biomass potential (MINBIOGAS1), - # Sludge biomass (MINBIOSLU1), + df["group"] = df["E-Comm"].map(glossary.group) + df["commodity"] = df["E-Comm"].map(glossary.description) - df = df.loc[year, scenario, :] + to_rename = { + "NUTS2 Potential available by Bio Commodity": "potential", + "NUST2": "NUTS2", + } + df.rename(columns=to_rename, inplace=True) + + # fill up with NUTS0 if NUTS2 is not given + df.NUTS2 = df.apply(lambda x: x.NUTS0 if x.NUTS2 == '-' else x.NUTS2, axis=1) - grouper = {v: k for k, vv in config["classes"].items() for v in vv} - df = df.groupby(grouper, axis=1).sum() + # convert PJ to TWh + df.potential /= 3.6 + df.Unit = "TWh/a" - df.index.name = "MWh/a" + dff = df.query("Year == @year and Scenario == @scenario") - df.to_csv(snakemake.output.biomass_potentials) + bio = dff.groupby(["NUTS2", "commodity"]).potential.sum().unstack() + + # currently Serbia and Kosovo not split, so aggregate + bio.loc["RS"] += bio.loc["XK"] + bio.drop("XK", inplace=True) + + return bio + + +def disaggregate_nuts0(bio): + """ + Some commodities are only given on NUTS0 level. + These are disaggregated here using the NUTS2 + population as distribution key. + + Parameters + ---------- + bio : pd.DataFrame + from enspreso_biomass_potentials() + + Returns + ------- + pd.DataFrame + """ + + pop = build_nuts_population_data() + + # get population in nuts2 + pop_nuts2 = pop.loc[pop.index.str.len() == 4] + by_country = pop_nuts2.total.groupby(pop_nuts2.ct).sum() + pop_nuts2["fraction"] = pop_nuts2.total / pop_nuts2.ct.map(by_country) + + # distribute nuts0 data to nuts2 by population + bio_nodal = bio.loc[pop_nuts2.ct] + bio_nodal.index = pop_nuts2.index + bio_nodal = bio_nodal.mul(pop_nuts2.fraction, axis=0) + + # update inplace + bio.update(bio_nodal) + + return bio + + +def build_nuts2_shapes(): + """ + - load NUTS2 geometries + - add RS, AL, BA country shapes (not covered in NUTS 2013) + - consistently name ME, MK + """ + + nuts2 = gpd.GeoDataFrame(gpd.read_file(snakemake.input.nuts2).set_index('id').geometry) + + countries = gpd.read_file(snakemake.input.country_shapes).set_index('name') + missing = countries.loc[["AL", "RS", "BA"]] + nuts2.rename(index={"ME00": "ME", "MK00": "MK"}, inplace=True) + + return nuts2.append(missing) + + +def area(gdf): + """Returns area of GeoDataFrame geometries in square kilometers.""" + return gdf.to_crs(epsg=3035).area.div(1e6) + + +def convert_nuts2_to_regions(bio_nuts2, regions): + """ + Converts biomass potentials given in NUTS2 to PyPSA-Eur regions based on the + overlay of both GeoDataFrames in proportion to the area. + + Parameters + ---------- + bio_nuts2 : gpd.GeoDataFrame + JRC ENSPRESO biomass potentials indexed by NUTS2 shapes. + regions : gpd.GeoDataFrame + PyPSA-Eur clustered onshore regions + + Returns + ------- + gpd.GeoDataFrame + """ + + # calculate area of nuts2 regions + bio_nuts2["area_nuts2"] = area(bio_nuts2) + + overlay = gpd.overlay(regions, bio_nuts2) + + # calculate share of nuts2 area inside region + overlay["share"] = area(overlay) / overlay["area_nuts2"] + + # multiply all nuts2-level values with share of nuts2 inside region + adjust_cols = overlay.columns.difference({"name", "area_nuts2", "geometry", "share"}) + overlay[adjust_cols] = overlay[adjust_cols].multiply(overlay["share"], axis=0) + + bio_regions = overlay.groupby("name").sum() + + bio_regions.drop(["area_nuts2", "share"], axis=1, inplace=True) + + return bio_regions if __name__ == "__main__": @@ -57,12 +196,28 @@ if __name__ == "__main__": from helper import mock_snakemake snakemake = mock_snakemake('build_biomass_potentials') + config = snakemake.config['biomass'] + year = config["year"] + scenario = config["scenario"] - # This is a hack, to be replaced once snakemake is unicode-conform + enspreso = enspreso_biomass_potentials(year, scenario) - solid_biomass = snakemake.config['biomass']['classes']['solid biomass'] - if 'Secondary Forestry residues sawdust' in solid_biomass: - solid_biomass.remove('Secondary Forestry residues sawdust') - solid_biomass.append('Secondary Forestry residues – sawdust') + enspreso = disaggregate_nuts0(enspreso) - build_biomass_potentials() + nuts2 = build_nuts2_shapes() + + df_nuts2 = gpd.GeoDataFrame(nuts2.geometry).join(enspreso) + + regions = gpd.read_file(snakemake.input.regions_onshore) + + df = convert_nuts2_to_regions(df_nuts2, regions) + + df.to_csv(snakemake.output.biomass_potentials_all) + + grouper = {v: k for k, vv in config["classes"].items() for v in vv} + df = df.groupby(grouper, axis=1).sum() + + df *= 1e6 # TWh/a to MWh/a + df.index.name = "MWh/a" + + df.to_csv(snakemake.output.biomass_potentials) diff --git a/scripts/build_biomass_transport_costs.py b/scripts/build_biomass_transport_costs.py new file mode 100644 index 00000000..aaec215b --- /dev/null +++ b/scripts/build_biomass_transport_costs.py @@ -0,0 +1,90 @@ +""" +Reads biomass transport costs for different countries of the JRC report + + "The JRC-EU-TIMES model. + Bioenergy potentials + for EU and neighbouring countries." + (2015) + +converts them from units 'EUR per km/ton' -> 'EUR/ (km MWh)' + +assuming as an approximation energy content of wood pellets + +@author: bw0928 +""" + +import pandas as pd +import tabula as tbl + +ENERGY_CONTENT = 4.8 # unit MWh/t (wood pellets) + +def get_countries(): + + pandas_options = dict( + skiprows=range(6), + header=None, + index_col=0 + ) + + return tbl.read_pdf( + str(snakemake.input.transport_cost_data), + pages="145", + multiple_tables=False, + pandas_options=pandas_options + )[0].index + + +def get_cost_per_tkm(page, countries): + + pandas_options = dict( + skiprows=range(6), + header=0, + sep=' |,', + engine='python', + index_col=False, + ) + + sc = tbl.read_pdf( + str(snakemake.input.transport_cost_data), + pages=page, + multiple_tables=False, + pandas_options=pandas_options + )[0] + sc.index = countries + sc.columns = sc.columns.str.replace("€", "EUR") + + return sc + + +def build_biomass_transport_costs(): + + countries = get_countries() + + sc1 = get_cost_per_tkm(146, countries) + sc2 = get_cost_per_tkm(147, countries) + + # take mean of both supply chains + to_concat = [sc1["EUR/km/ton"], sc2["EUR/km/ton"]] + transport_costs = pd.concat(to_concat, axis=1).mean(axis=1) + + # convert tonnes to MWh + transport_costs /= ENERGY_CONTENT + transport_costs.name = "EUR/km/MWh" + + # rename country names + to_rename = { + "UK": "GB", + "XK": "KO", + "EL": "GR" + } + transport_costs.rename(to_rename, inplace=True) + + # add missing Norway with data from Sweden + transport_costs["NO"] = transport_costs["SE"] + + transport_costs.to_csv(snakemake.output[0]) + + +if __name__ == "__main__": + + build_biomass_transport_costs() diff --git a/scripts/build_energy_totals.py b/scripts/build_energy_totals.py index 4b4dc2bc..b21200b5 100644 --- a/scripts/build_energy_totals.py +++ b/scripts/build_energy_totals.py @@ -213,6 +213,12 @@ def idees_per_country(ct, year): assert df.index[47] == "Electricity" ct_totals["electricity residential"] = df[47] + assert df.index[46] == "Derived heat" + ct_totals["Derived heat residential"] = df[46] + + assert df.index[50] == 'Thermal uses' + ct_totals["thermal uses residential"] = df[50] + # services df = pd.read_excel(fn_tertiary, "SER_hh_fec", index_col=0)[year] @@ -240,6 +246,13 @@ def idees_per_country(ct, year): assert df.index[50] == "Electricity" ct_totals["electricity services"] = df[50] + assert df.index[49] == "Derived heat" + ct_totals["derived heat services"] = df[49] + + assert df.index[53] == 'Thermal uses' + ct_totals["thermal uses services"] = df[53] + + # agriculture, forestry and fishing start = "Detailed split of energy consumption (ktoe)" @@ -371,6 +384,7 @@ def build_idees(countries, year): with mp.Pool(processes=nprocesses) as pool: totals_list = list(tqdm(pool.imap(func, countries), **tqdm_kwargs)) + totals = pd.concat(totals_list, axis=1) # convert ktoe to TWh @@ -380,6 +394,13 @@ def build_idees(countries, year): # convert TWh/100km to kWh/km totals.loc["passenger car efficiency"] *= 10 + # district heating share + district_heat = totals.loc[["derived heat residential", + "derived heat services"]].sum() + total_heat = totals.loc[["thermal uses residential", + "thermal uses services"]].sum() + totals.loc["district heat share"] = district_heat.div(total_heat) + return totals.T @@ -522,7 +543,7 @@ def build_energy_totals(countries, eurostat, swiss, idees): for purpose in ["passenger", "freight"]: attrs = [f"total domestic aviation {purpose}", f"total international aviation {purpose}"] - df.loc[missing, f"total aviation {purpose}"] = df.loc[missing, attrs].sum(axis=1) + df.loc[missing, f"total aviation {purpose}"] = df.loc[missing, attrs].sum(axis=1) if "BA" in df.index: # fill missing data for BA (services and road energy data) @@ -531,6 +552,14 @@ def build_energy_totals(countries, eurostat, swiss, idees): ratio = df.at["BA", "total residential"] / df.at["RS", "total residential"] df.loc['BA', missing] = ratio * df.loc["RS", missing] + # Missing district heating share + dh_share = pd.read_csv(snakemake.input.district_heat_share, + index_col=0, usecols=[0, 1]) + # make conservative assumption and take minimum from both data sets + df["district heat share"] = (pd.concat([df["district heat share"], + dh_share.reindex(index=df.index)/100], + axis=1).min(axis=1)) + return df diff --git a/scripts/build_industrial_energy_demand_per_country_today.py b/scripts/build_industrial_energy_demand_per_country_today.py index 1d906b24..0adf84e7 100644 --- a/scripts/build_industrial_energy_demand_per_country_today.py +++ b/scripts/build_industrial_energy_demand_per_country_today.py @@ -103,6 +103,7 @@ def add_ammonia_energy_demand(demand): demand['Basic chemicals (without ammonia)'] = demand["Basic chemicals"] - demand["Ammonia"] demand['Basic chemicals (without ammonia)'].clip(lower=0, inplace=True) + demand.drop(columns='Basic chemicals', inplace=True) return demand @@ -114,6 +115,11 @@ def add_non_eu28_industrial_energy_demand(demand): fn = snakemake.input.industrial_production_per_country production = pd.read_csv(fn, index_col=0) / 1e3 + #recombine HVC, Chlorine and Methanol to Basic chemicals (without ammonia) + chemicals = ["HVC", "Chlorine", "Methanol"] + production["Basic chemicals (without ammonia)"] = production[chemicals].sum(axis=1) + production.drop(columns=chemicals, inplace=True) + eu28_production = production.loc[eu28].sum() eu28_energy = demand.groupby(level=1).sum() eu28_averages = eu28_energy / eu28_production diff --git a/scripts/build_industrial_production_per_country.py b/scripts/build_industrial_production_per_country.py index 1754752a..eadfb224 100644 --- a/scripts/build_industrial_production_per_country.py +++ b/scripts/build_industrial_production_per_country.py @@ -179,8 +179,8 @@ def industry_production(countries): return demand -def add_ammonia_demand_separately(demand): - """Include ammonia demand separately and remove ammonia from basic chemicals.""" +def separate_basic_chemicals(demand): + """Separate basic chemicals into ammonia, chlorine, methanol and HVC.""" ammonia = pd.read_csv(snakemake.input.ammonia_production, index_col=0) @@ -189,7 +189,7 @@ def add_ammonia_demand_separately(demand): print("Following countries have no ammonia demand:", missing) - demand.insert(2, "Ammonia", 0.) + demand["Ammonia"] = 0. demand.loc[there, "Ammonia"] = ammonia.loc[there, str(year)] @@ -198,9 +198,13 @@ def add_ammonia_demand_separately(demand): # EE, HR and LT got negative demand through subtraction - poor data demand['Basic chemicals'].clip(lower=0., inplace=True) - to_rename = {"Basic chemicals": "Basic chemicals (without ammonia)"} - demand.rename(columns=to_rename, inplace=True) + # assume HVC, methanol, chlorine production proportional to non-ammonia basic chemicals + distribution_key = demand["Basic chemicals"] / demand["Basic chemicals"].sum() + demand["HVC"] = config["HVC_production_today"] * 1e3 * distribution_key + demand["Chlorine"] = config["chlorine_production_today"] * 1e3 * distribution_key + demand["Methanol"] = config["methanol_production_today"] * 1e3 * distribution_key + demand.drop(columns=["Basic chemicals"], inplace=True) if __name__ == '__main__': if 'snakemake' not in globals(): @@ -211,12 +215,14 @@ if __name__ == '__main__': year = snakemake.config['industry']['reference_year'] + config = snakemake.config["industry"] + jrc_dir = snakemake.input.jrc eurostat_dir = snakemake.input.eurostat demand = industry_production(countries) - add_ammonia_demand_separately(demand) + separate_basic_chemicals(demand) fn = snakemake.output.industrial_production_per_country demand.to_csv(fn, float_format='%.2f') diff --git a/scripts/build_industrial_production_per_country_tomorrow.py b/scripts/build_industrial_production_per_country_tomorrow.py index ba69e0a6..ccf31839 100644 --- a/scripts/build_industrial_production_per_country_tomorrow.py +++ b/scripts/build_industrial_production_per_country_tomorrow.py @@ -39,11 +39,14 @@ if __name__ == '__main__': al_primary_fraction = get(config["Al_primary_fraction"], investment_year) fraction_persistent_primary = al_primary_fraction * total_aluminium.sum() / production[key_pri].sum() - + production[key_pri] = fraction_persistent_primary * production[key_pri] production[key_sec] = total_aluminium - production[key_pri] - production["Basic chemicals (without ammonia)"] *= config['HVC_primary_fraction'] + production["HVC (mechanical recycling)"] = get(config["HVC_mechanical_recycling_fraction"], investment_year) * production["HVC"] + production["HVC (chemical recycling)"] = get(config["HVC_chemical_recycling_fraction"], investment_year) * production["HVC"] + + production["HVC"] *= get(config['HVC_primary_fraction'], investment_year) fn = snakemake.output.industrial_production_per_country_tomorrow production.to_csv(fn, float_format='%.2f') diff --git a/scripts/build_industrial_production_per_node.py b/scripts/build_industrial_production_per_node.py index b5361e6b..4ceffee9 100644 --- a/scripts/build_industrial_production_per_node.py +++ b/scripts/build_industrial_production_per_node.py @@ -9,7 +9,11 @@ sector_mapping = { 'Integrated steelworks': 'Iron and steel', 'DRI + Electric arc': 'Iron and steel', 'Ammonia': 'Chemical industry', - 'Basic chemicals (without ammonia)': 'Chemical industry', + 'HVC': 'Chemical industry', + 'HVC (mechanical recycling)': 'Chemical industry', + 'HVC (chemical recycling)': 'Chemical industry', + 'Methanol': 'Chemical industry', + 'Chlorine': 'Chemical industry', 'Other chemicals': 'Chemical industry', 'Pharmaceutical products etc.': 'Chemical industry', 'Cement': 'Cement', @@ -40,12 +44,12 @@ def build_nodal_industrial_production(): countries = keys.country.unique() sectors = industrial_production.columns - + for country, sector in product(countries, sectors): buses = keys.index[keys.country == country] mapping = sector_mapping.get(sector, "population") - + key = keys.loc[buses, mapping] nodal_production.loc[buses, sector] = industrial_production.at[country, sector] * key diff --git a/scripts/build_industry_sector_ratios.py b/scripts/build_industry_sector_ratios.py index adfb1d3c..49c82138 100644 --- a/scripts/build_industry_sector_ratios.py +++ b/scripts/build_industry_sector_ratios.py @@ -279,7 +279,7 @@ def chemicals_industry(): df = pd.DataFrame(index=index) - # Basid chemicals + # Basic chemicals sector = "Basic chemicals" @@ -374,52 +374,82 @@ def chemicals_industry(): # putting in ammonia demand for H2 and electricity separately s_emi = idees["emi"][3:57] - s_out = idees["out"][8:9] assert s_emi.index[0] == sector - assert sector in str(s_out.index) - ammonia = pd.read_csv(snakemake.input.ammonia_production, index_col=0) - - # ktNH3/a - ammonia_total = ammonia.loc[ammonia.index.intersection(eu28), str(year)].sum() - - s_out -= ammonia_total + # convert from MtHVC/a to ktHVC/a + s_out = config["HVC_production_today"] * 1e3 # tCO2/t material df.loc["process emission", sector] += ( s_emi["Process emissions"] - config["petrochemical_process_emissions"] * 1e3 - config["NH3_process_emissions"] * 1e3 - ) / s_out.values + ) / s_out # emissions originating from feedstock, could be non-fossil origin # tCO2/t material df.loc["process emission from feedstock", sector] += ( config["petrochemical_process_emissions"] * 1e3 - ) / s_out.values + ) / s_out # convert from ktoe/a to GWh/a sources = ["elec", "biomass", "methane", "hydrogen", "heat", "naphtha"] df.loc[sources, sector] *= toe_to_MWh + # subtract ammonia energy demand (in ktNH3/a) + ammonia = pd.read_csv(snakemake.input.ammonia_production, index_col=0) + ammonia_total = ammonia.loc[ammonia.index.intersection(eu28), str(year)].sum() df.loc["methane", sector] -= ammonia_total * config["MWh_CH4_per_tNH3_SMR"] df.loc["elec", sector] -= ammonia_total * config["MWh_elec_per_tNH3_SMR"] - # MWh/t material - df.loc[sources, sector] = df.loc[sources, sector] / s_out.values + # subtract chlorine demand + chlorine_total = config["chlorine_production_today"] + df.loc["hydrogen", sector] -= chlorine_total * config["MWh_H2_per_tCl"] + df.loc["elec", sector] -= chlorine_total * config["MWh_elec_per_tCl"] - to_rename = {sector: f"{sector} (without ammonia)"} - df.rename(columns=to_rename, inplace=True) + # subtract methanol demand + methanol_total = config["methanol_production_today"] + df.loc["methane", sector] -= methanol_total * config["MWh_CH4_per_tMeOH"] + df.loc["elec", sector] -= methanol_total * config["MWh_elec_per_tMeOH"] + + # MWh/t material + df.loc[sources, sector] = df.loc[sources, sector] / s_out + + df.rename(columns={sector: "HVC"}, inplace=True) + + # HVC mechanical recycling + + sector = "HVC (mechanical recycling)" + df[sector] = 0.0 + df.loc["elec", sector] = config["MWh_elec_per_tHVC_mechanical_recycling"] + + # HVC chemical recycling + + sector = "HVC (chemical recycling)" + df[sector] = 0.0 + df.loc["elec", sector] = config["MWh_elec_per_tHVC_chemical_recycling"] # Ammonia sector = "Ammonia" - df[sector] = 0.0 - df.loc["hydrogen", sector] = config["MWh_H2_per_tNH3_electrolysis"] df.loc["elec", sector] = config["MWh_elec_per_tNH3_electrolysis"] + # Chlorine + + sector = "Chlorine" + df[sector] = 0.0 + df.loc["hydrogen", sector] = config["MWh_H2_per_tCl"] + df.loc["elec", sector] = config["MWh_elec_per_tCl"] + + # Methanol + + sector = "Methanol" + df[sector] = 0.0 + df.loc["methane", sector] = config["MWh_CH4_per_tMeOH"] + df.loc["elec", sector] = config["MWh_elec_per_tMeOH"] + # Other chemicals sector = "Other chemicals" diff --git a/scripts/plot_network.py b/scripts/plot_network.py index cd74d3ea..5182123b 100644 --- a/scripts/plot_network.py +++ b/scripts/plot_network.py @@ -289,7 +289,7 @@ def plot_h2_map(network): title='Electrolyzer capacity', handler_map=make_handler_map_to_scale_circles_as_in(ax) ) - + ax.add_artist(l2) handles = [] @@ -398,7 +398,8 @@ def plot_series(network, carrier="AC", name="test"): supply = pd.DataFrame(index=n.snapshots) for c in n.iterate_components(n.branch_components): - for i in range(2): + n_port = 4 if c.name=='Link' else 2 + for i in range(n_port): supply = pd.concat((supply, (-1) * c.pnl["p" + str(i)].loc[:, c.df.index[c.df["bus" + str(i)].isin(buses)]].groupby(c.df.carrier, @@ -522,10 +523,11 @@ if __name__ == "__main__": snakemake = mock_snakemake( 'plot_network', simpl='', - clusters=48, - lv=1.0, - sector_opts='Co2L0-168H-T-H-B-I-solar3-dist1', - planning_horizons=2050, + clusters=45, + lv=1.5, + opts='', + sector_opts='Co2L0-168H-T-H-B-I-solar+p3-dist1', + planning_horizons=2030, ) overrides = override_component_attrs(snakemake.input.overrides) diff --git a/scripts/prepare_sector_network.py b/scripts/prepare_sector_network.py index c256c251..645dc9e1 100644 --- a/scripts/prepare_sector_network.py +++ b/scripts/prepare_sector_network.py @@ -19,6 +19,56 @@ from helper import override_component_attrs import logging logger = logging.getLogger(__name__) +from types import SimpleNamespace +spatial = SimpleNamespace() + + +def define_spatial(nodes): + """ + Namespace for spatial + + Parameters + ---------- + nodes : list-like + """ + + global spatial + global options + + spatial.nodes = nodes + + # biomass + + spatial.biomass = SimpleNamespace() + + if options["biomass_transport"]: + spatial.biomass.nodes = nodes + " solid biomass" + spatial.biomass.locations = nodes + spatial.biomass.industry = nodes + " solid biomass for industry" + spatial.biomass.industry_cc = nodes + " solid biomass for industry CC" + else: + spatial.biomass.nodes = ["EU solid biomass"] + spatial.biomass.locations = ["EU"] + spatial.biomass.industry = ["solid biomass for industry"] + spatial.biomass.industry_cc = ["solid biomass for industry CC"] + + spatial.biomass.df = pd.DataFrame(vars(spatial.biomass), index=nodes) + + # co2 + + spatial.co2 = SimpleNamespace() + + if options["co2_network"]: + spatial.co2.nodes = nodes + " co2 stored" + spatial.co2.locations = nodes + spatial.co2.vents = nodes + " co2 vent" + else: + spatial.co2.nodes = ["co2 stored"] + spatial.co2.locations = ["EU"] + spatial.co2.vents = ["co2 vent"] + + spatial.co2.df = pd.DataFrame(vars(spatial.co2), index=nodes) + def emission_sectors_from_opts(opts): @@ -58,6 +108,40 @@ def get(item, investment_year=None): return item +def create_network_topology(n, prefix, connector=" -> "): + """ + Create a network topology like the power transmission network. + + Parameters + ---------- + n : pypsa.Network + prefix : str + connector : str + + Returns + ------- + pd.DataFrame with columns bus0, bus1 and length + """ + + ln_attrs = ["bus0", "bus1", "length"] + lk_attrs = ["bus0", "bus1", "length", "underwater_fraction"] + + candidates = pd.concat([ + n.lines[ln_attrs], + n.links.loc[n.links.carrier == "DC", lk_attrs] + ]).fillna(0) + + positive_order = candidates.bus0 < candidates.bus1 + candidates_p = candidates[positive_order] + swap_buses = {"bus0": "bus1", "bus1": "bus0"} + candidates_n = candidates[~positive_order].rename(columns=swap_buses) + candidates = pd.concat([candidates_p, candidates_n]) + + topo = candidates.groupby(["bus0", "bus1"], as_index=False).mean() + topo.index = topo.apply(lambda c: prefix + c.bus0 + connector + c.bus1, axis=1) + return topo + + def co2_emissions_year(countries, opts, year): """ Calculate CO2 emissions in one specific year (e.g. 1990 or 2018). @@ -79,7 +163,7 @@ def co2_emissions_year(countries, opts, year): co2_emissions = co2_totals.loc[countries, sectors].sum().sum() # convert MtCO2 to GtCO2 - co2_emissions *= 0.001 + co2_emissions *= 0.001 return co2_emissions @@ -106,17 +190,14 @@ def build_carbon_budget(o, fn): #emissions at the beginning of the path (last year available 2018) e_0 = co2_emissions_year(countries, opts, year=2018) - - #emissions in 2019 and 2020 assumed equal to 2018 and substracted - carbon_budget -= 2 * e_0 - + planning_horizons = snakemake.config['scenario']['planning_horizons'] t_0 = planning_horizons[0] if "be" in o: # final year in the path - t_f = t_0 + (2 * carbon_budget / e_0).round(0) + t_f = t_0 + (2 * carbon_budget / e_0).round(0) def beta_decay(t): cdf_term = (t - t_0) / (t_f - t_0) @@ -148,6 +229,53 @@ def add_lifetime_wind_solar(n, costs): n.generators.loc[gen_i, "lifetime"] = costs.at[carrier, 'lifetime'] +def create_network_topology(n, prefix, connector=" -> ", bidirectional=True): + """ + Create a network topology like the power transmission network. + + Parameters + ---------- + n : pypsa.Network + prefix : str + connector : str + bidirectional : bool, default True + True: one link for each connection + False: one link for each connection and direction (back and forth) + + Returns + ------- + pd.DataFrame with columns bus0, bus1 and length + """ + + ln_attrs = ["bus0", "bus1", "length"] + lk_attrs = ["bus0", "bus1", "length", "underwater_fraction"] + + candidates = pd.concat([ + n.lines[ln_attrs], + n.links.loc[n.links.carrier == "DC", lk_attrs] + ]).fillna(0) + + positive_order = candidates.bus0 < candidates.bus1 + candidates_p = candidates[positive_order] + swap_buses = {"bus0": "bus1", "bus1": "bus0"} + candidates_n = candidates[~positive_order].rename(columns=swap_buses) + candidates = pd.concat([candidates_p, candidates_n]) + + def make_index(c): + return prefix + c.bus0 + connector + c.bus1 + + topo = candidates.groupby(["bus0", "bus1"], as_index=False).mean() + topo.index = topo.apply(make_index, axis=1) + + if not bidirectional: + topo_reverse = topo.copy() + topo_reverse.rename(columns=swap_buses, inplace=True) + topo_reverse.index = topo_reverse.apply(make_index, axis=1) + topo = topo.append(topo_reverse) + + return topo + + # TODO merge issue with PyPSA-Eur def update_wind_solar_costs(n, costs): """ @@ -277,6 +405,9 @@ def patch_electricity_network(n): update_wind_solar_costs(n, costs) n.loads["carrier"] = "electricity" n.buses["location"] = n.buses.index + # remove trailing white space of load index until new PyPSA version after v0.18. + n.loads.rename(lambda x: x.strip(), inplace=True) + n.loads_t.p_set.rename(lambda x: x.strip(), axis=1, inplace=True) def add_co2_tracking(n, options): @@ -303,26 +434,26 @@ def add_co2_tracking(n, options): ) # this tracks CO2 stored, e.g. underground - n.add("Bus", - "co2 stored", - location="EU", + n.madd("Bus", + spatial.co2.nodes, + location=spatial.co2.locations, carrier="co2 stored" ) - n.add("Store", - "co2 stored", + n.madd("Store", + spatial.co2.nodes, e_nom_extendable=True, - e_nom_max=options['co2_sequestration_potential'] * 1e6, + e_nom_max=np.inf, capital_cost=options['co2_sequestration_cost'], carrier="co2 stored", - bus="co2 stored" + bus=spatial.co2.nodes ) if options['co2_vent']: - n.add("Link", - "co2 vent", - bus0="co2 stored", + n.madd("Link", + spatial.co2.vents, + bus0=spatial.co2.nodes, bus1="co2 atmosphere", carrier="co2 vent", efficiency=1., @@ -330,6 +461,28 @@ def add_co2_tracking(n, options): ) +def add_co2_network(n, costs): + + logger.info("Adding CO2 network.") + co2_links = create_network_topology(n, "CO2 pipeline ") + + cost_onshore = (1 - co2_links.underwater_fraction) * costs.at['CO2 pipeline', 'fixed'] * co2_links.length + cost_submarine = co2_links.underwater_fraction * costs.at['CO2 submarine pipeline', 'fixed'] * co2_links.length + capital_cost = cost_onshore + cost_submarine + + n.madd("Link", + co2_links.index, + bus0=co2_links.bus0.values + " co2 stored", + bus1=co2_links.bus1.values + " co2 stored", + p_min_pu=-1, + p_nom_extendable=True, + length=co2_links.length.values, + capital_cost=capital_cost.values, + carrier="CO2 pipeline", + lifetime=costs.at['CO2 pipeline', 'lifetime'] + ) + + def add_dac(n, costs): heat_carriers = ["urban central heat", "services urban decentral heat"] @@ -340,10 +493,9 @@ def add_dac(n, costs): efficiency3 = -(costs.at['direct air capture', 'heat-input'] - costs.at['direct air capture', 'compression-heat-output']) n.madd("Link", - locations, - suffix=" DAC", + heat_buses.str.replace(" heat", " DAC"), bus0="co2 atmosphere", - bus1="co2 stored", + bus1=spatial.co2.df.loc[locations, "nodes"].values, bus2=locations.values, bus3=heat_buses, carrier="DAC", @@ -487,6 +639,8 @@ def prepare_data(n): nodal_energy_totals = energy_totals.loc[pop_layout.ct].fillna(0.) nodal_energy_totals.index = pop_layout.index + # district heat share not weighted by population + district_heat_share = nodal_energy_totals["district heat share"].round(2) nodal_energy_totals = nodal_energy_totals.multiply(pop_layout.fraction, axis=0) # copy forward the daily average heat demand into each hour, so it can be multipled by the intraday profile @@ -609,7 +763,7 @@ def prepare_data(n): ) - return nodal_energy_totals, heat_demand, ashp_cop, gshp_cop, solar_thermal, transport, avail_profile, dsm_profile, nodal_transport_data + return nodal_energy_totals, heat_demand, ashp_cop, gshp_cop, solar_thermal, transport, avail_profile, dsm_profile, nodal_transport_data, district_heat_share # TODO checkout PyPSA-Eur script @@ -775,7 +929,8 @@ def insert_electricity_distribution_grid(n, costs): marginal_cost=n.generators.loc[solar, 'marginal_cost'], capital_cost=costs.at['solar-rooftop', 'fixed'], efficiency=n.generators.loc[solar, 'efficiency'], - p_max_pu=n.generators_t.p_max_pu[solar] + p_max_pu=n.generators_t.p_max_pu[solar], + lifetime=costs.at['solar-rooftop', 'lifetime'] ) n.add("Carrier", "home battery") @@ -823,7 +978,7 @@ def insert_gas_distribution_costs(n, costs): # TODO options? f_costs = options['gas_distribution_grid_cost_factor'] - + print("Inserting gas distribution grid with investment cost factor of", f_costs) capital_cost = costs.loc['electricity distribution grid']["fixed"] * f_costs @@ -832,7 +987,7 @@ def insert_gas_distribution_costs(n, costs): gas_b = n.links.index[n.links.carrier.str.contains("gas boiler") & (~n.links.carrier.str.contains("urban central"))] n.links.loc[gas_b, "capital_cost"] += capital_cost - + # micro CHPs mchp = n.links.index[n.links.carrier.str.contains("micro gas")] n.links.loc[mchp, "capital_cost"] += capital_cost @@ -994,10 +1149,11 @@ def add_storage(n, costs): if options['methanation']: n.madd("Link", - nodes + " Sabatier", + spatial.nodes, + suffix=" Sabatier", bus0=nodes + " H2", bus1="EU gas", - bus2="co2 stored", + bus2=spatial.co2.nodes, p_nom_extendable=True, carrier="Sabatier", efficiency=costs.at["methanation", "efficiency"], @@ -1009,10 +1165,11 @@ def add_storage(n, costs): if options['helmeth']: n.madd("Link", - nodes + " helmeth", + spatial.nodes, + suffix=" helmeth", bus0=nodes, bus1="EU gas", - bus2="co2 stored", + bus2=spatial.co2.nodes, carrier="helmeth", p_nom_extendable=True, efficiency=costs.at["helmeth", "efficiency"], @@ -1025,11 +1182,12 @@ def add_storage(n, costs): if options['SMR']: n.madd("Link", - nodes + " SMR CC", + spatial.nodes, + suffix=" SMR CC", bus0="EU gas", bus1=nodes + " H2", bus2="co2 atmosphere", - bus3="co2 stored", + bus3=spatial.co2.nodes, p_nom_extendable=True, carrier="SMR CC", efficiency=costs.at["SMR CC", "efficiency"], @@ -1080,7 +1238,7 @@ def add_land_transport(n, costs): suffix=" EV battery", carrier="Li ion" ) - + p_set = electric_share * (transport[nodes] + cycling_shift(transport[nodes], 1) + cycling_shift(transport[nodes], 2)) / 3 n.madd("Load", @@ -1091,8 +1249,8 @@ def add_land_transport(n, costs): p_set=p_set ) - - p_nom = nodal_transport_data["number cars"] * options.get("bev_charge_rate", 0.011) * electric_share + + p_nom = nodal_transport_data["number cars"] * options.get("bev_charge_rate", 0.011) * electric_share n.madd("Link", nodes, @@ -1124,7 +1282,7 @@ def add_land_transport(n, costs): if electric_share > 0 and options["bev_dsm"]: - e_nom = nodal_transport_data["number cars"] * options.get("bev_energy", 0.05) * options["bev_availability"] * electric_share + e_nom = nodal_transport_data["number cars"] * options.get("bev_energy", 0.05) * options["bev_availability"] * electric_share n.madd("Store", nodes, @@ -1184,12 +1342,11 @@ def add_heat(n, costs): sectors = ["residential", "services"] - nodes = create_nodes_for_heat_sector() + + nodes, dist_fraction, urban_fraction = create_nodes_for_heat_sector() #NB: must add costs of central heating afterwards (EUR 400 / kWpeak, 50a, 1% FOM from Fraunhofer ISE) - urban_fraction = options['central_fraction'] * pop_layout["urban"] / pop_layout[["urban", "rural"]].sum(axis=1) - # exogenously reduce space heat demand if options["reduce_space_heat_exogenously"]: dE = get(options["reduce_space_heat_exogenously_factor"], investment_year) @@ -1204,7 +1361,7 @@ def add_heat(n, costs): "services urban decentral", "urban central" ] - + for name in heat_systems: name_type = "central" if name == "urban central" else "decentral" @@ -1220,15 +1377,22 @@ def add_heat(n, costs): ## Add heat load for sector in sectors: + # heat demand weighting if "rural" in name: factor = 1 - urban_fraction[nodes[name]] - elif "urban" in name: - factor = urban_fraction[nodes[name]] + elif "urban central" in name: + factor = dist_fraction[nodes[name]] + elif "urban decentral" in name: + factor = urban_fraction[nodes[name]] - \ + dist_fraction[nodes[name]] + else: + raise NotImplementedError(f" {name} not in " f"heat systems: {heat_systems}") + if sector in name: heat_load = heat_demand[[sector + " water",sector + " space"]].groupby(level=1,axis=1).sum()[nodes[name]].multiply(factor) if name == "urban central": - heat_load = heat_demand.groupby(level=1,axis=1).sum()[nodes[name]].multiply(urban_fraction[nodes[name]] * (1 + options['district_heating_loss'])) + heat_load = heat_demand.groupby(level=1,axis=1).sum()[nodes[name]].multiply(factor * (1 + options['district_heating']['district_heating_loss'])) n.madd("Load", nodes[name], @@ -1286,16 +1450,16 @@ def add_heat(n, costs): p_nom_extendable=True ) - + if isinstance(options["tes_tau"], dict): tes_time_constant_days = options["tes_tau"][name_type] else: logger.warning("Deprecated: a future version will require you to specify 'tes_tau' ", "for 'decentral' and 'central' separately.") tes_time_constant_days = options["tes_tau"] if name_type == "decentral" else 180. - + # conversion from EUR/m^3 to EUR/MWh for 40 K diff and 1.17 kWh/m^3/K - capital_cost = costs.at[name_type + ' water tank storage', 'fixed'] / 0.00117 / 40 + capital_cost = costs.at[name_type + ' water tank storage', 'fixed'] / 0.00117 / 40 n.madd("Store", nodes[name] + f" {name} water tanks", @@ -1378,7 +1542,7 @@ def add_heat(n, costs): bus1=nodes[name], bus2=nodes[name] + " urban central heat", bus3="co2 atmosphere", - bus4="co2 stored", + bus4=spatial.co2.df.loc[nodes[name], "nodes"].values, carrier="urban central gas CHP CC", p_nom_extendable=True, capital_cost=costs.at['central gas CHP', 'fixed']*costs.at['central gas CHP', 'efficiency'] + costs.at['biomass CHP capture', 'fixed']*costs.at['gas', 'CO2 intensity'], @@ -1508,37 +1672,54 @@ def create_nodes_for_heat_sector(): # rural are areas with low heating density and individual heating # urban are areas with high heating density # urban can be split into district heating (central) and individual heating (decentral) - + + ct_urban = pop_layout.urban.groupby(pop_layout.ct).sum() + # distribution of urban population within a country + pop_layout["urban_ct_fraction"] = pop_layout.urban / pop_layout.ct.map(ct_urban.get) + sectors = ["residential", "services"] - + nodes = {} + urban_fraction = pop_layout.urban / pop_layout[["rural", "urban"]].sum(axis=1) + for sector in sectors: nodes[sector + " rural"] = pop_layout.index + nodes[sector + " urban decentral"] = pop_layout.index - if options["central"]: - # TODO: this looks hardcoded, move to config - urban_decentral_ct = pd.Index(["ES", "GR", "PT", "IT", "BG"]) - nodes[sector + " urban decentral"] = pop_layout.index[pop_layout.ct.isin(urban_decentral_ct)] - else: - nodes[sector + " urban decentral"] = pop_layout.index - - # for central nodes, residential and services are aggregated - nodes["urban central"] = pop_layout.index.symmetric_difference(nodes["residential urban decentral"]) - - return nodes + # maximum potential of urban demand covered by district heating + central_fraction = options['district_heating']["potential"] + # district heating share at each node + dist_fraction_node = district_heat_share * pop_layout["urban_ct_fraction"] / pop_layout["fraction"] + nodes["urban central"] = dist_fraction_node.index + # if district heating share larger than urban fraction -> set urban + # fraction to district heating share + urban_fraction = pd.concat([urban_fraction, dist_fraction_node], + axis=1).max(axis=1) + # difference of max potential and today's share of district heating + diff = (urban_fraction * central_fraction) - dist_fraction_node + progress = get(options["district_heating"]["potential"], investment_year) + dist_fraction_node += diff * progress + print( + "The current district heating share compared to the maximum", + f"possible is increased by a progress factor of\n{progress}", + f"resulting in a district heating share of\n{dist_fraction_node}" + ) + + return nodes, dist_fraction_node, urban_fraction def add_biomass(n, costs): print("adding biomass") - # biomass distributed at country level - i.e. transport within country allowed - countries = n.buses.country.dropna().unique() - biomass_potentials = pd.read_csv(snakemake.input.biomass_potentials, index_col=0) - n.add("Carrier", "biogas") + if options["biomass_transport"]: + biomass_potentials_spatial = biomass_potentials.rename(index=lambda x: x + " solid biomass") + else: + biomass_potentials_spatial = biomass_potentials.sum() + n.add("Carrier", "biogas") n.add("Carrier", "solid biomass") n.add("Bus", @@ -1547,9 +1728,9 @@ def add_biomass(n, costs): carrier="biogas" ) - n.add("Bus", - "EU solid biomass", - location="EU", + n.madd("Bus", + spatial.biomass.nodes, + location=spatial.biomass.locations, carrier="solid biomass" ) @@ -1557,18 +1738,18 @@ def add_biomass(n, costs): "EU biogas", bus="EU biogas", carrier="biogas", - e_nom=biomass_potentials.loc[countries, "biogas"].sum(), + e_nom=biomass_potentials["biogas"].sum(), marginal_cost=costs.at['biogas', 'fuel'], - e_initial=biomass_potentials.loc[countries, "biogas"].sum() + e_initial=biomass_potentials["biogas"].sum() ) - n.add("Store", - "EU solid biomass", - bus="EU solid biomass", + n.madd("Store", + spatial.biomass.nodes, + bus=spatial.biomass.nodes, carrier="solid biomass", - e_nom=biomass_potentials.loc[countries, "solid biomass"].sum(), + e_nom=biomass_potentials_spatial["solid biomass"], marginal_cost=costs.at['solid biomass', 'fuel'], - e_initial=biomass_potentials.loc[countries, "solid biomass"].sum() + e_initial=biomass_potentials_spatial["solid biomass"] ) n.add("Link", @@ -1583,6 +1764,32 @@ def add_biomass(n, costs): p_nom_extendable=True ) + if options["biomass_transport"]: + + transport_costs = pd.read_csv( + snakemake.input.biomass_transport_costs, + index_col=0, + squeeze=True + ) + + # add biomass transport + biomass_transport = create_network_topology(n, "biomass transport ", bidirectional=False) + + # costs + bus0_costs = biomass_transport.bus0.apply(lambda x: transport_costs[x[:2]]) + bus1_costs = biomass_transport.bus1.apply(lambda x: transport_costs[x[:2]]) + biomass_transport["costs"] = pd.concat([bus0_costs, bus1_costs], axis=1).mean(axis=1) + + n.madd("Link", + biomass_transport.index, + bus0=biomass_transport.bus0 + " solid biomass", + bus1=biomass_transport.bus1 + " solid biomass", + p_nom_extendable=True, + length=biomass_transport.length.values, + marginal_cost=biomass_transport.costs * biomass_transport.length.values, + capital_cost=1, + carrier="solid biomass transport" + ) #AC buses with district heating urban_central = n.buses.index[n.buses.carrier == "urban central heat"] @@ -1593,7 +1800,7 @@ def add_biomass(n, costs): n.madd("Link", urban_central + " urban central solid biomass CHP", - bus0="EU solid biomass", + bus0=spatial.biomass.df.loc[urban_central, "nodes"].values, bus1=urban_central, bus2=urban_central + " urban central heat", carrier="urban central solid biomass CHP", @@ -1607,11 +1814,11 @@ def add_biomass(n, costs): n.madd("Link", urban_central + " urban central solid biomass CHP CC", - bus0="EU solid biomass", + bus0=spatial.biomass.df.loc[urban_central, "nodes"].values, bus1=urban_central, bus2=urban_central + " urban central heat", bus3="co2 atmosphere", - bus4="co2 stored", + bus4=spatial.co2.df.loc[urban_central, "nodes"].values, carrier="urban central solid biomass CHP CC", p_nom_extendable=True, capital_cost=costs.at[key, 'fixed'] * costs.at[key, 'efficiency'] + costs.at['biomass CHP capture', 'fixed'] * costs.at['solid biomass', 'CO2 intensity'], @@ -1633,36 +1840,39 @@ def add_industry(n, costs): # 1e6 to convert TWh to MWh industrial_demand = pd.read_csv(snakemake.input.industrial_demand, index_col=0) * 1e6 - solid_biomass_by_country = industrial_demand["solid biomass"].groupby(pop_layout.ct).sum() - - n.add("Bus", - "solid biomass for industry", - location="EU", + n.madd("Bus", + spatial.biomass.industry, + location=spatial.biomass.locations, carrier="solid biomass for industry" ) - n.add("Load", - "solid biomass for industry", - bus="solid biomass for industry", + if options["biomass_transport"]: + p_set = industrial_demand.loc[spatial.biomass.locations, "solid biomass"].rename(index=lambda x: x + " solid biomass for industry") / 8760 + else: + p_set = industrial_demand["solid biomass"].sum() / 8760 + + n.madd("Load", + spatial.biomass.industry, + bus=spatial.biomass.industry, carrier="solid biomass for industry", - p_set=solid_biomass_by_country.sum() / 8760 + p_set=p_set ) - n.add("Link", - "solid biomass for industry", - bus0="EU solid biomass", - bus1="solid biomass for industry", + n.madd("Link", + spatial.biomass.industry, + bus0=spatial.biomass.nodes, + bus1=spatial.biomass.industry, carrier="solid biomass for industry", p_nom_extendable=True, efficiency=1. ) - n.add("Link", - "solid biomass for industry CC", - bus0="EU solid biomass", - bus1="solid biomass for industry", + n.madd("Link", + spatial.biomass.industry_cc, + bus0=spatial.biomass.nodes, + bus1=spatial.biomass.industry, bus2="co2 atmosphere", - bus3="co2 stored", + bus3=spatial.co2.nodes, carrier="solid biomass for industry CC", p_nom_extendable=True, capital_cost=costs.at["cement capture", "fixed"] * costs.at['solid biomass', 'CO2 intensity'], @@ -1695,12 +1905,13 @@ def add_industry(n, costs): efficiency2=costs.at['gas', 'CO2 intensity'] ) - n.add("Link", - "gas for industry CC", + n.madd("Link", + spatial.co2.locations, + suffix=" gas for industry CC", bus0="EU gas", bus1="gas for industry", bus2="co2 atmosphere", - bus3="co2 stored", + bus3=spatial.co2.nodes, carrier="gas for industry CC", p_nom_extendable=True, capital_cost=costs.at["cement capture", "fixed"] * costs.at['gas', 'CO2 intensity'], @@ -1759,9 +1970,9 @@ def add_industry(n, costs): if shipping_hydrogen_share < 1: shipping_oil_share = 1 - shipping_hydrogen_share - + p_set = shipping_oil_share * nodal_energy_totals.loc[nodes, all_navigation].sum(axis=1) * 1e6 / 8760. - + n.madd("Load", nodes, suffix=" shipping oil", @@ -1769,7 +1980,7 @@ def add_industry(n, costs): carrier="shipping oil", p_set=p_set ) - + co2 = shipping_oil_share * nodal_energy_totals.loc[nodes, all_navigation].sum().sum() * 1e6 / 8760 * costs.at["oil", "CO2 intensity"] n.add("Load", @@ -1788,7 +1999,7 @@ def add_industry(n, costs): ) if "EU oil Store" not in n.stores.index: - + #could correct to e.g. 0.001 EUR/kWh * annuity and O&M n.add("Store", "EU oil Store", @@ -1810,7 +2021,7 @@ def add_industry(n, costs): if options["oil_boilers"]: - nodes_heat = create_nodes_for_heat_sector() + nodes_heat = create_nodes_for_heat_sector()[0] for name in ["residential rural", "services rural", "residential urban decentral", "services urban decentral"]: @@ -1831,7 +2042,7 @@ def add_industry(n, costs): nodes + " Fischer-Tropsch", bus0=nodes + " H2", bus1="EU oil", - bus2="co2 stored", + bus2=spatial.co2.nodes, carrier="Fischer-Tropsch", efficiency=costs.at["Fischer-Tropsch", 'efficiency'], capital_cost=costs.at["Fischer-Tropsch", 'fixed'], @@ -1920,11 +2131,12 @@ def add_industry(n, costs): ) #assume enough local waste heat for CC - n.add("Link", - "process emissions CC", + n.madd("Link", + spatial.co2.locations, + suffix=" process emissions CC", bus0="process emissions", bus1="co2 atmosphere", - bus2="co2 stored", + bus2=spatial.co2.nodes, carrier="process emissions CC", p_nom_extendable=True, capital_cost=costs.at["cement capture", "fixed"], @@ -2020,7 +2232,7 @@ def add_agriculture(n, costs): def decentral(n): - """Removes the electricity transmission system.""" + """Removes the electricity transmission system.""" n.lines.drop(n.lines.index, inplace=True) n.links.drop(n.links.index[n.links.carrier.isin(["DC", "B2B"])], inplace=True) @@ -2053,7 +2265,7 @@ def maybe_adjust_costs_and_potentials(n, opts): if attr == 'p_nom_max': comps = {"Generator", "Link", "StorageUnit"} elif attr == 'e_nom_max': - comps = {"Store"} + comps = {"Store"} else: comps = {"Generator", "Link", "StorageUnit", "Store"} for c in n.iterate_components(comps): @@ -2072,17 +2284,18 @@ def limit_individual_line_extension(n, maxext): hvdc = n.links.index[n.links.carrier == 'DC'] n.links.loc[hvdc, 'p_nom_max'] = n.links.loc[hvdc, 'p_nom'] + maxext - +#%% if __name__ == "__main__": if 'snakemake' not in globals(): from helper import mock_snakemake snakemake = mock_snakemake( 'prepare_sector_network', simpl='', - clusters=48, + opts="", + clusters="37", lv=1.0, sector_opts='Co2L0-168H-T-H-B-I-solar3-dist1', - planning_horizons=2020, + planning_horizons="2020", ) logging.basicConfig(level=snakemake.config['logging_level']) @@ -2107,8 +2320,10 @@ if __name__ == "__main__": patch_electricity_network(n) + define_spatial(pop_layout.index) + if snakemake.config["foresight"] == 'myopic': - + add_lifetime_wind_solar(n, costs) conventional = snakemake.config['existing_capacities']['conventional_carriers'] @@ -2129,11 +2344,13 @@ if __name__ == "__main__": if o[:4] == "dist": options['electricity_distribution_grid'] = True options['electricity_distribution_grid_cost_factor'] = float(o[4:].replace("p", ".").replace("m", "-")) + if o == "biomasstransport": + options["biomass_transport"] = True - nodal_energy_totals, heat_demand, ashp_cop, gshp_cop, solar_thermal, transport, avail_profile, dsm_profile, nodal_transport_data = prepare_data(n) + nodal_energy_totals, heat_demand, ashp_cop, gshp_cop, solar_thermal, transport, avail_profile, dsm_profile, nodal_transport_data, district_heat_share = prepare_data(n) if "nodistrict" in opts: - options["central"] = False + options["district_heating"]["progress"] = 0.0 if "T" in opts: add_land_transport(n, costs) @@ -2162,6 +2379,9 @@ if __name__ == "__main__": if "noH2network" in opts: remove_h2_network(n) + if options["co2_network"]: + add_co2_network(n, costs) + for o in opts: m = re.match(r'^\d+h$', o, re.IGNORECASE) if m is not None: diff --git a/scripts/solve_network.py b/scripts/solve_network.py index a46acc30..eac1f581 100644 --- a/scripts/solve_network.py +++ b/scripts/solve_network.py @@ -3,6 +3,7 @@ import pypsa import numpy as np +import pandas as pd from pypsa.linopt import get_var, linexpr, define_constraints @@ -19,12 +20,47 @@ pypsa.pf.logger.setLevel(logging.WARNING) def add_land_use_constraint(n): - #warning: this will miss existing offwind which is not classed AC-DC and has carrier 'offwind' - for carrier in ['solar', 'onwind', 'offwind-ac', 'offwind-dc']: - existing = n.generators.loc[n.generators.carrier == carrier, "p_nom"].groupby(n.generators.bus.map(n.buses.location)).sum() - existing.index += " " + carrier + "-" + snakemake.wildcards.planning_horizons - n.generators.loc[existing.index, "p_nom_max"] -= existing + if 'm' in snakemake.wildcards.clusters: + _add_land_use_constraint_m(n) + else: + _add_land_use_constraint(n) + +def _add_land_use_constraint(n): + #warning: this will miss existing offwind which is not classed AC-DC and has carrier 'offwind' + + for carrier in ['solar', 'onwind', 'offwind-ac', 'offwind-dc']: + existing = n.generators.loc[n.generators.carrier==carrier,"p_nom"].groupby(n.generators.bus.map(n.buses.location)).sum() + existing.index += " " + carrier + "-" + snakemake.wildcards.planning_horizons + n.generators.loc[existing.index,"p_nom_max"] -= existing + + n.generators.p_nom_max.clip(lower=0, inplace=True) + + +def _add_land_use_constraint_m(n): + # if generators clustering is lower than network clustering, land_use accounting is at generators clusters + + planning_horizons = snakemake.config["scenario"]["planning_horizons"] + grouping_years = snakemake.config["existing_capacities"]["grouping_years"] + current_horizon = snakemake.wildcards.planning_horizons + + for carrier in ['solar', 'onwind', 'offwind-ac', 'offwind-dc']: + + existing = n.generators.loc[n.generators.carrier==carrier,"p_nom"] + ind = list(set([i.split(sep=" ")[0] + ' ' + i.split(sep=" ")[1] for i in existing.index])) + + previous_years = [ + str(y) for y in + planning_horizons + grouping_years + if y < int(snakemake.wildcards.planning_horizons) + ] + + for p_year in previous_years: + ind2 = [i for i in ind if i + " " + carrier + "-" + p_year in existing.index] + sel_current = [i + " " + carrier + "-" + current_horizon for i in ind2] + sel_p_year = [i + " " + carrier + "-" + p_year for i in ind2] + n.generators.loc[sel_current, "p_nom_max"] -= existing.loc[sel_p_year].rename(lambda x: x[:-4] + current_horizon) + n.generators.p_nom_max.clip(lower=0, inplace=True) @@ -150,8 +186,26 @@ def add_chp_constraints(n): define_constraints(n, lhs, "<=", 0, 'chplink', 'backpressure') +def add_co2_sequestration_limit(n, sns): + + co2_stores = n.stores.loc[n.stores.carrier=='co2 stored'].index + + if co2_stores.empty or ('Store', 'e') not in n.variables.index: + return + + vars_final_co2_stored = get_var(n, 'Store', 'e').loc[sns[-1], co2_stores] + + lhs = linexpr((1, vars_final_co2_stored)).sum() + rhs = n.config["sector"].get("co2_sequestration_potential", 200) * 1e6 + + name = 'co2_sequestration_limit' + define_constraints(n, lhs, "<=", rhs, 'GlobalConstraint', + 'mu', axes=pd.Index([name]), spec=name) + + def extra_functionality(n, snapshots): add_battery_constraints(n) + add_co2_sequestration_limit(n, snapshots) def solve_network(n, config, opts='', **kwargs):