Subgradient hub

mpisppy/cylinders/hub.py

@@ -10,6 +10,7 @@
 import abc
 import logging
 import mpisppy.log
+from mpisppy.opt.subgradient import Subgradient
 from mpisppy.opt.aph import APH
 
 from mpisppy import MPI
@@ -538,6 +539,7 @@ def is_converged(self):
         ## might as well get a bound, in this case
         if self.opt._PHIter == 1 and self.use_trivial_bound:
             self.BestOuterBound = self.OuterBoundUpdate(self.opt.trivial_bound)
+        self.BestOuterBound = self.OuterBoundUpdate(self.opt.best_bound_obj_val)
 
         if not self.has_innerbound_spokes:
             if self.opt._PHIter == 1:
@@ -553,7 +555,7 @@ def is_converged(self):
             return False
 
         if not self.has_outerbound_spokes:
-            if self.opt._PHIter == 1:
+            if self.opt._PHIter == 1 and not isinstance(self.opt, Subgradient):
                 global_toc(
                     "Without outer bound spokes, no progress "
                     "will be made on the Best Bound")

mpisppy/generic_cylinders.py

@@ -62,6 +62,7 @@ def _parse_args(m):
     cfg.two_sided_args()
     cfg.ph_args()
     cfg.aph_args()
+    cfg.sub_args()
     cfg.fixer_args()    
     cfg.integer_relax_then_enforce_args()
     cfg.gapper_args()    
@@ -143,6 +144,16 @@ def _do_decomp(module, cfg, scenario_creator, scenario_creator_kwargs, scenario_
                                    rho_setter = rho_setter,
                                    all_nodenames = all_nodenames,
                                    )
+    elif cfg.SUB:
+        # Vanilla Subgradient hub
+        hub_dict = vanilla.subgradient_hub(
+                       *beans,
+                       scenario_creator_kwargs=scenario_creator_kwargs,
+                       ph_extensions=None,
+                       ph_converger=ph_converger,
+                       rho_setter = rho_setter,
+                       all_nodenames = all_nodenames,
+                   )
     else:
         # Vanilla PH hub
         hub_dict = vanilla.ph_hub(*beans,

mpisppy/opt/aph.py

@@ -1067,6 +1067,7 @@ def APH_main(self, spcomm=None, finalize=True):
         # End APH-specific Prep
 
         trivial_bound = self.Iter0()
+        self.best_bound_obj_val = trivial_bound
 
         self.setup_Lens()
         self.setup_dispatchrecord()

mpisppy/opt/ph.py

@@ -61,6 +61,7 @@ def ph_main(self, finalize=True):
         if (verbose):
             print('Calling PH Iter0 on global rank {}'.format(global_rank))
         trivial_bound = self.Iter0()
+        self.best_bound_obj_val = trivial_bound
         if (verbose):
             print ('Completed PH Iter0 on global rank {}'.format(global_rank))
         if ('asynchronousPH' in self.options) and (self.options['asynchronousPH']):

mpisppy/opt/subgradient.py

@@ -0,0 +1,146 @@
+###############################################################################
+# 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.
+###############################################################################
+
+import mpisppy.phbase
+import mpisppy.MPI as mpi
+
+from pyomo.common.collections import ComponentSet
+
+_global_rank = mpi.COMM_WORLD.Get_rank()
+
+class Subgradient(mpisppy.phbase.PHBase):
+    """ Subgradient Algorithm """
+
+    def subgradient_main(self, finalize=True):
+        """ Execute the subgradient algorithm.
+
+        Args:
+            finalize (bool, optional, default=True):
+                If True, call `PH.post_loops()`, if False, do not,
+                and return None for Eobj
+
+        Returns:
+            tuple:
+                Tuple containing
+
+                conv (float):
+                    The convergence value (not easily interpretable).
+                Eobj (float or `None`):
+                    If `finalize=True`, this is the expected, weighted
+                    objective value. This value is not directly useful.
+                    If `finalize=False`, this value is `None`.
+                trivial_bound (float):
+                    The "trivial bound", computed by solving the model with no
+                    nonanticipativity constraints (immediately after iter 0).
+        """
+        verbose = self.options['verbose']
+        smoothed = self.options['smoothed']
+        if smoothed != 0:
+            raise RuntimeError("Cannnot use smoothing with Subgradient algorithm")
+        self.create_fixed_nonant_cache()
+        self.PH_Prep(attach_prox=False, attach_smooth=smoothed)
+
+        if (verbose):
+            print('Calling Subgradient Iter0 on global rank {}'.format(_global_rank))
+        trivial_bound = self.Iter0()
+        self.best_bound_obj_val = trivial_bound
+        if (verbose):
+            print('Completed Subgradient Iter0 on global rank {}'.format(_global_rank))
+
+        self.iterk_loop()
+
+        if finalize:
+            Eobj = self.post_loops(self.extensions)
+        else:
+            Eobj = None
+
+        return self.conv, Eobj, trivial_bound
+
+    def ph_main(self, finalize=True):
+        # for working with a PHHub
+        return self.subgradient_main(finalize=finalize)
+
+    def solve_loop(self,
+                   solver_options=None,
+                   use_scenarios_not_subproblems=False,
+                   dtiming=False,
+                   dis_W=False,
+                   dis_prox=False,
+                   gripe=False,
+                   disable_pyomo_signal_handling=False,
+                   tee=False,
+                   verbose=False,
+                   need_solution=True,
+                   ):
+        """ Loop over `local_subproblems` and solve them in a manner
+        dicated by the arguments.
+
+        In addition to changing the Var values in the scenarios, this function
+        also updates the `_PySP_feas_indictor` to indicate which scenarios were
+        feasible/infeasible.
+
+        Args:
+            solver_options (dict, optional):
+                The scenario solver options.
+            use_scenarios_not_subproblems (boolean, optional):
+                If True, solves individual scenario problems, not subproblems.
+                This distinction matters when using bundling. Default is False.
+            dtiming (boolean, optional):
+                If True, reports solve timing information. Default is False.
+            dis_W (boolean, optional):
+                If True, duals weights (Ws) are disabled before solve, then
+                re-enabled after solve. Default is False.
+            dis_prox (boolean, optional):
+                If True, prox terms are disabled before solve, then
+                re-enabled after solve. Default is False.
+            gripe (boolean, optional):
+                If True, output a message when a solve fails. Default is False.
+            disable_pyomo_signal_handling (boolean, optional):
+                True for asynchronous PH; ignored for persistent solvers.
+                Default False.
+            tee (boolean, optional):
+                If True, displays solver output. Default False.
+            verbose (boolean, optional):
+                If True, displays verbose output. Default False.
+            need_solution (boolean, optional):
+                If True, raises an exception if a solution is not available.
+                Default True
+        """
+        super().solve_loop(
+            solver_options=solver_options,
+            use_scenarios_not_subproblems=use_scenarios_not_subproblems,
+            dtiming=dtiming,
+            dis_W=dis_W,
+            dis_prox=dis_prox,
+            gripe=gripe,
+            disable_pyomo_signal_handling=disable_pyomo_signal_handling,
+            tee=tee,
+            verbose=verbose,
+            need_solution=need_solution,
+        )
+
+        # set self.best_bound_obj_val if we don't have any additional fixed variables
+        if self.can_update_best_bound():
+            self.best_bound_obj_val = self.Ebound(verbose)
+
+
+    def create_fixed_nonant_cache(self):
+        self._initial_fixed_varibles = ComponentSet()
+        for s in self.local_scenarios.values():
+            for v in s._mpisppy_data.nonant_indices.values():
+                if v.fixed:
+                    self._initial_fixed_varibles.add(v)
+
+    def can_update_best_bound(self):
+        for s in self.local_scenarios.values():
+            for v in s._mpisppy_data.nonant_indices.values():
+                if v.fixed:
+                    if v not in self._initial_fixed_varibles:
+                        return False
+        return True

mpisppy/spbase.py

@@ -97,6 +97,9 @@ def __init__(
         self.first_stage_solution_available = False
         self.best_solution_obj_val = None
 
+        # sometimes we know a best bound
+        self.best_bound_obj_val = None
+
         if options.get("toc", True):
             global_toc("Initializing SPBase")
 

mpisppy/utils/cfg_vanilla.py

@@ -18,6 +18,7 @@
 from mpisppy.opt.ph import PH
 from mpisppy.opt.aph import APH
 from mpisppy.opt.lshaped import LShapedMethod
+from mpisppy.opt.subgradient import Subgradient
 from mpisppy.fwph.fwph import FWPH
 from mpisppy.utils.xhat_eval import Xhat_Eval
 import mpisppy.utils.sputils as sputils
@@ -175,6 +176,47 @@ def aph_hub(cfg,
 
     return hub_dict
 
+def subgradient_hub(cfg,
+    scenario_creator,
+    scenario_denouement,
+    all_scenario_names,
+    scenario_creator_kwargs=None,
+    ph_extensions=None,
+    extension_kwargs=None,
+    ph_converger=None,
+    rho_setter=None,
+    variable_probability=None,
+    all_nodenames=None,
+):
+    shoptions = shared_options(cfg)
+    options = copy.deepcopy(shoptions)
+    options["convthresh"] = cfg.intra_hub_conv_thresh
+    options["bundles_per_rank"] = cfg.bundles_per_rank
+    options["smoothed"] = 0
+
+    hub_dict = {
+        "hub_class": PHHub,
+        "hub_kwargs": {"options": {"rel_gap": cfg.rel_gap,
+                                   "abs_gap": cfg.abs_gap,
+                                   "max_stalled_iters": cfg.max_stalled_iters}},
+        "opt_class": Subgradient,
+        "opt_kwargs": {
+            "options": options,
+            "all_scenario_names": all_scenario_names,
+            "scenario_creator": scenario_creator,
+            "scenario_creator_kwargs": scenario_creator_kwargs,
+            "scenario_denouement": scenario_denouement,
+            "rho_setter": rho_setter,
+            "variable_probability": variable_probability,
+            "extensions": ph_extensions,
+            "extension_kwargs": extension_kwargs,
+            "ph_converger": ph_converger,
+            "all_nodenames": all_nodenames
+        }
+    }
+    add_wxbar_read_write(hub_dict, cfg)
+    add_ph_tracking(hub_dict, cfg)
+    return hub_dict
 
 def extension_adder(hub_dict,ext_class):
     # TBD March 2023: this is not really good enough

mpisppy/utils/config.py

@@ -420,6 +420,12 @@ def aph_args(self):
                             domain=float,
                             default=0.01)
 
+    def sub_args(self):
+
+        self.add_to_config(name="SUB",
+                           description="Use subgradient hub instead of PH (default False)",
+                           domain=bool,
+                           default=False)
 
     def fixer_args(self):