Merge remote-tracking branch 'rael/wecc' into plots

Author: staadecker

docs/Performance.md

@@ -0,0 +1,74 @@
+# Performance
+
+Memory use and solve time are two important factors that we try to keep to a minimum in our models. There are multiple
+things one can do to improve performance.
+
+## Solving methods
+
+By far the biggest factor that impacts performance is the method used by Gurobi. The fastest method is barrier solve
+without crossover (use `--recommended-fast`)
+however this method often returns a suboptimal solution. The next fastest is barrier solve followed by crossover and
+simplex (use `--recommended`) which almost always works. In some cases barrier solve encounters numerical issues (
+see [`Numerical Issues.md`](./Numerical%20Issues.md))
+in which case the slower Simplex method must be used (`--solver-options-string method=1`).
+
+## Solver interface
+
+Solver interfaces are how Pyomo communicates with Gurobi (or any solver).
+
+There are two solver interfaces that you should know about: `gurobi` and `gurobi_direct`.
+
+- When using `gurobi`, Pyomo will write the entire model to a temporary text file and then start a *separate Gurobi
+  process* that will read the file, solve the model and write the results to another temporary text file. Once Gurobi
+  finishes writing the results Pyomo will read the results text file and load the results back into the Python program
+  before running post_solve (e.g. generate csv files, create graphs, etc). Note that these temporary text files are
+  stored in `/tmp` but if you use `--recommended-debug` Pyomo and Gurobi will instead use a `temp` folder in your model.
+
+- `gurobi_direct` uses Gurobi's Python library to create and solve the model directly in Python without the use of
+  intermediate text files.
+
+In theory `gurobi_direct` should be faster and more efficient however in practice we find that that's not the case. As
+such we recommend using `gurobi` and all our defaults do so. If someone has the time they could profile `gurobi_direct`
+to improve performance at which point we could make `gurobi_direct` the default (and enable `--save-warm-start` by default, see below).
+
+The `gurobi` interface has the added advantage of separating Gurobi and Pyomo into separate threads. This means that
+while Gurobi is solving and Pyomo is idle, the operating system can automatically move Pyomo's memory usage
+to [virtual memory](https://serverfault.com/questions/48486/what-is-swap-memory)
+which will free up more memory for Gurobi.
+
+## Warm starting
+
+Warm starting is the act of using a solution from a previous similar model to start the solver closer to your expected
+solution. Theoretically this can help performance however in practice there are several limitations. For this section, *
+previous solution* refers to the results from an already solved model that you are using to warm start the solver. *
+Current solution* refers to the solution you are trying to find while using the warm start feature.
+
+- To warm start a model use `switch solve --warm-start <path_to_previous_solution>`.
+
+- Warm starting only works if the previous solution does not break any constraints of the current solution. This usually
+  only happens if a) the model has the exact same set of variables b)
+  the previous solution was "harder" (e.g. it had more constraints to satisfy).
+
+- Warm starting always uses the slower Simplex method. This means unless you expect the previous solution and current
+  solution to be very similar, it may be faster to solve without warm start using the barrier method.
+
+- If your previous solution didn't use crossover (e.g. you used `--recommended-fast`) then warm starting will be even
+  slower since the solver will need to first run crossover before warm starting.
+
+- Our implementation of warm starting only works if your previous solution has an `outputs/warm_start.pickle`
+  file. This file is only generated when you use `--save-warm-start`.
+
+- `--save-warm-start` and `--warm-start` both use an extension of the `gurobi_direct` solver interface which is
+  generally slower than the `gurobi` solver interface (see section above).
+  
+## Tools for improving performance
+
+- [Memory profiler](https://pypi.org/project/memory-profiler/) for generating plots of the memory
+use over time. Use `mprof run --interval 60 --multiprocess switch solve ...` and once solving is done
+  run `mprof plot -o profile.png` to make the plot.
+  
+- [Fil Profiler](https://pypi.org/project/filprofiler/) is an amazing tool for seeing which parts of the code are
+using up memory during peak memory usage.
+  
+- Using `switch_model.utilities.StepTimer` to measure how long certain code blocks take to run. See examples
+throughout the code.

switch_model/solve.py

@@ -35,6 +35,7 @@
 from switch_model.tools.graph.cli_graph import main as graph_main
 from switch_model.utilities.patches import patch_pyomo
 from switch_model.utilities.results_info import save_info, add_info, ResultsInfoSection
+import switch_model.utilities.gurobi_aug  # We keep this line here to ensure that 'gurobi_aug' gets registered as a solver
 
 
 def main(
@@ -155,9 +156,6 @@ def debug(type, value, tb):
         # create an instance (also reports time spent reading data and loading into model)
         instance = model.load_inputs(attach_data_portal=attach_data_portal)
 
-        if instance.options.warm_start:
-            warm_start(instance)
-
         #### Below here, we refer to instance instead of model ####
 
         instance.pre_solve()
@@ -190,6 +188,13 @@ def debug(type, value, tb):
         # We no longer need model (only using instance) so we can garbage collect it to optimize our memory usage
         del model
 
+        if instance.options.warm_start_mip:
+            if instance.options.verbose:
+                timer.step_time()
+            warm_start_mip(instance)
+            if instance.options.verbose:
+                print(f"Loaded warm start inputs in {timer.step_time_as_str()}.")
+
         if instance.options.reload_prior_solution:
             print("Loading prior solution...")
             reload_prior_solution_from_pickle(instance, instance.options.outputs_dir)
@@ -276,14 +281,14 @@ def debug(type, value, tb):
         )
 
 
-def warm_start(instance):
+def warm_start_mip(instance):
     """
-    This function loads in the variables from a previous run
-    and starts out our model at these variables to make it reach
-    a solution faster.
+    This function loads the results from a previous run into the Pyomo variables.
+    This allows Gurobi's Mixed Integer Programming algorithm to "warm start" (start closer to the solution).
+    Warm starting only works in Gurobi if the initial values don't violate any constraints
+    (i.e. valid but not optimal solution).
     """
-    warm_start_timer = StepTimer()
-    warm_start_dir = os.path.join(instance.options.warm_start, "outputs")
+    warm_start_dir = os.path.join(instance.options.warm_start_mip, "outputs")
     if not os.path.isdir(warm_start_dir):
         warnings.warn(
             f"Path {warm_start_dir} does not exist and cannot be used to warm start solver. Warm start skipped."
@@ -310,8 +315,6 @@ def warm_start(instance):
                 # If the index isn't valid that's ok, just don't warm start that variable
                 pass
 
-    print(f"Loaded warm start inputs in {warm_start_timer.step_time_as_str()}.")
-
 
 def reload_prior_solution_from_pickle(instance, outdir):
     with open(os.path.join(outdir, "results.pickle"), "rb") as fh:
@@ -524,10 +527,16 @@ def define_arguments(argparser):
     # whether that does the same thing as --solver-options-string so we don't reuse the same name.
     argparser.add_argument(
         "--solver-options-string",
-        default=None,
+        default="",
         help="A quoted string of options to pass to the model solver. Each option must be of the form option=value. "
         "(e.g., --solver-options-string \"mipgap=0.001 primalopt='' advance=2 threads=1\")",
     )
+    argparser.add_argument(
+        "--solver-method",
+        default=None,
+        type=int,
+        help="Specify the solver method to use.",
+    )
     argparser.add_argument(
         "--keepfiles",
         action="store_true",
@@ -639,6 +648,12 @@ def define_arguments(argparser):
         action="store_true",
         help="Save the solution to a pickle file after model is solved to allow for later inspection via --reload-prior-solution.",
     )
+    argparser.add_argument(
+        "--save-warm-start",
+        default=False,
+        action="store_true",
+        help="Save warm_start.pickle to the outputs which allows future runs to warm start from this one.",
+    )
     argparser.add_argument(
         "--interact",
         default=False,
@@ -685,10 +700,27 @@ def define_arguments(argparser):
         help="Number of threads to be used while solving. Currently only supported for Gurobi",
     )
 
+    argparser.add_argument(
+        "--warm-start-mip",
+        default=None,
+        help="Enables warm start for a Mixed Integer problem by specifying the "
+        "path to a previous scenario. Warm starting only works if the solution to the previous solution"
+        "is also a feasible (but not necessarily optimal) solution to the current scenario.",
+    )
+
     argparser.add_argument(
         "--warm-start",
         default=None,
-        help="Path to folder of directory to use for warm start",
+        help="Enables warm start for a LP Problem by specifying the path to the previous scenario. Note"
+        " that all variables must be the same between the previous and current scenario.",
+    )
+
+    argparser.add_argument(
+        "--gurobi-make-mps",
+        default=False,
+        action="store_true",
+        help="Instead of solving just output a Gurobi .mps file that can be used for debugging numerical properties."
+        " See https://github.com/staadecker/lp-analyzer/ for details.",
     )
 
 
@@ -719,7 +751,7 @@ def add_recommended_args(argparser):
         "--recommended-debug",
         default=False,
         action="store_true",
-        help='Equivalent to running with all of the following options: --solver gurobi -v --sorted-output --keepfiles --tempdir temp --stream-output --symbolic-solver-labels --log-run --debug --solver-options-string "method=2 BarHomogeneous=1 FeasibilityTol=1e-5"',
+        help="Same as --recommended but adds the flags --keepfiles --tempdir temp --symbolic-solver-labels which are useful when debugging Gurobi.",
     )
 
 
@@ -748,8 +780,10 @@ def parse_recommended_args(args):
         "--log-run",
         "--debug",
         "--graph",
+        "--solver-method",
+        "2",  # Method 2 is barrier solve which is much faster than default
     ] + args
-    solver_options_string = "BarHomogeneous=1 FeasibilityTol=1e-5 method=2"
+    solver_options_string = "BarHomogeneous=1 FeasibilityTol=1e-5"
     if options.recommended_fast:
         solver_options_string += " crossover=0"
     args = ["--solver-options-string", solver_options_string] + args
@@ -851,6 +885,18 @@ def solve(model):
         solver = model.solver
         solver_manager = model.solver_manager
     else:
+        # If we need warm start switch the solver to our augmented version that supports warm starting
+        if model.options.warm_start is not None or model.options.save_warm_start:
+            if model.options.solver not in ("gurobi", "gurobi_direct"):
+                raise NotImplementedError(
+                    "Warm start functionality requires --solver gurobi"
+                )
+            model.options.solver = "gurobi_aug"
+
+        if model.options.warm_start is not None:
+            # Method 1 (dual simplex) is required since it supports warm starting.
+            model.options.solver_method = 1
+
         # Create a solver object the first time in. We don't do this until a solve is
         # requested, because sometimes a different solve function may be used,
         # with its own solver object (e.g., with runph or a parallel solver server).
@@ -860,28 +906,35 @@ def solve(model):
         # Note previously solver was saved in model however this is very memory inefficient.
         solver = SolverFactory(model.options.solver, solver_io=model.options.solver_io)
 
-        # If this option is enabled, gurobi will output an IIS to outputs\iis.ilp.
-        if model.options.gurobi_find_iis:
-            # Enable symbolic labels since otherwise we can't debug the .ilp file.
-            model.options.symbolic_solver_labels = True
+        if model.options.gurobi_find_iis and model.options.gurobi_make_mps:
+            raise Exception("Can't use --gurobi-find-iis with --gurobi-make-mps.")
 
-            # If no string is passed make the string empty so we can add to it
-            if model.options.solver_options_string is None:
-                model.options.solver_options_string = ""
+        if model.options.gurobi_find_iis or model.options.gurobi_make_mps:
+            # If we are outputting a file we want to enable symbolic labels to help debugging
+            model.options.symbolic_solver_labels = True
 
+        # If this option is enabled, gurobi will output an IIS to outputs\iis.ilp.
+        if model.options.gurobi_find_iis:
             # Add to the solver options 'ResultFile=iis.ilp'
             # https://stackoverflow.com/a/51994135/5864903
-            iis_file_path = os.path.join(model.options.outputs_dir, "iis.ilp")
-            model.options.solver_options_string += " ResultFile={}".format(
-                iis_file_path
+            model.options.solver_options_string += " ResultFile=iis.ilp"
+        if model.options.gurobi_make_mps:
+            # Output the input file and set time limit to zero to ensure it doesn't actually solve
+            model.options.solver_options_string += (
+                f" ResultFile=problem.mps TimeLimit=0"
             )
 
         if model.options.threads:
+            model.options.solver_options_string += f" Threads={model.options.threads}"
+
+        if model.options.solver_method is not None:
             # If no string is passed make the string empty so we can add to it
             if model.options.solver_options_string is None:
                 model.options.solver_options_string = ""
 
-            model.options.solver_options_string += f" Threads={model.options.threads}"
+            model.options.solver_options_string += (
+                f" method={model.options.solver_method}"
+            )
 
         solver_manager = SolverManagerFactory(model.options.solver_manager)
 
@@ -896,9 +949,17 @@ def solve(model):
         else None,
     )
 
-    if model.options.warm_start is not None:
+    if model.options.warm_start_mip is not None or model.options.warm_start is not None:
         solver_args["warmstart"] = True
 
+    if model.options.warm_start is not None:
+        solver_args["read_warm_start"] = os.path.join(
+            model.options.warm_start, "outputs", "warm_start.pickle"
+        )
+
+    if model.options.save_warm_start:
+        solver_args["write_warm_start"] = os.path.join("outputs", "warm_start.pickle")
+
     # drop all the unspecified options
     solver_args = {k: v for (k, v) in solver_args.items() if v is not None}