Skip to content

Surface Plant time series integrations fix (Electricity Provided + other profiles final year value correction) #50

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.

Sign up for GitHub

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

Closed
wants to merge 11 commits into from
Closed

Surface Plant time series integrations fix (Electricity Provided + other profiles final year value correction) #50

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
47e5a4b
Adjust dx passed to trapz if last slice has fewer timesteps in heat k…
softwareengineerprogrammer Mar 5, 2025
49886ba
Adjust dx for slice size for all integrations in SurfacePlant.annual_…
softwareengineerprogrammer Mar 5, 2025
d0ac2d7
SurfacePlant.integrate_time_series_slice utility method; incorporated…
softwareengineerprogrammer Mar 5, 2025
e4eff98
SurfacePlantAbsorptionChiller fix
softwareengineerprogrammer Mar 5, 2025
3179a86
fix trapz/slice integration in AC & HP surface plants (unit tests sti…
softwareengineerprogrammer Mar 5, 2025
aabadb8
Temporarily revert slice integration behavior to original and thus al…
softwareengineerprogrammer Mar 5, 2025
dace656
Use correct slice integration behavior and update unit tests.
softwareengineerprogrammer Mar 5, 2025
6af58c7
fix/update s-dac profile parsing test
softwareengineerprogrammer Mar 5, 2025
591949b
Don't wrap dx_steps as single-member tuple
softwareengineerprogrammer Mar 5, 2025
3abcfc8
Run SBT example tests last to reduce wait time under some development…
softwareengineerprogrammer Mar 5, 2025
7351f3c
unit test time series integration
softwareengineerprogrammer Mar 6, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
53 changes: 36 additions & 17 deletions src/geophires_x/SurfacePlant.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,16 +1,35 @@
import sys
import os
import numpy as np

from .GeoPHIRESUtils import quantity
from .OptionList import EndUseOptions, PlantType
from .Parameter import floatParameter, intParameter, strParameter, OutputParameter, ReadParameter, \
from .Parameter import floatParameter, intParameter, OutputParameter, ReadParameter, \
coerce_int_params_to_enum_values
from .Units import *
import geophires_x.Model as Model
import pandas as pd


class SurfacePlant:
@staticmethod
def integrate_time_series_slice(
series: np.ndarray,
_i: int,
time_steps_per_year: int,
utilization_factor: float
) -> np.float64:
_slice = series[(0 + _i * time_steps_per_year):((_i + 1) * time_steps_per_year) + 1]

# Note that len(_slice) - 1 may be less than time_steps_per_year for the last slice.

if len(_slice) == 1:
return _slice[0]

dx_steps = len(_slice) - 1

return np.trapz(
_slice,
dx=1. / dx_steps * 365. * 24.
) * 1000. * utilization_factor

Comment on lines +12 to +32
Copy link
Owner Author

@softwareengineerprogrammer softwareengineerprogrammer Mar 5, 2025 •
edited
Loading

Choose a reason for hiding this comment

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

The crux of this fix is that the previous behavior was doing the equivalent of erroneously setting dx_steps to time_steps_per_year. The last slices of the time series in question always have one less data point than time_steps_per_year, which was causing the trapezoidal rule to behave as if the slice had a 'phantom' entry with value=0 appended to it, thus lowering the estimated integral value of the slice.

def remaining_reservoir_heat_content(self, InitialReservoirHeatContent: np.ndarray, HeatkWhExtracted: np.ndarray) -> np.ndarray:
"""
Calculate reservoir heat content
Expand Down Expand Up @@ -152,15 +171,15 @@ def electricity_heat_production(self, enduse_option: EndUseOptions, availability

return ElectricityProduced, HeatExtracted, HeatProduced, HeatExtractedTowardsElectricity

def annual_electricity_pumping_power(self, plant_lifetime: int, enduse_option: EndUseOptions, HeatExtracted: np.ndarray,
timestepsperyear: np.ndarray, utilization_factor: float, PumpingPower: np.ndarray,
def annual_electricity_pumping_power(self, plant_lifetime: int,enduse_option: EndUseOptions, HeatExtracted: np.ndarray,
time_steps_per_year: int, utilization_factor: float, PumpingPower: np.ndarray,
ElectricityProduced: np.ndarray, NetElectricityProduced: np.ndarray, HeatProduced: np.ndarray) -> tuple:
"""
Calculate annual electricity/heat production
:param plant_lifetime: plant lifetime
:param enduse_option: end-use option
:param HeatExtracted: heat extracted
:param timestepsperyear: timesteps per year
:param time_steps_per_year: time steps per year
:param utilization_factor: utilization factor
:param PumpingPower: pumping power
:param ElectricityProduced: electricity produced
Expand All @@ -176,11 +195,13 @@ def annual_electricity_pumping_power(self, plant_lifetime: int, enduse_option: E
NetkWhProduced = np.zeros(plant_lifetime)
HeatkWhProduced = np.zeros(plant_lifetime)

def _integrate_slice(series: np.ndarray, _i: int) -> np.float64:
return SurfacePlant.integrate_time_series_slice(series, _i, time_steps_per_year, utilization_factor)

for i in range(0, plant_lifetime):
HeatkWhExtracted[i] = np.trapz(HeatExtracted[(0 + i * timestepsperyear):((i + 1) * timestepsperyear) + 1],
dx = 1. / timestepsperyear * 365. * 24.) * 1000. * utilization_factor
PumpingkWh[i] = np.trapz(PumpingPower[(0 + i * timestepsperyear):((i + 1) * timestepsperyear) + 1],
dx = 1. / timestepsperyear * 365. * 24.) * 1000. * utilization_factor
HeatkWhExtracted[i] = _integrate_slice(HeatExtracted, i)
PumpingkWh[i] = _integrate_slice(PumpingPower, i)


if enduse_option in [EndUseOptions.ELECTRICITY, EndUseOptions.COGENERATION_TOPPING_EXTRA_HEAT,
EndUseOptions.COGENERATION_TOPPING_EXTRA_ELECTRICITY,
Expand All @@ -192,16 +213,14 @@ def annual_electricity_pumping_power(self, plant_lifetime: int, enduse_option: E
TotalkWhProduced = np.zeros(plant_lifetime)
NetkWhProduced = np.zeros(plant_lifetime)
for i in range(0, plant_lifetime):
TotalkWhProduced[i] = np.trapz(ElectricityProduced[(0 + i * timestepsperyear):((i + 1) * timestepsperyear) + 1],
dx=1. / timestepsperyear * 365. * 24.) * 1000. * utilization_factor
NetkWhProduced[i] = np.trapz(NetElectricityProduced[(0 + i * timestepsperyear):((i + 1) * timestepsperyear) + 1],
dx=1. / timestepsperyear * 365. * 24.) * 1000. * utilization_factor
TotalkWhProduced[i] = _integrate_slice(ElectricityProduced, i)
NetkWhProduced[i] = _integrate_slice(NetElectricityProduced, i)

if enduse_option is not EndUseOptions.ELECTRICITY:
# all those end-use options have a direct-use component
HeatkWhProduced = np.zeros(plant_lifetime)
for i in range(0, plant_lifetime):
HeatkWhProduced[i] = np.trapz(HeatProduced[(0 + i * timestepsperyear):((i + 1) * timestepsperyear) + 1],
dx=1. / timestepsperyear * 365. * 24.) * 1000. * utilization_factor
HeatkWhProduced[i] = _integrate_slice(HeatProduced, i)

return HeatkWhExtracted, PumpingkWh, TotalkWhProduced, NetkWhProduced, HeatkWhProduced

Expand Down
21 changes: 9 additions & 12 deletions src/geophires_x/SurfacePlantAGS.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -738,26 +738,23 @@ def Calculate(self, model: Model) -> None:
self.HeatExtracted.value = self.HeatExtracted.value / 1000.0
# useful direct-use heat provided to application [MWth]
self.HeatProduced.value = self.HeatExtracted.value * self.enduseefficiencyfactor.value

def _integrate_slice(series: np.ndarray, _i: int) -> np.float64:
return SurfacePlant.integrate_time_series_slice(
series, _i, model.economics.timestepsperyear.value, self.utilization_factor.value
)

for i in range(0, self.plant_lifetime.value):
self.HeatkWhExtracted.value[i] = np.trapz(self.HeatExtracted.value[
(i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.PumpingkWh.value[i] = np.trapz(model.wellbores.PumpingPower.value[
(i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.HeatkWhExtracted.value[i] = _integrate_slice(self.HeatExtracted.value, i)
self.PumpingkWh.value[i] = _integrate_slice(model.wellbores.PumpingPower.value, i)

self.RemainingReservoirHeatContent.value = model.reserv.InitialReservoirHeatContent.value - np.cumsum(
self.HeatkWhExtracted.value) * 3600 * 1E3 / 1E15

if self.End_use is not EndUseOptions.ELECTRICITY:
self.HeatkWhProduced.value = np.zeros(self.plant_lifetime.value)
for i in range(0, self.plant_lifetime.value):
self.HeatkWhProduced.value[i] = np.trapz(self.HeatProduced.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.HeatkWhProduced.value[i] = _integrate_slice(self.HeatProduced.value, i)
else:
# copy some arrays so we have a GEOPHIRES equivalent
self.TotalkWhProduced.value = self.Annual_electricity_production.copy()
Expand Down
25 changes: 9 additions & 16 deletions src/geophires_x/SurfacePlantAbsorptionChiller.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -114,29 +114,22 @@ def Calculate(self, model: Model) -> None:
self.HeatkWhExtracted.value = np.zeros(self.plant_lifetime.value)
self.PumpingkWh.value = np.zeros(self.plant_lifetime.value)

def _integrate_slice(series, _i):
return SurfacePlant.integrate_time_series_slice(
series, _i, model.economics.timestepsperyear.value, self.utilization_factor.value
)

for i in range(0, self.plant_lifetime.value):
self.HeatkWhExtracted.value[i] = np.trapz(self.HeatExtracted.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.PumpingkWh.value[i] = np.trapz(model.wellbores.PumpingPower.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.HeatkWhExtracted.value[i] = _integrate_slice(self.HeatExtracted.value, i)
self.PumpingkWh.value[i] = _integrate_slice(model.wellbores.PumpingPower.value, i)

self.HeatkWhProduced.value = np.zeros(self.plant_lifetime.value)
for i in range(0, self.plant_lifetime.value):
self.HeatkWhProduced.value[i] = np.trapz(self.HeatProduced.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.HeatkWhProduced.value[i] = _integrate_slice(self.HeatProduced.value, i)

self.cooling_kWh_Produced.value = np.zeros(self.plant_lifetime.value)
for i in range(0, self.plant_lifetime.value):
self.cooling_kWh_Produced.value[i] = np.trapz(self.cooling_produced.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.cooling_kWh_Produced.value[i] = _integrate_slice(self.cooling_produced.value, i)

# calculate reservoir heat content
self.RemainingReservoirHeatContent.value = SurfacePlant.remaining_reservoir_heat_content(
Expand Down
50 changes: 26 additions & 24 deletions src/geophires_x/SurfacePlantDistrictHeating.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -215,35 +215,37 @@ def Calculate(self, model: Model) -> None:
self.HeatkWhExtracted.value = np.zeros(self.plant_lifetime.value)
self.PumpingkWh.value = np.zeros(self.plant_lifetime.value)

def _integrate_slice(series, _i, util_factor):
return SurfacePlant.integrate_time_series_slice(
series, _i, model.economics.timestepsperyear.value, util_factor
)


for i in range(0, self.plant_lifetime.value):
if self.plant_type.value == PlantType.DISTRICT_HEATING: # for district heating, we have a util_factor_array
self.HeatkWhExtracted.value[i] = np.trapz(self.HeatExtracted.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * \
self.util_factor_array.value[i]
self.PumpingkWh.value[i] = np.trapz(model.wellbores.PumpingPower.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * \
self.util_factor_array.value[i]
if self.plant_type.value == PlantType.DISTRICT_HEATING:
self.HeatkWhExtracted.value[i] = _integrate_slice(
self.HeatExtracted.value,
i,
self.util_factor_array.value[i]
)

self.PumpingkWh.value[i] = _integrate_slice(
model.wellbores.PumpingPower.value,
i,
# for district heating, we have a util_factor_array
self.util_factor_array.value[i]
)
else:
self.HeatkWhExtracted.value[i] = np.trapz(self.HeatExtracted.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.PumpingkWh.value[i] = np.trapz(model.wellbores.PumpingPower.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * self.utilization_factor.value
self.HeatkWhExtracted.value[i] = _integrate_slice(self.HeatExtracted.value, i, self.utilization_factor.value)
self.PumpingkWh.value[i] = _integrate_slice(model.wellbores.PumpingPower.value, i, self.utilization_factor.value)

self.HeatkWhProduced.value = np.zeros(self.plant_lifetime.value)
for i in range(0, self.plant_lifetime.value):
self.HeatkWhProduced.value[i] = np.trapz(self.HeatProduced.value[
(0 + i * model.economics.timestepsperyear.value):((
i + 1) * model.economics.timestepsperyear.value) + 1],
dx=1. / model.economics.timestepsperyear.value * 365. * 24.) * 1000. * \
self.util_factor_array.value[i]
self.HeatkWhProduced.value[i] = _integrate_slice(
self.HeatProduced.value,
i,
self.util_factor_array.value[i]
)

# calculate reservoir heat content
self.RemainingReservoirHeatContent.value = SurfacePlant.remaining_reservoir_heat_content(
Expand Down
Loading