PHDualSpoke

doc/src/hubs.rst

@@ -47,11 +47,11 @@ this value (default 1e-1).
 If only the binary terms should be 
 approximated, the option `linearize_binary_proximal_terms` can be used. 
 
-PHNonant
+PHPrimal
 --------
 This is exactly like the PH hub, except it does not send W values.
 It can be useful when some other cylinder is providing W values
-(e.g., RelaxedPHSpoke).
+(e.g., PHDualSpoke, RelaxedPHSpoke).
 
 
 lshaped

doc/src/spokes.rst

@@ -174,15 +174,23 @@ cross scenario
 
 Computes and passes cross scenario cuts.
 
-relaxed_ph
+Relaxed PH
 ^^^^^^^^^^
 
-For S-MIPs, runs progressive hedging on the linear programming relaxation
-of the scenario subproblems. Provides Ws and "relaxed nonants", the former
-of which can be utilized for Lagrangian to compute lower bounds, and the later
-of which can be utilized to inform fixings for the hub via the relaxed_ph_fixer
-extention. This method can be effective when the subproblem linear programming
-relaxation is "strong".
-
-The option ``relaxed_ph_rescale_rho_factor``, which defaults to 1, rescales
-rho for this spoke between iteration 0 and iteration 1.
+For S-MIPs, executes PH on the continuous relaxation of the problem. Provides W values
+to the Lagrangian or Reduced Costs spokes, as well as relaxed
+nonants which can be used with the relaxed PH fixer for the PH Primal Hub. This spoke
+can be useful, along with the PH Primal Hub, when the continuous relaxation allows for
+good dual (but not primal) convergence of the original problem. This typically happens
+for problems with a "strong" continuous relaxation. The option ``relaxed_ph_rescale_rho_factor``
+allows the user to adjust provided rho values by a constant multiplier across all variables,
+which occurs between PH iteration 0 and PH iteration 1.
+
+PH Dual
+^^^^^^^
+
+Executes another copy of PH, providing W values to the Lagrangian and/or Reduceds Costs
+spokes. This spoke can be useful, along with the PH Primal Hub, for problems which need
+different rho strategies for primal and dual convergence. The option ``ph_dual_rescale_rho_factor``
+allows the user to adjust provided rho values by a constant multiplier across all variables,
+which occurs between PH iteration 0 and PH iteration 1.

examples/run_all.py

@@ -158,6 +158,7 @@ def do_one_mmw(dirname, runefstring, npyfile, mmwargstring):
 
         os.remove(npyfile)
     os.chdir("..")
+
 do_one("farmer", "farmer_ef.py", 1,
        "1 3 {}".format(solver_name))
 # for farmer_cylinders, the first arg is num_scens and is required
@@ -191,6 +192,11 @@ def do_one_mmw(dirname, runefstring, npyfile, mmwargstring):
        "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
        "--fwph-stop-check-tol 0.1 "
        "--default-rho=1 --solver-name={} --lagrangian --xhatshuffle --fwph".format(solver_name))
+do_one("farmer", "../../mpisppy/generic_cylinders.py", 4,
+       "--module-name farmer "
+       "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
+       "--ph-primal-hub --ph-dual --ph-dual-rescale-rho-factor=0.1 "
+       "--default-rho=1 --solver-name={} --lagrangian --xhatshuffle".format(solver_name))
 do_one("farmer", "farmer_cylinders.py", 2,
        "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
        "--default-rho=1 "
@@ -229,6 +235,18 @@ def do_one_mmw(dirname, runefstring, npyfile, mmwargstring):
        "farmer_cylinders.py",
        2,
        f"--num-scens 30 --max-iterations=10 --default-rho=1.0 --display-progress  --bundles-per-rank=0 --xhatshuffle --aph-gamma=1.0 --aph-nu=1.0 --aph-frac-needed=1.0 --aph-dispatch-frac=0.5 --abs-gap=1 --aph-sleep-seconds=0.01 --run-async --bundles-per-rank=5 --solver-name={solver_name}")
+do_one("farmer",
+       "../../mpisppy/generic_cylinders.py",
+       4,
+       "--module-name farmer --farmer-with-integer "
+       "--num-scens=3 "
+       "--lagrangian --ph-primal-hub "
+       "--max-iterations=10 --default-rho=0.1 "
+       "--relaxed-ph-rescale-rho-factor=10 "
+       "--relaxed-ph --relaxed-ph-fixer --xhatshuffle "
+       "--linearize-proximal-terms "
+       "--rel-gap=0.0 "
+       "--solver-name={}".format(solver_name))
 
 do_one("farmer",
        "farmer_ama.py",
@@ -275,10 +293,10 @@ def do_one_mmw(dirname, runefstring, npyfile, mmwargstring):
        4,
        "--instance-name=sslp_15_45_10 --bundles-per-rank=0 "
        "--integer-relax-then-enforce "
-       "--integer-relax-then-enforce-ratio=0.95 "
+       "--integer-relax-then-enforce-ratio=0.8 "
        "--lagrangian "
        "--reduced-costs-rho "
-       "--max-iterations=100 --default-rho=1e-6 "
+       "--max-iterations=20 --default-rho=1e-6 "
        "--reduced-costs --rc-fixer --xhatshuffle "
        "--linearize-proximal-terms "
        "--rel-gap=0.0 --surrogate-nonant "

mpisppy/cylinders/hub.py

@@ -202,7 +202,7 @@ def sync_bounds(self):
         self.send_boundsout()
 
 
-class PHNonantHub(Hub):
+class PHPrimalHub(Hub):
     """
     Like PHHub, but only sends nonants and omits Ws. To be used
     when another cylinder is supplying Ws (like RelaxedPHSpoke).
@@ -320,9 +320,9 @@ def send_ws(self):
         pass
 
 
-class PHHub(PHNonantHub):
-    send_fields = (*PHNonantHub.send_fields, Field.DUALS, )
-    receive_fields = (*PHNonantHub.receive_fields,)
+class PHHub(PHPrimalHub):
+    send_fields = (*PHPrimalHub.send_fields, Field.DUALS, )
+    receive_fields = (*PHPrimalHub.receive_fields,)
 
     def send_ws(self):
         """ Send dual weights to the appropriate spokes

mpisppy/cylinders/ph_dual_spoke.py

@@ -0,0 +1,58 @@
+###############################################################################
+# mpi-sppy: MPI-based Stochastic Programming in PYthon
+#
+# Copyright (c) 2024, Lawrence Livermore National Security, LLC, Alliance for
+# Sustainable Energy, LLC, The Regents of the University of California, et al.
+# All rights reserved. Please see the files COPYRIGHT.md and LICENSE.md for
+# full copyright and license information.
+###############################################################################
+
+from mpisppy.cylinders.spwindow import Field
+from mpisppy.cylinders.spoke import Spoke
+from mpisppy.cylinders.hub import PHHub
+
+
+class _PHDualSpokeBase(Spoke, PHHub):
+
+    send_fields = (*Spoke.send_fields, Field.DUALS,) 
+    receive_fields = (*Spoke.receive_fields, )
+
+    def send_boundsout(self):
+        # overwrite PHHub.send_boundsout (not a hub)
+        return
+
+    def update_rho(self):
+        rho_factor = self.opt.options.get("rho_factor", 1.0)
+        if rho_factor == 1.0:
+            return
+        for s in self.opt.local_scenarios.values():
+            for rho in s._mpisppy_model.rho.values():
+                rho._value = rho_factor * rho._value
+
+    def main(self):
+        # setup, PH Iter0
+        smoothed = self.options.get('smoothed', 0)
+        attach_prox = True
+        self.opt.PH_Prep(attach_prox=attach_prox, attach_smooth = smoothed)
+        trivial_bound = self.opt.Iter0()
+        if self.opt._can_update_best_bound():
+            self.opt.best_bound_obj_val = trivial_bound
+
+        # update the rho
+        self.update_rho()
+
+        # rest of PH
+        self.opt.iterk_loop()
+
+        return self.opt.conv, None, trivial_bound
+
+
+class PHDualSpoke(_PHDualSpokeBase):
+
+    @property
+    def nonant_field(self):
+        return None
+
+    def send_nonants(self):
+        # overwrite PHHub.send_nonants (don't send anything)
+        return

mpisppy/cylinders/relaxed_ph_spoke.py

@@ -8,50 +8,23 @@
 ###############################################################################
 
 from mpisppy.cylinders.spwindow import Field
-from mpisppy.cylinders.spoke import Spoke
-from mpisppy.cylinders.hub import PHHub
+from mpisppy.cylinders.ph_dual_spoke import _PHDualSpokeBase
 
 import pyomo.environ as pyo
 
-class RelaxedPHSpoke(Spoke, PHHub):
+class RelaxedPHSpoke(_PHDualSpokeBase):
 
-    send_fields = (*Spoke.send_fields, Field.DUALS, Field.RELAXED_NONANT, )
-    receive_fields = (*Spoke.receive_fields, )
+    send_fields = (*_PHDualSpokeBase.send_fields, Field.RELAXED_NONANT, )
+    receive_fields = (*_PHDualSpokeBase.receive_fields, )
 
     @property
     def nonant_field(self):
         return Field.RELAXED_NONANT
 
-    def send_boundsout(self):
-        # overwrite PHHub.send_boundsout (not a hub)
-        return
-
-    def update_rho(self):
-        rho_factor = self.opt.options.get("relaxed_ph_rho_factor", 1.0)
-        if rho_factor == 1.0:
-            return
-        for s in self.opt.local_scenarios.values():
-            for rho in s._mpisppy_model.rho.values():
-                rho._value = rho_factor * rho._value
-
     def main(self):
         # relax the integers
         integer_relaxer = pyo.TransformationFactory("core.relax_integer_vars")
         for s in self.opt.local_scenarios.values():
             integer_relaxer.apply_to(s)
 
-        # setup, PH Iter0
-        smoothed = self.options.get('smoothed', 0)
-        attach_prox = True
-        self.opt.PH_Prep(attach_prox=attach_prox, attach_smooth = smoothed)
-        trivial_bound = self.opt.Iter0()
-        if self.opt._can_update_best_bound():
-            self.opt.best_bound_obj_val = trivial_bound
-
-        # update the rho
-        self.update_rho()
-
-        # rest of PH
-        self.opt.iterk_loop()
-
-        return self.opt.conv, None, trivial_bound
+        return super().main()

mpisppy/generic_cylinders.py

@@ -82,7 +82,8 @@ def _parse_args(m):
     cfg.integer_relax_then_enforce_args()
     cfg.gapper_args()    
     cfg.gapper_args(name="lagrangian")
-    cfg.ph_nonant_args()
+    cfg.ph_primal_args()
+    cfg.ph_dual_args()
     cfg.relaxed_ph_args()
     cfg.fwph_args()
     cfg.lagrangian_args()
@@ -192,8 +193,8 @@ def _do_decomp(module, cfg, scenario_creator, scenario_creator_kwargs, scenario_
                        rho_setter = rho_setter,
                        all_nodenames = all_nodenames,
                    )
-    elif cfg.ph_nonant_hub:
-        hub_dict = vanilla.ph_nonant_hub(*beans,
+    elif cfg.ph_primal_hub:
+        hub_dict = vanilla.ph_primal_hub(*beans,
                                   scenario_creator_kwargs=scenario_creator_kwargs,
                                   ph_extensions=None,
                                   ph_converger=ph_converger,
@@ -343,6 +344,20 @@ def _do_decomp(module, cfg, scenario_creator, scenario_creator_kwargs, scenario_
         if cfg.sensi_rho:
             vanilla.add_sensi_rho(ph_ob_spoke, cfg)
 
+    # dual ph spoke
+    if cfg.ph_dual:
+        ph_dual_spoke = vanilla.ph_dual_spoke(*beans,
+                                          scenario_creator_kwargs=scenario_creator_kwargs,
+                                          rho_setter = rho_setter,
+                                          all_nodenames = all_nodenames,
+                                          )
+        if cfg.sep_rho:
+            vanilla.add_sep_rho(ph_dual_spoke, cfg)
+        if cfg.coeff_rho:
+            vanilla.add_coeff_rho(ph_dual_spoke, cfg)
+        if cfg.sensi_rho:
+            vanilla.add_sensi_rho(ph_dual_spoke, cfg)
+
     # relaxed ph spoke
     if cfg.relaxed_ph:
         relaxed_ph_spoke = vanilla.relaxed_ph_spoke(*beans,
@@ -404,6 +419,8 @@ def _do_decomp(module, cfg, scenario_creator, scenario_creator_kwargs, scenario_
         list_of_spoke_dict.append(lagrangian_spoke)
     if cfg.ph_ob:
         list_of_spoke_dict.append(ph_ob_spoke)
+    if cfg.ph_dual:
+        list_of_spoke_dict.append(ph_dual_spoke)
     if cfg.relaxed_ph:
         list_of_spoke_dict.append(relaxed_ph_spoke)
     if cfg.subgradient:

mpisppy/utils/cfg_vanilla.py

@@ -37,7 +37,8 @@
 from mpisppy.cylinders.cross_scen_spoke import CrossScenarioCutSpoke
 from mpisppy.cylinders.reduced_costs_spoke import ReducedCostsSpoke
 from mpisppy.cylinders.relaxed_ph_spoke import RelaxedPHSpoke
-from mpisppy.cylinders.hub import PHNonantHub, PHHub, SubgradientHub, APHHub, FWPHHub
+from mpisppy.cylinders.ph_dual_spoke import PHDualSpoke
+from mpisppy.cylinders.hub import PHPrimalHub, PHHub, SubgradientHub, APHHub, FWPHHub
 from mpisppy.extensions.extension import MultiExtension
 from mpisppy.extensions.fixer import Fixer
 from mpisppy.extensions.mipgapper import Gapper
@@ -154,7 +155,7 @@ def ph_hub(
     add_timed_mipgap(hub_dict, cfg)
     return hub_dict
 
-def ph_nonant_hub(
+def ph_primal_hub(
         cfg,
         scenario_creator,
         scenario_denouement,
@@ -181,7 +182,7 @@ def ph_nonant_hub(
         all_nodenames=all_nodenames,
     )
     # use PHNonantHub instead of PHHub
-    hub_dict["hub_class"] = PHNonantHub
+    hub_dict["hub_class"] = PHPrimalHub
     return hub_dict
 
 def aph_hub(cfg,
@@ -765,6 +766,44 @@ def subgradient_spoke(
     return subgradient_spoke
 
 
+def ph_dual_spoke(
+    cfg,
+    scenario_creator,
+    scenario_denouement,
+    all_scenario_names,
+    scenario_creator_kwargs=None,
+    rho_setter=None,
+    all_nodenames=None,
+    ph_extensions=None,
+    extension_kwargs=None,
+):
+    ph_dual_spoke = _PHBase_spoke_foundation(
+        PHDualSpoke,
+        cfg,
+        scenario_creator,
+        scenario_denouement,
+        all_scenario_names,
+        scenario_creator_kwargs=scenario_creator_kwargs,
+        rho_setter=rho_setter,
+        all_nodenames=all_nodenames,
+        ph_extensions=ph_extensions,
+        extension_kwargs=extension_kwargs,
+    )
+    options = ph_dual_spoke["opt_kwargs"]["options"]
+    if cfg.ph_dual_rescale_rho_factor is not None:
+        options["rho_factor"] = cfg.ph_dual_rescale_rho_factor
+
+    # make sure this spoke doesn't hit the time or iteration limit
+    options["time_limit"] = None
+    options["PHIterLimit"] = cfg.max_iterations * 1_000_000
+    options["display_progress"] = False
+    options["display_convergence_detail"] = False
+
+    add_ph_tracking(ph_dual_spoke, cfg, spoke=True)
+
+    return ph_dual_spoke
+
+
 def relaxed_ph_spoke(
     cfg,
     scenario_creator,
@@ -790,7 +829,7 @@ def relaxed_ph_spoke(
     )
     options = relaxed_ph_spoke["opt_kwargs"]["options"]
     if cfg.relaxed_ph_rescale_rho_factor is not None:
-        options["relaxed_ph_rho_factor"] = cfg.relaxed_ph_rescale_rho_factor
+        options["rho_factor"] = cfg.relaxed_ph_rescale_rho_factor
 
     # make sure this spoke doesn't hit the time or iteration limit
     options["time_limit"] = None