Merge branch 'main' into scaling-uninitialized

.github/workflows/test_branches.yml

@@ -519,7 +519,7 @@ jobs:
         $BARON_DIR = "${env:TPL_DIR}/baron"
         echo "$BARON_DIR" | `
             Out-File -FilePath $env:GITHUB_PATH -Encoding utf8 -Append
-        $URL = "https://www.minlp.com/downloads/xecs/baron/current/"
+        $URL = "https://minlp.com/downloads/xecs/baron/current/"
         if ( "${{matrix.TARGET}}" -eq "win" ) {
             $INSTALLER = "${env:DOWNLOAD_DIR}/baron_install.exe"
             $URL += "baron-win64.exe"

.github/workflows/test_pr_and_main.yml

@@ -562,7 +562,7 @@ jobs:
         $BARON_DIR = "${env:TPL_DIR}/baron"
         echo "$BARON_DIR" | `
             Out-File -FilePath $env:GITHUB_PATH -Encoding utf8 -Append
-        $URL = "https://www.minlp.com/downloads/xecs/baron/current/"
+        $URL = "https://minlp.com/downloads/xecs/baron/current/"
         if ( "${{matrix.TARGET}}" -eq "win" ) {
             $INSTALLER = "${env:DOWNLOAD_DIR}/baron_install.exe"
             $URL += "baron-win64.exe"

pyomo/contrib/appsi/solvers/maingo.py

@@ -57,33 +57,13 @@
 from pyomo.repn.util import valid_expr_ctypes_minlp
 
 
-def _import_SolverModel():
-    try:
-        from . import maingo_solvermodel
-    except ImportError:
-        raise
-    return maingo_solvermodel
-
-
-maingo_solvermodel, solvermodel_available = attempt_import(
-    "maingo_solvermodel", importer=_import_SolverModel
-)
-
-MaingoVar = namedtuple("MaingoVar", "type name lb ub init")
-
 logger = logging.getLogger(__name__)
-
-
-def _import_maingopy():
-    try:
-        import maingopy
-    except ImportError:
-        MAiNGO._available = MAiNGO.Availability.NotFound
-        raise
-    return maingopy
-
-
-maingopy, maingopy_available = attempt_import("maingopy", importer=_import_maingopy)
+MaingoVar = namedtuple("MaingoVar", "type name lb ub init")
+maingopy, maingopy_available = attempt_import("maingopy")
+# Note that importing maingo_solvermodel will trigger the import of
+# maingopy, so we defer that import using attempt_import (which will
+# always succeed, even if maingopy is not available)
+maingo_solvermodel = attempt_import("pyomo.contrib.appsi.solvers.maingo_solvermodel")[0]
 
 
 class MAiNGOConfig(MIPSolverConfig):
@@ -185,9 +165,11 @@ def __init__(self, only_child_vars=False):
         self._last_results_object: Optional[MAiNGOResults] = None
 
     def available(self):
-        if not maingopy_available:
-            return self.Availability.NotFound
-        self._available = True
+        if self._available is None:
+            if maingopy_available:
+                MAiNGO._available = True
+            else:
+                MAiNGO._available = MAiNGO.Availability.NotFound
         return self._available
 
     def version(self):

pyomo/contrib/appsi/solvers/maingo_solvermodel.py

@@ -28,15 +28,7 @@
 from pyomo.repn.util import valid_expr_ctypes_minlp
 
 
-def _import_maingopy():
-    try:
-        import maingopy
-    except ImportError:
-        raise
-    return maingopy
-
-
-maingopy, maingopy_available = attempt_import("maingopy", importer=_import_maingopy)
+maingopy, maingopy_available = attempt_import("maingopy")
 
 _plusMinusOne = {1, -1}
 
@@ -219,7 +211,7 @@ def _linear_to_maingo(self, node):
         return sum(values)
 
 
-class SolverModel(maingopy.MAiNGOmodel):
+class SolverModel(maingopy.MAiNGOmodel if maingopy_available else object):
     def __init__(self, var_list, objective, con_list, idmap, logger):
         maingopy.MAiNGOmodel.__init__(self)
         self._var_list = var_list

pyomo/contrib/incidence_analysis/scc_solver.py

@@ -66,7 +66,14 @@ def generate_strongly_connected_components(
             )
         )
 
-    assert len(variables) == len(constraints)
+    if len(variables) != len(constraints):
+        nvar = len(variables)
+        ncon = len(constraints)
+        raise RuntimeError(
+            "generate_strongly_connected_components only supports systems with the"
+            f" same numbers of variables and equality constraints. Got {nvar}"
+            f" variables and {ncon} constraints."
+        )
     if igraph is None:
         igraph = IncidenceGraphInterface()
 
@@ -78,6 +85,8 @@ def generate_strongly_connected_components(
         subsets, include_fixed=include_fixed
     ):
         # TODO: How does len scale for reference-to-list?
+        # If this assert fails, it may be due to a bug in block_triangularize
+        # or generate_subsystem_block.
         assert len(block.vars) == len(block.cons)
         yield (block, inputs)
 

pyomo/contrib/incidence_analysis/tests/test_scc_solver.py

@@ -501,5 +501,22 @@ def test_with_inequalities(self):
         self.assertEqual(m.x[3].value, 1.0)
 
 
+@unittest.skipUnless(scipy_available, "SciPy is not available")
+@unittest.skipUnless(networkx_available, "NetworkX is not available")
+class TestExceptions(unittest.TestCase):
+    def test_nonsquare_system(self):
+        m = pyo.ConcreteModel()
+        m.x = pyo.Var([1, 2], initialize=1)
+        m.eq = pyo.Constraint(expr=m.x[1] + m.x[2] == 1)
+
+        msg = "Got 2 variables and 1 constraints"
+        with self.assertRaisesRegex(RuntimeError, msg):
+            list(
+                generate_strongly_connected_components(
+                    constraints=[m.eq], variables=[m.x[1], m.x[2]]
+                )
+            )
+
+
 if __name__ == "__main__":
     unittest.main()

pyomo/core/kernel/conic.py

@@ -632,7 +632,7 @@ def as_domain(cls, r, x):
         b = block()
         b.r = variable_tuple([variable(lb=0) for i in range(len(r))])
         b.x = variable()
-        b.c = _build_linking_constraints(list(r) + [x], list(b.r) + [x])
+        b.c = _build_linking_constraints(list(r) + [x], list(b.r) + [b.x])
         b.q = cls(r=b.r, x=b.x)
         return b
 
@@ -934,7 +934,7 @@ def as_domain(cls, r, x):
         b = block()
         b.r = variable_tuple([variable(lb=0) for i in range(len(r))])
         b.x = variable()
-        b.c = _build_linking_constraints(list(r) + [x], list(b.r) + [x])
+        b.c = _build_linking_constraints(list(r) + [x], list(b.r) + [b.x])
         b.q = cls(r=b.r, x=b.x)
         return b
 

pyomo/core/tests/unit/kernel/test_conic.py

@@ -35,6 +35,8 @@
     primal_power,
     dual_exponential,
     dual_power,
+    primal_geomean,
+    dual_geomean,
 )
 
 
@@ -784,6 +786,40 @@ def test_as_domain(self):
         x[1].value = None
 
 
+# These mosek 10 constraints can't be evaluated, pprinted, checked for convexity,
+# pickled, etc., so I won't use the _conic_tester_base for them
+class Test_primal_geomean(unittest.TestCase):
+    def test_as_domain(self):
+        b = primal_geomean.as_domain(r=[2, 3], x=6)
+        self.assertIs(type(b), block)
+        self.assertIs(type(b.q), primal_geomean)
+        self.assertIs(type(b.r), variable_tuple)
+        self.assertIs(type(b.x), variable)
+        self.assertIs(type(b.c), constraint_tuple)
+        self.assertExpressionsEqual(b.c[0].body, b.r[0])
+        self.assertExpressionsEqual(b.c[0].rhs, 2)
+        self.assertExpressionsEqual(b.c[1].body, b.r[1])
+        self.assertExpressionsEqual(b.c[1].rhs, 3)
+        self.assertExpressionsEqual(b.c[2].body, b.x)
+        self.assertExpressionsEqual(b.c[2].rhs, 6)
+
+
+class Test_dual_geomean(unittest.TestCase):
+    def test_as_domain(self):
+        b = dual_geomean.as_domain(r=[2, 3], x=6)
+        self.assertIs(type(b), block)
+        self.assertIs(type(b.q), dual_geomean)
+        self.assertIs(type(b.r), variable_tuple)
+        self.assertIs(type(b.x), variable)
+        self.assertIs(type(b.c), constraint_tuple)
+        self.assertExpressionsEqual(b.c[0].body, b.r[0])
+        self.assertExpressionsEqual(b.c[0].rhs, 2)
+        self.assertExpressionsEqual(b.c[1].body, b.r[1])
+        self.assertExpressionsEqual(b.c[1].rhs, 3)
+        self.assertExpressionsEqual(b.c[2].body, b.x)
+        self.assertExpressionsEqual(b.c[2].rhs, 6)
+
+
 class TestMisc(unittest.TestCase):
     def test_build_linking_constraints(self):
         c = _build_linking_constraints([], [])

pyomo/neos/__init__.py

@@ -30,7 +30,7 @@
     'minos': 'SLC NLP solver',
     'minto': 'MILP solver',
     'mosek': 'Interior point NLP solver',
-    'octeract': 'Deterministic global MINLP solver',
+    #'octeract': 'Deterministic global MINLP solver',
     'ooqp': 'Convex QP solver',
     'path': 'Nonlinear MCP solver',
     'snopt': 'SQP NLP solver',

pyomo/neos/tests/test_neos.py

@@ -79,6 +79,9 @@ def test_doc(self):
         doc = pyomo.neos.doc
         dockeys = set(doc.keys())
 
+        # Octeract interface is disabled, see #3321
+        amplsolvers.remove('octeract')
+
         self.assertEqual(amplsolvers, dockeys)
 
         # gamssolvers = set(v[0].lower() for v in tmp if v[1]=='GAMS')
@@ -149,8 +152,11 @@ def test_minto(self):
     def test_mosek(self):
         self._run('mosek')
 
-    def test_octeract(self):
-        self._run('octeract')
+    # [16 Jul 24] Octeract is erroring.  We will disable the interface
+    # (and testing) until we have time to resolve #3321
+    #
+    # def test_octeract(self):
+    #     self._run('octeract')
 
     def test_ooqp(self):
         if self.sense == pyo.maximize:

pyomo/repn/plugins/nl_writer.py

@@ -1988,6 +1988,18 @@ def _record_named_expression_usage(self, named_exprs, src, comp_type):
             elif info[comp_type] != src:
                 info[comp_type] = 0
 
+    def _resolve_subexpression_args(self, nl, args):
+        final_args = []
+        for arg in args:
+            if arg in self.var_id_to_nl_map:
+                final_args.append(self.var_id_to_nl_map[arg])
+            else:
+                _nl, _ids, _ = self.subexpression_cache[arg][1].compile_repn(
+                    self.visitor
+                )
+                final_args.append(self._resolve_subexpression_args(_nl, _ids))
+        return nl % tuple(final_args)
+
     def _write_nl_expression(self, repn, include_const):
         # Note that repn.mult should always be 1 (the AMPLRepn was
         # compiled before this point).  Omitting the assertion for
@@ -2007,18 +2019,7 @@ def _write_nl_expression(self, repn, include_const):
                     nl % tuple(map(self.var_id_to_nl_map.__getitem__, args))
                 )
             except KeyError:
-                final_args = []
-                for arg in args:
-                    if arg in self.var_id_to_nl_map:
-                        final_args.append(self.var_id_to_nl_map[arg])
-                    else:
-                        _nl, _ids, _ = self.subexpression_cache[arg][1].compile_repn(
-                            self.visitor
-                        )
-                        final_args.append(
-                            _nl % tuple(map(self.var_id_to_nl_map.__getitem__, _ids))
-                        )
-                self.ostream.write(nl % tuple(final_args))
+                self.ostream.write(self._resolve_subexpression_args(nl, args))
 
         elif include_const:
             self.ostream.write(self.template.const % repn.const)

pyomo/repn/tests/ampl/test_nlv2.py

@@ -2703,3 +2703,72 @@ def test_presolve_check_invalid_monomial_constraints(self):
             r"\(fixed body value 5.0 outside bounds \[10, None\]\)\.",
         ):
             nl_writer.NLWriter().write(m, OUT, linear_presolve=True)
+
+    def test_nested_external_expressions(self):
+        # This tests nested external functions in a single expression
+        DLL = find_GSL()
+        if not DLL:
+            self.skipTest("Could not find the amplgsl.dll library")
+
+        m = ConcreteModel()
+        m.hypot = ExternalFunction(library=DLL, function="gsl_hypot")
+        m.p = Param(initialize=1, mutable=True)
+        m.x = Var(bounds=(None, 3))
+        m.y = Var(bounds=(3, None))
+        m.z = Var(initialize=1)
+        m.o = Objective(expr=m.z**2 * m.hypot(m.z, m.hypot(m.x, m.y)) ** 2)
+        m.c = Constraint(expr=m.x == m.y)
+
+        OUT = io.StringIO()
+        nl_writer.NLWriter().write(
+            m, OUT, symbolic_solver_labels=True, linear_presolve=False
+        )
+        self.assertEqual(
+            *nl_diff(
+                """g3 1 1 0	#problem unknown
+ 3 1 1 0 1	#vars, constraints, objectives, ranges, eqns
+ 0 1 0 0 0 0	#nonlinear constrs, objs; ccons: lin, nonlin, nd, nzlb
+ 0 0	#network constraints: nonlinear, linear
+ 0 3 0	#nonlinear vars in constraints, objectives, both
+ 0 1 0 1	#linear network variables; functions; arith, flags
+ 0 0 0 0 0	#discrete variables: binary, integer, nonlinear (b,c,o)
+ 2 3	#nonzeros in Jacobian, obj. gradient
+ 1 1	#max name lengths: constraints, variables
+ 0 0 0 0 0	#common exprs: b,c,o,c1,o1
+F0 1 -1 gsl_hypot
+C0	#c
+n0
+O0 0	#o
+o2	#*
+o5	#^
+v0	#z
+n2
+o5	#^
+f0 2	#hypot
+v0	#z
+f0 2	#hypot
+v1	#x
+v2	#y
+n2
+x1	#initial guess
+0 1	#z
+r	#1 ranges (rhs's)
+4 0	#c
+b	#3 bounds (on variables)
+3	#z
+1 3	#x
+2 3	#y
+k2	#intermediate Jacobian column lengths
+0
+1
+J0 2	#c
+1 1
+2 -1
+G0 3	#o
+0 0
+1 0
+2 0
+""",
+                OUT.getvalue(),
+            )
+        )

pyomo/solvers/plugins/solvers/xpress_direct.py

@@ -1036,10 +1036,8 @@ def _load_slacks(self, cons_to_load=None):
             if xpress_con in self._range_constraints:
                 ## for xpress, the slack on a range constraint
                 ## is based on the upper bound
-                ## FIXME: This looks like a bug - there is no variable named
-                ## `con` - there is, however, `xpress_con` and `pyomo_con`
-                lb = con.lb
-                ub = con.ub
+                lb = xpress_con.lb
+                ub = xpress_con.ub
                 ub_s = val
                 expr_val = ub - ub_s
                 lb_s = lb - expr_val