Merge pull request NREL#185 from malcolm-dsider/main

Adding overpressure calculation

Author: malcolm-dsider

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.InputDepth.quantity().to('m').magnitude)
             )
 
             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.InputDepth.quantity().to('m').magnitude)
             )
             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.InputDepth.quantity().to('m').magnitude))
                 # 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.InputDepth.quantity().to('m').magnitude)
                 ) / 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.InputDepth.quantity().to('m').magnitude)
                 ) * 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.InputDepth.quantity().to('m').magnitude)
             )
             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.InputDepth.quantity().to('m').magnitude)
                 ) * 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.InputDepth.quantity().to('m').magnitude)
                 ) * 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.InputDepth.quantity().to('m').magnitude),
             )
 
 
             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.InputDepth.quantity().to('m').magnitude),
             )  # Need this for surface plant output calculation
 
             # set pumping power to zero for all times, assuming that the thermosphere wil always

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.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=model.reserv.lithostatic_pressure()
+            pressure=self.lithostatic_pressure(model.reserv.rhorock.value, model.reserv.depth.quantity().to('m').magnitude)
         )
 
         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_kg_per_m3: float, depth_m: float) -> PlainQuantity:
+            return quantity(static_pressure_MPa(rho_rock_kg_per_m3, depth_m), 'MPa')