Ph dual grad rho

doc/src/extensions.rst

@@ -194,7 +194,7 @@ CoeffRho
 ^^^^^^^^
 
 Set per variable rho values proportional to the cost coefficient on each non-anticipative variable,
-with an optional multiplier (default = 1.0). If the coefficient is 0, the default rho value is used instead.
+with an optional multiplier (default = 1.0) that is applied to the computed value. If the coefficient is 0, the default rho value is used instead.
 
 primal_dual_rho
 ^^^^^^^^^^^^^^^
@@ -228,7 +228,7 @@ There are options in ``cfg`` to control dynamic updates.
 mult_rho_updater
 ^^^^^^^^^^^^^^^^
 
-This extension does a simple multiplicative update of rho.
+This extension does a simple multiplicative update of rho; consequently, the update is cumulative.
 
 cross-scenario cuts
 ^^^^^^^^^^^^^^^^^^^

examples/run_all.py

@@ -158,6 +158,7 @@ def do_one_mmw(dirname, modname, runefstring, npyfile, mmwargstring):
 
         os.remove(npyfile)
     os.chdir("..")
+    os.chdir("..")  # moved to CI directory
 
 do_one("farmer/CI", "farmer_ef.py", 1,
        "1 3 {}".format(solver_name))
@@ -190,13 +191,25 @@ def do_one_mmw(dirname, modname, runefstring, npyfile, mmwargstring):
        "--solver-name={}".format(solver_name))
 do_one("farmer/archive", "farmer_cylinders.py", 4,
        "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
-       "--fwph-stop-check-tol 0.1 "
+       "--fwph-stop-check-tol 0.1 --rel-gap 0.001 "
        "--default-rho=1 --solver-name={} --lagrangian --xhatshuffle --fwph".format(solver_name))
+do_one("farmer", "../../mpisppy/generic_cylinders.py", 3,
+       "--module-name farmer --num-scens 6 "
+       "--rel-gap 0.001 --max-iterations=50 "
+       "--grad-rho --grad-order-stat 0.5 "
+       "--default-rho=2 --solver-name={} --lagrangian --xhatshuffle".format(solver_name))
 do_one("farmer", "../../mpisppy/generic_cylinders.py", 4,
-       "--module-name farmer "
-       "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
+       "--module-name farmer --num-scens 6 "
+       "--rel-gap 0.001 --max-iterations=50 "
        "--ph-primal-hub --ph-dual --ph-dual-rescale-rho-factor=0.1 "
+       "--default-rho=2 --solver-name={} --lagrangian --xhatshuffle".format(solver_name))
+do_one("farmer", "../../mpisppy/generic_cylinders.py", 4,
+       "--module-name farmer "
+       "--num-scens 6 --max-iterations=50 --grad-rho --grad-order-stat 0.5 "
+       "--ph-dual-grad-order-stat 0.3 "
+       "--ph-primal-hub --ph-dual --ph-dual-rescale-rho-factor=0.1 --ph-dual-rho-multiplier 0.2 "
        "--default-rho=1 --solver-name={} --lagrangian --xhatshuffle".format(solver_name))
+
 do_one("farmer/archive", "farmer_cylinders.py", 2,
        "--num-scens 6 --bundles-per-rank=2 --max-iterations=50 "
        "--default-rho=1 "

mpisppy/extensions/grad_rho.py

@@ -312,4 +312,4 @@ def register_receive_fields(self):
             self.best_xhat_spoke_index,
         )
 
-        return
+        return

mpisppy/generic_cylinders.py

@@ -11,6 +11,7 @@
 import sys
 import os
 import json
+import copy
 import shutil
 import numpy as np
 import pyomo.environ as pyo
@@ -340,14 +341,21 @@ def _do_decomp(module, cfg, scenario_creator, scenario_creator_kwargs, scenario_
                                           rho_setter = rho_setter,
                                           all_nodenames = all_nodenames,
                                           )
+        if cfg.sep_rho or cfg.coeff_rho or cfg.sensi_rho or cfg.grad_rho:
+            # Note that this deepcopy might be expensive if certain wrappers were used.
+            # (Could we do the modification to cfg in ph_dual to obviate the need?)
+            modified_cfg = copy.deepcopy(cfg)
+            modified_cfg["grad_rho_multiplier"] = cfg.ph_dual_rho_multiplier            
         if cfg.sep_rho:
-            vanilla.add_sep_rho(ph_dual_spoke, cfg)
+            vanilla.add_sep_rho(ph_dual_spoke, modified_cfg)
         if cfg.coeff_rho:
-            vanilla.add_coeff_rho(ph_dual_spoke, cfg)
+            vanilla.add_coeff_rho(ph_dual_spoke, modified_cfg)
         if cfg.sensi_rho:
-            vanilla.add_sensi_rho(ph_dual_spoke, cfg)
-        # TBD xxxxx add grad rho asap after PR#559 is merged
-        # Note: according to DLW, it is non-trivial to deal with the cfg args
+            vanilla.add_sensi_rho(ph_dual_spoke, modified_cfg)
+        if cfg.grad_rho:
+            modified_cfg["grad_order_stat"] = cfg.ph_dual_grad_order_stat
+            vanilla.add_grad_rho(ph_dual_spoke, modified_cfg)
+
 
     # relaxed ph spoke
     if cfg.relaxed_ph:

mpisppy/utils/config.py

@@ -143,21 +143,30 @@ def get(self, name, ifmissing=None):
     def checker(self):
         """Verify that options *selected* make sense with respect to each other
         """
-        def _bad_rho_setters(msg):
-            raise ValueError("Rho setter options do not make sense together:\n"
+        def _bad_options(msg):
+            raise ValueError("Options do not make sense together:\n"
                              f"{msg}")
-
-        if self.get("grad_rho") and self.get("sensi_rho"):
-            _bad_rho_setters("Only one rho setter can be active.")
+
+        # remember that True is 1 and False is 0
+        if (self.get("grad_rho") + self.get("sensi_rho") + self.get("coeff_rho") + self.get("reduced_costs_rho") + self.get("sep_rho")) > 1:
+            _bad_options("Only one rho setter can be active.")
         if not (self.get("grad_rho")
                 or self.get("sensi_rho")
                 or self.get("sep_rho")
                 or self.get("reduced_costs_rho")):
             if self.get("dynamic_rho_primal_crit") or self.get("dynamic_rho_dual_crit"):
-                _bad_rho_setters("dynamic rho only works with grad-, sensi-, and sep-rho")
+                _bad_options("dynamic rho only works with an automated rho setter")
+        if self.get("grad_rho") and self.get("bundles_per_rank") != 0:
+            _bad_options("Grad rho does not work with loose bundling (--bundles-per-rank).\n "
+                         "Also note that loose bundling is being deprecated in favor of proper bundles.")
+
+        if self.get("ph_primal_hub")\
+           and not (self.get("ph_dual") or self.get("relaxed_ph")):
+            _bad_options("--ph-primal-hub is used only when there is a cylinder that provideds Ws "
+                         "such as --ph-dual or --relaxed-ph")
+
         if self.get("rc_fixer") and not self.get("reduced_costs"):
-            _bad_rho_setters("--rc-fixer requires --reduced-costs")
-
+            _bad_options("--rc-fixer requires --reduced-costs")
 
     def add_solver_specs(self, prefix=""):
         sstr = f"{prefix}_solver" if prefix != "" else "solver"
@@ -732,7 +741,8 @@ def subgradient_bounder_args(self):
 
 
     def ph_ob_args(self):
-        raise RuntimeError("ph_ob (the --ph-ob option) and ph_ob_args were deprecated and replaced with ph_dual August 2025")
+        raise RuntimeError("ph_ob (the --ph-ob option) and ph_ob_args were deprecated and replaced with ph_dual August 2025\n"
+                           "To get the same effect as ph_ob, use --ph-dual with --ph-dual-grad-rho")
 
     def relaxed_ph_args(self):
 
@@ -753,9 +763,19 @@ def ph_dual_args(self):
                             domain=bool,
                             default=False)
         self.add_to_config("ph_dual_rescale_rho_factor",
-                            description="Used to rescale rho initially (default=1.0)",
+                            description="Used to rescale rho initially (default=0.1)",
+                            domain=float,
+                            default=0.1)
+        self.add_to_config("ph_dual_rho_multiplier",
+                            description="Rescale factor for dynamic updates in ph_dual if ph_dual and a rho setter are chosen;"
+                            " note that it is not cummulative (default=1.0)",
                             domain=float,
                             default=1.0)
+        self.add_to_config("ph_dual_grad_order_stat",
+                            description="Order stat for selecting rho if ph_dual and ph_dual_grad_rho are chosen;"
+                            " note that this is impacted by the multiplier (default=0.0)",
+                            domain=float,
+                            default=0.0)
 
 
     def xhatlooper_args(self):