Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Linopy integration for Power and Sector-Coupled Modelling #1172

Open
wants to merge 62 commits into
base: main
Choose a base branch
from

Conversation

ekatef
Copy link
Member

@ekatef ekatef commented Nov 5, 2024

Closes #494

(and hopefully #968)

Changes proposed in this Pull Request

Integrated changes needed to update PyPSA version and use linopy.

To have a better overview, the PR has been implemented in three steps:

  1. Update PyPSA version (follows Update PyPSA version #1065)
  2. Integrate linopy into power part of the model (follows Linopy transition #796, the diffs can be checked in this local PR)
  3. Integrate linopy into cross-sectoral part of the model (implemented currently as a dedicated PR to the Sec- repo)

Checklist

  • I consent to the release of this PR's code under the AGPLv3 license and non-code contributions under CC0-1.0 and CC-BY-4.0.
  • I tested my contribution locally and it seems to work fine.
  • Code and workflow changes are sufficiently documented.
  • Newly introduced dependencies are added to envs/environment.yaml and doc/requirements.txt.
  • Changes in configuration options are added in all of config.default.yaml and config.tutorial.yaml.
  • Add a test config or line additions to test/ (note tests are changing the config.tutorial.yaml)
  • Changes in configuration options are also documented in doc/configtables/*.csv and line references are adjusted in doc/configuration.rst and doc/tutorial.rst.
  • A note for the release notes doc/release_notes.rst is amended in the format of previous release notes, including reference to the requested PR.

github-actions bot and others added 3 commits October 30, 2024 14:47
* Add a zenodo link to natura.tiff

* Update environment

* Revise structure definition for lines

* Remove get_aggregation_strategies

* Fix typo aggregation_strategies

* Replace aggregategenerators with aggregateoneport

* Add aggregation strategies as a parameter

* Re-define aggregation strategies

* Update aggregation strategies

* Update aggregation strategies for lines

* Update aggregation strategies for buses

* Fix typo

* Put aggregation strategies into a variable

* Parametrize the aggregation strategies

* Refactor update of the aggregation strategies

* Clean-up the code

* Revert "Add a zenodo link to natura.tiff"

This reverts commit 7700759.

* Define an explicit clustering strategy for v_nom

* Add a release note

* Get glpk back

* Specify v_nom for buses explicitly

* Revert "Specify v_nom for buses explicitly"

This reverts commit 20192e6.

* Add a version restriction to the environment specification

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Adjust naming

* Move the variable definition

* Move the variable

* Upgrade PyPSA version

---------

Co-authored-by: Davide Fioriti <fioritidavidesubs@gmail.com>
Co-authored-by: Davide Fioriti <67809479+davide-f@users.noreply.github.com>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
* Add a zenodo link to natura.tiff

* Update environment

* Revise structure definition for lines

* Remove get_aggregation_strategies

* Fix typo aggregation_strategies

* Replace aggregategenerators with aggregateoneport

* Add aggregation strategies as a parameter

* Re-define aggregation strategies

* Update aggregation strategies

* Update aggregation strategies for lines

* Update aggregation strategies for buses

* Fix typo

* Put aggregation strategies into a variable

* Parametrize the aggregation strategies

* Refactor update of the aggregation strategies

* Clean-up the code

* Revert "Add a zenodo link to natura.tiff"

This reverts commit 7700759.

* Define an explicit clustering strategy for v_nom

* Add a release note

* Get glpk back

* Specify v_nom for buses explicitly

* Revert "Specify v_nom for buses explicitly"

This reverts commit 20192e6.

* Add a version restriction to the environment specification

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Adjust naming

* Move the variable definition

* Move the variable

* Upgrade PyPSA version

* Update docstring

* Fix imports duplication

* Update imports

* Update the carrier-capacity constraint

* Add docstring

* Update the equity constraint

* Add docstring

* Update BAU constraint

* Update SAFE constraint

* Add docstring

* Update operational reserve margin constraint

* Add docstring

* Add an new argument to the RM constraint

* Update the update of capacity constraints

* Update adding an operational reserve margin constraint

* Update docstring

* Update battery constraint

* Add docstring

* Update a constraint related to a RES share

* Fix usage of add_ERS_constraints

* Update solving script

* Update a solving run

* Fix typos

---------

Co-authored-by: Davide Fioriti <fioritidavidesubs@gmail.com>
Co-authored-by: Davide Fioriti <67809479+davide-f@users.noreply.github.com>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
@ekatef
Copy link
Member Author

ekatef commented Nov 5, 2024

The PR has been tested against main, and the change in the objective function has been found to be less than 1%.

Codewise, the implementation is same as in #1065. Manual reimplementation has been needed because of #1170

@ekatef
Copy link
Member Author

ekatef commented Nov 5, 2024

@GbotemiB As discussed, we can use a branch of this PR (pypsa_linopy_update located in upstream directly) to accommodate your linopy PR from PyPSA-Earth-Sec. What do you think?

@ekatef ekatef mentioned this pull request Nov 5, 2024
8 tasks
@GbotemiB
Copy link
Contributor

GbotemiB commented Nov 7, 2024

@GbotemiB As discussed, we can use a branch of this PR (pypsa_linopy_update located in upstream directly) to accommodate your linopy PR from PyPSA-Earth-Sec. What do you think?

Thanks @ekatef. I can start working on the linopy PR for the sec aspect.

Copy link
Member

@davide-f davide-f left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great :)
Very minor things, basically ready to merge :)

nodes = n.buses.index[n.buses.carrier == "battery"]
if nodes.empty or ("Link", "p_nom") not in n.variables.index:
# TODO Check if the second part of the condition can make sense
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this still apply?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

My understanding is that this second part addresses some special case when the network doesn't contain any batteries, despite of the fact that there are nodes with battery carrier. Not sure when exactly can it be the case and just leaving it as is for now, as it seems not to relate to interfacing with linopy. Does it make sense for you or am I missing something?

lgrouper = n.loads.bus.map(n.buses.country)
# TODO drop load
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this for a later stage?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, my feeling is that the enclosing add_RES_constraints( ) will need some additional revisions which should go beyond this PR.

In our current upstream implementation, the function throws a warning stating that it's only a preliminary implementation. I think, it may be a good idea to open an issue to track further improvements

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In the updates version of add_RES_constraints( ), this TODO is not needed anymore

# Generators
# TODO restore grouping by countries un-commenting calls of groupby()
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this not work? for later?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, linopy.groupby( ) has been not in the best shape when I have been debugging this part. So, I have left it out for future improvements of add_RES_constraints( )

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe we should add a warning&issue then?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Have fixed the functionality. It has been more straight-forward ;)

n.lines.loc[n.lines.index.str.contains("new"), "s_nom_min"] = (
snakemake.params.augmented_line_connection.get("min_expansion")
)
# TODO Double-check handling the augmented case
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why comment this? Does it trigger an issue?
It seems look it generates issues, what do you think?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A great catch! That is just an old comment needed to debug the main functionality. Now it must be lifted

Testing it

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Have get this part back, and testing worked well

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can we drop this then?

scripts/solve_network.py Outdated Show resolved Hide resolved
scripts/solve_network.py Show resolved Hide resolved
yerbol-akhmetov and others added 13 commits November 8, 2024 20:47
Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>
* Add a zenodo link to natura.tiff

* Update environment

* Revise structure definition for lines

* Remove get_aggregation_strategies

* Fix typo aggregation_strategies

* Replace aggregategenerators with aggregateoneport

* Add aggregation strategies as a parameter

* Re-define aggregation strategies

* Update aggregation strategies

* Update aggregation strategies for lines

* Update aggregation strategies for buses

* Fix typo

* Put aggregation strategies into a variable

* Parametrize the aggregation strategies

* Refactor update of the aggregation strategies

* Clean-up the code

* Revert "Add a zenodo link to natura.tiff"

This reverts commit 7700759.

* Define an explicit clustering strategy for v_nom

* Add a release note

* Get glpk back

* Specify v_nom for buses explicitly

* Revert "Specify v_nom for buses explicitly"

This reverts commit 20192e6.

* Add a version restriction to the environment specification

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Adjust naming

* Move the variable definition

* Move the variable

* Upgrade PyPSA version

---------

Co-authored-by: Davide Fioriti <fioritidavidesubs@gmail.com>
Co-authored-by: Davide Fioriti <67809479+davide-f@users.noreply.github.com>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
* Add a zenodo link to natura.tiff

* Update environment

* Revise structure definition for lines

* Remove get_aggregation_strategies

* Fix typo aggregation_strategies

* Replace aggregategenerators with aggregateoneport

* Add aggregation strategies as a parameter

* Re-define aggregation strategies

* Update aggregation strategies

* Update aggregation strategies for lines

* Update aggregation strategies for buses

* Fix typo

* Put aggregation strategies into a variable

* Parametrize the aggregation strategies

* Refactor update of the aggregation strategies

* Clean-up the code

* Revert "Add a zenodo link to natura.tiff"

This reverts commit 7700759.

* Define an explicit clustering strategy for v_nom

* Add a release note

* Get glpk back

* Specify v_nom for buses explicitly

* Revert "Specify v_nom for buses explicitly"

This reverts commit 20192e6.

* Add a version restriction to the environment specification

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Adjust naming

* Move the variable definition

* Move the variable

* Upgrade PyPSA version

* Update docstring

* Fix imports duplication

* Update imports

* Update the carrier-capacity constraint

* Add docstring

* Update the equity constraint

* Add docstring

* Update BAU constraint

* Update SAFE constraint

* Add docstring

* Update operational reserve margin constraint

* Add docstring

* Add an new argument to the RM constraint

* Update the update of capacity constraints

* Update adding an operational reserve margin constraint

* Update docstring

* Update battery constraint

* Add docstring

* Update a constraint related to a RES share

* Fix usage of add_ERS_constraints

* Update solving script

* Update a solving run

* Fix typos

---------

Co-authored-by: Davide Fioriti <fioritidavidesubs@gmail.com>
Co-authored-by: Davide Fioriti <67809479+davide-f@users.noreply.github.com>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
@ekatef ekatef changed the title Integrate PyPSA and linopy update for power and cross-sectoral parts of the model Linopy integration for power Nov 14, 2024
@ekatef ekatef changed the title Linopy integration for power Linopy integration for Power Nov 14, 2024
@ekatef ekatef force-pushed the pypsa_linopy_update branch from 7ce0894 to 8172112 Compare December 20, 2024 23:28
Co-authored-by: Davide Fioriti <67809479+davide-f@users.noreply.github.com>
@davide-f davide-f mentioned this pull request Dec 20, 2024
5 tasks
@ekatef
Copy link
Member Author

ekatef commented Dec 21, 2024

Contribution looks quite interesting and nearly finalized.

It is advisable to update the branch too and check the new results of the run. I've compared the results of the objective values of the CI and they seem a bit different:

  • the first power test: 3809010949.0 [this PR] vs 2540477210.0 main
  • the sec test: 7.5107178e+11 [this PR] vs 7.270497e+11[main]

I remember you crosschecked them already; do you recall?

@davide-f thanks for the review! 😄

I confirm that me and Emmanuel have checked that the PR doesn't change the objective function. Though, the may be some pitfalls with testing extra-functionality functions. Have to check it.

Currently, CI is failing as it's using pinned environments which must be re-build using the updated environment.yaml. Happy to look into it but if you have some experience on that, some hints would be much appreciated 🙂

@ekatef ekatef force-pushed the pypsa_linopy_update branch from 6efea61 to 75b60ca Compare December 21, 2024 14:38
ekatef and others added 6 commits December 21, 2024 15:40
* Update pinned environment files for all platforms

---------

Co-authored-by: ekatef <30229437+ekatef@users.noreply.github.com>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Copy link
Member

@davide-f davide-f left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Added few comments related to the reserve constraints, hoping to help

lhs = merge(
dispatch * 1,
reserve * 1,
)

if not ext_i.empty:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In pypsa-eur, there is no need for this if case, though I don't think is the reason of the issue

Copy link
Member Author

@ekatef ekatef Dec 24, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

True. Probably p_max_pu[ext_i] is safe for empty ext_i. Have added a TODO to revise this part

.loc[vres_i.intersection(ext_i)]
.rename({"Generator-ext": "Generator"})
)
lhs = merge(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see you have made extensive use of the merge function.
I'm dropping a comment hoping it may help.
The documentation of merge informs that when a list of expressions is merged, it assumes they share the same index basically: in case of a list, the first index overrides the others.
https://linopy.readthedocs.io/en/latest/generated/linopy.expressions.merge.html

Not sure if it may have unintended consequences.
In pypsa-eur, I've seen a lower usage of merge but it may not be the cause of this issue

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agree and revised, except for add_CCL_constraints where it's not straght-forward (the current design implies forming an expression as a list, while pypsa-eur implementation has changed significantly)

)
lhs = merge(
lhs,
(renewable_capacity_variables * (-EPSILON_VRES * capacity_factor)).sum(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

if we rename p_nom_vres the renewable capacity variables, maybe we shorten a bit the text, but minor comment

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In pypsa-eur, they also include a xr.DataArray(capacity_factor), not sure if it has implications into the calculations

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agree, aligned with pypsa-eur implementation

@ekatef
Copy link
Member Author

ekatef commented Dec 23, 2024

Update: re-testing of the objective function is working well in some testing cases, but leads to the differences in others. In particular, the tutorial-based test is not working properly. I'm able to reproduce the issue locally, and it appears to be linked with the changes introduced before solve_network script.

Basically, the value of the objective function does not depend on the solving script but on the procedure before: solving the network prepared by main workflow with the solve_network of this PR leads to the same objective value as in this PR. Visa versa is true as well: solving the prepared network of this PR with solve_network of main gives the same objective value as in this PR.

A similar issue has been observed before for power-only model, and has been resolved by setting v_nom needed to resolve the negative objective issue #961. It looks like a similar pitfall has appeared due to the changes introduced by some cross-sectoral related changes.

@davide-f
Copy link
Member

Update: re-testing of the objective function is working well in some testing cases, but leads to the differences in others. In particular, the tutorial-based test is not working properly. I'm able to reproduce the issue locally, and it appears to be linked with the changes introduced before solve_network script.

Basically, the value of the objective function does not depend on the solving script but on the procedure before: solving the network prepared by main workflow with the solve_network of this PR leads to the same objective value as in this PR. Visa versa is true as well: solving the prepared network of this PR with solve_network of main gives the same objective value as in this PR.

A similar issue has been observed before for power-only model, and has been resolved by setting v_nom needed to resolve the negative objective issue #961. It looks like a similar pitfall has appeared due to the changes introduced by some cross-sectoral related changes.

Many thanks for investigating. What is your insight on this issue and status of the PR therefore?

can we track what has been merged that may have led to the issue?

@davide-f
Copy link
Member

Let's remind to change the dependencies also in the file doc/requirements.txt.
The API reference was not compiling because of that

@davide-f
Copy link
Member

davide-f commented Dec 24, 2024

I've been testing the model and to me, it seems like the change in the config.tutorial.yaml are not due to changes in solve_network but something else.
In particular, when debugging the solve_network in this PR and in the main, only 3 functions are actually activated add_battery_constraints, add_chp_constraints (but irrelevant and skipped) and add_co2_sequestration_limit.
Only add_battery_constraint actually changes the model but they seem fine

The checks on the network seem to suggest that the input networks match but something changes for some reason.

At deeper analysis, it seems like the onwind (and marginally offwind) are different.
In particular, the new formulation has half onwind than the base case, which may lead to different outputs

Update: the onwind p_max_pu is expected to be the same after elec.nc [add_electricity]; probably there is something going on in the clustering steps that is different

@ekatef
Copy link
Member Author

ekatef commented Dec 24, 2024

@davide-f many thanks for testing (and providing the fantastic infrastructure for inter-comparing the networks!). The issues indeed are quite weird and it has been precious to have your independent workflow to be sure that they are reproducible. Writing down the results of the testing below to possibly improve further both implementations 🙂

It has been found that the network structure and parameters are exactly the same for both implementations. The reason of the differences in the objective function are linked with the difference in p_max_pu time series for wind generators (a quick way to reproduce: n_gen_onwind = n.generators[n.generators.carrier == "onwind"].index; n.generators_t.p_max_pu[n_gen_onwind].sum(axis=1).sum()). This difference stems from the fact that aggregation_strategies are defined in the different way for the main and the updated implementation.

In the older version of PyPSA, a default strategy have has been used for p_max_pu time-series, while now we must use custom definitions for aggregation strategies. If we change aggregation_strategies:generators:p_max_pu between mean and max, the objective value changes between ~4.5 [financial units] and ~1.7 [financial units] for the updated implementation with ~2.5 [financial units] for the implementation in main.

Points of action
Given a dramatic change in the underlaying computational approaches, it doesn't make any sense trying to keep the objective. What must be done instead is tracing across the code values affected by the aggregation in simplify_network and cluster_network and check that conservation laws are respected. Potentially, it means that checks or improvements can be also needed to be added into the current main implementation which will remain as a stable backward compatible version after this PR will be merged.

@davide-f great that we have been able to discover the reason, as it appears to be a potentially major issue also for the current implementation 🎉 🎉 🎉 Thank you so much for that. I think it's a great result, and we are very close to both finalise this PR and make the current implementation really stable 😄

@ekatef
Copy link
Member Author

ekatef commented Dec 28, 2024

It appears, the trick was to use an explicit definition of an aggregation strategy as p_max_pu: weighted_average. That brings the objective function to 2512675751.0 which is quite close to 2512675455.0 from main.

The full testing results looks like follows:

main

INFO:main:Objective function: 2512675455.0
INFO:main:Objective function: 2341159599.0
INFO:main:Objective function: 2321096575.0
INFO:main:Objective function: 336330443.3
INFO:main:Objective function: 4248565071.0
INFO:main:Objective function: 4602716893.0
INFO:main:Objective function: 5723516039.0
INFO:main:Objective function: 476366499.4
INFO:main:Objective function: 10926432310.0
INFO:main:Objective function: 1461452593.0
INFO:main:Objective function: 0.211751626
INFO:main:Objective function: 4376884296.0
INFO:main:Objective function: -0.006962630898
INFO:main:Objective function: 722220262793.1165

this PR

INFO:main:Objective function: 2512675751.0
INFO:main:Objective function: 2341214199.0
INFO:main:Objective function: 2321096682.0
INFO:main:Objective function: 476365695.9
INFO:main:Objective function: 10747113970.0
INFO:main:Objective function: 1461452491.0
INFO:main:Objective function: 4324387954.0
INFO:main:Objective function: 0.001945018768
INFO:main:Objective function: 5723515621.0
INFO:main:Objective function: 336330259.2
INFO:main:Objective function: 4503925925.0
INFO:main:Objective function: 4186016571.0
INFO:main:Objective function: 5.44154318e-07
INFO:main:Objective function: 722265045983.7013

That gives the following differences:

config.custom.yaml

1.178027e-07

config.landlock.yaml

2.332123e-05

config.sector.yaml

4.609890e-08

config.test_myopic.yaml

6.200382e-05

For Monte-Carlo test, the differences may be very high (2.939659e-01 6.046785e-01 -2.149413e+00 -3.235436e-01 -2.449161e+11 -9.090421e-01 [10] -3.345290e+00 1.000000e+00 -4.559650e-02) but probably it's intentional.

@ekatef
Copy link
Member Author

ekatef commented Dec 28, 2024

@davide-f I think we have done it 😄

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Integrate Linopy
4 participants