import numpy as np
from scipy.optimize import minimize
initial_guess = [50, 4, 0.2, 2.22]
def modified_cost(n, a, b, c, d):
return a + (b * n + c) ** d
def original_cost(n):
return modified_cost(n, *initial_guess)
# minimize deviation at n < 7 and keep ~3 policies ahead at n = 12 to ~25
def objective(params):
a, b, c, d = params
ns = np.arange(0, 7)
original = original_cost(ns)
modified = modified_cost(ns, a, b, c, d)
low_n_error = np.sum((modified - original) ** 2) / len(ns) # per point error
# Maintain ~3 ahead in cost: modified(n) ≈ original(n + 3)
ns_mid = np.arange(12, 26)
target = original_cost(ns_mid + 3)
modified_mid = modified_cost(ns_mid, a, b, c, d)
offset_error = np.sum((modified_mid - target) ** 2) / len(ns_mid)
# Scale goes up by at least 10x, so weight early more to preserve behavior for small n
return 20 * low_n_error + offset_error
# note c cannot go negative or we will get a complex number
result = minimize(objective, initial_guess, bounds=[(40, 60), (3, 5), (0, 1), (2, 2.5)])
print(result.x)