[CyIpopt] Automatically use Hessian approximation when grey-box exact Hessian isn't available

[Robbybp]

Summary

When testing some simple grey-box models for #3870, I noticed that we don't automatically set hessian_approximation="limited-memory" when using a grey-box model that doesn't support Hessians. This should be simple enough to do. Without manually setting this option, we can get some very bad behavior.

Steps to reproduce

Here's an example:

import pyomo.environ as pyo
from pyomo.contrib.pynumero.interfaces.external_grey_box import ExternalGreyBoxBlock
from pyomo.contrib.pynumero.examples.external_grey_box.react_example.reactor_model_outputs import (
    ReactorConcentrationsOutputModel,
)
m = pyo.ConcreteModel()
m.reactor = ExternalGreyBoxBlock(external_model=ReactorConcentrationsOutputModel())
m.k1con = pyo.Constraint(expr=m.reactor.inputs["k1"] == 5 / 6)
m.k2con = pyo.Constraint(expr=m.reactor.inputs["k2"] == 5 / 3)
m.k3con = pyo.Constraint(expr=m.reactor.inputs["k3"] == 1 / 6000)
m.cafcon = pyo.Constraint(expr=m.reactor.inputs["caf"] == 10000)
m.obj = pyo.Objective(expr=m.reactor.outputs["cb"], sense=pyo.maximize)
solver = pyo.SolverFactory("cyipopt")
solver.config.options["print_user_options"] = "yes"
results = solver.solve(m, tee=True)
pyo.assert_optimal_termination(results)

output

List of user-set options:

                                    Name   Value                used
                      print_user_options = yes                   yes

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

This is Ipopt version 3.14.11, running with linear solver ma27.

Number of nonzeros in equality constraint Jacobian...:       28
Number of nonzeros in inequality constraint Jacobian.:        0
Number of nonzeros in Lagrangian Hessian.............:        0

Total number of variables............................:        9
                     variables with only lower bounds:        5
                variables with lower and upper bounds:        0
                     variables with only upper bounds:        0
Total number of equality constraints.................:        8
Total number of inequality constraints...............:        0
        inequality constraints with only lower bounds:        0
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:        0

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0 -1.0000000e+00 6.30e+03 2.00e+00  -1.0 0.00e+00    -  0.00e+00 0.00e+00   0
   1 -6.5266413e+02 7.28e+03 4.83e+08  -1.0 2.78e+03    -  3.57e-04 3.75e-01f  1
   2 -1.0559804e+03 3.59e+03 2.59e+06  -1.0 1.94e+03    -  2.38e-01 1.00e+00f  1
   3 -1.3120593e+03 1.82e+03 1.75e+05  -1.0 1.22e+03    -  6.44e-01 1.00e+00f  1
   4 -1.2420075e+03 9.04e+02 1.87e+05  -1.0 1.18e+03    -  1.00e+00 1.00e+00h  1
   5 -1.0199862e+03 3.00e+02 2.03e+05  -1.0 7.59e+02    -  6.56e-01 1.00e+00h  1
   6 -1.1420440e+03 6.28e+02 1.49e+04  -1.0 1.97e+05    -  3.80e-03 4.72e-03f  2
   7 -1.1427067e+03 7.07e+02 2.19e+04  -1.0 6.68e+04    -  1.58e-02 6.47e-03f  3
   8 -1.1757736e+03 8.06e+02 1.48e+04  -1.0 1.44e+06    -  1.12e-03 3.65e-04f  2
   9 -1.1763624e+03 8.88e+02 2.30e+04  -1.0 3.40e+05    -  4.53e-02 1.26e-03f  3
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  10 -1.1902226e+03 9.07e+02 3.42e+03  -1.0 2.73e+06    -  5.93e-04 7.97e-05f  3
  11 -1.2150801e+03 1.07e+03 5.87e+03  -1.0 3.37e+06    -  7.72e-04 2.14e-04f  1
  12 -1.2934436e+03 7.91e+03 1.76e+06  -1.0 4.68e+05    -  2.58e-05 3.80e-03f  1
  13 -1.4090049e+03 8.08e+03 1.80e+06  -1.0 1.54e+07    -  4.00e-04 7.41e-05f  1
  14 -1.4102949e+03 8.08e+03 1.79e+06  -1.0 6.39e+06    -  2.84e-04 1.79e-06f  1
  15 -7.1974579e+02 1.91e+03 7.71e+05  -1.0 3.45e+03    -  1.00e+00 1.00e+00h  1
  16 -1.3094965e+03 2.80e+03 6.16e+05  -1.0 4.47e+03    -  2.70e-01 7.64e-01f  1
  17 -1.3299674e+03 1.44e+03 2.17e+04  -1.0 1.31e+03    -  1.00e+00 8.82e-01f  1
  18 -1.1019077e+03 6.59e+02 2.87e+05  -1.0 1.22e+03    -  6.70e-01 1.00e+00h  1
  19 -8.1471373e+02 3.99e+02 2.63e+05  -1.0 1.04e+03    -  1.76e-01 1.00e+00h  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  20 -9.8896790e+02 7.54e+02 1.42e+04  -1.0 4.25e+04    -  7.35e-04 2.19e-02f  2
  21 -1.0380726e+03 8.15e+02 1.49e+03  -1.0 6.25e+06    -  7.62e-05 7.39e-05f  3
  22 -1.0621291e+03 8.69e+02 1.13e+05  -1.0 5.25e+06    -  2.77e-04 8.52e-05f  3
  23 -1.0637442e+03 8.80e+02 5.00e+04  -1.0 5.17e+05    -  3.32e-03 4.37e-04f  3
  24 -1.0642945e+03 8.80e+02 5.68e+04  -1.0 4.77e+06    -  1.69e-04 5.71e-06f  7
  25 -1.0935404e+03 7.35e+03 2.08e+06  -1.0 3.88e+06    -  3.40e-03 4.48e-04f  1
  26 -1.2733620e+03 7.74e+03 2.82e+06  -1.0 5.48e+07    -  1.37e-05 3.28e-05f  1
  27 -1.2754702e+03 7.74e+03 2.83e+06  -1.0 2.59e+06    -  4.15e-03 6.93e-06f  1
  28 -1.1255570e+03 4.01e+03 9.28e+05  -1.0 1.13e+03    -  1.00e+00 1.00e+00h  1
  29 -1.3240320e+03 2.33e+03 4.16e+05  -1.0 1.37e+03    -  3.62e-01 7.55e-01f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  30 -1.2777896e+03 1.09e+03 2.03e+05  -1.0 1.32e+03    -  1.00e+00 1.00e+00h  1
  31 -9.9704897e+02 4.97e+02 2.92e+05  -1.0 1.13e+03    -  6.48e-01 1.00e+00h  1
  32 -1.1239214e+03 6.42e+02 4.02e+05  -1.0 1.65e+03    -  2.67e-03 1.00e+00f  1
  33 -1.1618303e+03 9.39e+02 3.23e+04  -1.0 1.96e+06    -  8.50e-04 5.27e-04f  1
  34 -1.1912870e+03 1.01e+03 2.07e+05  -1.0 6.30e+04    -  2.67e-04 7.25e-03f  3
  35 -1.1946390e+03 1.09e+03 5.03e+05  -1.0 2.40e+06    -  1.87e-04 1.80e-04f  3
  36 -1.2036953e+03 1.10e+03 3.94e+05  -1.0 1.89e+07    -  9.02e-05 8.36e-06f  3
  37 -1.2082689e+03 1.10e+03 3.16e+05  -1.0 1.90e+07    -  2.36e-04 6.74e-06f  3
  38 -1.2101901e+03 1.11e+03 2.60e+05  -1.0 1.82e+07    -  1.60e-04 5.60e-06f  3
  39 -1.2119471e+03 1.14e+03 1.24e+05  -1.0 8.06e+05    -  2.97e-01 4.00e-04f  1
...

Note the lack of hessian_approximation = limited-memory in the list of user options. This never converges, as far as I can tell. This is particularly bad as we have no indication that anything is wrong other than the solve not converging. By (a) not providing a grey-box Hessian and (b) not using a limited-memory approximation, we are telling IPOPT that these constraints have a Hessian of zero.

Preferred solution

• If a grey-box block is being used and the underlying model doesn't support Hessians, set hessian_approximation="limited-memory". Reset this option after the solve.
• If hessian_approximation="exact" and we have a grey-box block that doesn't support Hessians, log a warning before overriding the option.

[Robbybp]

This will require a has_hessian_support method on PyomoNLPWithGreyBoxBlocks to reuse the Hessian detection logic in _ExternalGreyBoxAsNLP.

[mrmundt]

@Robbybp - .... Use the template, Robby. Don't delete it. 😜