Merge remote-tracking branch 'upstream/main' into lshaped_fix

.github/workflows/pull_push_regression.yml

@@ -1,4 +1,4 @@
-name: mpi-sppy (Pyomo released)
+name: test_ef_ph (Pyomo released)
 
 on:
   push:

.github/workflows/schur_complement.yml

@@ -4,7 +4,9 @@ name: schur-complement
 
 on:
   push:
+    branches: [ main ]
   pull_request:
+    branches: [ main ]
 
 jobs:
   build:

.github/workflows/testmmw.yml

@@ -0,0 +1,40 @@
+# aph (pyomo released)
+
+name: mmw tests
+
+on:
+  push:
+    branches: [ main ]
+  pull_request:
+    branches: [ main ]
+
+defaults:
+    run:
+        shell: bash -l {0}
+
+jobs:
+  build:
+
+    runs-on: ubuntu-latest
+
+    steps:
+      - uses: actions/checkout@v2
+      - uses: conda-incubator/setup-miniconda@v2
+        with:
+          activate-environment: test_env
+          python-version: 3.8
+          auto-activate-base: false
+      - name: Install dependencies
+        run: |
+          conda install mpi4py numpy setuptools 
+          pip install pyomo xpress cplex scipy
+
+      - name: setup the program
+        run: |
+          python setup.py develop
+
+      - name: run tests
+        timeout-minutes: 10
+        run: |
+          cd mpisppy/tests
+          python test_mmw.py 

doc/src/extensions.rst

@@ -19,7 +19,7 @@ creation of more. Extensions can be found in ``mpisppy.extensions``.
 Multiple Extensions
 -------------------
 
-To employ multiple PH extensions, use ``mpisppy.extensions.extension import MultiPHExtension``
+To employ multiple PH extensions, use ``mpisppy.extensions.extension import MultiExtension``
 that allows you to give a list of extensions that will fire in order
 at each callout point. See, e.g. ``examples.sizes.sizes_demo.py`` for an
 example of use.

doc/src/index.rst

@@ -17,6 +17,7 @@ MPI is used.
    hubs.rst
    spokes.rst
    extensions.rst
+   mmw.rst
    pysp.rst
    ef.rst
    secretmenu.rst

doc/src/mmw.rst

@@ -0,0 +1,38 @@
+.. _MMW confidence interval:
+
+MMW confidence interval
+=======================
+
+If we want to assess the quality of a given candidate solution ``xhat``, we could
+try and evaluate the optimality gap, i.e. the gap between the value of the objective function
+at ``xhat`` and the value of the solution to our problem.
+The class ``MMWConfidenceIntervals`` compute an estimator of the optimality gap
+as described in [mmw1999]_ (Section 3.2) and an asymptotic confidence interval for
+this gap. 
+
+We will document two steps in the process : finding a candidate solution ``xhat``, 
+and evaluating it
+
+
+Finding a candidate solution
+----------------------------
+
+Computing this confidence interval means that we need to find a solution to 
+an approximate problem, and evaluate how good a solution to this approximate problem ``xhat`` is.
+In order to use MMW, ``xhat`` must be written using one of two functions 
+``ef_ROOT_nonants_npy_serializer`` or ``write_spin_the_wheel_first_stage_solution``.
+These functions write ``xhat`` to a file and can be read using ``read_xhat``.
+
+Computing a confidence interval
+-------------------------------
+
+The first step in computing a confidence interval is creating a ``MMWConfidenceIntervals`` object
+that takes as an argument an ``xhat`` and options.
+This object has a ``run`` method that returns a gap estimator and a confidence interval on the gap.
+
+Example
+-------
+
+An example of use, with the ``farmer`` problem, can be found in the main of ``mmwci.py``.
+
+

doc/src/refs.rst

@@ -8,3 +8,5 @@ References
 .. [gade2016] Gade, D., Hackebeil, G., Ryan, S., Watson, J.P., Wets, R., Woodruff, D.: Obtaining lower bounds from the progressive hedging algorithm for stochastic mixed-integer programming. Mathematical Programming 157(1), 47–67 (2016)
 
 .. [birge2011] Birge, J.R., Louveaux, F.: Introduction to Stochastic Programming, Second Ed. Springer Verlag (2011)
+
+.. [mmw1999] Mak, W.K., Morton,D.P., Wood, R.K.: Monte Carlo bounding techniques for determining solution quality in stochastic programs. Operations Research Letters 24(1-2), 47-56 (1999)

examples/acopf3/ccopf2wood.py

@@ -12,7 +12,7 @@
 from mpisppy.cylinders.hub import PHHub
 # Make it all go
 from mpisppy.utils.sputils import spin_the_wheel
-from mpisppy.utils.xhat_tryer import XhatTryer
+from mpisppy.utils.xhat_eval import Xhat_Eval
 
 # the problem
 import ACtree as etree
@@ -112,34 +112,34 @@ def main():
                                              rootnode.ScenarioList)+1)]
     all_nodenames = scenario_creator_kwargs["etree"].All_Nonleaf_Nodenames()
 
-    PHoptions = dict()
+    options = dict()
     if scenario_creator_kwargs["convex_relaxation"]:
-        PHoptions["solvername"] = solvername
-        if "gurobi" in PHoptions["solvername"]:
-            PHoptions["iter0_solver_options"] = {"BarHomogeneous": 1}
-            PHoptions["iterk_solver_options"] = {"BarHomogeneous": 1}
+        options["solvername"] = solvername
+        if "gurobi" in options["solvername"]:
+            options["iter0_solver_options"] = {"BarHomogeneous": 1}
+            options["iterk_solver_options"] = {"BarHomogeneous": 1}
         else:
-            PHoptions["iter0_solver_options"] = None
-            PHoptions["iterk_solver_options"] = None
+            options["iter0_solver_options"] = None
+            options["iterk_solver_options"] = None
     else:
-        PHoptions["solvername"] = solvername  # needs to be ipopt
-        PHoptions["iter0_solver_options"] = None
-        PHoptions["iterk_solver_options"] = None
-    PHoptions["PHIterLimit"] = PHIterLimit
-    PHoptions["defaultPHrho"] = 1
-    PHoptions["convthresh"] = 0.001
-    PHoptions["subsolvedirectives"] = None
-    PHoptions["verbose"] = False
-    PHoptions["display_timing"] = False
-    PHoptions["display_progress"] = True
-    PHoptions["iter0_solver_options"] = None
-    PHoptions["iterk_solver_options"] = None
-    PHoptions["branching_factors"] = branching_factors
+        options["solvername"] = solvername  # needs to be ipopt
+        options["iter0_solver_options"] = None
+        options["iterk_solver_options"] = None
+    options["PHIterLimit"] = PHIterLimit
+    options["defaultPHrho"] = 1
+    options["convthresh"] = 0.001
+    options["subsolvedirectives"] = None
+    options["verbose"] = False
+    options["display_timing"] = False
+    options["display_progress"] = True
+    options["iter0_solver_options"] = None
+    options["iterk_solver_options"] = None
+    options["branching_factors"] = branching_factors
 
     # try to do something interesting for bundles per rank
     if scenperbun > 0:
         nscen = branching_factors[0] * branching_factors[1]
-        PHoptions["bundles_per_rank"] = int((nscen / n_proc) / scenperbun)
+        options["bundles_per_rank"] = int((nscen / n_proc) / scenperbun)
     if global_rank == 0:
         appfile = "acopf.app"
         if not os.path.isfile(appfile):
@@ -150,53 +150,53 @@ def main():
                 f.write(", PH_lastiter, PH_wallclock, aph_frac_needed")
                 f.write(", APH_IB, APH_OB")
                 f.write(", APH_lastiter, APH_wallclock")
-        if "bundles_per_rank" in PHoptions:
-            nbunstr = str(PHoptions["bundles_per_rank"])
+        if "bundles_per_rank" in options:
+            nbunstr = str(options["bundles_per_rank"])
         else:
             nbunstr = "0"
         oline = "\n"+ str(start_time)+","+socket.gethostname()
         oline += ","+str(branching_factors[0])+","+str(branching_factors[1])
-        oline += ", "+str(seed) + ", "+str(PHoptions["solvername"])
+        oline += ", "+str(seed) + ", "+str(options["solvername"])
         oline += ", "+str(n_proc) + ", "+ nbunstr
-        oline += ", "+str(PHoptions["PHIterLimit"])
-        oline += ", "+str(PHoptions["convthresh"])
+        oline += ", "+str(options["PHIterLimit"])
+        oline += ", "+str(options["convthresh"])
 
         with open(appfile, "a") as f:
             f.write(oline)
         print(oline)
 
 
     # PH hub
-    PHoptions["tee-rank0-solves"] = True
+    options["tee-rank0-solves"] = True
     hub_dict = {
         "hub_class": PHHub,
         "hub_kwargs": {"options": None},
         "opt_class": PH,
         "opt_kwargs": {
-            "PHoptions": PHoptions,
+            "options": options,
             "all_scenario_names": all_scenario_names,
             "scenario_creator": pysp2_callback,
             'scenario_denouement': scenario_denouement,
             "scenario_creator_kwargs": scenario_creator_kwargs,
             "rho_setter": rho_setter.ph_rhosetter_callback,
-            "PH_extensions": None,
+            "extensions": None,
             "all_nodenames":all_nodenames,
         }
     }
 
-    xhat_options = PHoptions.copy()
+    xhat_options = options.copy()
     xhat_options['bundles_per_rank'] = 0 #  no bundles for xhat
     xhat_options["xhat_specific_options"] = {"xhat_solver_options":
-                                          PHoptions["iterk_solver_options"],
+                                          options["iterk_solver_options"],
                                           "xhat_scenario_dict": xhat_dict,
                                           "csvname": "specific.csv"}
 
     ub2 = {
         'spoke_class': XhatSpecificInnerBound,
         "spoke_kwargs": dict(),
-        "opt_class": XhatTryer,
+        "opt_class": Xhat_Eval,
         'opt_kwargs': {
-            'PHoptions': xhat_options,
+            'options': xhat_options,
             'all_scenario_names': all_scenario_names,
             'scenario_creator': pysp2_callback,
             'scenario_denouement': scenario_denouement,

examples/acopf3/fourstage.py

@@ -12,7 +12,7 @@
 from mpisppy.cylinders.hub import PHHub
 # Make it all go
 from mpisppy.utils.sputils import spin_the_wheel
-from mpisppy.utils.xhat_tryer import XhatTryer
+from mpisppy.utils.xhat_eval import Xhat_Eval
 
 # the problem
 import ACtree as etree
@@ -119,34 +119,34 @@ def main():
                                              rootnode.ScenarioList)+1)]
     all_nodenames = scenario_creator_kwargs["etree"].All_Nonleaf_Nodenames()
 
-    PHoptions = dict()
+    options = dict()
     if convex_relaxation:
-        PHoptions["solvername"] = solvername
-        if "gurobi" in PHoptions["solvername"]:
-            PHoptions["iter0_solver_options"] = {"BarHomogeneous": 1}
-            PHoptions["iterk_solver_options"] = {"BarHomogeneous": 1}
+        options["solvername"] = solvername
+        if "gurobi" in options["solvername"]:
+            options["iter0_solver_options"] = {"BarHomogeneous": 1}
+            options["iterk_solver_options"] = {"BarHomogeneous": 1}
         else:
-            PHoptions["iter0_solver_options"] = None
-            PHoptions["iterk_solver_options"] = None
+            options["iter0_solver_options"] = None
+            options["iterk_solver_options"] = None
     else:
-        PHoptions["solvername"] = "ipopt"
-        PHoptions["iter0_solver_options"] = None
-        PHoptions["iterk_solver_options"] = None
-    PHoptions["PHIterLimit"] = PHIterLimit
-    PHoptions["defaultPHrho"] = 1
-    PHoptions["convthresh"] = 0.001
-    PHoptions["subsolvedirectives"] = None
-    PHoptions["verbose"] = False
-    PHoptions["display_timing"] = False
-    PHoptions["display_progress"] = True
-    PHoptions["iter0_solver_options"] = None
-    PHoptions["iterk_solver_options"] = None
-    PHoptions["branching_factors"] = branching_factors
+        options["solvername"] = "ipopt"
+        options["iter0_solver_options"] = None
+        options["iterk_solver_options"] = None
+    options["PHIterLimit"] = PHIterLimit
+    options["defaultPHrho"] = 1
+    options["convthresh"] = 0.001
+    options["subsolvedirectives"] = None
+    options["verbose"] = False
+    options["display_timing"] = False
+    options["display_progress"] = True
+    options["iter0_solver_options"] = None
+    options["iterk_solver_options"] = None
+    options["branching_factors"] = branching_factors
 
     # try to do something interesting for bundles per rank
     if scenperbun > 0:
         nscen = branching_factors[0] * branching_factors[1]
-        PHoptions["bundles_per_rank"] = int((nscen / n_proc) / scenperbun)
+        options["bundles_per_rank"] = int((nscen / n_proc) / scenperbun)
     if global_rank == 0:
         appfile = "acopf.app"
         if not os.path.isfile(appfile):
@@ -157,53 +157,53 @@ def main():
                 f.write(", PH_lastiter, PH_wallclock, aph_frac_needed")
                 f.write(", APH_IB, APH_OB")
                 f.write(", APH_lastiter, APH_wallclock")
-        if "bundles_per_rank" in PHoptions:
-            nbunstr = str(PHoptions["bundles_per_rank"])
+        if "bundles_per_rank" in options:
+            nbunstr = str(options["bundles_per_rank"])
         else:
             nbunstr = "0"
         oline = "\n"+ str(start_time)+","+socket.gethostname()
         oline += ","+str(branching_factors[0])+","+str(branching_factors[1])
-        oline += ", "+str(seed) + ", "+str(PHoptions["solvername"])
+        oline += ", "+str(seed) + ", "+str(options["solvername"])
         oline += ", "+str(n_proc) + ", "+ nbunstr
-        oline += ", "+str(PHoptions["PHIterLimit"])
-        oline += ", "+str(PHoptions["convthresh"])
+        oline += ", "+str(options["PHIterLimit"])
+        oline += ", "+str(options["convthresh"])
 
         with open(appfile, "a") as f:
             f.write(oline)
         print(oline)
 
 
     # PH hub
-    PHoptions["tee-rank0-solves"] = True
+    options["tee-rank0-solves"] = True
     hub_dict = {
         "hub_class": PHHub,
         "hub_kwargs": {"options": None},
         "opt_class": PH,
         "opt_kwargs": {
-            "PHoptions": PHoptions,
+            "options": options,
             "all_scenario_names": all_scenario_names,
             "scenario_creator": pysp2_callback,
             'scenario_denouement': scenario_denouement,
             "scenario_creator_kwargs": scenario_creator_kwargs,
             "rho_setter": rho_setter.ph_rhosetter_callback,
-            "PH_extensions": None,
+            "extensions": None,
             "all_nodenames":all_nodenames,
         }
     }
 
-    xhat_options = PHoptions.copy()
+    xhat_options = options.copy()
     xhat_options['bundles_per_rank'] = 0 #  no bundles for xhat
     xhat_options["xhat_specific_options"] = {"xhat_solver_options":
-                                          PHoptions["iterk_solver_options"],
+                                          options["iterk_solver_options"],
                                           "xhat_scenario_dict": xhat_dict,
                                           "csvname": "specific.csv"}
 
     ub2 = {
         'spoke_class': XhatSpecificInnerBound,
         "spoke_kwargs": dict(),
-        "opt_class": XhatTryer,
+        "opt_class": Xhat_Eval,
         'opt_kwargs': {
-            'PHoptions': xhat_options,
+            'options': xhat_options,
             'all_scenario_names': all_scenario_names,
             'scenario_creator': pysp2_callback,
             'scenario_denouement': scenario_denouement,