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

.bumpversion.cfg

@@ -1,5 +1,5 @@
 [bumpversion]
-current_version = 3.7.21
+current_version = 3.7.23
 commit = True
 tag = True
 

.cookiecutterrc

@@ -54,7 +54,7 @@ default_context:
     sphinx_doctest: "no"
     sphinx_theme: "sphinx-py3doc-enhanced-theme"
     test_matrix_separate_coverage: "no"
-    version: 3.7.21
+    version: 3.7.23
     version_manager: "bump2version"
     website: "https://github.com/NREL"
     year_from: "2023"

README.rst

@@ -56,9 +56,9 @@ Free software: `MIT license <LICENSE>`__
     :alt: Supported implementations
     :target: https://pypi.org/project/geophires-x
 
-.. |commits-since| image:: https://img.shields.io/github/commits-since/softwareengineerprogrammer/GEOPHIRES-X/v3.7.21.svg
+.. |commits-since| image:: https://img.shields.io/github/commits-since/softwareengineerprogrammer/GEOPHIRES-X/v3.7.23.svg
     :alt: Commits since latest release
-    :target: https://github.com/softwareengineerprogrammer/GEOPHIRES-X/compare/v3.7.21...main
+    :target: https://github.com/softwareengineerprogrammer/GEOPHIRES-X/compare/v3.7.23...main
 
 .. |docs| image:: https://readthedocs.org/projects/GEOPHIRES-X/badge/?style=flat
     :target: https://nrel.github.io/GEOPHIRES-X

docs/conf.py

@@ -18,7 +18,7 @@
 year = '2025'
 author = 'NREL'
 copyright = f'{year}, {author}'
-version = release = '3.7.21'
+version = release = '3.7.23'
 
 pygments_style = 'trac'
 templates_path = ['./templates']

setup.py

@@ -13,7 +13,7 @@ def read(*names, **kwargs):
 
 setup(
     name='geophires-x',
-    version='3.7.21',
+    version='3.7.23',
     license='MIT',
     description='GEOPHIRES is a free and open-source geothermal techno-economic simulator.',
     long_description='{}\n{}'.format(

src/geophires_x/SurfacePlant.py

@@ -1,16 +1,43 @@
-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_start_index = _i * time_steps_per_year
+        slice_end_index = ((_i + 1) * time_steps_per_year) + 1
+        _slice = list(series[slice_start_index:slice_end_index])
+
+        # Note that len(_slice) - 1 may be less than time_steps_per_year for the last slice.
+
+        if len(_slice) == 1:
+            extrapolated_future_datapoint = _slice[0]
+            if slice_start_index - 1 > 0:
+                delta = series[slice_start_index] - series[slice_start_index - 1]
+                extrapolated_future_datapoint = _slice[0] + delta
+            _slice.append(extrapolated_future_datapoint)
+
+        dx_steps = len(_slice) - 1
+
+        integral = np.trapz(
+            _slice,
+            dx=1. / dx_steps * 365. * 24.
+        )
+
+        return integral * 1000. * utilization_factor
+
     def remaining_reservoir_heat_content(self, InitialReservoirHeatContent: np.ndarray, HeatkWhExtracted:  np.ndarray) -> np.ndarray:
         """
         Calculate reservoir heat content
@@ -152,15 +179,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
@@ -176,11 +203,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,
@@ -192,16 +221,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
 

src/geophires_x/SurfacePlantAGS.py

@@ -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()

src/geophires_x/SurfacePlantAbsorptionChiller.py

@@ -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(

src/geophires_x/SurfacePlantDistrictHeating.py

@@ -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(