diff --git a/README.md b/README.md index e617b61f8..a4459da4d 100644 --- a/README.md +++ b/README.md @@ -18,4 +18,6 @@ SWITCH. This will build HTML documentation files from python doc strings which will include descriptions of each module, their intentions, model components they define, and what input files they expect. -To learn about **numerical solvers** read [`docs/Numerical Solvers.md`](/docs/Numerical%20Solvers.md) \ No newline at end of file +To learn about **numerical solvers** read [`docs/Numerical Solvers.md`](/docs/Numerical%20Solvers.md) + +To learn about **numerical issues** read [`docs/Numerical Issues.md`](/docs/Numerical%20Issues.md) diff --git a/docs/Numerical Issues.md b/docs/Numerical Issues.md new file mode 100644 index 000000000..6844366d7 --- /dev/null +++ b/docs/Numerical Issues.md @@ -0,0 +1,235 @@ +## Numerical Issues while Using Solvers + +### What are numerical issues and why do they occur? + +All computers store numbers in binary, that is, all numbers are represented as a finite sequence of 0s and 1s. +Therefore, not all numbers can be perfectly stored. For example, the fraction +`1/3`, when stored on a computer, might actually be stored as `0.33333334...` since representing an infinite number of 3s +is not possible with only a finite number of 0s and 1s. + +When solving linear programs, these small errors can accumulate and cause significant deviations from the +'true' result. When this occurs, we say that we've encountered *numerical issues*. + +Numerical issues most often arise when our linear program contains numerical values of very small or very large +magnitudes (e.g. 10-10 or 1010). This is because—due to how numbers are stored in binary—very large or +very small values are stored less accurately (and therefore with a greater error). + +For more details on why numerical issues occur, the curious can read +about [floating-point arithmetic](https://en.wikipedia.org/wiki/Floating-point_arithmetic) +(how arithmetic is done on computers) and the [IEEE 754 standard](https://en.wikipedia.org/wiki/IEEE_754) +(the standard used by almost all computers today). For a deep dive into the topic, +read [What Every Computer Scientist Should Know About Floating-Point Arithmetic](https://www.itu.dk/~sestoft/bachelor/IEEE754_article.pdf). + +### How to detect numerical issues in Gurobi? + +Most solvers, including Gurobi, will have tests and thresholds to detect when the numerical errors become +significant. If the solver detects that we've exceeded its thresholds, it will display warnings or errors. Based +on [Gurobi's documentation](https://www.gurobi.com/documentation/9.1/refman/does_my_model_have_numeric.html), here are +some warnings or errors that Gurobi might display. + +``` +Warning: Model contains large matrix coefficient range +Consider reformulating model or setting NumericFocus parameter +to avoid numerical issues. +Warning: Markowitz tolerance tightened to 0.5 +Warning: switch to quad precision +Numeric error +Numerical trouble encountered +Restart crossover... +Sub-optimal termination +Warning: ... variables dropped from basis +Warning: unscaled primal violation = ... and residual = ... +Warning: unscaled dual violation = ... and residual = ... +``` + +Many of these warnings indicate that Gurobi is trying to workaround the numerical issues. The following list gives +examples of what Gurobi does internally to workaround numerical issues. + +- If Gurobi's normal barrier method fails due to numerical issues, Gurobi will switch to the slower but more reliable + dual simplex method (often indicated by the message: `Numerical trouble encountered`). + + +- If Gurobi's dual simplex method encounters numerical issues, Gurobi might switch to quadruple precision + (indicated by the warning: `Warning: switch to quad precision`). This is 20 to + 80 times slower (on my computer) but will represent numbers using 128-bits instead of the normal 64-bits, allowing + much greater precision and avoiding numerical issues. + +Sometimes Gurobi's internal mechanisms to avoid numerical issues are insufficient or not desirable +(e.g. too slow). In this case, we need to resolve the numerical issues ourselves. One way to do this is by scaling our +model. + +### Scaling our model to resolve numerical issues + +#### Introduction, an example of scaling + +As mentioned, numerical issues occur when our linear program contains numerical values of very small or very large +magnitude (e.g. 10-10 or 1010). We can get rid of these very large or small values by scaling our model. Consider +the following example of a linear program. + +``` +Maximize +3E17 * x + 2E10 * y +With constraints +500 * x + 1E-5 * y < 1E-5 +``` + +This program contains many large and small coefficients that we wish to get rid of. If we multiply our objective +function by 10-10, and the constraint by 105 we get: + +``` +Maximize +3E7 * x + 2 * y +With constraints +5E7 * x + y < 0 +``` + +Then if we define a new variable `x'` as 107 times the value of `x` we get: + +``` +Maximize +3 * x' + 2 * y +With constraints +5 * x' + y < 0 +``` + +This last model can be solved without numerical issues since the coefficients are neither too +small or too large. Once we solve the model, +all that's left to do is dividing `x'` by 107 to get `x`. + +This example, although not very realistic, gives an example +of what it means to scale a model. Scaling is often the best solution to resolving numerical issues +and can be easily done with Pyomo (see below). In some cases scaling is insufficient and other +changes need to be made, this is explained in the section *Other techniques to resolve numerical issues*. + +#### Gurobi Guidelines for ranges of values + +An obvious question is, what is considered too small or too large? +Gurobi provides some guidelines on what is a reasonable range +for numerical values ([here](https://www.gurobi.com/documentation/9.1/refman/recommended_ranges_for_var.html) and [here](https://www.gurobi.com/documentation/9.1/refman/advanced_user_scaling.html)). +Here's a summary: + +- Right-hand sides of inequalities and variable domains (i.e. variable bounds) should + be on the order of 104 or less. + +- The objective function's optimal value (i.e. the solution) should be between 1 and 105. + +- The matrix coefficients should span a range of six orders of magnitude + or less and ideally between 10-3 and 106. + + + + +#### Scaling our model with Pyomo + +Scaling an objective function or constraint is easy. +Simply multiply the expression by the scaling factor. For example, + +```python +# Objective function without scaling +model.TotalCost = Objective(..., rule=lambda m, ...: some_expression) +# Objective function ith scaling +scaling_factor = 1e-7 +model.TotalCost = Objective(..., rule=lambda m, ...: (some_expression) * scaling_factor) + +# Constraint without scaling +model.SomeConstraint = Constraint(..., rule=lambda m, ...: left_hand_side < right_hand_side) +# Constraint with scaling +scaling_factor = 1e-2 +model.SomeConstraint = Constraint( + ..., + rule=lambda m, ...: (left_hand_side) * scaling_factor < (right_hand_side) * scaling_factor +) +``` + +Scaling a variable is more of a challenge since the variable +might be used in multiple places, and we don't want to need +to change our code in multiple places. The trick is to define an expression that equals our unscaled variable. +We can use this expression throughout our model, and Pyomo will still use the underlying +scaled variable when solving. Here's what I mean. + +```python +from pyomo.environ import Var, Expression +... +scaling_factor = 1e5 +# Define the variable +model.ScaledVariable = Var(...) +# Define an expression that equals the variable but unscaled. This is what we use elsewhere in our model. +model.UnscaledExpression = Expression(..., rule=lambda m, *args: m.ScaledVariable[args] / scaling_factor) +... +``` + +Thankfully, I've written the `ScaledVariable` class that will do this for us. +When we add a `ScaledVariable` to our model, the actual scaled +variable is created behind the scenes and what's returned is the unscaled expression that +we can use elsewhere in our code. Here's how to use `ScaledVariable` in practice. + + +```python +# Without scaling +from pyomo.environ import Var +model.SomeVariable = Var(...) + +# With scaling +from switch_model.utilities.scaling import ScaledVariable +model.SomeVariable = ScaledVariable(..., scaling_factor=5) +``` + +Here, we can use `SomeVariable` throughout our code just as before. +Internally, Pyomo (and the solver) will be using a scaled version of `SomeVariable`. +In this case the solver will use a variable that is 5 times greater +than the value we reference in our code. This means the +coefficients in front of `SomeVariable` will be 1/5th of the usual. + +#### How much to scale by ? + +Earlier we shared Gurobi's guidelines on the ideal range for our matrix coefficients, +right-hand side values, variable bounds and objective function. Yet how do +we know where our values currently lie to determine how much to scale them by? + +For large models, determining which variables to scale and by how much can be a challenging task. + +First, when solving with the flag `--stream-solver` and `-v`, +Gurobi will print to console helpful information for preliminary analysis. +Here's an example of what Gurobi's output might look like. It might also +include warnings about ranges that are too wide. + +``` +Coefficient statistics: + Matrix range [2e-05, 1e+06] + Objective range [2e-05, 6e+04] + Bounds range [3e-04, 4e+04] + RHS range [1e-16, 3e+05] +``` + +Second, if we solved our model with the flags `--keepfiles`, `--tempdir` and `--symbolic-solver-labels`, then +we can read the `.lp` file in the temporary folder that contains the entire model including the coefficients. +This is the file Gurobi reads before solving and has all the information we might need. +However, this file is very hard to analyze manually due its size. + +Third, I'm working on a tool to automatically analyze the `.lp` file and return information +useful that would help resolve numerical issues. The tool is available [here](https://github.com/staadecker/lp-analyzer). + + +### Other techniques to resolve numerical issues + +In some cases, scaling might not be sufficient to resolve numerical issues. +For example, if variables within the same set have values that span too wide of a range, +scaling will not be able to reduce the range since scaling affects all variables +within a set equally. + +Other than scaling, some techniques to resolve numerical issues are: + +- Reformulating the model + +- Avoiding unnecessarily large penalty terms + +- Changing the solver's method + +- Loosening tolerances (at the risk of getting less accurate, or inaccurate results) + +One can learn more about these methods +by reading [Gurobi's guidelines on Numerical Issues](https://www.gurobi.com/documentation/9.1/refman/guidelines_for_numerical_i.html) +or a [shorter PDF](http://files.gurobi.com/Numerics.pdf) that Gurobi has released. + + + diff --git a/docs/Numerical Solvers.md b/docs/Numerical Solvers.md index 324259f44..346de74bb 100644 --- a/docs/Numerical Solvers.md +++ b/docs/Numerical Solvers.md @@ -1,250 +1,55 @@ -# A Guide on Numerical Solvers +# Numerical Solvers by Martin Staadecker -## Content +## Introduction -Numerical solvers, such as Gurobi, are tools used to solve [linear programs](https://en.wikipedia.org/wiki/Linear_programming). +Numerical solvers, such as [Gurobi](https://gurobi.com), are tools used to solve [linear programs](https://en.wikipedia.org/wiki/Linear_programming). This document is a record of what I've learnt about using numerical solvers -while working with Switch, Pyomo and Gurobi. It includes conceptual explanations, -solutions to problems one might encounter, and techniques that I've found to be useful. +while working with Switch, Pyomo and Gurobi. -Currently, the topics covered are: +## Gurobi resources -- Numerical issues while using solvers +The best resource when working with Gurobi is [Gurobi's reference manual](https://www.gurobi.com/documentation/9.1/refman/index.html). -## Numerical Issues while Using Solvers - -### What are numerical issues and why do they occur? +A few especially useful pages are: -All computers store numbers in binary, that is, all numbers are represented as a finite sequence of 0s and 1s. -Therefore, not all numbers can be perfectly stored. For example, the fraction -`1/3`, when stored on a computer, might actually be stored as `0.33333334...` since representing an infinite number of 3s -is not possible with only a finite number of 0s and 1s. - -When solving linear programs, these small errors can accumulate and cause significant deviations from the -'true' result. When this occurs, we say that we've encountered *numerical issues*. +- Gurobi's [parameter guidelines for continuous models](https://www.gurobi.com/documentation/9.1/refman/continuous_models.html) -Numerical issues most often arise when our linear program contains numerical values of very small or very large -magnitudes (e.g. 10-10 or 1010). This is because—due to how numbers are stored in binary—very large or -very small values are stored less accurately (and therefore with a greater error). +- Gurobi's [parameter list](https://www.gurobi.com/documentation/9.1/refman/parameters.html#sec:Parameters) especially the [`Method` parameter](https://www.gurobi.com/documentation/9.1/refman/method.html). -For more details on why numerical issues occur, the curious can read -about [floating-point arithmetic](https://en.wikipedia.org/wiki/Floating-point_arithmetic) -(how arithmetic is done on computers) and the [IEEE 754 standard](https://en.wikipedia.org/wiki/IEEE_754) -(the standard used by almost all computers today). For a deep dive into the topic, -read [What Every Computer Scientist Should Know About Floating-Point Arithmetic](https://www.itu.dk/~sestoft/bachelor/IEEE754_article.pdf). +- Gurobi's [guidelines for numerical issues](https://www.gurobi.com/documentation/9.1/refman/guidelines_for_numerical_i.html) (see `docs/Numerical Issues.md`). -### How to detect numerical issues in Gurobi? +## Specifying parameters -Most solvers, including Gurobi, will have tests and thresholds to detect when the numerical errors become -significant. If the solver detects that we've exceeded its thresholds, it will display warnings or errors. Based -on [Gurobi's documentation](https://www.gurobi.com/documentation/9.1/refman/does_my_model_have_numeric.html), here are -some warnings or errors that Gurobi might display. +To specify a Gurobi parameter use the following format: -``` -Warning: Model contains large matrix coefficient range -Consider reformulating model or setting NumericFocus parameter -to avoid numerical issues. -Warning: Markowitz tolerance tightened to 0.5 -Warning: switch to quad precision -Numeric error -Numerical trouble encountered -Restart crossover... -Sub-optimal termination -Warning: ... variables dropped from basis -Warning: unscaled primal violation = ... and residual = ... -Warning: unscaled dual violation = ... and residual = ... -``` +`switch solve --solver gurobi --solver-options-string "Parameter1=Value Parameter2=Value"`. -Many of these warnings indicate that Gurobi is trying to workaround the numerical issues. The following list gives -examples of what Gurobi does internally to workaround numerical issues. +We recommend always using `"method=2 BarHomogeneous=1 FeasibilityTol=1e-5 crossover=0"` +as your base parameters (this is what `switch solve --recommended` does). -- If Gurobi's normal barrier method fails due to numerical issues, Gurobi will switch to the slower but more reliable - dual simplex method (often indicated by the message: `Numerical trouble encountered`). +## Solving Methods +There are two methods that exist when solving a linear program. +The first is the Simplex Method and the second is the Barrier +solve method also known as interior-point method (IPM). -- If Gurobi's dual simplex method encounters numerical issues, Gurobi might switch to quadruple precision - (indicated by the warning: `Warning: switch to quad precision`). This is 20 to - 80 times slower (on my computer) but will represent numbers using 128-bits instead of the normal 64-bits, allowing - much greater precision and avoiding numerical issues. - -Sometimes Gurobi's internal mechanisms to avoid numerical issues are insufficient or not desirable -(e.g. too slow). In this case, we need to resolve the numerical issues ourselves. One way to do this is by scaling our -model. - -### Scaling our model to resolve numerical issues - -#### Introduction, an example of scaling - -As mentioned, numerical issues occur when our linear program contains numerical values of very small or very large -magnitude (e.g. 10-10 or 1010). We can get rid of these very large or small values by scaling our model. Consider -the following example of a linear program. - -``` -Maximize -3E17 * x + 2E10 * y -With constraints -500 * x + 1E-5 * y < 1E-5 -``` - -This program contains many large and small coefficients that we wish to get rid of. If we multiply our objective -function by 10-10, and the constraint by 105 we get: - -``` -Maximize -3E7 * x + 2 * y -With constraints -5E7 * x + y < 0 -``` - -Then if we define a new variable `x'` as 107 times the value of `x` we get: - -``` -Maximize -3 * x' + 2 * y -With constraints -5 * x' + y < 0 -``` - -This last model can be solved without numerical issues since the coefficients are neither too -small or too large. Once we solve the model, -all that's left to do is dividing `x'` by 107 to get `x`. - -This example, although not very realistic, gives an example -of what it means to scale a model. Scaling is often the best solution to resolving numerical issues -and can be easily done with Pyomo (see below). In some cases scaling is insufficient and other -changes need to be made, this is explained in the section *Other techniques to resolve numerical issues*. - -#### Gurobi Guidelines for ranges of values - -An obvious question is, what is considered too small or too large? -Gurobi provides some guidelines on what is a reasonable range -for numerical values ([here](https://www.gurobi.com/documentation/9.1/refman/recommended_ranges_for_var.html) and [here](https://www.gurobi.com/documentation/9.1/refman/advanced_user_scaling.html)). -Here's a summary: - -- Right-hand sides of inequalities and variable domains (i.e. variable bounds) should - be on the order of 104 or less. - -- The objective function's optimal value (i.e. the solution) should be between 1 and 105. - -- The matrix coefficients should span a range of six orders of magnitude - or less and ideally between 10-3 and 106. - - - - -#### Scaling our model with Pyomo - -Scaling an objective function or constraint is easy. -Simply multiply the expression by the scaling factor. For example, - -```python -# Objective function without scaling -model.TotalCost = Objective(..., rule=lambda m, ...: some_expression) -# Objective function ith scaling -scaling_factor = 1e-7 -model.TotalCost = Objective(..., rule=lambda m, ...: (some_expression) * scaling_factor) - -# Constraint without scaling -model.SomeConstraint = Constraint(..., rule=lambda m, ...: left_hand_side < right_hand_side) -# Constraint with scaling -scaling_factor = 1e-2 -model.SomeConstraint = Constraint( - ..., - rule=lambda m, ...: (left_hand_side) * scaling_factor < (right_hand_side) * scaling_factor -) -``` - -Scaling a variable is more of a challenge since the variable -might be used in multiple places, and we don't want to need -to change our code in multiple places. The trick is to define an expression that equals our unscaled variable. -We can use this expression throughout our model, and Pyomo will still use the underlying -scaled variable when solving. Here's what I mean. - -```python -from pyomo.environ import Var, Expression -... -scaling_factor = 1e5 -# Define the variable -model.ScaledVariable = Var(...) -# Define an expression that equals the variable but unscaled. This is what we use elsewhere in our model. -model.UnscaledExpression = Expression(..., rule=lambda m, *args: m.ScaledVariable[args] / scaling_factor) -... -``` - -Thankfully, I've written the `ScaledVariable` class that will do this for us. -When we add a `ScaledVariable` to our model, the actual scaled -variable is created behind the scenes and what's returned is the unscaled expression that -we can use elsewhere in our code. Here's how to use `ScaledVariable` in practice. - - -```python -# Without scaling -from pyomo.environ import Var -model.SomeVariable = Var(...) - -# With scaling -from switch_model.utilities.scaling import ScaledVariable -model.SomeVariable = ScaledVariable(..., scaling_factor=5) -``` - -Here, we can use `SomeVariable` throughout our code just as before. -Internally, Pyomo (and the solver) will be using a scaled version of `SomeVariable`. -In this case the solver will use a variable that is 5 times greater -than the value we reference in our code. This means the -coefficients in front of `SomeVariable` will be 1/5th of the usual. - -#### How much to scale by ? - -Earlier we shared Gurobi's guidelines on the ideal range for our matrix coefficients, -right-hand side values, variable bounds and objective function. Yet how do -we know where our values currently lie to determine how much to scale them by? - -For large models, determining which variables to scale and by how much can be a challenging task. - -First, when solving with the flag `--stream-solver` and `-v`, -Gurobi will print to console helpful information for preliminary analysis. -Here's an example of what Gurobi's output might look like. It might also -include warnings about ranges that are too wide. - -``` -Coefficient statistics: - Matrix range [2e-05, 1e+06] - Objective range [2e-05, 6e+04] - Bounds range [3e-04, 4e+04] - RHS range [1e-16, 3e+05] -``` - -Second, if we solved our model with the flags `--keepfiles`, `--tempdir` and `--symbolic-solver-labels`, then -we can read the `.lp` file in the temporary folder that contains the entire model including the coefficients. -This is the file Gurobi reads before solving and has all the information we might need. -However, this file is very hard to analyze manually due its size. - -Third, I'm working on a tool to automatically analyze the `.lp` file and return information -useful that would help resolve numerical issues. The tool is available [here](https://github.com/staadecker/lp-analyzer). - - -### Other techniques to resolve numerical issues - -In some cases, scaling might not be sufficient to resolve numerical issues. -For example, if variables within the same set have values that span too wide of a range, -scaling will not be able to reduce the range since scaling affects all variables -within a set equally. - -Other than scaling, some techniques to resolve numerical issues are: - -- Reformulating the model +- The Simplex Method is more robust than IPM (it can find + an optimal solution even with numerically challenging problems). -- Avoiding unnecessarily large penalty terms - -- Changing the solver's method - -- Loosening tolerances (at the risk of getting less accurate, or inaccurate results) - -One can learn more about these methods -by reading [Gurobi's guidelines on Numerical Issues](https://www.gurobi.com/documentation/9.1/refman/guidelines_for_numerical_i.html) -or a [shorter PDF](http://files.gurobi.com/Numerics.pdf) that Gurobi has released. - +- The Simplex Method uses only 1 thread while the Barrier method can +leverage multiple threads. + +- The Barrier method is significantly faster for our model sizes. +- Running `switch solve --recommended` selects the Barrier method. +- By default, when the Barrier method finishes it converts its solution +to a simplex solution in what is called the crossover phase (see [details](https://www.gurobi.com/documentation/9.1/refman/barrier_logging.html)). + This crossover phase takes the most time and is not necessary. Therefore is gets + disabled by the `--recommended` flag. + +- The crossover phase *is* important if the barrier method produces a sub-optimal solution. + In this case use `--recommended-robust` to enable the crossover. + diff --git a/switch_model/reporting/__init__.py b/switch_model/reporting/__init__.py index eecc9603a..56cac92c0 100644 --- a/switch_model/reporting/__init__.py +++ b/switch_model/reporting/__init__.py @@ -52,13 +52,13 @@ def define_arguments(argparser): help="List of expressions to save in addition to variables; can also be 'all' or 'none'." ) -def get_cell_formatter(sig_digits): +def get_cell_formatter(sig_digits, zero_cutoff): sig_digits_formatter = "{0:." + str(sig_digits) + "g}" def format_cell(c): if not isinstance(c, float): return c - if abs(c) < 1e-8: + if abs(c) < zero_cutoff: return 0 else: return sig_digits_formatter.format(c) @@ -75,7 +75,7 @@ def write_table(instance, *indexes, output_file=None, **kwargs): # don't know what that is. if output_file is None: raise Exception("Must specify output_file in write_table()") - cell_formatter = get_cell_formatter(instance.options.sig_figs_output) + cell_formatter = get_cell_formatter(instance.options.sig_figs_output, instance.options.zero_cutoff_output) if 'df' in kwargs: df = kwargs.pop('df') @@ -133,7 +133,7 @@ def post_solve(instance, outdir): def save_generic_results(instance, outdir, sorted_output): - cell_formatter = get_cell_formatter(instance.options.sig_figs_output) + cell_formatter = get_cell_formatter(instance.options.sig_figs_output, instance.options.zero_cutoff_output) components = list(instance.component_objects(Var)) # add Expression objects that should be saved, if any diff --git a/switch_model/solve.py b/switch_model/solve.py index e2bc2004c..3abc7fb88 100755 --- a/switch_model/solve.py +++ b/switch_model/solve.py @@ -4,7 +4,7 @@ from __future__ import print_function from pyomo.environ import * -from pyomo.opt import SolverFactory, SolverStatus, TerminationCondition +from pyomo.opt import SolverFactory, SolverStatus, TerminationCondition, SolutionStatus import pyomo.version import pandas as pd @@ -580,6 +580,10 @@ def define_arguments(argparser): "--sig-figs-output", default=5, type=int, help='The number of significant digits to include in the output by default' ) + argparser.add_argument( + "--zero-cutoff-output", default=1e-5, type=float, + help="If the magnitude of an output value is less than this value, it is rounded to 0." + ) argparser.add_argument( "--sorted-output", default=False, action='store_true', @@ -609,7 +613,15 @@ def add_recommended_args(argparser): """ argparser.add_argument( "--recommended", default=False, action='store_true', - help='Equivalent to running with all of the following options: --solver gurobi -v --sorted-output --stream-output --log-run --solver-io python --solver-options-string "method=2 BarHomogeneous=1 FeasibilityTol=1e-5"' + help='Includes several flags that are recommended including --solver gurobi --verbose --stream-output and more. ' + 'See parse_recommended_args() in solve.py for the full list of recommended flags.' + ) + + argparser.add_argument( + "--recommended-robust", default=False, action='store_true', + help='Equivalent to --recommended however enables crossover during solving. Crossover is useful' + ' if the solution found by the barrier method is suboptimal. If you find that the solver returns' + ' a suboptimal solution use this flag.' ) argparser.add_argument( @@ -623,15 +635,13 @@ def parse_recommended_args(args): add_recommended_args(argparser) options = argparser.parse_known_args(args)[0] - if options.recommended and options.recommended_debug: - raise Exception("Can't use both --recommended and --recommended-debug") - if not options.recommended and not options.recommended_debug: + flags_used = options.recommended + options.recommended_robust + options.recommended_debug + if flags_used > 1: + raise Exception("Must pick between --recommended-debug, --recommended-fast or --recommended.") + if flags_used == 0: return args - # If you change the flags, make sure to update the help text in - # the function above. - # Note we don't append but rather prepend so that flags can override the --recommend or - # --recommend-debug flags. + # Note we don't append but rather prepend so that flags can override the --recommend flags. args = [ '--solver', 'gurobi', '-v', @@ -640,12 +650,11 @@ def parse_recommended_args(args): '--log-run', '--debug', '--graph', - '--solver-options-string', - 'method=2 BarHomogeneous=1 FeasibilityTol=1e-5' ] + args - # We use to use solver-io=python however the seperation of processes is useful in helping the OS - # if options.recommended: - # args = ['--solver-io', 'python'] + args + solver_options_string = "BarHomogeneous=1 FeasibilityTol=1e-5" + if not options.recommended_robust: + solver_options_string += " crossover=0 method=2" + args = ['--solver-options-string', solver_options_string] + args if options.recommended_debug: args = ['--keepfiles', '--tempdir', 'temp', '--symbolic-solver-labels'] + args @@ -816,42 +825,30 @@ def solve(model): "https://docs.python.org/3/library/pdb.html for instructions on using pdb.") breakpoint() - # Treat infeasibility as an error, rather than trying to load and save the results - # (note: in this case, results.solver.status may be SolverStatus.warning instead of - # SolverStatus.error) - if (results.solver.termination_condition == TerminationCondition.infeasible): - if hasattr(model, "iis"): - print("Model was infeasible; irreducibly inconsistent set (IIS) returned by solver:") - print("\n".join(sorted(c.name for c in model.iis))) - else: - print("Model was infeasible; if the solver can generate an irreducibly inconsistent set (IIS),") - print("more information may be available by setting the appropriate flags in the ") - print('solver_options_string and calling this script with "--suffixes iis" or "--gurobi-find-iis".') - raise RuntimeError("Infeasible model") - - # Report any warnings; these are written to stderr so users can find them in - # error logs (e.g. on HPC systems). These can occur, e.g., if solver reaches - # time limit or iteration limit but still returns a valid solution - if results.solver.status == SolverStatus.warning: - warn( - "Solver terminated with warning.\n" - + " Solution Status: {}\n".format(model.solutions[-1].status) - + " Termination Condition: {}".format(results.solver.termination_condition) - ) + solver_status = results.solver.status + solver_message = results.solver.message + termination_condition = results.solver.termination_condition + solution_status = model.solutions[-1].status if len(model.solutions) != 0 else None - if(results.solver.status != SolverStatus.ok or - results.solver.termination_condition != TerminationCondition.optimal): + if solver_status != SolverStatus.ok or termination_condition != TerminationCondition.optimal: warn( - f"Solver terminated with status '{results.solver.status}' and termination condition" - f" {results.solver.termination_condition}" + f"Solver termination status is not 'ok' or not 'optimal':\n" + f"\t- Termination condition: {termination_condition}\n" + f"\t- Solver status: {solver_status}\n" + f"\t- Solver message: {solver_message}\n" + f"\t- Solution status: {solution_status}" ) + if solution_status == SolutionStatus.feasible and solver_status == SolverStatus.warning: + print("If you used --recommended, you might want to try --recommended-robust.") + + if query_yes_no("Do you want to abort and exit?", default=None): + raise SystemExit() + if model.options.verbose: - print("") - print("Optimization termination condition was {}.".format( - results.solver.termination_condition)) - if str(results.solver.message) != '': - print('Solver message: {}'.format(results.solver.message)) + print(f"\nOptimization termination condition was {termination_condition}.") + if str(solver_message) != '': + print(f'Solver message: {solver_message}') print("") if model.options.save_solution: