dws-city-res-solve/main.py
2026-06-08 11:35:22 -04:00

725 lines
28 KiB
Python

"""
City Resource Optimization -> CP-SAT (Google OR-Tools)
Maximise Electrum * Brass * Steel after the gains of step 5.
The objective is a product of three variables, so this is a constraint-
programming / nonlinear problem, hence CP-SAT (cp_model) rather than the
LP/MIP solver. Read the comments at:
- PARAMETERS (initial pools + arrival schedule + bonus mode)
- "OBJECTIVE IS SET HERE" (the product being maximised -- tweak freely)
"""
from ortools.sat.python import cp_model
import printer
# ======================================================================
# PARAMETERS -- edit these
# ======================================================================
# Starting resource pools at the start of step 1: (Electrum, Brass, Steel, Capital)
INITIAL = (3, 3, 3, 3)
# Arrival schedule. Key = step (1..5), value = list of city types that arrive
# at the START of that step. Types: 'H' Hub, 'F' Foundry, 'M' Metropolis, 'N' Monument.
# Total cities across all steps must be <= 7. Arriving cities act that same step.
ARRIVALS = {
1: ["H", "F", "H", "H", "N"],
2: [],
3: [],
4: [],
5: [],
}
# Collect Bonus (b) adds +1 to whatever a Collect gives. On a foundry that is
# +1 of the resource collected (i.e. the chosen vat's value + 1). This is the
# same uniform rule as Hub (+1 Capital) and Metropolis (+1 resource pick).
NUM_STEPS = 5
MAX_RES = (
200 # upper bound on any resource pool (raise if you expect more; affects speed)
)
MAX_VAT = 15 # upper bound on a foundry vat value (1 + 2*nsteps is plenty)
# ======================================================================
# MODEL
# ======================================================================
def solve(
initial=INITIAL,
arrivals=ARRIVALS,
max_res=MAX_RES,
max_vat=MAX_VAT,
time_limit=60.0,
num_workers=8,
verbose=True,
):
# ---- build the city list -----------------------------------------
cities = [] # list of (arrival_step, arrival_type)
for s in range(1, NUM_STEPS + 1):
for typ in arrivals.get(s, []):
cities.append((s, typ))
N = len(cities)
# assert N <= 7, f"At most 7 cities allowed, got {N}"
m = cp_model.CpModel()
def AND(a, b):
"""Boolean AND of two 0/1 vars, returned as a new 0/1 var."""
c = m.NewBoolVar("")
m.AddMultiplicationEquality(c, [a, b])
return c
# ---- state variables ---------------------------------------------
# Indexed [i][t]; t in 1..NUM_STEPS+1 for state (t=NUM_STEPS+1 == "after step 5").
isH, isF, isM, isMon, present = {}, {}, {}, {}, {}
hasA, hasB, hasD = {}, {}, {}
vE, vB, vS = {}, {}, {} # foundry vats at start of step t
# action variables, t in 1..NUM_STEPS
col, ua, ub, ud = {}, {}, {}, {} # collect / upgrade a,b,d
ow = {} # overwork (global limit: 1 per step)
rH, rF, rM = {}, {}, {} # renovate -> Hub/Foundry/Metro
cvE, cvB, cvS = {}, {}, {} # foundry: which vat collected (normal collect)
owcvE, owcvB, owcvS = {}, {}, {} # foundry: which vat collected (overwork)
mE, mB, mS, mC = {}, {}, {}, {} # metropolis: resources picked (+1 each)
owmE, owmB, owmS, owmC = {}, {}, {}, {} # metropolis: resources picked (overwork)
for i in range(N):
a_step, a_type = cities[i]
for t in range(1, NUM_STEPS + 2):
isH[i, t] = m.NewBoolVar(f"isH_{i}_{t}")
isF[i, t] = m.NewBoolVar(f"isF_{i}_{t}")
isM[i, t] = m.NewBoolVar(f"isM_{i}_{t}")
isMon[i, t] = m.NewBoolVar(f"isMon_{i}_{t}")
present[i, t] = m.NewBoolVar(f"present_{i}_{t}")
hasA[i, t] = m.NewBoolVar(f"hasA_{i}_{t}")
hasB[i, t] = m.NewBoolVar(f"hasB_{i}_{t}")
hasD[i, t] = m.NewBoolVar(f"hasD_{i}_{t}")
vE[i, t] = m.NewIntVar(0, max_vat, f"vE_{i}_{t}")
vB[i, t] = m.NewIntVar(0, max_vat, f"vB_{i}_{t}")
vS[i, t] = m.NewIntVar(0, max_vat, f"vS_{i}_{t}")
# presence: present from arrival step onward, persists
m.Add(present[i, t] == (1 if t >= a_step else 0))
# exactly one type iff present
m.Add(isH[i, t] + isF[i, t] + isM[i, t] + isMon[i, t] == present[i, t])
for t in range(1, NUM_STEPS + 1):
col[i, t] = m.NewBoolVar(f"col_{i}_{t}")
ua[i, t] = m.NewBoolVar(f"ua_{i}_{t}")
ub[i, t] = m.NewBoolVar(f"ub_{i}_{t}")
ud[i, t] = m.NewBoolVar(f"ud_{i}_{t}")
ow[i, t] = m.NewBoolVar(f"ow_{i}_{t}")
rH[i, t] = m.NewBoolVar(f"rH_{i}_{t}")
rF[i, t] = m.NewBoolVar(f"rF_{i}_{t}")
rM[i, t] = m.NewBoolVar(f"rM_{i}_{t}")
cvE[i, t] = m.NewBoolVar(f"cvE_{i}_{t}")
cvB[i, t] = m.NewBoolVar(f"cvB_{i}_{t}")
cvS[i, t] = m.NewBoolVar(f"cvS_{i}_{t}")
owcvE[i, t] = m.NewBoolVar(f"owcvE_{i}_{t}")
owcvB[i, t] = m.NewBoolVar(f"owcvB_{i}_{t}")
owcvS[i, t] = m.NewBoolVar(f"owcvS_{i}_{t}")
mE[i, t] = m.NewIntVar(0, 3, f"mE_{i}_{t}")
mB[i, t] = m.NewIntVar(0, 3, f"mB_{i}_{t}")
mS[i, t] = m.NewIntVar(0, 3, f"mS_{i}_{t}")
mC[i, t] = m.NewIntVar(0, 3, f"mC_{i}_{t}")
owmE[i, t] = m.NewIntVar(0, 6, f"owmE_{i}_{t}")
owmB[i, t] = m.NewIntVar(0, 6, f"owmB_{i}_{t}")
owmS[i, t] = m.NewIntVar(0, 6, f"owmS_{i}_{t}")
owmC[i, t] = m.NewIntVar(0, 6, f"owmC_{i}_{t}")
# ---- per-step gain/cost accumulators (linear expressions) ---------
gain_E = {t: [] for t in range(1, NUM_STEPS + 1)}
gain_B = {t: [] for t in range(1, NUM_STEPS + 1)}
gain_S = {t: [] for t in range(1, NUM_STEPS + 1)}
gain_C = {t: [] for t in range(1, NUM_STEPS + 1)}
cost_C = {t: [] for t in range(1, NUM_STEPS + 1)} # capital cost (collects)
cost_S = {t: [] for t in range(1, NUM_STEPS + 1)} # steel cost (upgrades b,d)
# ---- per-city logic -----------------------------------------------
for i in range(N):
a_step, a_type = cities[i]
# initial type at arrival step
init = {"H": isH, "F": isF, "M": isM, "N": isMon}
m.Add(init[a_type][i, a_step] == 1)
# no upgrades at arrival
m.Add(hasA[i, a_step] == 0)
m.Add(hasB[i, a_step] == 0)
m.Add(hasD[i, a_step] == 0)
# vats at arrival
m.Add(vE[i, a_step] == 1)
m.Add(vB[i, a_step] == 1)
m.Add(vS[i, a_step] == 1)
# before arrival: everything zero
for t in range(1, a_step):
for v in (isH, isF, isM, isMon, hasA, hasB, hasD, vE, vB, vS):
m.Add(v[i, t] == 0)
if t <= NUM_STEPS:
for v in (
col,
ua,
ub,
ud,
ow,
rH,
rF,
rM,
cvE,
cvB,
cvS,
owcvE,
owcvB,
owcvS,
mE,
mB,
mS,
mC,
owmE,
owmB,
owmS,
owmC,
):
m.Add(v[i, t] == 0)
# action + transition logic for active steps
for t in range(a_step, NUM_STEPS + 1):
P = present[i, t]
ren = m.NewBoolVar("") # any renovation this step
m.Add(ren == rH[i, t] + rF[i, t] + rM[i, t])
no_ren = m.NewBoolVar("")
m.Add(no_ren == 1 - ren)
# exactly one action while present
m.Add(
col[i, t]
+ ua[i, t]
+ ub[i, t]
+ ud[i, t]
+ ow[i, t]
+ rH[i, t]
+ rF[i, t]
+ rM[i, t]
== P
)
# renovation must change type (renovate-to-X requires not-X now)
m.Add(isH[i, t] == 0).OnlyEnforceIf(rH[i, t])
m.Add(isF[i, t] == 0).OnlyEnforceIf(rF[i, t])
m.Add(isM[i, t] == 0).OnlyEnforceIf(rM[i, t])
# upgrade legality
m.Add(hasA[i, t] == 0).OnlyEnforceIf(ua[i, t]) # can't re-acquire
m.Add(hasB[i, t] == 0).OnlyEnforceIf(ub[i, t])
m.Add(hasD[i, t] == 0).OnlyEnforceIf(ud[i, t])
m.Add(isF[i, t] == 1).OnlyEnforceIf(ud[i, t]) # d is foundry-only
# LLM allowed renovation into metropolis, need to prevent that now
m.Add(rM[i, t] == 0)
# Monument constraints due to maximizing resources and not wanting enemy to get free upgrades
m.Add(col[i, t] == 0).OnlyEnforceIf(isMon[i, t])
m.Add(ua[i, t] == 0).OnlyEnforceIf(isMon[i, t])
m.Add(ub[i, t] == 0).OnlyEnforceIf(isMon[i, t])
m.Add(ud[i, t] == 0).OnlyEnforceIf(isMon[i, t])
m.Add(ow[i, t] == 0).OnlyEnforceIf(isMon[i, t])
# overwork cooldown: city that overworked last step can't collect or overwork
if t > a_step:
m.Add(col[i, t] == 0).OnlyEnforceIf(ow[i, t - 1])
m.Add(ow[i, t] == 0).OnlyEnforceIf(ow[i, t - 1])
# ---- type transition t -> t+1 ----
m.Add(isH[i, t + 1] == isH[i, t]).OnlyEnforceIf(no_ren)
m.Add(isF[i, t + 1] == isF[i, t]).OnlyEnforceIf(no_ren)
m.Add(isM[i, t + 1] == isM[i, t]).OnlyEnforceIf(no_ren)
for r, target in ((rH[i, t], isH), (rF[i, t], isF), (rM[i, t], isM)):
m.Add(target[i, t + 1] == 1).OnlyEnforceIf(r)
m.Add(isF[i, t + 1] == 0).OnlyEnforceIf(rH[i, t])
m.Add(isM[i, t + 1] == 0).OnlyEnforceIf(rH[i, t])
m.Add(isH[i, t + 1] == 0).OnlyEnforceIf(rF[i, t])
m.Add(isM[i, t + 1] == 0).OnlyEnforceIf(rF[i, t])
m.Add(isH[i, t + 1] == 0).OnlyEnforceIf(rM[i, t])
m.Add(isF[i, t + 1] == 0).OnlyEnforceIf(rM[i, t])
# ---- upgrade transition t -> t+1 ----
# a, b survive renovation and are monotone
m.AddMaxEquality(hasA[i, t + 1], [hasA[i, t], ua[i, t]])
m.AddMaxEquality(hasB[i, t + 1], [hasB[i, t], ub[i, t]])
# d is stripped on renovation, else monotone
m.Add(hasD[i, t + 1] == 0).OnlyEnforceIf(ren)
d_keep = m.NewBoolVar("")
m.AddMaxEquality(d_keep, [hasD[i, t], ud[i, t]])
m.Add(hasD[i, t + 1] == d_keep).OnlyEnforceIf(no_ren)
# ---- collect sub-choices ----
fcol = AND(isF[i, t], col[i, t])
mcol = AND(isM[i, t], col[i, t])
hcol = AND(isH[i, t], col[i, t])
# foundry: exactly one vat chosen iff foundry collects
m.Add(cvE[i, t] + cvB[i, t] + cvS[i, t] == fcol)
# metropolis: pick (2 + hasB) resources (+1 each); else nothing
npick = m.NewIntVar(0, 3, "")
m.Add(npick == 2 + hasB[i, t]).OnlyEnforceIf(mcol)
m.Add(npick == 0).OnlyEnforceIf(mcol.Not())
m.Add(mE[i, t] + mB[i, t] + mS[i, t] + mC[i, t] == npick)
# ================= COSTS =================
# foundry & metropolis collect each cost 1 Capital
cost_C[t].append(fcol)
cost_C[t].append(mcol)
# upgrades b,d cost 2 Steel, reduced by 1 if cost-reduction (a) already held
for upg in (ub[i, t], ud[i, t]):
c = m.NewIntVar(0, 2, "")
both = AND(upg, hasA[i, t])
m.Add(c == 2).OnlyEnforceIf(upg, hasA[i, t].Not())
m.Add(c == 1).OnlyEnforceIf(both)
m.Add(c == 0).OnlyEnforceIf(upg.Not())
cost_S[t].append(c)
# ================= GAINS =================
# Hub collect: +2 Capital (+1 more if collect-bonus b)
hub_gain = m.NewIntVar(0, 3, "")
m.Add(hub_gain == 2 + hasB[i, t]).OnlyEnforceIf(hcol)
m.Add(hub_gain == 0).OnlyEnforceIf(hcol.Not())
gain_C[t].append(hub_gain)
# Metropolis collect: +1 per pick
gain_E[t].append(mE[i, t])
gain_B[t].append(mB[i, t])
gain_S[t].append(mS[i, t])
gain_C[t].append(mC[i, t])
# Foundry collect: gain chosen vat's value as that resource
gEf = m.NewIntVar(0, max_vat, "")
gBf = m.NewIntVar(0, max_vat, "")
gSf = m.NewIntVar(0, max_vat, "")
m.AddMultiplicationEquality(gEf, [cvE[i, t], vE[i, t]])
m.AddMultiplicationEquality(gBf, [cvB[i, t], vB[i, t]])
m.AddMultiplicationEquality(gSf, [cvS[i, t], vS[i, t]])
gain_E[t].append(gEf)
gain_B[t].append(gBf)
gain_S[t].append(gSf)
# Collect Bonus (b): adds +1 to the amount a Collect gives. For a
# foundry that means +1 of the collected resource (vat value + 1),
# the same uniform "+1 to what Collect gives" rule as Hub/Metro.
gain_E[t].append(AND(cvE[i, t], hasB[i, t]))
gain_B[t].append(AND(cvB[i, t], hasB[i, t]))
gain_S[t].append(AND(cvS[i, t], hasB[i, t]))
# ---- overwork sub-choices (double-collect; no Capital cost) ----
fow = AND(isF[i, t], ow[i, t])
mow = AND(isM[i, t], ow[i, t])
how = AND(isH[i, t], ow[i, t])
# foundry overwork: exactly one vat chosen iff foundry overworks
m.Add(owcvE[i, t] + owcvB[i, t] + owcvS[i, t] == fow)
# metropolis overwork: pick 2*(2 + hasB) resources; else nothing
ow_npick = m.NewIntVar(0, 6, "")
m.Add(ow_npick == 2 * (2 + hasB[i, t])).OnlyEnforceIf(mow)
m.Add(ow_npick == 0).OnlyEnforceIf(mow.Not())
m.Add(owmE[i, t] + owmB[i, t] + owmS[i, t] + owmC[i, t] == ow_npick)
# ================= OVERWORK GAINS (2x normal collect) =================
# Hub overwork: 2*(2 + hasB) Capital
hub_ow_gain = m.NewIntVar(0, 6, "")
m.Add(hub_ow_gain == 2 * (2 + hasB[i, t])).OnlyEnforceIf(how)
m.Add(hub_ow_gain == 0).OnlyEnforceIf(how.Not())
gain_C[t].append(hub_ow_gain)
# Metropolis overwork: +1 per pick (picks are already doubled via ow_npick)
gain_E[t].append(owmE[i, t])
gain_B[t].append(owmB[i, t])
gain_S[t].append(owmS[i, t])
gain_C[t].append(owmC[i, t])
# Foundry overwork: 2 * chosen vat's value
owgEf = m.NewIntVar(0, 2 * max_vat, "")
owgBf = m.NewIntVar(0, 2 * max_vat, "")
owgSf = m.NewIntVar(0, 2 * max_vat, "")
_owE = m.NewIntVar(0, max_vat, "")
_owB = m.NewIntVar(0, max_vat, "")
_owS = m.NewIntVar(0, max_vat, "")
m.AddMultiplicationEquality(_owE, [owcvE[i, t], vE[i, t]])
m.AddMultiplicationEquality(_owB, [owcvB[i, t], vB[i, t]])
m.AddMultiplicationEquality(_owS, [owcvS[i, t], vS[i, t]])
m.Add(owgEf == 2 * _owE)
m.Add(owgBf == 2 * _owB)
m.Add(owgSf == 2 * _owS)
gain_E[t].append(owgEf)
gain_B[t].append(owgBf)
gain_S[t].append(owgSf)
# Collect Bonus (b) for foundry overwork: 2 * (+1) = +2 of that resource
ow_be = AND(owcvE[i, t], hasB[i, t])
ow_bb = AND(owcvB[i, t], hasB[i, t])
ow_bs = AND(owcvS[i, t], hasB[i, t])
gain_E[t].append(ow_be)
gain_E[t].append(ow_be) # added twice == *2
gain_B[t].append(ow_bb)
gain_B[t].append(ow_bb)
gain_S[t].append(ow_bs)
gain_S[t].append(ow_bs)
# ---- vat update producing vat[i, t+1] ----
# increment added to the two non-collected vats (1, or 2 with upgrade d)
inc = m.NewIntVar(1, 2, "")
m.Add(inc == 1 + hasD[i, t])
# vat_next = result of this step's action (only meaningful if foundry)
vEn = m.NewIntVar(0, max_vat, "")
vBn = m.NewIntVar(0, max_vat, "")
vSn = m.NewIntVar(0, max_vat, "")
# collect E: E->0, B,S += inc
m.Add(vEn == 0).OnlyEnforceIf(cvE[i, t])
m.Add(vBn == vB[i, t] + inc).OnlyEnforceIf(cvE[i, t])
m.Add(vSn == vS[i, t] + inc).OnlyEnforceIf(cvE[i, t])
# collect B
m.Add(vBn == 0).OnlyEnforceIf(cvB[i, t])
m.Add(vEn == vE[i, t] + inc).OnlyEnforceIf(cvB[i, t])
m.Add(vSn == vS[i, t] + inc).OnlyEnforceIf(cvB[i, t])
# collect S
m.Add(vSn == 0).OnlyEnforceIf(cvS[i, t])
m.Add(vEn == vE[i, t] + inc).OnlyEnforceIf(cvS[i, t])
m.Add(vBn == vB[i, t] + inc).OnlyEnforceIf(cvS[i, t])
# overwork vat transitions: same reset/increment as normal collect
m.Add(vEn == 0).OnlyEnforceIf(owcvE[i, t])
m.Add(vBn == vB[i, t] + inc).OnlyEnforceIf(owcvE[i, t])
m.Add(vSn == vS[i, t] + inc).OnlyEnforceIf(owcvE[i, t])
m.Add(vBn == 0).OnlyEnforceIf(owcvB[i, t])
m.Add(vEn == vE[i, t] + inc).OnlyEnforceIf(owcvB[i, t])
m.Add(vSn == vS[i, t] + inc).OnlyEnforceIf(owcvB[i, t])
m.Add(vSn == 0).OnlyEnforceIf(owcvS[i, t])
m.Add(vEn == vE[i, t] + inc).OnlyEnforceIf(owcvS[i, t])
m.Add(vBn == vB[i, t] + inc).OnlyEnforceIf(owcvS[i, t])
# foundry but neither collecting nor overworking: vats unchanged
# f_noncollect_noow = isF AND NOT fcol AND NOT fow = isF - fcol - fow
f_noncollect_noow = m.NewBoolVar("")
m.Add(f_noncollect_noow == isF[i, t] - fcol - fow)
m.Add(vEn == vE[i, t]).OnlyEnforceIf(f_noncollect_noow)
m.Add(vBn == vB[i, t]).OnlyEnforceIf(f_noncollect_noow)
m.Add(vSn == vS[i, t]).OnlyEnforceIf(f_noncollect_noow)
# assign vat[i, t+1]:
# renovate-to-foundry -> reset to 1
# continuing foundry -> vat_next
# otherwise (not foundry next) -> 0
cont_F = AND(isF[i, t], no_ren) # stays a foundry next step
for vnext, vn in (
(vE[i, t + 1], vEn),
(vB[i, t + 1], vBn),
(vS[i, t + 1], vSn),
):
m.Add(vnext == 1).OnlyEnforceIf(rF[i, t])
m.Add(vnext == vn).OnlyEnforceIf(cont_F)
m.Add(isF[i, t + 1] == 0).OnlyEnforceIf(
rF[i, t].Not(), cont_F.Not()
) # tautology guard
not_F_next = isF[i, t + 1].Not()
m.Add(vE[i, t + 1] == 0).OnlyEnforceIf(not_F_next)
m.Add(vB[i, t + 1] == 0).OnlyEnforceIf(not_F_next)
m.Add(vS[i, t + 1] == 0).OnlyEnforceIf(not_F_next)
# ---- global overwork limit: at most 1 city overworks per step ----
for t in range(1, NUM_STEPS + 1):
m.Add(sum(ow[i, t] for i in range(N)) <= 1)
# ---- resource pool recursion --------------------------------------
E = {1: m.NewIntVar(initial[0], initial[0], "E1")}
B = {1: m.NewIntVar(initial[1], initial[1], "B1")}
S = {1: m.NewIntVar(initial[2], initial[2], "S1")}
C = {1: m.NewIntVar(initial[3], initial[3], "C1")}
for t in range(1, NUM_STEPS + 1):
E[t + 1] = m.NewIntVar(0, max_res, f"E_{t + 1}")
B[t + 1] = m.NewIntVar(0, max_res, f"B_{t + 1}")
S[t + 1] = m.NewIntVar(0, max_res, f"S_{t + 1}")
C[t + 1] = m.NewIntVar(0, max_res, f"C_{t + 1}")
# costs are paid from the START-of-step pool (must work out before gains)
m.Add(C[t] - sum(cost_C[t]) >= 0)
m.Add(S[t] - sum(cost_S[t]) >= 0)
# next pool = start - costs + gains
m.Add(E[t + 1] == E[t] + sum(gain_E[t]))
m.Add(B[t + 1] == B[t] + sum(gain_B[t]))
m.Add(S[t + 1] == S[t] - sum(cost_S[t]) + sum(gain_S[t]))
m.Add(C[t + 1] == C[t] - sum(cost_C[t]) + sum(gain_C[t]))
finalE, finalB, finalS = E[NUM_STEPS + 1], B[NUM_STEPS + 1], S[NUM_STEPS + 1]
# ---- Phase 1: real per-resource ceilings (fast LINEAR maximisations) ----
# The triple product E*B*S has a very loose relaxation, so proving optimality
# directly is slow. We first find the true max each resource can reach on its
# own (a linear objective, solved to optimality quickly), then clamp the final
# pools to those ceilings. This tightens the product's propagation enough to
# prove optimality, and is valid because no feasible solution can exceed an
# individual resource's standalone maximum.
def _ceiling(var):
s = cp_model.CpSolver()
s.parameters.max_time_in_seconds = 20.0
s.parameters.num_search_workers = num_workers
m.Maximize(var)
st = s.Solve(m)
return (
int(s.ObjectiveValue())
if st in (cp_model.OPTIMAL, cp_model.FEASIBLE)
else max_res
)
capE = _ceiling(finalE)
capB = _ceiling(finalB)
capS = _ceiling(finalS)
m.Add(finalE <= capE)
m.Add(finalB <= capB)
m.Add(finalS <= capS)
# ======================================================================
# OBJECTIVE IS SET HERE -- maximise Electrum * Brass * Steel (post step 5)
# To change the objective, edit the three "final" pools and/or the product
# below. (finalE/finalB/finalS are the pools after step 5's gains.)
# ======================================================================
## NOTE: product can be changed here
# prodEB = m.NewIntVar(0, capE * capB, "prodEB")
# m.AddMultiplicationEquality(prodEB, [finalE, finalB])
# obj = m.NewIntVar(0, capE * capB * capS, "obj")
# m.AddMultiplicationEquality(obj, [prodEB, finalS])
# m.Maximize(obj)
# Linear objective instead (it sucks)
# m.Maximize(finalE + finalB + finalS)
# New Product
def Eprod(v):
return v * v
def Bprod(v):
return v
def Sprod(v):
return v * v
prodEE = m.NewIntVar(0, Eprod(capE), "prodEE")
m.AddMultiplicationEquality(prodEE, [finalE])
prodSS = m.NewIntVar(0, Sprod(capS), "prodSS")
m.AddMultiplicationEquality(prodSS, [finalS])
prodBB = m.NewIntVar(0, Bprod(capB), "prodBB")
m.AddMultiplicationEquality(prodBB, [finalB])
prodEB = m.NewIntVar(0, Eprod(capE) * Bprod(capB), "prodEB")
m.AddMultiplicationEquality(prodEB, [prodEE, prodBB])
obj = m.NewIntVar(0, Eprod(capE) * Bprod(capB) * Sprod(capS), "obj")
m.AddMultiplicationEquality(obj, [prodEB, prodSS])
m.Maximize(obj)
# ---- Phase 2: solve the product to optimality ----
solver = cp_model.CpSolver()
solver.parameters.max_time_in_seconds = time_limit
solver.parameters.num_search_workers = num_workers
status = solver.Solve(
m,
printer.IntermediateSolutionPrinter(
{"electrum": finalE, "brass": finalB, "steel": finalS}
),
)
if verbose:
print(f"(resource ceilings used: E<={capE} B<={capB} S<={capS})")
_report(
solver,
status,
cities,
N,
isH,
isF,
isM,
isMon,
col,
ua,
ub,
ud,
ow,
rH,
rF,
rM,
cvE,
cvB,
cvS,
owcvE,
owcvB,
owcvS,
mE,
mB,
mS,
mC,
owmE,
owmB,
owmS,
owmC,
hasA,
hasB,
hasD,
vE,
vB,
vS,
E,
B,
S,
C,
finalE,
finalB,
finalS,
)
return solver, status
def _report(
solver,
status,
cities,
N,
isH,
isF,
isM,
isMon,
col,
ua,
ub,
ud,
ow,
rH,
rF,
rM,
cvE,
cvB,
cvS,
owcvE,
owcvB,
owcvS,
mE,
mB,
mS,
mC,
owmE,
owmB,
owmS,
owmC,
hasA,
hasB,
hasD,
vE,
vB,
vS,
E,
B,
S,
C,
finalE,
finalB,
finalS,
):
print("status:", solver.StatusName(status))
if status not in (cp_model.OPTIMAL, cp_model.FEASIBLE):
return
typ_name = {}
for i in range(N):
for t in range(1, NUM_STEPS + 1):
if solver.Value(isH[i, t]):
typ_name[i, t] = "Hub"
elif solver.Value(isF[i, t]):
typ_name[i, t] = "Foundry"
elif solver.Value(isM[i, t]):
typ_name[i, t] = "Metro"
elif solver.Value(isMon[i, t]):
typ_name[i, t] = "Monument"
else:
typ_name[i, t] = "-"
def action_str(i, t):
if solver.Value(col[i, t]):
if typ_name[i, t] == "Foundry":
v = (
"E"
if solver.Value(cvE[i, t])
else "B"
if solver.Value(cvB[i, t])
else "S"
)
return f"Collect vat {v}"
if typ_name[i, t] == "Metro":
picks = []
for nm, var in (("E", mE), ("B", mB), ("S", mS), ("C", mC)):
picks += [nm] * solver.Value(var[i, t])
return "Collect {" + ",".join(picks) + "}"
return "Collect (+Capital)"
if solver.Value(ow[i, t]):
if typ_name[i, t] == "Foundry":
v = (
"E"
if solver.Value(owcvE[i, t])
else "B"
if solver.Value(owcvB[i, t])
else "S"
)
return f"Overwork vat {v} (2x)"
if typ_name[i, t] == "Metro":
picks = []
for nm, var in (("E", owmE), ("B", owmB), ("S", owmS), ("C", owmC)):
picks += [nm] * solver.Value(var[i, t])
return "Overwork {" + ",".join(picks) + "} (2x)"
return "Overwork (+Capital 2x)"
if solver.Value(ua[i, t]):
return "Upgrade a (cost-reduction)"
if solver.Value(ub[i, t]):
return "Upgrade b (collect-bonus)"
if solver.Value(ud[i, t]):
return "Upgrade d (vat-increment)"
if solver.Value(rH[i, t]):
return "Renovate -> Hub"
if solver.Value(rF[i, t]):
return "Renovate -> Foundry"
if solver.Value(rM[i, t]):
return "Renovate -> Metro"
return "(inactive)"
print("\nPer-step pools (start of step):")
print(" step: E B S C")
for t in range(1, NUM_STEPS + 2):
label = f"after5" if t == NUM_STEPS + 1 else f"start {t}"
print(
f" {label:>8}: {solver.Value(E[t]):3} {solver.Value(B[t]):3} "
f"{solver.Value(S[t]):3} {solver.Value(C[t]):3}"
)
print("\nActions:")
for i in range(N):
a_step, a_type = cities[i]
print(f" City {i} (arrives step {a_step} as {a_type}):")
for t in range(a_step, NUM_STEPS + 1):
ups = "".join(
n
for n, h in (("a", hasA), ("b", hasB), ("d", hasD))
if solver.Value(h[i, t])
)
extra = (
f" vats(E{solver.Value(vE[i, t])},B{solver.Value(vB[i, t])},S{solver.Value(vS[i, t])})"
if typ_name[i, t] == "Foundry"
else ""
)
print(
f" step {t}: [{typ_name[i, t]:7}] {action_str(i, t):28}"
f" upg[{ups}]{extra}"
)
fe, fb, fs = solver.Value(finalE), solver.Value(finalB), solver.Value(finalS)
print(
f"\nFINAL E={fe} B={fb} S={fs} product = {fe * fb * fs} sum = {fe + fb + fs}"
)
if __name__ == "__main__":
solve()