Petroleum starter - Decline forecasting

MODEL "Petroleum starter - Decline forecasting"
DESCRIPTION "Synthetic Arps screening with exponential, hyperbolic and harmonic limits. Rate in STB/day, nominal decline per year, 365-day years; constant operating conditions assumed."
VERSION "1.0.0"
TAGS "petroleum", "starter", "decline", "scenarios"

strict

# Inputs: the interval 0 <= b <= 1 is a deliberate teaching scope.
input InitialRate_stbd = 1000
  VALIDATE BETWEEN 0.001 AND 1000000 MESSAGE "Initial rate must be positive (0.001 to 1,000,000 STB/day)."
input InitialDecline_per_year = 0.60
  VALIDATE BETWEEN 0 AND 5 MESSAGE "Nominal decline must be between 0 and 5 per year. Zero means constant rate."
input HyperbolicExponent = 0.50
  VALIDATE BETWEEN 0 AND 1 MESSAGE "Use b=0 for exponential, b=1 for harmonic, or 0<b<1 for hyperbolic decline."
input EconomicLimit_stbd = 100
  VALIDATE BETWEEN 0.001 AND 1000000 MESSAGE "Economic-limit rate must be positive. At or above the initial rate, the limit is already reached."
input Horizon_years = 10
  VALIDATE BETWEEN 0 AND 50 MESSAGE "Forecast horizon must be between 0 and 50 years."

# Stable elementary differences near the exponential and harmonic limits.
log_one_plus = (value) => ABS(value) < 0.00001 THEN value - value^2 / 2 + value^3 / 3 - value^4 / 4 + value^5 / 5 ELSE LN(1 + value)
exponential_minus_one = (value) => ABS(value) < 0.00001 THEN value + value^2 / 2 + value^3 / 6 + value^4 / 24 + value^5 / 120 ELSE EXP(value) - 1
rate_at = (year) => HyperbolicExponent = 0 THEN InitialRate_stbd * EXP(-InitialDecline_per_year * year) ELSE InitialRate_stbd * EXP(-log_one_plus(HyperbolicExponent * InitialDecline_per_year * year) / HyperbolicExponent)
cumulative_at = (year) => InitialDecline_per_year = 0 THEN 365 * InitialRate_stbd * year
  ELSE HyperbolicExponent = 0 THEN -365 * InitialRate_stbd / InitialDecline_per_year * exponential_minus_one(-InitialDecline_per_year * year)
  ELSE HyperbolicExponent = 1 THEN 365 * InitialRate_stbd / InitialDecline_per_year * log_one_plus(InitialDecline_per_year * year)
  ELSE -365 * InitialRate_stbd / (InitialDecline_per_year * (1 - HyperbolicExponent)) * exponential_minus_one((HyperbolicExponent - 1) / HyperbolicExponent * log_one_plus(HyperbolicExponent * InitialDecline_per_year * year))

# A zero-decline well above the limit has no finite economic-limit time.
output LimitStatus = EconomicLimit_stbd >= InitialRate_stbd THEN "Limit already reached at year zero."
  ELSE InitialDecline_per_year = 0 THEN "Constant rate: no finite time to the lower rate limit."
  ELSE "Finite time to the selected rate limit; this is a rate cutoff, not a discounted economic assessment."
output LimitTime_years = EconomicLimit_stbd >= InitialRate_stbd THEN 0
  ELSE InitialDecline_per_year = 0 THEN #N/A
  ELSE HyperbolicExponent = 0 THEN LN(InitialRate_stbd / EconomicLimit_stbd) / InitialDecline_per_year
  ELSE exponential_minus_one(HyperbolicExponent * LN(InitialRate_stbd / EconomicLimit_stbd)) / (HyperbolicExponent * InitialDecline_per_year)
output RateAtHorizon_stbd IS number = rate_at(Horizon_years)
output CumulativeAtHorizon_stb IS number = cumulative_at(Horizon_years)
output Assumptions = "Continuous production; 365-day years; Arps under constant conditions and boundary-dominated behavior. No transient-flow fit, terminal decline transition, taxes, prices or reserves classification. Horizon totals are uncurtailed by the rate cutoff."
output ScenarioGuide = "Save Baseline; set b to 0 and save Exponential, then b to 1 and save Harmonic. Set decline to 0 to check constant rate; N/A limit time then means no finite lower-rate crossing."

Years = MAP(SEQUENCE(61, 1, 0, 1), (index) => Horizon_years * index / 60)
Rates_stbd = MAP(Years, (year) => rate_at(year))
Forecast = MAP(Years, (year) => {year: year, rate: rate_at(year)})
A1 = "Rate at horizon (STB/day)"
B1 = RateAtHorizon_stbd
A2 = "Uncurtailed cumulative production (STB)"
B2 = CumulativeAtHorizon_stb
A3 = LimitStatus
A4 = Assumptions

Charts_1!type = "charts"
Charts_1!config = <toml>
  version = 1
  kind = "charts"
  title = "Decline forecasting"

  [layout]
  title = "Uncurtailed Arps rate forecast"
  x_axis = "Time (365-day years)"
  y_axis = "Oil production rate (STB/day)"

  [[charts]]
  type = "line"
  bind = "Forecast"
  x = "year"
  y = "rate"
  label = "Oil rate"
</toml>

END MODEL