diff --git a/biosteam/_tea.py b/biosteam/_tea.py index 34cb4c6f..4c43d78a 100644 --- a/biosteam/_tea.py +++ b/biosteam/_tea.py @@ -182,6 +182,7 @@ def taxable_and_nontaxable_cashflows( start, years, lang_factor, accumulate_interest_during_construction, + f_inflation, ): # Cash flow data and parameters # C_FC: Fixed capital @@ -202,6 +203,13 @@ def taxable_and_nontaxable_cashflows( ) for i in unit_capital_costs: add_all_replacement_costs_to_cashflow_array(i, C_FC, years, start, lang_factor) + + # Multiply for inflation factors, if inflation None the factors are 1. + C *= f_inflation + S *= f_inflation + C_FC *= f_inflation + C_WC *= f_inflation + if finance_interest: interest = finance_interest years = finance_years @@ -220,6 +228,11 @@ def taxable_and_nontaxable_cashflows( nontaxable_cashflow = D - C_FC - C_WC return taxable_cashflow, nontaxable_cashflow +def build_nominal_factors(factors, rate): + operating_index = np.arange(factors.size) + factors[:] = (1.0 + rate)**operating_index + return factors + def NPV_with_sales( sales, taxable_cashflow, @@ -274,7 +287,10 @@ class TEA: system : Should contain feed and product streams. IRR : - Internal rate of return (fraction). + Real internal rate of return (fraction). If `inflation_rate` is given, + cashflows are escalated to nominal values and this real IRR is internally + converted to a nominal discount rate for NPV calculations. This is done by + using Fisher equation: nom_IRR = (1 + IRR)*(1 + inflation_rate) - 1. duration : Start and end year of venture (e.g. (2018, 2038)). depreciation : @@ -306,11 +322,21 @@ class TEA: WC_over_FCI : Working capital as a fraction of fixed capital investment. finance_interest : - Yearly interest of capital cost financing as a fraction. + Yearly interest of capital cost financing as a fraction. If `inflation_rate` + is provided, nominal yearly interest of capiltal cost financing must be + given. finance_years : Number of years the loan is paid for. finance_fraction : Fraction of capital cost that needs to be financed. + inflation_rate : + Annual constant inflation rate as a fraction. If provided, + operating costs, sales, capital expenses, replacement costs + and working capital flows are escalated to nominal dollars + using this rate. NPV calculations then use a nominal + discount rate computed from `IRR` and `inflation_rate`. + If `None`, no inflation is applied and cashflows are treated + as base-year dollars. Warning ------- @@ -333,7 +359,8 @@ class TEA: '_startup_schedule', '_operating_days', '_duration', '_depreciation_key', '_depreciation', '_years', '_duration', '_start', 'IRR', '_IRR', '_sales', - '_duration_array_cache', 'accumulate_interest_during_construction') + '_duration_array_cache', 'accumulate_interest_during_construction', + '_inflation_rate', '_inflation_factors') #: Available depreciation schedules. Defaults include modified #: accelerated cost recovery system from U.S. IRS publication 946 (MACRS), @@ -446,16 +473,24 @@ def __init__(self, system: System, IRR: float, duration: tuple[int, int], startup_months: float, startup_FOCfrac: float, startup_VOCfrac: float, startup_salesfrac: float, WC_over_FCI: float, finance_interest: float, finance_years: int, finance_fraction: float, - accumulate_interest_during_construction: bool=False): + accumulate_interest_during_construction: bool=False, + inflation_rate: float|None = None): #: System being evaluated. self.system: System = system + # Inflation + self._inflation_rate = inflation_rate + self._inflation_factors = None + + # Time periods self.duration = duration self.depreciation = depreciation self.construction_schedule = construction_schedule self.startup_months = startup_months self.operating_days = operating_days + self.update_inflation_factors() + #: Internal rate of return (fraction). self.IRR: float = IRR @@ -493,7 +528,7 @@ def __init__(self, system: System, IRR: float, duration: tuple[int, int], #: Whether to immediately pay interest before operation or to accumulate interest during construction self.accumulate_interest_during_construction = accumulate_interest_during_construction - + #: For convenience, set a TEA attribute for the system system._TEA = self @@ -567,6 +602,7 @@ def duration(self, duration): start, end = [int(i) for i in duration] self._duration = (start, end) self._years = end - start + self.update_inflation_factors() @property def depreciation(self) -> str|NDArray[float]: @@ -644,6 +680,7 @@ def construction_schedule(self) -> Sequence[float]: def construction_schedule(self, schedule): self._construction_schedule = np.array(schedule, dtype=float) self._start = len(schedule) + self.update_inflation_factors() @property def startup_months(self) -> float: @@ -728,6 +765,51 @@ def PBP(self) -> float: net_earnings = self.net_earnings return FCI/net_earnings + @property + def inflation_rate(self): + """Annual constant inflation rate used to escalate cashflows to nominal dollars""" + return self._inflation_rate + + @inflation_rate.setter + def inflation_rate(self, value): + self._inflation_rate = value + self.update_inflation_factors() + + @property + def inflation_factors(self): + """Multiplicative nominal escalation factors aligned with the cashflow array""" + return self._inflation_factors + + def update_inflation_factors(self): + """Update inflation factors""" + if hasattr(self,"_start") and hasattr(self, "_years"): + self._inflation_factors = self._build_inflation_factors() + + def _build_inflation_factors(self): + """ + Build nominal escalation factors for all construction and operating years. + + Returns an array of length `self._start + self._years`. If `inflation_rate` + is None, all factors are 1. Otherwise, factors are compounded as: + factor[t] = (1 + inflation_rate)**t + """ + length = self._start + self._years + inflation = self.inflation_rate + + if inflation is None: + return np.ones(length, dtype=float) + + if isinstance(inflation, (int, float)): + return build_nominal_factors(np.ones(length, dtype=float), inflation) + + raise TypeError("inflation_rate must be None or float annual rate.") + + def _get_discount_rate(self): + """Return the discount rate consistent with the cashflow basis.""" + if self.inflation_rate is None: + return self.IRR + return (1.0 + self.IRR) * (1.0 + self.inflation_rate) - 1.0 + def _get_duration_array(self): key = start, years = (self._start, self._years) if key in _duration_array_cache: @@ -756,7 +838,17 @@ def _fill_depreciation_array(self, D, start, years, TDC): D[start:start + N_depreciation_years] = TDC * depreciation_array def get_cashflow_table(self): - """Return DataFrame of the cash flow analysis.""" + """ + Return DataFrame of the cash flow analysis. + + If `inflation_rate` is provided annual, annual costs, sales, capital expenses, + replacement costs and working capital are reported in nominal dollars. + Discount factors are computed using the nominal discount rate derived from + the real `IRR` and `inflation_rate`. + + If `inflation_rate` is None, values are reported in real or base-year dollars and + discounted directly with `IRR`. + """ # Cash flow data and parameters # index: Year since construction until end of venture # C_D: Depreciable capital @@ -778,16 +870,16 @@ def get_cashflow_table(self): # DF: Discount factor # NPV: Net present value # CNPV: Cumulative NPV - TDC = self.TDC - FCI = self._FCI(TDC) + TDC0 = self.TDC + FCI = self._FCI(TDC0) start = self._start years = self._years + f_inflation = self.inflation_factors FOC = self._FOC(FCI) VOC = self.VOC sales = self.sales length = start + years C_D, C_FC, C_WC, D, L, LI, LP, LPl, C, S, T, I, TE, FL, NE, CF, DF, NPV, CNPV = data = np.zeros((19, length)) - self._fill_depreciation_array(D, start, years, TDC) w0 = self._startup_time % 1 w1 = 1. - w0 end_start = start + int(self._startup_time) @@ -800,15 +892,21 @@ def get_cashflow_table(self): start1 = end_start + 1 C[start1:] = VOC + FOC S[start1:] = sales + C *= f_inflation + S *= f_inflation WC = self.WC_over_FCI * FCI - C_D[:start] = TDC*self._construction_schedule + C_D[:start] = TDC0*self._construction_schedule + C_D *= f_inflation + self._fill_depreciation_array(D, start, years, C_D[:start].sum()) C_FC[:start] = FCI*self._construction_schedule C_WC[start-1] = WC C_WC[-1] = -WC + C_WC *= f_inflation system = self.system lang_factor = system.lang_factor unit_capital_costs = system.unit_capital_costs.values() if isinstance(system, bst.AgileSystem) else system.cost_units for i in unit_capital_costs: add_all_replacement_costs_to_cashflow_array(i, C_FC, years, start, lang_factor) + C_FC *= f_inflation if self.finance_interest: interest = self.finance_interest years = self.finance_years @@ -848,7 +946,7 @@ def get_cashflow_table(self): ) NE[:] = taxable_cashflow + I - T CF[:] = NE + nontaxable_cashflow - DF[:] = 1/(1.+self.IRR)**self._get_duration_array() + DF[:] = 1/(1.+self._get_discount_rate())**self._get_duration_array() NPV[:] = CF * DF CNPV[:] = NPV.cumsum() DF *= 1e6 @@ -858,7 +956,13 @@ def get_cashflow_table(self): columns=cashflow_columns) @property def NPV(self) -> float: - """Net present value.""" + """ + Net present value. + + Uses cashflows consistent with the inflation setting. With inflation, cashflows are + nominal and discounted with the internally computed nominal discount rate. without + inflation, cashflows are treated as base-year values and discounted with `IRR`. + """ taxable_cashflow, nontaxable_cashflow, depreciation = self._taxable_nontaxable_depreciation_cashflows() tax = np.zeros_like(taxable_cashflow) incentives = tax.copy() @@ -868,7 +972,7 @@ def NPV(self) -> float: nontaxable_cashflow, tax, depreciation ) cashflow = nontaxable_cashflow + taxable_cashflow + incentives - tax - return NPV_at_IRR(self.IRR, cashflow, self._get_duration_array()) + return NPV_at_IRR(self._get_discount_rate(), cashflow, self._get_duration_array()) def _AOC(self, FCI): """Return AOC at given FCI""" @@ -886,21 +990,23 @@ def _taxable_nontaxable_depreciation_cashflows(self): # S: Sales # NE: Net earnings # CF: Cash flow - TDC = self.TDC - FCI = self._FCI(TDC) + TDC0 = self.TDC + FCI = self._FCI(TDC0) start = self._start years = self._years + inflation_factors = self.inflation_factors + TDC_nom = (TDC0 * self.construction_schedule * inflation_factors[:start]).sum() FOC = self._FOC(FCI) VOC = self.VOC D, C_FC, C_WC, Loan, LP, C, S = np.zeros((7, start + years)) - self._fill_depreciation_array(D, start, years, TDC) + self._fill_depreciation_array(D, start, years, TDC_nom) WC = self.WC_over_FCI * FCI system = self.system return ( *taxable_and_nontaxable_cashflows( system.unit_capital_costs if isinstance(system, bst.AgileSystem) else system.cost_units, D, C, S, C_FC, C_WC, Loan, LP, - FCI, WC, TDC, VOC, FOC, self.sales, + FCI, WC, TDC0, VOC, FOC, self.sales, self._startup_time, self.startup_VOCfrac, self.startup_FOCfrac, @@ -912,6 +1018,7 @@ def _taxable_nontaxable_depreciation_cashflows(self): start, years, self.lang_factor, self.accumulate_interest_during_construction, + inflation_factors, ), D ) @@ -990,8 +1097,23 @@ def total_production_cost(self, products: Collection[bst.Stream], with_annual_de return self.AOC - coproduct_sales def solve_IRR(self, financing=True, bounds=None): - """Return the IRR at the break even point (NPV = 0) through cash flow analysis.""" + """ + Return the real internal rate of return at the break-even point. + + If `inflation_rate` is provided, this method solves IRR for inflated + cashflows and the resulting nominal IRR is converted to real IRR applying + Fisher equation: real_IRR = (1 + nominal_IRR)/(1 + inflation_rate) - 1. + + If bounds are provided and `inflation_rate` is not None, bounds must be given + as nominal IRR values because the solver operates on nominal cashflows. + + """ IRR = self._IRR + + # Use nominal IRR as initial guess to solve nominal cashflows + if self.inflation_rate is not None: + IRR = (1 + IRR)*(1 + self.inflation_rate) - 1 + if not IRR or np.isnan(IRR) or IRR < 0.: IRR = 0.01 if financing: args = (self.cashflow_array, self._get_duration_array()) @@ -1022,6 +1144,11 @@ def solve_IRR(self, financing=True, bounds=None): ) finally: self.finance_fraction, self.finance_interest = financing_values + + # convert nominal IRR to real IRR + if self.inflation_rate is not None: + IRR = (1 + IRR)/(1 + self.inflation_rate) - 1 + self._IRR = IRR return IRR @@ -1070,10 +1197,14 @@ def FOC_table(self): def solve_sales(self): """ Return the required additional sales [USD] to reach the breakeven - point (NPV = 0) through cash flow analysis. + point (NPV = 0) through cash flow analysis. + + The returned value is expressed in base-year dollars. If `inflation_rate` + is provided, this additional sales value is escalated through time using + `inflation_factors` before calculating NPV. """ - discount_factors = (1 + self.IRR)**self._get_duration_array() + discount_factors = (1 + self._get_discount_rate())**self._get_duration_array() sales_coefficients = np.ones_like(discount_factors, dtype=float) start = self._start sales_coefficients[:start] = 0 @@ -1081,6 +1212,7 @@ def solve_sales(self): end_start = start + int(self._startup_time) sales_coefficients[end_start] = w0 * self.startup_salesfrac + (1. - w0) sales_coefficients[start:end_start] = self.startup_salesfrac + sales_coefficients *= self.inflation_factors taxable_cashflow, nontaxable_cashflow, depreciation = self._taxable_nontaxable_depreciation_cashflows() if np.isnan(taxable_cashflow).any(): warn('nan encountered in cashflow array; resimulating system', category=RuntimeWarning) @@ -1112,7 +1244,7 @@ def __repr__(self): def _info(self): return (f'{type(self).__name__}: {self.system}\n' - f'NPV: {self.NPV:,.0f} USD at {self.IRR:.1%} IRR') + f'NPV: {self.NPV:,.0f} USD at {self.IRR:.1%} real IRR') def show(self): """Prints information on unit.""" diff --git a/tests/test_tea.py b/tests/test_tea.py index c3448b0b..f82c6e9f 100644 --- a/tests/test_tea.py +++ b/tests/test_tea.py @@ -337,6 +337,118 @@ def _FCI(self, TDC): # fixed capital investment table['Loan principal [MM$]'].iloc[0]) assert_allclose(total_interest_payment1, total_interest_payment2, atol=1e-4) +def test_tea_with_inflation(): + cost = bst.decorators.cost + bst.CE = CE = bst.design_tools.CEPCI_by_year[2013] + + @cost('Fake scaler', 'Lumped cost', CE=CE, cost=1e6, S=1, n=1, BM=1) + class LumpedCost(bst.Unit): + '''Does nothing but adding given costs.''' + _units = {'Fake scaler': ''} + + def _design(self): + self.design_results['Fake scaler'] = 1 + + class TEA(bst.TEA): + def __init__( + self, + system, + FOC_over_installed=0.5, + DPI_over_installed=2, + TDC_over_DPI=1.2*1.4, + FCI_over_TDC=1, + **kwargs, + ): + self.FOC_over_installed = FOC_over_installed + self.DPI_over_installed = DPI_over_installed + self.TDC_over_DPI = TDC_over_DPI + self.FCI_over_TDC = FCI_over_TDC + bst.TEA.__init__(self, system, **kwargs) + + def _FOC(self, installed_equipment_cost): + return installed_equipment_cost * self.FOC_over_installed + + def _DPI(self, installed_equipment_cost): + return installed_equipment_cost * self.DPI_over_installed + + def _TDC(self, DPI): + return DPI * self.TDC_over_DPI + + def _FCI(self, TDC): + return TDC * self.FCI_over_TDC + + bst.settings.set_thermo([bst.Chemical('Water')]) + reactant = bst.Stream('reactant', Water=1, units='kg/hr') + product = bst.Stream('product', Water=1, price=2.5e6/365/24, units='kg/hr') + + U101 = LumpedCost('U101', ins=reactant, outs=product) + sys = bst.System('sys_inflation', path=(U101,)) + sys.simulate() + + IRR = 0.10 + inflation_rate = 0.03 + nominal_IRR = (1 + IRR) * (1 + inflation_rate) - 1 + + tea = TEA( + system=sys, + IRR=IRR, + inflation_rate=inflation_rate, + duration=(2013, 2013+15), + income_tax=0.31, + construction_schedule=(1,), + depreciation='MACRS7', + operating_days=365, + startup_months=0, + startup_FOCfrac=1, + startup_VOCfrac=1, + startup_salesfrac=1, + lang_factor=None, + WC_over_FCI=0.05, + finance_interest=0.08, + finance_years=10, + finance_fraction=0.6, + accumulate_interest_during_construction=False, + ) + + table = tea.get_cashflow_table() + factors = (1 + inflation_rate) ** np.arange(tea._start + tea._years) + + assert_allclose(tea.inflation_factors, factors) + assert_allclose(tea._get_discount_rate(), nominal_IRR) + + # Check nominal escalation of representative cashflows. + assert_allclose(table['Sales [MM$]'].values[0], 0.0) + assert_allclose(table['Sales [MM$]'].values[1:], 2.5 * factors[1:]) + + assert_allclose(table['Annual operating cost (excluding depreciation) [MM$]'].values[0], 0.0) + assert_allclose( + table['Annual operating cost (excluding depreciation) [MM$]'].values[1:], + 1.68 * factors[1:], + ) + + # Check construction-year capital and final working capital recovery. + assert_allclose(table['Fixed capital investment [MM$]'].iloc[0], 3.36) + assert_allclose(table['Working capital [MM$]'].iloc[0], 0.168) + assert_allclose(table['Working capital [MM$]'].iloc[-1], -0.168 * factors[-1]) + + # Check nominal discounting. + assert_allclose( + table['Discount factor'], + 1 / (1 + nominal_IRR) ** tea._get_duration_array(), + ) + + # Table and NPV property should be consistent. + assert_allclose( + tea.NPV, + table['Cumulative NPV [MM$]'].iloc[-1] * 1e6, + atol=1e-4, + ) + + # Check solve_price with inflation; indirectly checks solve_sales. + price = tea.solve_price(product) + product.price = price + assert_allclose(tea.NPV, 0, atol=100) + def test_add_replacement_cost(): cashflow_array = np.zeros(12) @@ -358,5 +470,6 @@ def test_add_replacement_cost(): test_depreciation_schedule() test_cashflow_consistency() test_tea() + test_tea_with_inflation() test_tea_startup_months() - test_add_replacement_cost() + test_add_replacement_cost() \ No newline at end of file