Merge pull request NREL#284 from malcolm-dsider/sbt-integration

SBT integration

Author: malcolm-dsider

src/geophires_x/Economics.py

@@ -1,56 +1,13 @@
 import math
 import sys
-import os
 import numpy as np
 import numpy_financial as npf
 import geophires_x.Model as Model
 from geophires_x.OptionList import Configuration, WellDrillingCostCorrelation, EconomicModel, EndUseOptions, PlantType
 from geophires_x.Parameter import intParameter, floatParameter, OutputParameter, ReadParameter, boolParameter, \
     coerce_int_params_to_enum_values
 from geophires_x.Units import *
-
-
-def calculate_total_drilling_lengths_m(Configuration, numnonverticalsections: int, nonvertical_length_km: float,
-                                       InputDepth_km: float, OutputDepth_km: float, nprod:int, ninj:int) -> tuple:
-    """
-    returns the total length, vertical length, and non-vertical lengths, depending on the configuration
-    :param Configuration: Configuration of the well
-    :type Configuration: :class:`~geophires
-    :param numnonverticalsections: number of non-vertical sections
-    :type numnonverticalsections: int
-    :param nonvertical_length_km: length of non-vertical sections in km
-    :type nonvertical_length_km: float
-    :param InputDepth_km: depth of the well in km
-    :type InputDepth_km: float
-    :param OutputDepth_km: depth of the output end of the well in km, if U shaped, and not horizontal
-    :type OutputDepth_km: float
-    :param nprod: number of production wells
-    :type nprod: int
-    :param ninj: number of injection wells
-    :return: total length, vertical length, and horizontal lengths in meters
-    :rtype: tuple
-    """
-    if Configuration == Configuration.ULOOP:
-        # Total drilling depth of both wells and laterals in U-loop [m]
-        vertical_pipe_length_m = (nprod * InputDepth_km * 1000.0) + (ninj * OutputDepth_km * 1000.0)
-        nonvertical_pipe_length_m = numnonverticalsections * nonvertical_length_km * 1000.0
-    elif Configuration == Configuration.COAXIAL:
-        # Total drilling depth of well and lateral in co-axial case [m]
-        vertical_pipe_length_m = (nprod + ninj) * InputDepth_km * 1000.0
-        nonvertical_pipe_length_m = numnonverticalsections * nonvertical_length_km * 1000.0
-    elif Configuration == Configuration.VERTICAL:
-        # Total drilling depth of well in vertical case [m]
-        vertical_pipe_length_m = (nprod + ninj) * InputDepth_km * 1000.0
-        nonvertical_pipe_length_m = 0.0
-    elif Configuration == Configuration.L:
-        # Total drilling depth of well in L case [m]
-        vertical_pipe_length_m = (nprod + ninj) * InputDepth_km * 1000.0
-        nonvertical_pipe_length_m = numnonverticalsections * nonvertical_length_km * 1000.0
-    else:
-        raise ValueError(f'Invalid Configuration: {Configuration}')
-
-    tot_pipe_length_m = vertical_pipe_length_m + nonvertical_pipe_length_m
-    return tot_pipe_length_m, vertical_pipe_length_m, nonvertical_pipe_length_m
+from geophires_x.WellBores import calculate_total_drilling_lengths_m
 
 
 def calculate_cost_of_one_vertical_well(model: Model, depth_m: float, well_correlation: int,
@@ -1809,8 +1766,14 @@ def __init__(self, model: Model):
             PreferredUnits=CurrencyUnit.MDOLLARS,
             CurrentUnits=CurrencyUnit.MDOLLARS
         )
-        self.cost_nonvertical_section = self.OutputParameterDict[self.cost_nonvertical_section.Name] = OutputParameter(
-            Name="Cost of the non-vertical section of a well",
+        self.cost_lateral_section = self.OutputParameterDict[self.cost_lateral_section.Name] = OutputParameter(
+            Name="Cost of the entire (multi-) lateral section of a well",
+            UnitType=Units.CURRENCY,
+            PreferredUnits=CurrencyUnit.MDOLLARS,
+            CurrentUnits=CurrencyUnit.MDOLLARS
+        )
+        self.cost_to_junction_section = self.OutputParameterDict[self.cost_to_junction_section.Name] = OutputParameter(
+            Name="Cost of the entire section of a well from bottom of vertical to junction with laterals",
             UnitType=Units.CURRENCY,
             PreferredUnits=CurrencyUnit.MDOLLARS,
             CurrentUnits=CurrencyUnit.MDOLLARS
@@ -2220,7 +2183,7 @@ def Calculate(self, model: Model) -> None:
                                 (self.cost_one_injection_well.value * model.wellbores.ninj.value))
         else:
             if hasattr(model.wellbores, 'numnonverticalsections') and model.wellbores.numnonverticalsections.Provided:
-                self.cost_nonvertical_section.value = 0.0
+                self.cost_lateral_section.value = 0.0
                 if not model.wellbores.IsAGS.value:
                     input_vert_depth_km = model.reserv.depth.quantity().to('km').magnitude
                     output_vert_depth_km = 0.0
@@ -2229,13 +2192,13 @@ def Calculate(self, model: Model) -> None:
                     output_vert_depth_km = model.reserv.OutputDepth.quantity().to('km').magnitude
                 model.wellbores.injection_reservoir_depth.value = input_vert_depth_km
 
-                tot_m, tot_vert_m, tot_horiz_m = calculate_total_drilling_lengths_m(model.wellbores.Configuration.value,
-                                                                      model.wellbores.numnonverticalsections.value,
-                                                                      model.wellbores.Nonvertical_length.value / 1000.0,
-                                                                      input_vert_depth_km,
-                                                                      output_vert_depth_km,
-                                                                      model.wellbores.nprod.value,
-                                                                      model.wellbores.ninj.value)
+                tot_m, tot_vert_m, tot_horiz_m, _ = calculate_total_drilling_lengths_m(model.wellbores.Configuration.value,
+                                                                                    model.wellbores.numnonverticalsections.value,
+                                                                                    model.wellbores.Nonvertical_length.value / 1000.0,
+                                                                                    input_vert_depth_km,
+                                                                                    output_vert_depth_km,
+                                                                                    model.wellbores.nprod.value,
+                                                                                    model.wellbores.ninj.value)
 
             else:
                 tot_m = tot_vert_m = model.reserv.depth.quantity().to('km').magnitude
@@ -2261,20 +2224,20 @@ def Calculate(self, model: Model) -> None:
                                                                                          self.injection_well_cost_adjustment_factor.value)
 
             if hasattr(model.wellbores, 'numnonverticalsections') and model.wellbores.numnonverticalsections.Provided:
-                self.cost_nonvertical_section.value = calculate_cost_of_non_vertical_section(model, tot_horiz_m,
+                self.cost_lateral_section.value = calculate_cost_of_non_vertical_section(model, tot_horiz_m,
                                             self.wellcorrelation.value,
                                             self.Nonvertical_drilling_cost_per_m.value,
                                             model.wellbores.numnonverticalsections.value,
                                             self.per_injection_well_cost.Name,
                                             model.wellbores.NonverticalsCased.value,
                                             self.production_well_cost_adjustment_factor.value)
             else:
-                self.cost_nonvertical_section.value = 0.0
+                self.cost_lateral_section.value = 0.0
             # cost of the well field
             # 1.05 for 5% indirect costs
             self.Cwell.value = 1.05 * ((self.cost_one_production_well.value * model.wellbores.nprod.value) +
                                           (self.cost_one_injection_well.value * model.wellbores.ninj.value) +
-                                          self.cost_nonvertical_section.value)
+                                          self.cost_lateral_section.value)
 
         # reservoir stimulation costs (M$/injection well). These are calculated whether totalcapcost.Valid = 1
         if self.ccstimfixed.Valid:

src/geophires_x/Model.py

@@ -11,6 +11,9 @@
 from geophires_x.TDPReservoir import TDPReservoir
 from geophires_x.WellBores import WellBores
 from geophires_x.SurfacePlant import SurfacePlant
+from geophires_x.SBTEconomics import SBTEconomics
+from geophires_x.SBTWellbores import SBTWellbores
+from geophires_x.SBTReservoir import SBTReservoir
 from geophires_x.SurfacePlantIndustrialHeat import SurfacePlantIndustrialHeat
 from geophires_x.SurfacePlantSubcriticalORC import SurfacePlantSubcriticalOrc
 from geophires_x.SurfacePlantSupercriticalORC import SurfacePlantSupercriticalOrc
@@ -72,68 +75,107 @@ def __init__(self, enable_geophires_logging_config=True, input_file=None):
         if input_file is None and len(sys.argv) > 1:
             input_file = sys.argv[1]
 
+        # Key step - read the entire provided input file
         read_input_file(self.InputParameters, logger=self.logger, input_file_name=input_file)
 
+        # initiate the outputs object
+        output_file = 'HDR.out'
+        if len(sys.argv) > 2:
+            output_file = sys.argv[2]
+        self.outputs = Outputs(self, output_file=output_file)
+
+        # Initiate the elements of the Model object
+        # this is where you can change what class get initiated - the superclass, or one of the subclasses
+        self.logger.info("Initiate the elements of the Model")
+
+        # Assume that SDAC and add-ons are not used
         self.sdacgtoutputs = None
         self.sdacgteconomics = None
         self.addoutputs = None
         self.addeconomics = None
 
-        # Initiate the elements of the Model
-        # this is where you can change what class get initiated - the superclass, or one of the subclasses
-        self.logger.info("Initiate the elements of the Model")
-        # we need to decide which reservoir to instantiate based on the user input (InputParameters),
-        # which we just read above for the first time
-        # Default is Thermal drawdown percentage model (GETEM)
-        self.reserv: TDPReservoir = TDPReservoir(self)
-        if 'Reservoir Model' in self.InputParameters:
-            if self.InputParameters['Reservoir Model'].sValue == '0':
-                self.reserv: CylindricalReservoir = CylindricalReservoir(self)  # Simple Cylindrical Reservoir
-            elif self.InputParameters['Reservoir Model'].sValue == '1':
-                self.reserv: MPFReservoir = MPFReservoir(self)  # Multiple parallel fractures model (LANL)
-            elif self.InputParameters['Reservoir Model'].sValue == '2':
-                self.reserv: LHSReservoir = LHSReservoir(self)  # Multiple parallel fractures model (LANL)
-            elif self.InputParameters['Reservoir Model'].sValue == '3':
-                self.reserv: SFReservoir = SFReservoir(self)  # Drawdown parameter model (Tester)
-            elif self.InputParameters['Reservoir Model'].sValue == '5':
-                self.reserv: UPPReservoir = UPPReservoir(self)  # Generic user-provided temperature profile
-            elif self.InputParameters['Reservoir Model'].sValue == '6':
-                self.reserv: TOUGH2Reservoir = TOUGH2Reservoir(self)  # Tough2 is called
-            elif self.InputParameters['Reservoir Model'].sValue == '7':
-                self.reserv: SUTRAReservoir = SUTRAReservoir(self)  # SUTRA output is created
-
         # initialize the default objects
+        self.reserv: TDPReservoir = TDPReservoir(self)
         self.wellbores: WellBores = WellBores(self)
-        self.surfaceplant: SurfacePlant = SurfacePlant(self)
+        self.surfaceplant = SurfacePlantIndustrialHeat(self)  # default is Industrial Heat
         self.economics: Economics = Economics(self)
 
-        output_file = 'HDR.out'
-        if len(sys.argv) > 2:
-            output_file = sys.argv[2]
-
-        self.outputs = Outputs(self, output_file=output_file)
+        # Now we need to handle the creation of all the special case objects based on the user settings
+        # We have to access the user setting from the overall master table because the read_parameters functions
+        # have not been called on the specific objects yet, so the instantiated objects only contain default values
+        # not user set values. The user set value can only come from the master dictionary "InputParameters"
+        # based on the user input, which we just read above for the first time
 
+        # First, we need to decide which reservoir to instantiate
+        # For reservoirs, the default is Thermal drawdown percentage model (GETEM); see above where it is initialized.
+        # The user can change this in the input file, so check the values and initiate the appropriate reservoir object
         if 'Reservoir Model' in self.InputParameters:
-            if self.InputParameters['Reservoir Model'].sValue == '7':
-                # if we use SUTRA output for simulating reservoir thermal energy storage, we use a special wellbore object that can handle SUTRA data
+            if self.InputParameters['Reservoir Model'].sValue in ['0', 'Simple cylindrical']:
+                self.reserv: CylindricalReservoir = CylindricalReservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['1', 'Multiple Parallel Fractures']:
+                self.reserv: MPFReservoir = MPFReservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['2', '1-D Linear Heat Sweep']:
+                self.reserv: LHSReservoir = LHSReservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['3', 'Single Fracture m/A Thermal Drawdown']:
+                self.reserv: SFReservoir = SFReservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['5', 'User-Provided Temperature Profile']:
+                self.reserv: UPPReservoir = UPPReservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['6', 'TOUGH2 Simulator']:
+                self.reserv: TOUGH2Reservoir = TOUGH2Reservoir(self)
+            elif self.InputParameters['Reservoir Model'].sValue in ['7', 'SUTRA']:
+                # if we use SUTRA output for simulating reservoir thermal energy storage,
+                # we use a special wellbore object that handles SUTRA data, and special Economics and Outputs objects
+                self.logger.info('Setup the SUTRA elements of the Model and instantiate new attributes as needed')
+                self.reserv: SUTRAReservoir = SUTRAReservoir(self)
                 self.wellbores: WellBores = SUTRAWellBores(self)
                 self.surfaceplant: SurfacePlantSUTRA = SurfacePlantSUTRA(self)
                 self.economics: SUTRAEconomics = SUTRAEconomics(self)
                 self.outputs: SUTRAOutputs = SUTRAOutputs(self, output_file=output_file)
+            elif self.InputParameters['Reservoir Model'].sValue in ['8', 'SBT']:
+                self.logger.info('Setup the SBT elements of the Model and instantiate new attributes as needed')
+                self.reserv: SBTReservoir = SBTReservoir(self)
+                self.wellbores: SBTWellBores = SBTWellbores(self)
+                self.economics: SBTEconomics = SBTEconomics(self)
 
+        # Now handle the special cases for all AGS cases (CLGS, SBT, or CLGS)
         if 'Is AGS' in self.InputParameters:
             if self.InputParameters['Is AGS'].sValue in ['True', 'true', 'TRUE', 'T', '1']:
-                self.logger.info("Setup the AGS elements of the Model and instantiate new attributes as needed")
-                # If we are doing AGS, we need to replace the various objects we with versions of the objects
-                # that have AGS functionality.
-                # that means importing them, initializing them, then reading their parameters
-                # use the simple cylindrical reservoir for all AGS systems.
-                self.reserv: CylindricalReservoir = CylindricalReservoir(self)
-                self.wellbores: WellBores = AGSWellBores(self)
-                self.surfaceplant: SurfacePlantAGS = SurfacePlantAGS(self)
-                self.economics: AGSEconomics = AGSEconomics(self)
-                self.outputs: AGSOutputs = AGSOutputs(self, output_file=output_file)
+                self.logger.info('Setup the AGS elements of the Model and instantiate new attributes as needed')
                 self.wellbores.IsAGS.value = True
+                if not isinstance(self.reserv, SBTReservoir):
+                    if self.InputParameters['Economic Model'].sValue not in ['4', 'Simple (CLGS)']:
+                        # must be doing wangju approach, # so go back to using  default objects
+                        self.surfaceplant = SurfacePlant(self)
+                        self.economics = Economics(self)
+                # Must be doing CLGS, so we need to instantiate the right objects
+                    self.reserv: CylindricalReservoir = CylindricalReservoir(self)
+                    self.wellbores: WellBores = AGSWellBores(self)
+                    self.surfaceplant: SurfacePlantAGS = SurfacePlantAGS(self)
+                    self.economics: AGSEconomics = AGSEconomics(self)
+                    self.outputs: AGSOutputs = AGSOutputs(self, output_file=output_file)
+
+        # initialize the right Power Plant Type
+        if 'Power Plant Type' in self.InputParameters:
+            # electricity
+            if self.InputParameters['Power Plant Type'].sValue in ['1', 'Subcritical ORC']:
+                self.surfaceplant = SurfacePlantSubcriticalOrc(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['2', 'Supercritical ORC']:
+                self.surfaceplant = SurfacePlantSupercriticalOrc(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['3', 'Single-Flash']:
+                self.surfaceplant = SurfacePlantSingleFlash(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['4', 'Double-Flash']:
+                self.surfaceplant = SurfacePlantDoubleFlash(self)
+            # Heat applications
+            elif self.InputParameters['Power Plant Type'].sValue in ['5', 'Absorption Chiller']:
+                self.surfaceplant = SurfacePlantAbsorptionChiller(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['6', 'Heat Pump']:
+                self.surfaceplant = SurfacePlantHeatPump(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['7', 'District Heating']:
+                self.surfaceplant = SurfacePlantDistrictHeating(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['8', 'Reservoir Thermal Energy Storage']:
+                self.surfaceplant = SurfacePlantSUTRA(self)
+            elif self.InputParameters['Power Plant Type'].sValue in ['9', 'Industrial']:
+                self.surfaceplant = SurfacePlantIndustrialHeat(self)
 
         # if we find out we have an add-ons, we need to instantiate it, then read for the parameters
         if 'AddOn Nickname 1' in self.InputParameters:
@@ -150,6 +192,7 @@ def __init__(self, enable_geophires_logging_config=True, input_file=None):
 
         self.logger.info(f'Complete {__class__}: {__name__}')
 
+
     def __str__(self):
         return "Model"