Fix unit conversion issue that was causing incorrect royalty calculations

Author: softwareengineerprogrammer

src/geophires_x/Economics.py

@@ -3363,18 +3363,31 @@ def _calculate_sam_economics(self, model: Model) -> None:
                 self.sam_economics_calculations.royalties_opex.value[pre_revenue_years_slice_index:]
             )
 
-            self.royalties_average_annual_cost.value = average_annual_royalties  # TODO unit conversion
-            self.Coam.value += self.royalties_average_annual_cost.quantity().to(self.Coam.CurrentUnits.value).magnitude
+            self.royalties_average_annual_cost.value = (quantity(
+                average_annual_royalties,
+                self.sam_economics_calculations.royalties_opex.CurrentUnits
+            ).to(self.royalties_average_annual_cost.CurrentUnits).magnitude)
+
+            self.Coam.value += (self.royalties_average_annual_cost.quantity()
+                                .to(self.Coam.CurrentUnits.value).magnitude)
+
+            self.royalty_holder_npv.value = quantity(
+                calculate_npv(
+                    self.royalty_holder_discount_rate.value,
+                    self.sam_economics_calculations.royalties_opex.value,
+                    self.discount_initial_year_cashflow.value
+                ),
+                self.sam_economics_calculations.royalties_opex.CurrentUnits.get_currency_unit_str()
+            ).to(self.royalty_holder_npv.CurrentUnits).magnitude
 
-            self.royalty_holder_npv.value = calculate_npv(
-                self.royalty_holder_discount_rate.value,
-                self.sam_economics_calculations.royalties_opex.value,
-                self.discount_initial_year_cashflow.value
-            )
             self.royalty_holder_annual_revenue.value = self.royalties_average_annual_cost.value
-            self.royalty_holder_total_revenue.value = np.sum(  # TODO unit conversion
-                self.sam_economics_calculations.royalties_opex.value[pre_revenue_years_slice_index:]
-            )
+
+            self.royalty_holder_total_revenue.value = quantity(
+                np.sum(
+                    self.sam_economics_calculations.royalties_opex.value[pre_revenue_years_slice_index:]
+                ),
+                self.sam_economics_calculations.royalties_opex.CurrentUnits.get_currency_unit_str()
+            ).to(self.royalty_holder_total_revenue.CurrentUnits).magnitude
 
 
         self.wacc.value = self.sam_economics_calculations.wacc.value

src/geophires_x/EconomicsSam.py

@@ -173,12 +173,13 @@ def sf(_v: float, num_sig_figs: int = 5) -> float:
     sam_economics.project_npv.value = sf(single_owner.Outputs.project_return_aftertax_npv * 1e-6)
     sam_economics.capex.value = single_owner.Outputs.adjusted_installed_cost * 1e-6
 
-    royalty_rate = model.economics.royalty_rate.quantity().to('dimensionless').magnitude
-    ppa_revenue_row = _cash_flow_profile_row(cash_flow, 'PPA revenue ($)')
-    royalties_unit = sam_economics.royalties_opex.CurrentUnits.value.replace('/yr', '')
-    sam_economics.royalties_opex.value = [
-        quantity(x * royalty_rate, 'USD').to(royalties_unit).magnitude for x in ppa_revenue_row
-    ]
+    if model.economics.royalty_rate.Provided:
+        # Assumes that royalties opex is the only possible O&M production-based expense - this logic will need to be
+        # updated if more O&M production-based expenses are added to SAM-EM
+        sam_economics.royalties_opex.value = [
+            quantity(it, 'USD / year').to(sam_economics.royalties_opex.CurrentUnits).magnitude
+            for it in _cash_flow_profile_row(cash_flow, 'O&M production-based expense ($)')
+        ]
 
     sam_economics.nominal_discount_rate.value, sam_economics.wacc.value = _calculate_nominal_discount_rate_and_wacc(
         model, single_owner
@@ -408,28 +409,11 @@ def _get_single_owner_parameters(model: Model) -> dict[str, Any]:
     geophires_ptr_tenths = Decimal(econ.PTR.value)
     ret['property_tax_rate'] = float(geophires_ptr_tenths * Decimal(100))
 
-    ppa_price_schedule_per_kWh = _ppa_pricing_model(
-        model.surfaceplant.plant_lifetime.value,
-        econ.ElecStartPrice.value,
-        econ.ElecEndPrice.value,
-        econ.ElecEscalationStart.value,
-        econ.ElecEscalationRate.value,
-    )
+    ppa_price_schedule_per_kWh = _get_ppa_price_schedule_per_kWh(model)
     ret['ppa_price_input'] = ppa_price_schedule_per_kWh
 
-    royalty_rate_schedule = _get_royalty_rate_schedule(model)
     if model.economics.royalty_rate.Provided:
-        # For each year, calculate the royalty as a $/MWh variable cost.
-        # The royalty is a percentage of revenue (MWh * $/MWh). By setting the
-        # variable O&M rate to (PPA Price * Royalty Rate), SAM's calculation
-        # (Rate * MWh) will correctly yield the total royalty payment.
-        variable_om_schedule_per_MWh = [
-            (price_kwh * 1000) * royalty_fraction  # TODO use pint unit conversion instead
-            for price_kwh, royalty_fraction in zip(ppa_price_schedule_per_kWh, royalty_rate_schedule)
-        ]
-
-        # The PySAM parameter for variable operating cost in $/MWh is 'om_production'.
-        ret['om_production'] = variable_om_schedule_per_MWh
+        ret['om_production'] = _get_royalties_variable_om_per_MWh_schedule(model)
 
     # Debt/equity ratio ('Fraction of Investment in Bonds' parameter)
     ret['debt_percent'] = _pct(econ.FIB)
@@ -452,6 +436,24 @@ def _get_single_owner_parameters(model: Model) -> dict[str, Any]:
     return ret
 
 
+def _get_royalties_variable_om_per_MWh_schedule(model: Model):
+    """TODO price unit in method name"""
+
+    royalty_rate_schedule = _get_royalty_rate_schedule(model)
+    ppa_price_schedule_per_kWh = _get_ppa_price_schedule_per_kWh(model)
+
+    # For each year, calculate the royalty as a $/MWh variable cost.
+    # The royalty is a percentage of revenue (MWh * $/MWh). By setting the
+    # variable O&M rate to (PPA Price * Royalty Rate), SAM's calculation
+    # (Rate * MWh) will correctly yield the total royalty payment.
+    variable_om_schedule_per_MWh = [
+        (price_kWh * 1000) * royalty_fraction  # TODO use pint unit conversion instead
+        for price_kWh, royalty_fraction in zip(ppa_price_schedule_per_kWh, royalty_rate_schedule)
+    ]
+
+    return variable_om_schedule_per_MWh
+
+
 def _get_fed_and_state_tax_rates(geophires_ctr_tenths: float) -> tuple[list[float]]:
     geophires_ctr_tenths = Decimal(geophires_ctr_tenths)
     max_fed_rate_tenths = Decimal(0.21)
@@ -469,9 +471,22 @@ def _pct(econ_value: Parameter) -> float:
     return econ_value.quantity().to(convertible_unit('%')).magnitude
 
 
+def _get_ppa_price_schedule_per_kWh(model: Model) -> list[float]:
+    """TODO price unit in method name"""
+
+    econ = model.economics
+    return _ppa_pricing_model(
+        model.surfaceplant.plant_lifetime.value,
+        econ.ElecStartPrice.value,
+        econ.ElecEndPrice.value,
+        econ.ElecEscalationStart.value,
+        econ.ElecEscalationRate.value,
+    )
+
+
 def _ppa_pricing_model(
     plant_lifetime: int, start_price: float, end_price: float, escalation_start_year: int, escalation_rate: float
-) -> list:
+) -> list[float]:
     # See relevant comment in geophires_x.EconomicsUtils.BuildPricingModel re:
     # https://github.com/NREL/GEOPHIRES-X/issues/340?title=Price+Escalation+Start+Year+seemingly+off+by+1.
     # We use the same utility method here for the sake of consistency despite technical incorrectness.

src/geophires_x/EconomicsUtils.py

@@ -150,8 +150,8 @@ def royalty_cost_output_parameter() -> OutputParameter:
     return OutputParameter(
             Name='Royalty Cost',
             UnitType=Units.CURRENCYFREQUENCY,
-            PreferredUnits=CurrencyFrequencyUnit.MDOLLARSPERYEAR,
-            CurrentUnits=CurrencyFrequencyUnit.MDOLLARSPERYEAR,
+            PreferredUnits=CurrencyFrequencyUnit.DOLLARSPERYEAR,
+            CurrentUnits=CurrencyFrequencyUnit.DOLLARSPERYEAR,
             ToolTipText='The annual costs paid to a royalty holder, calculated as a percentage of the '
                         'project\'s gross annual revenue. This is modeled as a variable operating expense.'
         )

src/geophires_x/Units.py

@@ -201,6 +201,9 @@ class CurrencyFrequencyUnit(str, Enum):
     KMXNPERYEAR = "KMXN/yr"
     MXNPERYEAR = "MXN/yr"
 
+    def get_currency_unit_str(self) -> str:
+        return self.value.split('/')[0]
+
 
 class EnergyCostUnit(str, Enum):
     DOLLARSPERKWH = "USD/kWh"

tests/examples/example_SAM-single-owner-PPA-4.out

@@ -4,10 +4,10 @@
 
 Simulation Metadata
 ----------------------
- GEOPHIRES Version: 3.9.56
- Simulation Date: 2025-09-16
- Simulation Time:  09:04
- Calculation Time:      1.182 sec
+ GEOPHIRES Version: 3.9.59
+ Simulation Date: 2025-09-24
+ Simulation Time:  08:47
+ Calculation Time:      1.162 sec
 
                            ***SUMMARY OF RESULTS***
 
@@ -115,8 +115,8 @@ Simulation Metadata
          Wellfield maintenance costs:                         1.13 MUSD/yr
          Power plant maintenance costs:                       3.90 MUSD/yr
          Water costs:                                         1.58 MUSD/yr
-         Average Annual Royalty Cost:                         2.50 MUSD/yr
-      Total operating and maintenance costs:                  9.10 MUSD/yr
+         Average Annual Royalty Cost:                         4.71 MUSD/yr
+      Total operating and maintenance costs:                 11.31 MUSD/yr
 
 
                            ***SURFACE EQUIPMENT SIMULATION RESULTS***
@@ -417,6 +417,6 @@ Interest earned on reserves ($)                  	0              	0          	0
 
                                 ***EXTENDED ECONOMICS***
 
-      Royalty Holder NPV:                                   29.63 MUSD
-      Royalty Holder Average Annual Revenue:                 2.50 MUSD/yr
-      Royalty Holder Total Revenue:                         49.93 MUSD
+      Royalty Holder NPV:                                   54.18 MUSD
+      Royalty Holder Average Annual Revenue:                 4.71 MUSD/yr
+      Royalty Holder Total Revenue:                         94.17 MUSD

tests/geophires_x_tests/test_units.py

@@ -0,0 +1,8 @@
+from base_test_case import BaseTestCase
+from geophires_x.Units import CurrencyFrequencyUnit
+
+
+class UnitsTestCase(BaseTestCase):
+
+    def test_get_currency_frequency_unit_currency_unit_str(self):
+        self.assertEqual('USD', CurrencyFrequencyUnit.DOLLARSPERYEAR.get_currency_unit_str())