Persistent Solvers should support IndexVars

[leroum]

Summary

It would be very helpful if the persistent solver would support updating IndexedVars in the model.

Rationale

When using pyomo through an application like oemof, it would be awsome to use the persitant solver. But when the original model containts IndexedVars it is not possible to update those vars or single entrys of them. This makes the persitant solver functionality not useable with those applications.

[michaelbynum]

Can you be more specific? Which persistent solver are you trying to use? It should be possible to update indexed variables with:

for k, v in m.x.items():
    opt.update_var(v)

[mrmundt]

@leroum - Can you please answer @michaelbynum 's question? We see that this is a popular request, so we'd like to know more.

[leroum]

Sorry for the late answer. I wanted to investigate a bit more after the information that it should work with the loop.

I am using gurobi_persistent and for indexed vars the update_var() function does not work even when using the loop. Where as the set_instance() functionality works fine.

[michaelbynum]

@leroum - do you get an error? If so, can you post the stack trace?

[leroum]

No i dont get an error message. The model is simply not updated in the persititant solver.

[quantresearch1]

Hello,

I think we would need a minimal working example to confirm what issue @leroum is encountering. I just tested it in on pyomo 6.9.4 with python 3.11 and the update_var seems to work as expected in my test case shared below. I have also tested the case with multi dimensional indexing of the form m2.z = pyo.Var(m2.I, m2.J, bounds=(0, 10)) which worked fine.

import pyomo.environ as pyo
from pyomo.solvers.plugins.solvers.gurobi_persistent import GurobiPersistent


print("=== IndexedVar Update Test Case ===")
opt = GurobiPersistent()

# Create a simple model with IndexedVar
print("\n1. Creating test model...")
m = pyo.ConcreteModel()

# Create indexed sets
m.I = pyo.Set(initialize=[1, 2, 3])
m.J = pyo.Set(initialize=['a', 'b'])

# Create IndexedVar with initial bounds
m.x = pyo.Var(m.I, bounds=(0, 10), initialize=5)
m.y = pyo.Var(m.J, bounds=(0, 5), initialize=2)

# Create some constraints to make the problem meaningful
m.constraint1 = pyo.Constraint(expr=sum(m.x[i] for i in m.I) <= 20)
m.constraint2 = pyo.Constraint(expr=sum(m.y[j] for j in m.J) >= 3)
m.constraint3 = pyo.Constraint(expr=m.x[1] + m.y['a'] >= 4)

# Objective: maximize sum of all variables
m.obj = pyo.Objective(expr=sum(m.x[i] for i in m.I) + sum(m.y[j] for j in m.J), 
                        sense=pyo.maximize)

print("Model created with:")
print(f"  - x[{list(m.I)}] with bounds (0, 10)")
print(f"  - y[{list(m.J)}] with bounds (0, 5)")
print("  - Objective: maximize sum of all variables")

# Set up persistent solver
print("\n2. Setting up persistent solver...")
opt.set_instance(m)

# Solve initial problem
print("\n3. Solving initial problem...")
results = opt.solve()
print(f"Initial solve status: {results.solver.termination_condition}")

if results.solver.termination_condition == pyo.TerminationCondition.optimal:
    print("Initial solution:")
    for i in m.I:
        print(f"  x[{i}] = {pyo.value(m.x[i]):.2f}")
    for j in m.J:
        print(f"  y[{j}] = {pyo.value(m.y[j]):.2f}")
    print(f"  Objective value = {pyo.value(m.obj):.2f}")
    initial_obj = pyo.value(m.obj)

# Now test the IndexedVar update functionality
print("\n4. Testing IndexedVar updates...")
print("Updating bounds: x variables to (0, 3), y variables to (0, 2)")

# Update bounds on the Pyomo model
for i in m.I:
    m.x[i].setlb(0)
    m.x[i].setub(3)  # Much tighter upper bound

for j in m.J:
    m.y[j].setlb(0)
    m.y[j].setub(2)  # Much tighter upper bound

print("\n5. Applying Michael's approach for x variables...")
for k, v in m.x.items():
    print(f"  Updating x[{k}] (bounds: {v.lb}, {v.ub})")
    opt.update_var(v)
for k, v in m.y.items():
    print(f"  Updating y[{k}] (bounds: {v.lb}, {v.ub})")
    opt.update_var(v)
print("✓ Successfully updated all y variables")

# Solve again with updated bounds
print("\n7. Solving with updated bounds...")
results = opt.solve()
print(f"Updated solve status: {results.solver.termination_condition}")

if results.solver.termination_condition == pyo.TerminationCondition.optimal:
    print("Solution with updated bounds:")
    for i in m.I:
        print(f"  x[{i}] = {pyo.value(m.x[i]):.2f}")
    for j in m.J:
        print(f"  y[{j}] = {pyo.value(m.y[j]):.2f}")
    updated_obj = pyo.value(m.obj)
    print(f"  Objective value = {updated_obj:.2f}")
    
    # Verify the bounds were actually applied
    print(f"\n8. Verification:")
    print(f"  Initial objective: {initial_obj:.2f}")
    print(f"  Updated objective: {updated_obj:.2f}")
    
    if updated_obj < initial_obj:
        print("✓ SUCCESS: Objective decreased as expected due to tighter bounds")