Is there any out of the box newton like solver that allows bounds (and a supplied jacobian and hessian)

Question

Is there any out of the box newton like solver that allows bounds on the input parameters (and a supplied expression for the jacobian and hessian).

I am basically exactly looking for "trust-constr" in SciPy, except that one does not seem to work at all. While it accepts bounds as an argument it seems to just ignore them. See also this question Scipy: How can I use Bounds with trust-constr? about this issue.

---

It seems that the failure of trust-constr to take into account bounds is not completely general. It would probably be nicest to show an instance with something standard like the rosen function but my initial try did not work. So here is a simple counter example

def fake_function(x):
    if x[0]==0.523 and x[1]==1.43:
        return -0.0318285
    print("value at x not available")

def fake_grad(x):
    return [9.21296, -1.98147]

def fake_hess(x):
    return [[-467.451, -98.9485], [-98.9485, 28.6649]]

scipy.optimize.minimize(fake_function, [0.523, 1.43], method='trust-constr', jac=fake_grad, hess=fake_hess,
                                  bounds=[(0.5,0.55),(0.5,2)],
                                  options={'gtol': 1e-5, 'disp': True})

Obviously this optimization doesn't make sense after the first step because I defined the function in by its value at the initial point. However, it is enough to show that at the very first step trust-constr will ignore the bounds and ask for a point outside them. (Now clarified to show that the fake_function is a well defined twice differentiable function that is too big to include in the minimal example. We only need the values at the initial point however to see the algorithm misbehave.)

---

Also changing the bounds a bit I realize that it is not completely ignoring them. It chooses different points based on the bounds however it does not mind stepping outside the bounds which seems to be a terrible feature.

Answer 1

Despite the (imho) terrible setup for demonstrating your problem, consider explicitly enforcing feasiblity (docs):

Code

import numpy as np
from scipy.optimize import minimize, Bounds

def fake_function(x):
    print('EVAL @: ', x)

    if x[0]==0.523 and x[1]==1.43:
        return -0.0318285
    print("value at x not available")

def fake_grad(x):
    return [9.21296, -1.98147]

def fake_hess(x):
    return [[-467.451, -98.9485], [-98.9485, 28.6649]]

# !!!
bounds_explicit_feas = Bounds(np.array([0.5, 0.5]), np.array([0.55, 2.0]), keep_feasible = True)
bounds_original = [(0.5,0.55),(0.5,2)]
# !!!    

minimize(fake_function, [0.523, 1.43], method='trust-constr', jac=fake_grad, hess=fake_hess,
                                  bounds=bounds_explicit_feas,
                                  options={'gtol': 1e-5, 'disp': True})

Original Bounds: Output

EVAL @:  [0.523 1.43 ]
EVAL @:  [0.19537876 1.56399192]
value at x not available

Explicit Feasible: Output

EVAL @:  [0.523 1.43 ]
EVAL @:  [0.52256323 1.49080832]
value at x not available