Implementing a way to handle overpressure

environments
malcolm-dsider · Apr 15, 2024 · 322b27a · 322b27a
1 parent 0d0ebeb
commit 322b27a
Show file tree

Hide file tree

Showing 40 changed files with 738 additions and 235 deletions.
diff --git a/src/geophires_x/AGSWellBores.py b/src/geophires_x/AGSWellBores.py
@@ -871,12 +871,12 @@ def CalculateNonverticalPressureDrop(self, model:Model, time_operation: float, t
             # nonvertical wellbore fluid conditions based on current temperature
             rhowater = density_water_kg_per_m3(
                 self.NonverticalProducedTemperature.value[year],
-                pressure=model.reserv.lithostatic_pressure()
+                pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
             )
 
             muwater = viscosity_water_Pa_sec(
                 self.NonverticalProducedTemperature.value[year],
-                pressure=model.reserv.lithostatic_pressure()
+                pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
             )
             vhoriz = self.q_circulation / rhowater / (math.pi / 4. * self.nonverticalwellborediameter.value ** 2)
 
@@ -957,15 +957,15 @@ def Calculate(self, model: Model) -> None:
                 # MIR figure out how to calculate year and extract Tini from reserv Tresoutput array
                 year = math.trunc(self.time_operation.value / self.al)
                 self.NonverticalProducedTemperature.value[year] = inverselaplace(
-                    self, 16, 0, model.reserv.lithostatic_pressure())
+                    self, 16, 0, model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value))
                 # update alpha_fluid value based on next temperature of reservoir
 
                 self.alpha_fluid = self.WaterThermalConductivity.value / density_water_kg_per_m3(
                     self.NonverticalProducedTemperature.value[year],
-                    pressure=model.reserv.lithostatic_pressure()
+                    pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
                 ) / heat_capacity_water_J_per_kg_per_K(
                     self.NonverticalProducedTemperature.value[year],
-                    pressure=model.reserv.lithostatic_pressure()
+                    pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
                 ) * 24.0 * 3600.0
                 self.time_operation.value += self.al
 
@@ -979,7 +979,7 @@ def Calculate(self, model: Model) -> None:
             self.ProdTempDrop.value = self.tempdropprod.value
             model.reserv.cpwater.value = heat_capacity_water_J_per_kg_per_K(
                 self.NonverticalProducedTemperature.value[0],
-                pressure=model.reserv.lithostatic_pressure()
+                pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
             )
             if self.rameyoptionprod.value:
                 self.ProdTempDrop.value = RameyCalc(model.reserv.krock.value,
@@ -1002,13 +1002,13 @@ def Calculate(self, model: Model) -> None:
             if self.productionwellpumping.value:
                 self.rhowaterinj = density_water_kg_per_m3(
                     model.reserv.Tsurf.value,
-                    pressure=model.reserv.lithostatic_pressure()
+                    pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
                 ) * np.linspace(1, 1,
                                 len(self.ProducedTemperature.value))
 
                 self.rhowaterprod = density_water_kg_per_m3(
                     model.reserv.Trock.value,
-                    pressure=model.reserv.lithostatic_pressure()
+                    pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
                 ) * np.linspace(1, 1, len(self.ProducedTemperature.value))
 
                 self.DPProdWell.value, f3, vprod, self.rhowaterprod = WellPressureDrop(model,
@@ -1113,13 +1113,13 @@ def Calculate(self, model: Model) -> None:
 
             rho_water = density_water_kg_per_m3(
                 self.Tout[0],
-                pressure = model.reserv.lithostatic_pressure(),
+                pressure = model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value),
             )
 
 
             model.reserv.cpwater.value = heat_capacity_water_J_per_kg_per_K(
                 self.Tout[0],
-                pressure=model.reserv.lithostatic_pressure(),
+                pressure=model.reserv.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value),
             )  # Need this for surface plant output calculation
 
             # set pumping power to zero for all times, assuming that the thermosphere wil always

diff --git a/src/geophires_x/CylindricalReservoir.py b/src/geophires_x/CylindricalReservoir.py
@@ -244,20 +244,20 @@ def Calculate(self, model: Model) -> None:
                                                  ) / 1e15  # 10^15 J
         self.cpwater.value = heat_capacity_water_J_per_kg_per_K(
             model.wellbores.Tinj.value * 0.5 + (self.Trock.value * 0.9 + model.wellbores.Tinj.value * 0.1) * 0.5,
-            pressure=model.reserv.lithostatic_pressure()
+            pressure=self.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
         )
         self.rhowater.value = density_water_kg_per_m3(
             model.wellbores.Tinj.value * 0.5 + (self.Trock.value * 0.9 + model.wellbores.Tinj.value * 0.1) * 0.5,
-            pressure=model.reserv.lithostatic_pressure()
+            pressure=self.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.Trock.value)
         )
 
         model.logger.info(f'complete {str(__class__)}: {sys._getframe().f_code.co_name}')
 
-    def lithostatic_pressure(self) -> PlainQuantity:
+    #def lithostatic_pressure(self) -> PlainQuantity:
         """
         @override
 
         Standard reservoir implementation uses depth but CylindricalReservoir sets depth to total drilled length
         """
-        return quantity(static_pressure_MPa(self.rhorock.quantity().to('kg/m**3').magnitude,
-                                            self.InputDepth.quantity().to('m').magnitude), 'MPa')
+        def lithostatic_pressure(self, rho_rock: float, depth: float) -> PlainQuantity:
+            return quantity(static_pressure_MPa(rho_rock, depth), 'MPa')
diff --git a/src/geophires_x/Outputs.py b/src/geophires_x/Outputs.py
@@ -1005,9 +1005,9 @@ def PrintOutputs(self, model: Model):
                 reservoir_parameters.append(OutputTableItem('Reservoir impedance', '{0:10.2f}'.format(model.wellbores.impedance.value / 1000),
                                     model.wellbores.impedance.CurrentUnits.value))
             else:
-                reservoir_parameters.append(OutputTableItem('Reservoir hydrostatic pressure',
-                                                            '{0:10.2f}'.format(model.wellbores.Phydrostaticcalc.value),
-                                                            model.wellbores.Phydrostaticcalc.CurrentUnits.value))
+                reservoir_parameters.append(OutputTableItem('Average reservoir pressure',
+                                                            '{0:10.2f}'.format(np.average(model.wellbores.production_reservoir_pressure.value)),
+                                                            model.wellbores.production_reservoir_pressure.CurrentUnits.value))
                 reservoir_parameters.append(OutputTableItem('Plant outlet pressure', '{0:10.2f}'.format(
                     model.surfaceplant.plant_outlet_pressure.value),
                                                             model.surfaceplant.plant_outlet_pressure.CurrentUnits.value))
@@ -1672,7 +1672,7 @@ def PrintOutputs(self, model: Model):
                     if model.wellbores.impedancemodelused.value:
                         f.write(f'      Reservoir impedance:                              {model.wellbores.impedance.value/1000:10.2f} ' + model.wellbores.impedance.CurrentUnits.value + NL)
                     else:
-                        f.write(f'      Reservoir hydrostatic pressure:                   {model.wellbores.Phydrostaticcalc.value:10.2f} ' + model.wellbores.Phydrostaticcalc.CurrentUnits.value + NL)
+                        f.write(f'      Average reservoir pressure:                       {np.average(model.wellbores.production_reservoir_pressure.value):10.2f} ' + model.wellbores.production_reservoir_pressure.CurrentUnits.value + NL)
                         f.write(f'      Plant outlet pressure:                            {model.surfaceplant.plant_outlet_pressure.value:10.2f} ' + model.surfaceplant.plant_outlet_pressure.CurrentUnits.value + NL)
                         if model.wellbores.productionwellpumping.value:
                             f.write(f'      Production wellhead pressure:                     {model.wellbores.Pprodwellhead.value:10.2f} ' + model.wellbores.Pprodwellhead.CurrentUnits.value + NL)

diff --git a/src/geophires_x/Reservoir.py b/src/geophires_x/Reservoir.py
@@ -749,12 +749,14 @@ def Calculate(self, model: Model) -> None:
             # calculate reservoir water properties
             self.cpwater.value = heat_capacity_water_J_per_kg_per_K(
                 model.wellbores.Tinj.value * 0.5 + (self.Trock.value * 0.9 + model.wellbores.Tinj.value * 0.1) * 0.5,
-                pressure=self.lithostatic_pressure()
+                pressure=self.lithostatic_pressure(self.rhorock.quantity().to('kg/m**3').magnitude,
+                                            self.depth.quantity().to('m').magnitude)
             )
 
             self.rhowater.value = density_water_kg_per_m3(
                 model.wellbores.Tinj.value * 0.5 + (self.Trock.value * 0.9 + model.wellbores.Tinj.value * 0.1) * 0.5,
-                pressure=self.lithostatic_pressure()
+                pressure=self.lithostatic_pressure(self.rhorock.quantity().to('kg/m**3').magnitude,
+                                            self.depth.quantity().to('m').magnitude)
             )
 
             # temperature gain in injection wells
@@ -766,8 +768,7 @@ def Calculate(self, model: Model) -> None:
 
         model.logger.info(f'complete {str(__class__)}: {sys._getframe().f_code.co_name}')
 
-    def lithostatic_pressure(self) -> PlainQuantity:
-        return quantity(static_pressure_MPa(self.rhorock.quantity().to('kg/m**3').magnitude,
-                                            self.depth.quantity().to('m').magnitude), 'MPa')
+    def lithostatic_pressure(self, rho_rock: float, depth: float) -> PlainQuantity:
+        return quantity(static_pressure_MPa(rho_rock, depth), 'MPa')
 
 
diff --git a/src/geophires_x/Units.py b/src/geophires_x/Units.py
@@ -58,6 +58,8 @@ class Units(IntEnum):
     POWERPERUNITAREA = auto()
     HEATPERUNITVOLUME = auto()
     POWERPERUNITVOLUME = auto()
+    DECAY_RATE=auto()
+    INFLATION_RATE=auto()
 
 
 class TemperatureUnit(str, Enum):
@@ -129,7 +131,7 @@ class EnergyUnit(str, Enum):
 
 
 class PowerUnit(str, Enum):
-    """Power (electrcity or heat) Units"""
+    """Power (electricity or heat) Units"""
     W = "W"
     KW = "kW"
     MW = "MW"
@@ -338,3 +340,13 @@ class HeatPerUnitVolumeUnit(str,Enum):
 class PowerPerUnitVolumeUnit(str,Enum):
     """Population Density Units"""
     MWPERCUBICKM = "MW/km**3"
+
+
+class Decay_RateUnit(str,Enum):
+    """Decay rate Units"""
+    PERCENTPERYEAR = "%/yr"
+
+
+class Inflation_RateUnit(str,Enum):
+    """Decay rate Units"""
+    KPASCALPERYEAR = "kPa/yr"