From c521df1f7545bde6e64da6605aceadd933d208a9 Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Sat, 21 Dec 2024 19:18:54 -0700 Subject: [PATCH 1/9] only add inequalities for bounds that are not dominated (redundant) --- var_elim/algorithms/replace.py | 45 +++++++++++++------ var_elim/algorithms/tests/test_replacement.py | 27 ++++++++++- 2 files changed, 58 insertions(+), 14 deletions(-) diff --git a/var_elim/algorithms/replace.py b/var_elim/algorithms/replace.py index 544b19d..2c70b38 100644 --- a/var_elim/algorithms/replace.py +++ b/var_elim/algorithms/replace.py @@ -32,6 +32,7 @@ ) from pyomo.contrib.incidence_analysis import IncidenceGraphInterface from pyomo.contrib.incidence_analysis.config import IncidenceMethod +import pyomo.contrib.fbbt.fbbt as fbbt from pyomo.common.modeling import unique_component_name from pyomo.common.timing import HierarchicalTimer @@ -132,21 +133,27 @@ def add_bounds_to_expr(var, var_expr): constraints on the expression if the variables replaced were bounded Each constraint added to the bound_cons list is indexed by var_name_ub or - var_name_lb depending upon whihc bound it adds to the expression + var_name_lb depending upon which bound it adds to the expression """ - if var.ub is None and var.lb is None: - lb_expr = None - ub_expr = None - elif var.lb is not None and var.ub is None: - lb_expr = var_expr >= var.lb - ub_expr = None - elif var.ub is not None and var.lb is None: - lb_expr = None - ub_expr = var_expr <= var.ub + lb, ub = fbbt.compute_bounds_on_expr(var_expr) + # TODO: Use a tolerance for bound equivalence here? + if var.lb is not None and (lb is None or lb < var.lb): + # We add a lower bound constraint if the variable has a lower bound + # and it is not dominated by the bounds on its defining expression. + add_lb_con = True else: - lb_expr = var_expr >= var.lb - ub_expr = var_expr <= var.ub - + add_lb_con = False + if var.ub is not None and (ub is None or ub > var.ub): + add_ub_con = True + else: + add_ub_con = False + + # TODO: If our expression is a linear (monotonic) function of a single + # variable, we can propagate its bound back to the independent variable. + + ub_expr = var_expr <= var.ub if add_ub_con else None + lb_expr = var_expr >= var.lb if add_lb_con else None + return lb_expr, ub_expr @@ -355,10 +362,19 @@ def eliminate_variables( # Update data structures for eliminated variables for var in var_order: var_expr = substitution_map[id(var)] + # Here, we should return None if adding inequalities was not necessary. + # In these cases, we either (a) do nothing (the bounds are dominated by + # existing bounds on dependent variables) or (b) translate the bounds + # to the (single, linear) defining variable. + # If we're updating bounds on defining variables in-place, I need to make + # sure these bounds are accounted for in later steps. lb_expr, ub_expr = add_bounds_to_expr(var, var_expr) lb_name = var.name + "_lb" ub_name = var.name + "_ub" if lb_expr is not None and type(lb_expr) is not bool: + # Checking for trivial constraints here seems unnecessary. But in this + # way we're masking an infeasibility. I guess we've determined that these + # don't happen in the problems we care about? if lb_expr is False: raise RuntimeError("Lower bound resolved to trivial infeasible constraint") bound_con_set.add(lb_name) @@ -382,6 +398,9 @@ def eliminate_variables( # iterate over variables-to-eliminate, and perform the substitution for each # (new) constraint that they're adjacent to. This way we don't check every # constraint if we're only eliminating a small number of variables. + # + # This loop replaces variables in equality and inequality constraints, but + # doesn't convert bounds on eliminated variables to inequalities. for con in igraph.constraints: if ( id(con) not in elim_con_set diff --git a/var_elim/algorithms/tests/test_replacement.py b/var_elim/algorithms/tests/test_replacement.py index 1fb2c70..bf2ae5f 100644 --- a/var_elim/algorithms/tests/test_replacement.py +++ b/var_elim/algorithms/tests/test_replacement.py @@ -187,6 +187,8 @@ def test_constraints_are_added(self): vars_to_elim = [m.x[1], m.x[2]] cons_to_elim = [m.eq1, m.eq2] + m.x[1].setlb(1) + var_order, con_order = define_elimination_order(vars_to_elim, cons_to_elim) eliminate_variables(m, var_order, con_order) @@ -229,6 +231,27 @@ def test_same_solution(self): pyo.value(x2_lb_con.lower), ) + def test_dominated_bounds(self): + m = pyo.ConcreteModel() + m.x = pyo.Var([1, 2, 3]) + m.eq = pyo.Constraint(pyo.PositiveIntegers) + m.eq[1] = m.x[1] == m.x[2] + m.x[3]**2 + m.obj = pyo.Objective(expr=m.x[1]**2 + m.x[2]**2 + m.x[3]**2) + e = m.x[2] + m.x[3]**2 + m.x[2].setlb(-2) + m.x[2].setub(2) + m.x[3].setlb(-1) + m.x[3].setub(4) + # Bounds on expression are (-2, 18) + m.x[1].setlb(-2) + m.x[1].setub(20) + + var_elim = [m.x[1]] + con_elim = [m.eq[1]] + var_exprs, var_lb_map, var_ub_map = eliminate_variables(m, var_elim, con_elim) + assert len(m.replaced_variable_bounds_set) == 0 + assert len(m.replaced_variable_bounds) == 0 + class TestReplacementInInequalities: def _make_simple_model(self): @@ -260,7 +283,9 @@ def test_simple_replacement(self): # Make sure new model has correct number of constraints new_igraph = IncidenceGraphInterface(m, include_inequality=True) - assert len(new_igraph.constraints) == 6 + # New constraints only get added for upper bounds, as lower bounds are + # redundant. + assert len(new_igraph.constraints) == 4 # Make sure proper replacement happened here assert ComponentSet(identify_variables(m.eq4.expr)) == ComponentSet(m.y[:]) From e17e22c65df64516d032a671cd25f881cb22196f Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Sun, 22 Dec 2024 16:32:55 -0700 Subject: [PATCH 2/9] update CLI --- var_elim/scripts/analyze_solvetime.py | 1 + var_elim/scripts/collect_results.py | 40 +++++++++++++++---------- var_elim/scripts/write_command_lines.py | 4 +++ 3 files changed, 30 insertions(+), 15 deletions(-) diff --git a/var_elim/scripts/analyze_solvetime.py b/var_elim/scripts/analyze_solvetime.py index 3520fef..03a001b 100644 --- a/var_elim/scripts/analyze_solvetime.py +++ b/var_elim/scripts/analyze_solvetime.py @@ -152,6 +152,7 @@ def main(args): # We need to re-set the callback each time we solve a model htimer.start("solver") solver.config.intermediate_callback = Callback() + # TODO: Option to set tee=True? res = solver.solve(model, tee=False, timer=htimer) htimer.stop("solver") diff --git a/var_elim/scripts/collect_results.py b/var_elim/scripts/collect_results.py index 69221f7..c9dd879 100644 --- a/var_elim/scripts/collect_results.py +++ b/var_elim/scripts/collect_results.py @@ -54,7 +54,8 @@ def main(args): fnames = [basename + "-" + s + ".csv" for s in suffixes] filedir = os.path.dirname(__file__) fpaths = [ - os.path.join(filedir, "results", args.result_type, fname) + #os.path.join(filedir, "results", args.result_type, fname) + os.path.join(args.results_dir, args.result_type, fname) for fname in fnames ] print() @@ -71,7 +72,7 @@ def main(args): # output-suffix overrides suffix suff_str = f"-{args.output_suffix}" output_fname = args.result_type + suff_str + ".csv" - output_fpath = os.path.join(config.get_results_dir(), output_fname) + output_fpath = os.path.join(args.results_dir, output_fname) dfs = [pd.read_csv(fpath) for fpath in fpaths] output_df = pd.concat(dfs, ignore_index=True, join="inner") @@ -101,25 +102,34 @@ def main(args): ), default=None, ) + # We now append results_type to results_dir, so this isn't necessary. + #argparser.add_argument( + # "--output-results-dir", + # help=( + # "Results dir to store collected output" + # " (different from results dir where we look for results to collect)" + # ), + # default=config.get_results_dir(), + #) # HACK: We change the default of the argparser so we can handle it specially # if --method or --model are used. # It's unclear whether this hack will be worth the convenience, but let's try it. - argparser.set_defaults(results_dir=None) + #argparser.set_defaults(results_dir=None) args = argparser.parse_args() - if args.results_dir is None: - if args.method is None and args.model is None: - # If neither method nor model is used (we are collecting all results) - # we put results in the top-level results directory. - args.results_dir = config.get_results_dir() - else: - # If either method or model is used, we put the results in the - # results/structure subdirectory. This is because we don't want the - # top-level results getting polluted with a bunch of files. - resdir = os.path.join(config.get_results_dir(), "solvetime") - config.validate_dir(resdir) - args.results_dir = resdir + #if args.results_dir is None: + # if args.method is None and args.model is None: + # # If neither method nor model is used (we are collecting all results) + # # we put results in the top-level results directory. + # args.results_dir = config.get_results_dir() + # else: + # # If either method or model is used, we put the results in the + # # results/structure subdirectory. This is because we don't want the + # # top-level results getting polluted with a bunch of files. + # resdir = os.path.join(config.get_results_dir(), args.results_dir) + # config.validate_dir(resdir) + # args.results_dir = resdir main(args) diff --git a/var_elim/scripts/write_command_lines.py b/var_elim/scripts/write_command_lines.py index cb834a6..53a2b24 100644 --- a/var_elim/scripts/write_command_lines.py +++ b/var_elim/scripts/write_command_lines.py @@ -75,6 +75,10 @@ def main(args): if args.suffix is not None: for cl in cl_lists: cl.append(f"--suffix={args.suffix}") + if args.results_dir != config.get_results_dir(): + for cl in cl_lists: + cl.append(f"--results-dir={args.results_dir}") + else: results_dir = os.path.join(os.path.dirname(__file__), "results", "sweep") suff_str = "" if args.suffix is None else f"-{args.suffix}" From 0b107c7f2f5bb0a634cbd679ceb7959577ccd81d Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Sun, 22 Dec 2024 23:40:58 -0700 Subject: [PATCH 3/9] use results-dir arg in sweep command lines --- var_elim/scripts/write_sweep_command_lines.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/var_elim/scripts/write_sweep_command_lines.py b/var_elim/scripts/write_sweep_command_lines.py index 74e8b25..4ad3791 100644 --- a/var_elim/scripts/write_sweep_command_lines.py +++ b/var_elim/scripts/write_sweep_command_lines.py @@ -66,9 +66,11 @@ def main(args): # Note that sample arg (and output filename) is base-1 for i in range(1, nsamples_total + 1) ] - if args.suffix is not None: - for cmd in sample_commands: + for cmd in sample_commands: + if args.suffix is not None: cmd.append(f"--suffix={args.suffix}") + if args.results_dir is not None: + cmd.append(f"--results-dir={args.results_dir}") sample_commands_str = [" ".join(cmd) + "\n" for cmd in sample_commands] @@ -95,6 +97,8 @@ def main(args): # suffix will be used to identify the right input files, while # output-suffix will be used for the output file. collect_cmd.append(f"--output_suffix={args.output_suffix}") + if args.results_dir is not None: + collect_cmd.append(f"--results-dir={args.results_dir}") collect_commands.append(collect_cmd) # TODO: Maybe send other arguments to the plotting command? From 65b3069e4feeb6e1acd7ea7e9f330f38a2246cd7 Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Mon, 23 Dec 2024 00:04:56 -0700 Subject: [PATCH 4/9] add tee arg --- var_elim/scripts/analyze_solvetime.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/var_elim/scripts/analyze_solvetime.py b/var_elim/scripts/analyze_solvetime.py index 03a001b..3320f82 100644 --- a/var_elim/scripts/analyze_solvetime.py +++ b/var_elim/scripts/analyze_solvetime.py @@ -152,8 +152,10 @@ def main(args): # We need to re-set the callback each time we solve a model htimer.start("solver") solver.config.intermediate_callback = Callback() - # TODO: Option to set tee=True? - res = solver.solve(model, tee=False, timer=htimer) + solver.config.options["linear_solver"] = "ma57" + solver.config.options["ma57_pivot_order"] = 4 + solver.config.options["print_user_options"] = "yes" + res = solver.solve(model, tee=args.tee, timer=htimer) htimer.stop("solver") timer.toc("Solve model") @@ -295,6 +297,7 @@ def main(args): default=None, help="Basename for file to write results to", ) + argparser.add_argument("--tee", action="store_true", help="Stream solver log to stdout") # HACK: We change the default of the argparser so we can handle it specially # if --method or --model are used. # It's unclear whether this hack will be worth the convenience, but let's try it. From 908e5d39fc8dfccf76f11f05a54cb73c8dd6121b Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Fri, 3 Jan 2025 13:27:16 -0700 Subject: [PATCH 5/9] propagate bounds when eliminating a linear, degree-2 constraint --- var_elim/algorithms/replace.py | 42 ++++++++++++++++++- var_elim/algorithms/tests/test_replacement.py | 28 ++++++------- 2 files changed, 54 insertions(+), 16 deletions(-) diff --git a/var_elim/algorithms/replace.py b/var_elim/algorithms/replace.py index 2c70b38..0e32724 100644 --- a/var_elim/algorithms/replace.py +++ b/var_elim/algorithms/replace.py @@ -135,7 +135,13 @@ def add_bounds_to_expr(var, var_expr): Each constraint added to the bound_cons list is indexed by var_name_ub or var_name_lb depending upon which bound it adds to the expression """ - lb, ub = fbbt.compute_bounds_on_expr(var_expr) + + # This seems to be fast and eliminate many bounds (for problems that have them), + # but leads to slowdowns in some instances. Basically, I don't know when the + # bounds were and weren't "useful for the algorithm" even when they're not + # necessary. + #lb, ub = fbbt.compute_bounds_on_expr(var_expr) + lb, ub = None, None # TODO: Use a tolerance for bound equivalence here? if var.lb is not None and (lb is None or lb < var.lb): # We add a lower bound constraint if the variable has a lower bound @@ -150,6 +156,40 @@ def add_bounds_to_expr(var, var_expr): # TODO: If our expression is a linear (monotonic) function of a single # variable, we can propagate its bound back to the independent variable. + # - Use standard-repn to check if defining expr is linear-degree-1 + # - Extract constant and coefficient + # - Set bounds on x as e.g. (y^L-b)/m + propagate_bounds = True + # Propagate bounds from eliminated variable to "defining variable" if the + # two bounds are equivalent. + if propagate_bounds: + # TODO: Only compute standard repn if defined variable has bounds? + # ... or just cache and reuse standard repn... + repn = generate_standard_repn(var_expr, compute_values=True, quadratic=False) + if ( + len(repn.nonlinear_vars) == 0 # Expression is affine + and len(repn.linear_vars) == 1 # and only contains one variable. + # ^ Can linear_vars contain duplicates? + ): + offset = repn.constant + coef = repn.linear_coefs[0] + defining_var = repn.linear_vars[0] + + # Bounds implied by bounds on defined variable + lb = None if var.lb is None else (var.lb - offset) / coef + ub = None if var.ub is None else (var.ub - offset) / coef + lbkey = lambda b: -float("inf") if b is None else b + ubkey = lambda b: float("inf") if b is None else b + # Take the more restrictive of the two bounds + lb = max(lb, defining_var.lb, key=lbkey) + ub = min(ub, defining_var.ub, key=ubkey) + lb = None if lb == -float("inf") else lb + ub = None if ub == float("inf") else ub + defining_var.setlb(lb) + defining_var.setub(ub) + + add_ub_con = False + add_lb_con = False ub_expr = var_expr <= var.ub if add_ub_con else None lb_expr = var_expr >= var.lb if add_lb_con else None diff --git a/var_elim/algorithms/tests/test_replacement.py b/var_elim/algorithms/tests/test_replacement.py index bf2ae5f..cc7aa0d 100644 --- a/var_elim/algorithms/tests/test_replacement.py +++ b/var_elim/algorithms/tests/test_replacement.py @@ -231,26 +231,26 @@ def test_same_solution(self): pyo.value(x2_lb_con.lower), ) - def test_dominated_bounds(self): + def test_propagate_bounds_linear(self): + # Test that we correctly propagate bounds when eliminating a linear constraint + # with two variables. m = pyo.ConcreteModel() m.x = pyo.Var([1, 2, 3]) + for i in range(1, 4): + m.x[i].setlb(-i) + m.x[i].setub(i) m.eq = pyo.Constraint(pyo.PositiveIntegers) - m.eq[1] = m.x[1] == m.x[2] + m.x[3]**2 + m.eq[1] = m.x[1] == 2 * m.x[2] + 3 + m.eq[2] = m.x[2] * m.x[2] == m.x[3] m.obj = pyo.Objective(expr=m.x[1]**2 + m.x[2]**2 + m.x[3]**2) - e = m.x[2] + m.x[3]**2 - m.x[2].setlb(-2) - m.x[2].setub(2) - m.x[3].setlb(-1) - m.x[3].setub(4) - # Bounds on expression are (-2, 18) - m.x[1].setlb(-2) - m.x[1].setub(20) var_elim = [m.x[1]] con_elim = [m.eq[1]] var_exprs, var_lb_map, var_ub_map = eliminate_variables(m, var_elim, con_elim) - assert len(m.replaced_variable_bounds_set) == 0 - assert len(m.replaced_variable_bounds) == 0 + assert m.x[1] not in var_lb_map + assert m.x[1] not in var_ub_map + assert math.isclose(m.x[2].lb, -2.0, abs_tol=1e-8) + assert math.isclose(m.x[2].ub, -1.0, abs_tol=1e-8) class TestReplacementInInequalities: @@ -283,9 +283,7 @@ def test_simple_replacement(self): # Make sure new model has correct number of constraints new_igraph = IncidenceGraphInterface(m, include_inequality=True) - # New constraints only get added for upper bounds, as lower bounds are - # redundant. - assert len(new_igraph.constraints) == 4 + assert len(new_igraph.constraints) == 6 # Make sure proper replacement happened here assert ComponentSet(identify_variables(m.eq4.expr)) == ComponentSet(m.y[:]) From e01b3f8277c9b45dc598478ccc5ac094677edcda Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Fri, 3 Jan 2025 13:45:33 -0700 Subject: [PATCH 6/9] objective and tee --- var_elim/scripts/analyze_solvetime.py | 14 ++++++++++++-- var_elim/scripts/write_command_lines.py | 4 ++++ 2 files changed, 16 insertions(+), 2 deletions(-) diff --git a/var_elim/scripts/analyze_solvetime.py b/var_elim/scripts/analyze_solvetime.py index 3320f82..9552ef6 100644 --- a/var_elim/scripts/analyze_solvetime.py +++ b/var_elim/scripts/analyze_solvetime.py @@ -79,6 +79,14 @@ def solve_reduced(m, tee=True, callback=matching_elim_callback): return m +def get_objective_value(m): + objs = list(m.component_data_objects(pyo.Objective, active=True)) + if len(objs) > 1: + raise RuntimeError(f"Model has {len(objs)} objectives") + else: + return pyo.value(objs[0].expr) + + def main(args): if args.model is not None: models = [(args.model, config.CONSTRUCTOR_LOOKUP[args.model])] @@ -114,6 +122,7 @@ def main(args): "success": [], "feasible": [], "max-infeasibility": [], + "objective-value": [], "elim-time": [], "solve-time": [], "init-time": [], @@ -152,8 +161,6 @@ def main(args): # We need to re-set the callback each time we solve a model htimer.start("solver") solver.config.intermediate_callback = Callback() - solver.config.options["linear_solver"] = "ma57" - solver.config.options["ma57_pivot_order"] = 4 solver.config.options["print_user_options"] = "yes" res = solver.solve(model, tee=args.tee, timer=htimer) htimer.stop("solver") @@ -251,12 +258,15 @@ def main(args): print(f"Time spent initializing solver: {init_time}") print() + objval = get_objective_value(model) + data["model"].append(mname) data["method"].append(elim_name) data["elim-time"].append(elim_time) data["success"].append(success) data["feasible"].append(valid) data["max-infeasibility"].append(max_infeas) + data["objective-value"].append(objval) data["solve-time"].append(solve_time) data["init-time"].append(init_time) # This is time to build the model, and has nothing to do with the elimination diff --git a/var_elim/scripts/write_command_lines.py b/var_elim/scripts/write_command_lines.py index 53a2b24..5cf4106 100644 --- a/var_elim/scripts/write_command_lines.py +++ b/var_elim/scripts/write_command_lines.py @@ -78,6 +78,9 @@ def main(args): if args.results_dir != config.get_results_dir(): for cl in cl_lists: cl.append(f"--results-dir={args.results_dir}") + if args.tee: + for cl in cl_lists: + cl.append("--tee") else: results_dir = os.path.join(os.path.dirname(__file__), "results", "sweep") @@ -117,5 +120,6 @@ def main(args): help="Parallelize by model, method, or both. Default=both", default="both", ) + argparser.add_argument("--tee", action="store_true", help="Solver log to stdout") args = argparser.parse_args() main(args) From 5493bc403dbc8af61f0bdaaf5111520fc08f9215 Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Fri, 3 Jan 2025 13:49:42 -0700 Subject: [PATCH 7/9] switch bounds if dividing by negative coefficient --- var_elim/algorithms/replace.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/var_elim/algorithms/replace.py b/var_elim/algorithms/replace.py index 0e32724..b03fe69 100644 --- a/var_elim/algorithms/replace.py +++ b/var_elim/algorithms/replace.py @@ -178,6 +178,8 @@ def add_bounds_to_expr(var, var_expr): # Bounds implied by bounds on defined variable lb = None if var.lb is None else (var.lb - offset) / coef ub = None if var.ub is None else (var.ub - offset) / coef + if coef < 0.0: + lb, ub = ub, lb lbkey = lambda b: -float("inf") if b is None else b ubkey = lambda b: float("inf") if b is None else b # Take the more restrictive of the two bounds From 527bfaace17c0ed4186bb4972ba86c60eb7d0517 Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Fri, 3 Jan 2025 15:49:57 -0700 Subject: [PATCH 8/9] dont add inequalities when eliminating degree-1 constraints --- var_elim/algorithms/replace.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/var_elim/algorithms/replace.py b/var_elim/algorithms/replace.py index b03fe69..79710f4 100644 --- a/var_elim/algorithms/replace.py +++ b/var_elim/algorithms/replace.py @@ -192,6 +192,15 @@ def add_bounds_to_expr(var, var_expr): add_ub_con = False add_lb_con = False + elif len(repn.nonlinear_vars) == 0 and len(repn.linear_vars) == 0: + # We are eliminating a fixed variable + if ( + (var.ub is not None and repn.constant > var.ub) + or (var.lb is not None and repn.constant < var.lb) + ): + raise ValueError("Attempting to fix a variable outside its bounds") + add_ub_con = False + add_lb_con = False ub_expr = var_expr <= var.ub if add_ub_con else None lb_expr = var_expr >= var.lb if add_lb_con else None From 4d3cba9f8cb0c8b2531bf1aec16255b8b654430f Mon Sep 17 00:00:00 2001 From: Robert Parker Date: Fri, 3 Jan 2025 15:50:32 -0700 Subject: [PATCH 9/9] remove redundant --- var_elim/algorithms/replace.py | 74 +++++++++++++++++----------------- 1 file changed, 36 insertions(+), 38 deletions(-) diff --git a/var_elim/algorithms/replace.py b/var_elim/algorithms/replace.py index 79710f4..2b933bd 100644 --- a/var_elim/algorithms/replace.py +++ b/var_elim/algorithms/replace.py @@ -159,48 +159,46 @@ def add_bounds_to_expr(var, var_expr): # - Use standard-repn to check if defining expr is linear-degree-1 # - Extract constant and coefficient # - Set bounds on x as e.g. (y^L-b)/m - propagate_bounds = True # Propagate bounds from eliminated variable to "defining variable" if the # two bounds are equivalent. - if propagate_bounds: - # TODO: Only compute standard repn if defined variable has bounds? - # ... or just cache and reuse standard repn... - repn = generate_standard_repn(var_expr, compute_values=True, quadratic=False) + # TODO: Only compute standard repn if defined variable has bounds? + # ... or just cache and reuse standard repn... + repn = generate_standard_repn(var_expr, compute_values=True, quadratic=False) + if ( + len(repn.nonlinear_vars) == 0 # Expression is affine + and len(repn.linear_vars) == 1 # and only contains one variable. + # ^ Can linear_vars contain duplicates? + ): + offset = repn.constant + coef = repn.linear_coefs[0] + defining_var = repn.linear_vars[0] + + # Bounds implied by bounds on defined variable + lb = None if var.lb is None else (var.lb - offset) / coef + ub = None if var.ub is None else (var.ub - offset) / coef + if coef < 0.0: + lb, ub = ub, lb + lbkey = lambda b: -float("inf") if b is None else b + ubkey = lambda b: float("inf") if b is None else b + # Take the more restrictive of the two bounds + lb = max(lb, defining_var.lb, key=lbkey) + ub = min(ub, defining_var.ub, key=ubkey) + lb = None if lb == -float("inf") else lb + ub = None if ub == float("inf") else ub + defining_var.setlb(lb) + defining_var.setub(ub) + + add_ub_con = False + add_lb_con = False + elif len(repn.nonlinear_vars) == 0 and len(repn.linear_vars) == 0: + # We are eliminating a fixed variable if ( - len(repn.nonlinear_vars) == 0 # Expression is affine - and len(repn.linear_vars) == 1 # and only contains one variable. - # ^ Can linear_vars contain duplicates? + (var.ub is not None and repn.constant > var.ub) + or (var.lb is not None and repn.constant < var.lb) ): - offset = repn.constant - coef = repn.linear_coefs[0] - defining_var = repn.linear_vars[0] - - # Bounds implied by bounds on defined variable - lb = None if var.lb is None else (var.lb - offset) / coef - ub = None if var.ub is None else (var.ub - offset) / coef - if coef < 0.0: - lb, ub = ub, lb - lbkey = lambda b: -float("inf") if b is None else b - ubkey = lambda b: float("inf") if b is None else b - # Take the more restrictive of the two bounds - lb = max(lb, defining_var.lb, key=lbkey) - ub = min(ub, defining_var.ub, key=ubkey) - lb = None if lb == -float("inf") else lb - ub = None if ub == float("inf") else ub - defining_var.setlb(lb) - defining_var.setub(ub) - - add_ub_con = False - add_lb_con = False - elif len(repn.nonlinear_vars) == 0 and len(repn.linear_vars) == 0: - # We are eliminating a fixed variable - if ( - (var.ub is not None and repn.constant > var.ub) - or (var.lb is not None and repn.constant < var.lb) - ): - raise ValueError("Attempting to fix a variable outside its bounds") - add_ub_con = False - add_lb_con = False + raise ValueError("Attempting to fix a variable outside its bounds") + add_ub_con = False + add_lb_con = False ub_expr = var_expr <= var.ub if add_ub_con else None lb_expr = var_expr >= var.lb if add_lb_con else None