Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

FIX: allow for print_level greater than 1 and ENH: allow warm start for dual variables #239

Closed
wants to merge 10 commits into from
38 changes: 34 additions & 4 deletions cyipopt/scipy_interface.py
Original file line number Diff line number Diff line change
Expand Up @@ -385,7 +385,10 @@ def minimize_ipopt(fun,
constraints=(),
tol=None,
callback=None,
options=None):
options=None,
mult_g=[],
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Note that using mutable default arguments in function signatures can cause bugs if the function is called more than once. See the "Mutable Default Arguments" section of this article for an explanation and suggestion on how to avoid. A simple way to avoid this foot gun is to use None as the default arguments in the function signature and then manually create a new instance of an empty list within the function body if the value of the variable in question is None.

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for your patience! I did neglect this.

mult_x_L=[],
mult_x_U=[],):
"""
Minimization using Ipopt with an interface like
:py:func:`scipy.optimize.minimize`.
Expand Down Expand Up @@ -487,6 +490,16 @@ def minimize_ipopt(fun,
This argument is ignored by the default `method` (Ipopt).
If `method` is one of the SciPy methods, this is a callable that is
called once per iteration. See [2]_ for details.
mult_g : list, optional
Initial guess for the Lagrange multipliers of the constraints. A list
of real elements of length ``m``, where ``m`` is the number of
constraints.
mult_x_L : list, optional
Initial guess for the Lagrange multipliers of the lower bounds on the
variables. A list of real elements of length ``n``.
mult_x_U : list, optional
Initial guess for the Lagrange multipliers of the upper bounds on the
variables. A list of real elements of length ``n``.

References
----------
Expand Down Expand Up @@ -596,8 +609,11 @@ def minimize_ipopt(fun,
# Rename some default scipy options
replace_option(options, b'disp', b'print_level')
replace_option(options, b'maxiter', b'max_iter')
if getattr(options, 'print_level', False) is True:
options[b'print_level'] = 1
if b'print_level' in options:
if options[b'print_level'] is True:
options[b'print_level'] = 1
elif options[b'print_level'] is False:
options[b'print_level'] = 0
Comment on lines -599 to +616
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This change doesn't appear to make any difference to the way the code functions, except increasing the number of lines of code to do the same thing. Am I missing something? If so, please can you explain why this change is needed?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ipopt print levels are integers between 1 and 12, with the default being 5: https://coin-or.github.io/Ipopt/OPTIONS.html#OPT_print_level

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we should use something like the following. This will ensure that if a print_level of True is used that the default print level is used, not breaking previous code. It will also allow users to manually select the exact desired print level as per the Ipopt options.

Suggested change
if getattr(options, 'print_level', False) is True:
options[b'print_level'] = 1
if b'print_level' in options:
if options[b'print_level'] is True:
options[b'print_level'] = 1
elif options[b'print_level'] is False:
options[b'print_level'] = 0
if getattr(options, 'print_level') is True:
options[b'print_level'] = 5
else:
print_level = getattr(options, 'print_level' , 5)
if not isinstance(print_level, int) or (print_level < 1) or (print_level > 12):
raise ValueError(f'`print_level` option must be an int between 1 and 12, not {print_level}')
options[b'print_level'] = print_level

Copy link
Author

@zhengzl18 zhengzl18 Nov 5, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought it this way: The original version would force the print_level to be 1 if the user has specified this parameter, no matter what value is chosen. I tried to make it work like this:

  1. 'disp': True -> 'print_level': 1
  2. 'disp': False -> 'print_level': 0
  3. no 'disp' or 'print_level' -> 'print_level': 0 (according to Ipopt's document, a print_level of value 0 is acceptable)
  4. 'disp'/'print_level': 0-12 -> 'print_level': 0-12

About the default print_level in the second case, I did have some doubts. I wondered before why Ipopt said that print_level was by default 5, but that was not how it actually behavioured. So when I was changing these codes, I hesitated about whether to set the default value to 1 or 5. I chose 1 eventually because it was the original setting.

else:
options[b'print_level'] = 0
if b'tol' not in options:
Expand All @@ -614,7 +630,11 @@ def minimize_ipopt(fun,
msg = 'Invalid option for IPOPT: {0}: {1} (Original message: "{2}")'
raise TypeError(msg.format(option, value, e))

x, info = nlp.solve(x0)
_dual_initial_guess_validation("mult_g", mult_g, len(cl))
_dual_initial_guess_validation("mult_x_L", mult_x_L, len(lb))
_dual_initial_guess_validation("mult_x_U", mult_x_U, len(ub))

x, info = nlp.solve(x0, lagrange=mult_g, zl=mult_x_L, zu=mult_x_U)

return OptimizeResult(x=x,
success=info['status'] == 0,
Expand Down Expand Up @@ -689,3 +709,13 @@ def _minimize_ipopt_iv(fun, x0, args, kwargs, method, jac, hess, hessp,

return (fun, x0, args, kwargs, method, jac, hess, hessp,
bounds, constraints, tol, callback, options)

def _dual_initial_guess_validation(dual_name: str, dual_value: list, length: int):
if not isinstance(dual_value, list):
raise TypeError(f'`{dual_name}` must be a list.')
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why does this have to be a list? Couldn't any iterable or sequence work?

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since Problem.solve() uses empty lists as the default values for lagrange, zl and zu, and create zero numpy arrays when they are [], I just followed it. Indeed any iterable or sequence would work.

if lagrange == []:

if len(dual_value) > 0:
assert all(isinstance(x, (int, float)) for x in dual_value), \
f'All elements of `{dual_name}` must be numeric.'
if len(dual_value) != length:
raise ValueError(f'`{dual_name}` must be empty or have length '
f'`{length}`.')
148 changes: 148 additions & 0 deletions cyipopt/tests/unit/test_dual_warm_start.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,148 @@
import pytest
import numpy as np
from numpy.testing import assert_, assert_allclose
from cyipopt import minimize_ipopt


class TestDualWarmStart:
atol = 1e-7

def setup_method(self):
self.opts = {'disp': False}

def fun(self, d, sign=1.0):
"""
Arguments:
d - A list of two elements, where d[0] represents x and d[1] represents y
in the following equation.
sign - A multiplier for f. Since we want to optimize it, and the SciPy
optimizers can only minimize functions, we need to multiply it by
-1 to achieve the desired solution
Returns:
2*x*y + 2*x - x**2 - 2*y**2

"""
x = d[0]
y = d[1]
return sign*(2*x*y + 2*x - x**2 - 2*y**2)

def jac(self, d, sign=1.0):
"""
This is the derivative of fun, returning a NumPy array
representing df/dx and df/dy.

"""
x = d[0]
y = d[1]
dfdx = sign*(-2*x + 2*y + 2)
dfdy = sign*(2*x - 4*y)
return np.array([dfdx, dfdy], float)

def f_eqcon(self, x, sign=1.0):
""" Equality constraint """
return np.array([x[0] - x[1]])

def f_ieqcon(self, x, sign=1.0):
""" Inequality constraint """
return np.array([x[0] - x[1] - 1.0])

def f_ieqcon2(self, x):
""" Vector inequality constraint """
return np.asarray(x)

def fprime_ieqcon2(self, x):
""" Vector inequality constraint, derivative """
return np.identity(x.shape[0])

# minimize
def test_dual_warm_start_unconstrained_without(self):
# unconstrained, without warm start.
res = minimize_ipopt(self.fun, [-1.0, 1.0], args=(-1.0, ),
jac=self.jac, method=None, options=self.opts)
assert_(res['success'], res['message'])
assert_allclose(res.x, [2, 1])

@pytest.mark.xfail(raises=(ValueError,), reason="Initial guesses for dual variables have wrong shape")
def test_dual_warm_start_unconstrained_with(self):
# unconstrained, with warm start.
res = minimize_ipopt(self.fun, [-1.0, 1.0], args=(-1.0, ),
jac=self.jac, method=None, options=self.opts, mult_g=[1, 1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [2, 1])

def test_dual_warm_start_equality_without(self):
# equality constraint, without warm start.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'eq', 'fun':self.f_eqcon,
'args': (-1.0, )},
options=self.opts)
assert_(res['success'], res['message'])
assert_allclose(res.x, [1, 1])

def test_dual_warm_start_equality_with_right(self):
# equality constraint, with right warm start.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'eq', 'fun':self.f_eqcon,
'args': (-1.0, )},
options=self.opts, mult_g=[1], mult_x_L=[1, 1], mult_x_U=[-1, -1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [1, 1])

@pytest.mark.xfail(raises=(ValueError,), reason="Initial guesses for dual variables have wrong shape")
def test_dual_warm_start_equality_with_wrong_shape(self):
# equality constraint, with wrong warm start shape.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'eq', 'fun':self.f_eqcon,
'args': (-1.0, )},
options=self.opts, mult_g=[1], mult_x_U=[1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [1, 1])

@pytest.mark.xfail(raises=(TypeError,), reason="Initial guesses for dual variables have wrong type")
def test_dual_warm_start_equality_with_wrong_type(self):
# equality constraint, with wrong warm start type.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'eq', 'fun':self.f_eqcon,
'args': (-1.0, )},
options=self.opts, mult_x_L=[1, 1], mult_x_U=np.array([1, 1]))
assert_(res['success'], res['message'])
assert_allclose(res.x, [1, 1])

def test_dual_warm_start_inequality_with_right(self):
# inequality constraint, with right warm start.
res = minimize_ipopt(self.fun, [-1.0, 1.0], method=None,
jac=self.jac, args=(-1.0, ),
constraints={'type': 'ineq',
'fun': self.f_ieqcon,
'args': (-1.0, )},
options=self.opts, mult_x_L=[-1, 1], mult_x_U=[1, 1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [2, 1], atol=1e-3)

@pytest.mark.xfail(raises=(ValueError,), reason="Initial guesses for dual variables have wrong shape")
def test_dual_warm_start_inequality_vec_with_wrong_shape(self):
# vector inequality constraint, with wrong warm start shape.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'ineq',
'fun': self.f_ieqcon2,
'jac': self.fprime_ieqcon2},
options=self.opts, mult_g=[1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [2, 1])

@pytest.mark.xfail(raises=(AssertionError,), reason="Initial guesses for dual variables have wrong type")
def test_dual_warm_start_inequality_vec_with_wrong_element_type(self):
# vector inequality constraint, with wrong warm start element type.
res = minimize_ipopt(self.fun, [-1.0, 1.0], jac=self.jac,
method=None, args=(-1.0,),
constraints={'type': 'ineq',
'fun': self.f_ieqcon2,
'jac': self.fprime_ieqcon2},
options=self.opts, mult_g=['1', 1])
assert_(res['success'], res['message'])
assert_allclose(res.x, [2, 1])