or-fusion / pao Goto Github PK

As a user, I would like to add an indexed variable to the fixed list in pao.duality.plugins.create_linear_dual

As it stands, if the indexed variable v is to be fixed, then I need to add v[1], v[2], ... to the fixed list rather than just v.

Dear all,

I am wondering if there is a plan to integrate a solver for this kind of problems:

• Linear Equation Structure
• Integer Upper-Level Variables
• Integer Lower-Level Variables
• and k-Bilevel Multilevel Representation

Thank you for your contribution. It is extremely useful :)
Best regards,

Ramon

PAO only supports a limited subset. We should walk the model and confirm that we don't see others.

Develop interfaces for

• MibS - https://github.com/coin-or/MibS
• Fischetti et al - https://msinnl.github.io/pages/bilevel.html
• Pineda et al. - https://github.com/salvapineda/bilevel

This feature was not re-implemented in PAO after reworking the logic in pyomo.bilevel. This requires us to revisit the solver logic to consider the role of bilinear terms.

I tried to run the simple example from the documentation with PAO version 1.0.2 and Pyomo 6.6.2 using MibS 1.2 as the solver using Python 3.9.10 on Ubuntu in WSL.

import pyomo.environ as pe
from pao.pyomo import *

M = pe.ConcreteModel()

M.x = pe.Var(bounds=(0,None), domain=pe.NonNegativeIntegers)

M.y = pe.Var(bounds=(0,None))

M.o = pe.Objective(expr=M.x - 4*M.y)

M.L = SubModel(fixed=M.x)

M.L.o = pe.Objective(expr=M.y)

M.L.c1 = pe.Constraint(expr=   -M.x -   M.y <= -3)
M.L.c2 = pe.Constraint(expr= -2*M.x +   M.y <=  0)
M.L.c3 = pe.Constraint(expr=  2*M.x +   M.y <= 12)
M.L.c4 = pe.Constraint(expr=  3*M.x - 2*M.y <=  4)

with Solver('pao.pyomo.MIBS') as solver:
    results = solver.solve(M)
    print(results)
    
print(M.x.value)
print(M.y.value)

There was no error, but also just a solution of (0,0). I debugged and everything initially seemed correct. I could actually solve the instance that was created with MibS on the command line from the MPS and AUX files that were output. But when I traced the execution within PAO, it finally ended up on this line and I couldn't really tell what happened.

FYI, the output format of the solution has changed a little in MibS 1.2 so I was originally thinking that might be the issue, but it seems that no output was obtained for some reason. I didn't try it on any other platforms. I didn't allocate a huge amount of time to this, just wanted to see if I could make it work :).

Edit: As an aside, I had to install pyutilib manually. It didn't get installed with pip install pao.

I'm using google collab and am trying to import pao but ran into the following error: "ImportError: cannot import name 'SimpleBlock' from 'pyomo.core'." Am I missing something or is something wrong with the package?

Screenshot of error

Thanks in advance

When a term of the force u*x exists in a constraint where u is binary and x is continuous, this term is overlooked when collecting dual data.

Idea:

1. Reformulate lowest layer using KKT conditions
2. Solve resulting bilevel problem with PCCG

Hello Pao devs! First of all thank you for your awesome work!

I wanted to ask a question about how to correctly build the following bilevel optimization program. More specifically, I am struggling with the first constraint, as it features variables from both the model and the submodel.

$$\begin{align} %r0 & \underset{z_{j}}{\text{maximize}}~\underset{y_{ij}}{\text{minimize}} ~~ \sum_{i \in \mathcal{I}} \sum_{j \in \mathcal{J}} c_{ij}y_{ij}\\ &\\ s.t. \quad & \sum_{i \in \mathcal{I}} y_{ij} \leq z_j~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\forall j \in \mathcal{J}\\ & \sum_{j \in \mathcal{J}} y_{ij} \leq 1 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ \forall i \in \mathcal{I}\\ %r5 & \sum_{i \in \mathcal{I}} \sum_{j \in \mathcal{J}} y_{cp} = min(|I|,|J|)\\ & y_{ij} \in \{0,1\} ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ \forall (i,j) \in \mathcal{I} \times \mathcal{J}\\ & z_{j} \in \{0,1\} ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ \forall (i,j) \in \mathcal{J}\\ \end{align}$$

If I initialize the constraint within the submodel e.g

model.submodel.each_taxi_serves_one_user_max = pe.Constraint(
            model.users, rule=each_taxi_serves_one_user_max
        )

File "/home/franciscodelio/francisco_python_projects/Dynamic_Assignment_Problem/static_optimization/solvers.py", line 311, in find_optimal_solution
    solver.solve(self.model, mip_solver="cbc")
  File "/home/franciscodelio/.cache/pypoetry/virtualenvs/dynamic-assignment-problem-ASJTY80e-py3.10/lib/python3.10/site-packages/pao/pyomo/solvers/mpr_solvers.py", line 36, in solve
    return super().solve(model, **options)
  File "/home/franciscodelio/.cache/pypoetry/virtualenvs/dynamic-assignment-problem-ASJTY80e-py3.10/lib/python3.10/site-packages/pao/pyomo/solver.py", line 55, in solve
    mp, soln_manager = convert_pyomo2MultilevelProblem(model, inequalities=True)
  File "/home/franciscodelio/.cache/pypoetry/virtualenvs/dynamic-assignment-problem-ASJTY80e-py3.10/lib/python3.10/site-packages/pao/pyomo/convert.py", line 554, in convert_pyomo2MultilevelProblem
    node.initialize_level(L, inequalities, var, vidmap, levelmap)
  File "/home/franciscodelio/.cache/pypoetry/virtualenvs/dynamic-assignment-problem-ASJTY80e-py3.10/lib/python3.10/site-packages/pao/pyomo/convert.py", line 204, in initialize_level
    t, nid, j = vidmap[vid]
KeyError: 140037957243440

then I can not access the variable in the model $z$ and I get the error above, and if I initialize the constraint from the model e.g

model.each_taxi_serves_one_user_max = pe.Constraint(
            model.users, rule=each_taxi_serves_one_user_max
        )

the program executes with no apparent errors, however the solutions provided violate the 1st constraint.

Hope you can help! Leaving a reproducible example below where the first constraint is initialized from the model (case 2).

from scipy.spatial import distance
import numpy as np
import pandas as pd
import pyomo.environ as pe
import pyomo.opt as po
import pao


"""Initializing dummy variables and useful functions"""
def generate_2d_coordinate():
    return tuple(np.random.uniform(low=0.0, high=100.0, size=2))

def all_unique(lst):
  return len(lst) == len(set(lst))

available_taxis_ids = [0,1,2,3,4]
users_ids = [0,1,2,3,4,6,7]

n_assignments = min(len(available_taxis_ids),len(users_ids))

taxi_coordinates = [
    generate_2d_coordinate() for _ in available_taxis_ids
]
destination_coordinates = [
    generate_2d_coordinate() for _ in users_ids
]
cost_matrix = distance.cdist(
    taxi_coordinates, destination_coordinates, "euclidean"
)

assignment_cost = (
    pd.DataFrame(
        cost_matrix, index=available_taxis_ids, columns=users_ids
    )
    .stack()
    .to_dict()
)

"""Build and solves the Pyomo model for the MILP problem."""
model = pe.ConcreteModel()

model.taxis = pe.Set(initialize=available_taxis_ids)
model.users = pe.Set(initialize=users_ids)
model.z = pe.Var(model.users, domain=pe.NonNegativeReals)

model.submodel = pao.pyomo.SubModel(fixed=[model.z])
model.submodel.taxis = pe.Set(initialize=available_taxis_ids)
model.submodel.users = pe.Set(initialize=users_ids)

model.submodel.taxis = pe.Set(initialize=available_taxis_ids)
model.submodel.users = pe.Set(initialize=users_ids)
model.submodel.assignment_cost = pe.Param(
    model.submodel.taxis, model.submodel.users, initialize=assignment_cost
)
model.submodel.y = pe.Var(model.submodel.taxis, model.submodel.users, domain=pe.NonNegativeReals)

expression = sum(
    model.submodel.assignment_cost[c, k] * model.submodel.y[c, k]
    for c in model.submodel.taxis
    for k in model.submodel.users
)

model.obj = pe.Objective(sense=pe.maximize, expr=expression)
model.submodel.obj = pe.Objective(sense=pe.minimize, expr=expression)

def each_taxi_serves_one_user_max(model, k):
    # assign at most one taxi c to every user.
    constraint = sum(model.submodel.y[c, k] for c in model.submodel.taxis) <= model.z[k]
    return constraint

def each_user_is_served_by_one_taxi_max(submodel, c):
    # a taxi c can only be assign to one given user k.
    constraint = sum(submodel.y[c, k] for k in submodel.users) <= 1
    return constraint

def maximal_injective_assignment(submodel):
    # a minimum of min(|Taxis|,|Users|) matches ought to be made.
    constraint = (
        sum(submodel.y[c, k] for k in submodel.users for c in submodel.taxis)
        == n_assignments
    )
    return constraint

model.each_taxi_serves_one_user_max = pe.Constraint(
    model.users, rule=each_taxi_serves_one_user_max
)

model.submodel.each_user_is_served_by_one_taxi_max = pe.Constraint(
    model.submodel.taxis, rule=each_user_is_served_by_one_taxi_max
)
model.submodel.maximal_injective_assignment = pe.Constraint(
    rule=maximal_injective_assignment
)

with pao.pyomo.Solver('pao.pyomo.FA') as solver:
    solver.solve(model, mip_solver="cbc")

"""Processing Results"""

output_df = (
            pd.Series(model.submodel.y.extract_values())
            .reset_index()
            .rename(columns={"level_0": "taxis", "level_1": "users", 0: "y"})
        )

solution_is_binary = set(output_df.y.unique()).issubset({0, 1})
assert solution_is_binary

solution_df = output_df.loc[lambda x: x.y == 1, ["users", "taxis"]].reset_index(
            drop=True
        )

print(solution_df)
assert all_unique(solution_df.users)

users  taxis
0      0      0
1      2      1
2      6      2
3      0      3
4      7      4

AssertionError                            

See #44 for a simple example that we should be able to address when neither variable in a bilinear term is integer.

Existing methods in PAO only apply to optimistic formulations, where the follower can select their best-possible alternatives when multiple solutions have the same value. We need to consider pessimistic formulations as well.

The example is solved by PCCG and CPLEX as follows.
https://pao.readthedocs.io/en/latest/examples.html

import pyomo.environ as pe
from pao.pyomo import *

# Create a model object
M = pe.ConcreteModel()

# Define decision variables
M.x = pe.Var(bounds=(0,None))
M.y = pe.Var(bounds=(0,None))

# Define the upper-level objective
M.o = pe.Objective(expr=M.x - 4*M.y)

# Create a SubModel component to declare a lower-level problem
# The variable M.x is fixed in this lower-level problem
M.L = SubModel(fixed=M.x)

# Define the lower-level objective
M.L.o = pe.Objective(expr=M.y)

# Define lower-level constraints
M.L.c1 = pe.Constraint(expr=   -M.x -   M.y <= -3)
M.L.c2 = pe.Constraint(expr= -2*M.x +   M.y <=  0)
M.L.c3 = pe.Constraint(expr=  2*M.x +   M.y <= 12)
M.L.c4 = pe.Constraint(expr=  3*M.x - 2*M.y <=  4)

# Create a solver and apply it
with Solver('pao.pyomo.PCCG',mip_solver="cplex") as solver:
    results = solver.solve(M)

# The final solution is loaded into the model
print(M.x.value)
print(M.y.value)

The solution is x=3 and y=6, which is different from the solution x=4 and y=4 obtained by FA. Th solution, x=3 and y=6, is inorrect since y should be 2.5 to reduce its objective function if x=3.
Then I check the status as follows.

print(results.solver.termination_condition)
print(results.check_optimal_termination())

The outputs are

maxIterations
False

How to set the the max iteration?
The version of the CPLEX for this test is 22.1.0.0.

Thanks.

I am having issues installing pao and using it in jupyter notebook. I have installed and uninstalled several times using pip and git clone.

If I try "from pao.pyomo import *". I get "AttributeError : partially initialized module 'pao' has no attribute 'bilevel' (most likely due to a circular import)".

If I try "from pao import *". I get "AttributeError: partially initialized module 'pao' has no attribute 'bilevel' (most likely due to a circular import)"

I was able to load "from pao.bilevel import *" on the second try after getting an error the first time. however, I am not sure if I am able to get full access to its features.

update:
I was able to fix the issue by uninstalling pyomo and pao, and just installing pao which installed the correct version of pyomo for it to work.

Hi,
I copy the code of Simple Examples and modify the domain of the lower level variable to type integer.
M.y = pe.Var(bounds=(0,None),within=pe.Integers)
The whole code is provided as follows.

import pyomo.environ as pe
from pao.pyomo import *

# Create a model object
M = pe.ConcreteModel()

# Define decision variables
M.x = pe.Var(bounds=(0,None))
M.y = pe.Var(bounds=(0,None),within=pe.Integers)

# Define the upper-level objective
M.o = pe.Objective(expr=M.x - 4*M.y)

# Create a SubModel component to declare a lower-level problem
# The variable M.x is fixed in this lower-level problem
M.L = SubModel(fixed=M.x)

# Define the lower-level objective
M.L.o = pe.Objective(expr=M.y)

# Define lower-level constraints
M.L.c1 = pe.Constraint(expr=   -M.x -   M.y <= -3)
M.L.c2 = pe.Constraint(expr= -2*M.x +   M.y <=  0)
M.L.c3 = pe.Constraint(expr=  2*M.x +   M.y <= 12)
M.L.c4 = pe.Constraint(expr=  3*M.x - 2*M.y <=  4)

# Create a solver and apply it
with Solver('pao.pyomo.PCCG',mip_solver="cplex") as solver:
    results = solver.solve(M)

# The final solution is loaded into the model
print(M.x.value)
print(M.y.value)

Then an error returns
RuntimeError: ERROR! ERROR! Master: Could not find optimal solution,
however, a solution x=4, y=4 can be found when the domain of variables is set to real.

Please go to my fork, branch "issue_20201903" and run the test script:
/pao/bilevel/tests/test_linear_bilevel.py

This script applies the solver1.py 'pao.bilevel.ld' to the following test cases:

1. t5.py
2. bqp_example2_indexed.py

The bqp_example2_indexed.py fails on line 106 of solver1.py with the following error:
File ".../sandbox/pao/pao/bilevel/solvers/solver1.py", line 106, in apply_solver
for data in vuid.find_component_on(self._instance).values():
AttributeError: '_GeneralVarData' object has no attribute 'values'

Whereas t5.py has the vuid in tdata.fixed such that vuid.find_component_on(self._instance) returns an object of type pyomo.core.base.var.SimpleVar, the bqp_example2_indexed.py returns an object of type pyomo.core.base.var._GeneralVarData, which does not have the attribute "values()" and throws an error.

Since both SimpleVar and _GeneralVarData are leaf nodes on the expression tree, why does one have the "values()" functionality and the other does not? Do we have to parse these data objects differently, or is there a better way to fix this issue?

Hello,

Does there exist a persistent solver interface? Perhaps analogous to GurobiPersistent?

Thanks,
Jake

Restrict the API of PAO solvers to only the options that are configured.

This will make the documentation more straightforward, and it will make the API more robust. Subsolver options can only be passed-in by defining a solver object that is configured.

hello , how can I express bilinear terms in the objective function

@anyacastillo raised a concern that the user may not know all of the variables that need to be fixed.

For example, if the user transforms upper-level constraints, such transformations may create variables that the user doesn't know about from a modeling perspective.

Does it make sense to change the internal transformation bookkeeping to track unfixed variables? I could see a similar issue arising in the lower-level constraints. Where a user specifies that these variables are unfixed, but doesn't know to specify variables that are created through subsequent transformations.

Hi,

I hope you are doing very well.

I have seen you have added some examples of bilevel optimization problems in the bilevel/examples folder. This is great to see these examples, however, my issue is that we cannot solve them and we have no idea what solver we should use. I have seen on a question posted online that we need to transform our problem using TransformatonFactory to be able to be solved by a particular solver. Would you be able to either share a more real-world example of this so we could learn more or perhaps explain this? Unfortunately, not many examples are available to learn from.

Thanks,

If the upper level model has binery varibles, continous varibles and interger varibles,it's a linering mixed integer model. And the lower level model has binery varibles, continous varibles and interger varibles, too. And it's also a linering mixed integer model. PAO.
with Solver('pao.pyomo.PCCG',mip_solver="gurobi") as solver:
results = solver.solve(M)
if I use the PCCG solver or other solvers, could the problem be solved?
Intermediate linking variables may be binery varibles or continous varibles

I applied the pao.pyomo.FA solver to the following problem:

from pao.pyomo import *

M = pe.ConcreteModel()
M.x = pe.Var(bounds=(0, None), domain=pe.Reals)

M.L = SubModel(fixed=M.x)
M.L.y = pe.Var(bounds=(0, 10), domain=pe.Reals) #This is the relevant code line

M.L.ymax = 3
M.xmin = 2

def ul_obj_rule(M):
    """upper-level objective"""
    mo = M.model()
    return(mo.x + mo.L.y)

def ll_obj_rule(M):
    mo = M.model()
    return(mo.x + mo.L.y)

def ul_constraint_rule(M):
    mod = M.model()
    return(M.x >= mod.xmin)

def ll_constraint_rule(M):
    mod = M.model()
    return(mod.L.y <= mod.x + M.ymax)

M.obj = pe.Objective(rule=ul_obj_rule, sense=pe.minimize)
M.L.obj = pe.Objective(rule=ll_obj_rule, sense=pe.maximize)
M.c1 = pe.Constraint(rule=ul_constraint_rule)
M.L.c1 = pe.Constraint(rule=ll_constraint_rule)

solver = Solver('pao.pyomo.FA')
results = solver.solve(M, tee=True)

Since x is the upper-level variable and y is the lower-level variable, the solution should be x=2 and y=5 according to the objectives and constraints (with an objective value of 7).

Though, the solver sets y to 1e-4 with an objective value close to 2, so the upper-level seems to determine not only x, but also y. If I set the upper-bound of y to "None" in the "relevant code line", the expected objective value of 7 is found. For an upper bound of 1e5, the value for y is set to 1 (instead of 5).

Running the same code with the pao.pyomo.PCCG solver or pao.pyomo.REG, this behaviour doesn't occur.

Since my model contains only linear equations, real upper- and lower-level variables and is a bilevel problem, I assume that all requirements to apply the pao.pyomo.FA solver are fulfilled. Can anyone please help me understand why the described (unexpected) results are found?

Hi all, first of all thank you very much for PAO. I think pyomo is in general a great tool, and the extension to bilevel problems looks wonderful.

I am trying to get my head around the PAO functioning compared to the pyomo.bilevel package, and though I am able to run simple tests I am getting unexpected results. I have a main which looks like

from pao.bilevel import solvers
import models

if __name__ == "__main__":
   # global quadratic bilevel
    solver = solvers.solver4.BilevelSolver4()    

    # BARD. Optimal values are x=2, y=1        
    print("Bard 5.1.1")
    model = models.bard_511()
    solver.options.solver = "cbc"        
    results = solver.solve(model)
    print('optimal x:\t', model.x.value)
    print('optimal y:\t', model.y.value)
    print('Is constraint 1 satisfied?\t', model.sub.c1.expr())
    print('Is constraint 2 satisfied?\t', model.sub.c2.expr())
    print('Is constraint 3 satisfied?\t', model.sub.c3.expr())
    print('Is constraint 4 satisfied?\t', model.sub.c4.expr())
    print('=====')

    # TESIUS BASIC MODEL. Optimal values are x = 20, y1=10, y2=0
    print("Tesius 1.0")
    model = models.tesius_10()
    solver.options.solver = "couenne"        
    results = solver.solve(model)
    print('optimal x:\t', model.x.value)
    print('optimal y1:\t', model.sub.y1.value)
    print('optimal y2:\t', model.sub.y2.value)
    print('Is constraint satisfied?\t', model.sub.c.expr())
    print('=====')

where models.py contains the definition for the "bard_511" and "tesius_10" models

from pyomo.environ import ConcreteModel, Var, Objective, Constraint, maximize, minimize
from pao.bilevel import SubModel


def bard_511():
    M = ConcreteModel()
    M.x = Var(bounds=(0,None))
    M.y = Var(bounds=(0,None))
    M.o = Objective(expr=M.x - 4*M.y)

    M.sub = SubModel(fixed=M.x)
    M.sub.o = Objective(expr=M.y)
    M.sub.c1 = Constraint(expr=-M.x - M.y <= -3)
    M.sub.c2 = Constraint(expr=-2*M.x +   M.y <=  0)
    M.sub.c3 = Constraint(expr= 2*M.x +   M.y <= 12)
    M.sub.c4 = Constraint(expr=-3*M.x + 2*M.y <= -4)

    return M

def tesius_10():
    M = ConcreteModel()
    M.x = Var(bounds=(0,25))
    M.sub = SubModel(fixed=M.x)    
    M.sub.y1 = Var(bounds=(0,11))
    M.sub.y2 = Var(bounds=(0,11))    

    M.o = Objective(expr=M.x*M.sub.y1, sense=maximize)
    M.sub.o = Objective(expr=M.x*M.sub.y1 + 20*M.sub.y2, sense=minimize)

    M.sub.c = Constraint(expr=M.sub.y1 + M.sub.y2 == 10)

    return M

The tesius_10 aims at defining a problem where the "leader" maximizes x*y1, while the follower minimizes x*y1 + 20*y2. A constraint links the y1 and y2 variables, as their sum should equal 10. The economically optimal choice should be that the leader sets x=20, in order to make the follower indifferent between "turning on" y1=10, y2=0 and y2=10-y1. Should the leader "attempts" to raise x above the threshold level of 20, y1 (and thus the upper objective function) would be unfavored against y2. This would be my expected result... but what I see from the execution of the main is

Tesius 1.0
WARNING: Unused arguments in the bigM map! These arguments were not used by
    the transformation:
        None
optimal x:       25.0
optimal y1:      11.0
optimal y2:      10.999999999999947
Is constraint satisfied?         False

Am I misunderstanding the way "bilevel" structure is represented by PAO? Or am I missing something else in the way PAO should be used?

Incidentally, the bard_511 test looks fine, i.e. I get the expected result (with the same solver4, using cbc or couenne if that matters)

Thanks a lot

Enrico

Hello,

Can you please look at this problem that I face with a bilevel problem? I reproduce it with an example problem found in the documentation.

import pyomo.environ as pyo
from pao.pyomo import *

# Create a model object
M = pyo.ConcreteModel()

# Define decision variables
M.x = pyo.Var(bounds=(0,None))
M.y = pyo.Var(bounds=(0,None))

# Define the upper-level objective
M.o = pyo.Objective(expr=M.x - 4*M.y)

# Create a SubModel component to declare a lower-level problem
# The variable M.x is fixed in this lower-level problem
M.L = SubModel(fixed=M.x)

# Define the lower-level objective
M.L.o = pyo.Objective(expr=M.y)

# Define lower-level constraints
M.L.c1 = pyo.Constraint(expr=   -M.x - M.y <= -3)
M.L.c2 = pyo.Constraint(expr= -2*M.x + M.y <=  0)
M.L.c3 = pyo.Constraint(expr=  2*M.x + M.y <= 12)
M.L.c4 = pyo.Constraint(expr=  3*M.x - 2*M.y <=  4)

pyo.SolverFactory('gurobi').solve(M)

The error that I get is:

WARNING: Empty constraint block written in LP format - solver may error
WARNING: Loading a SolverResults object with a warning status into
    model.name="unknown";
      - termination condition: infeasibleOrUnbounded
      - message from solver: Problem proven to be infeasible or unbounded.

This should be natural and appropriate.

This probably requires a generalization of the MPRs

The documentation illustrates (a) bilevel problems with multiple sub-problems, and (b) tri-level problems.

Do we need to include more general examples? Are there motivating applications?

Right now, PAO only pulls in variables/expressions/constraints from Block/Objective/Constraint expressions.

We should (a) throw errors when other components are encountered or (b) pull-out the modeling representations from those components.

I have written an extensive explanation of my problem at: https://stackoverflow.com/questions/75161332/pyomo-interfaces-dont-recognize-linear-multilevel-problems-in-pyomo-pao-or-mis

Thanks for your time & attention.

Dear pao-Team,

thank you for the development of this amazing package! :)
I would like to report a minor bug in the function check_termination in pccg_solver.py.
If execute_PCCG_solver reaches its iteration-limit and the terminition is checked, while the solver config.quiet-parameter is set True, check_termination is called and raises the RuntimeError(f'Error: Upper bound greater than lower bound after {k} iterations and {elapsed} seconds: Obj={LB} UB={UB}'). However, this is resulting in a NameError as k (and elapsed) is not provided to the check_termination-function.
If PRs are considered, I would be happy to solve this Issue.

Thank you in advance!
Best,
Kirill Kuroptev