Repository: scipy/scipy
Pull request: BUG: optimize.minimize: set `trust-constr` `success=False` when constraints are not satisfied
<!-- 
Thanks for contributing a pull request! Please ensure that
your PR satisfies the checklist before submitting:
http://scipy.github.io/devdocs/dev/contributor/development_workflow.html#checklist-before-submitting-a-pr

Also, please name and describe your PR as you would write a
commit message:
http://scipy.github.io/devdocs/dev/contributor/development_workflow.html#writing-the-commit-message.
However, please only include an issue number in the description, not the title,
and please ensure that any code names containing underscores are enclosed in backticks.

Depending on your changes, you can skip CI operations and save time and energy: 
http://scipy.github.io/devdocs/dev/contributor/continuous_integration.html#skipping

Note that we are a team of volunteers; we appreciate your
patience during the review process.

Again, thanks for contributing!
-->

#### Reference issue
Fixes #18882.

#### What does this implement/fix?
Previously trust-constr returned success=True (i.e. status equal to 1 or 2) even when the constraints were not satisfied.  Now trust-constr returns success=False and status=4 with message="Constraints are not satisfied."  A test is implemented in test_minimized_constrained.py to check success=False when constr_violation>gtol.

#### Additional information
This fix is based in part on the PR created by @ellieLitwack at https://github.com/scipy/scipy/pull/18890#issue-1806356775, but approaches the solution to the bug in a way that does not require a modification to the stop criteria.  The feedback from @andyfaff to this previous PR was also considered in the creation of this new fix.


Problem statement:
BUG: minimize does not satisfy nonlinear constraint even if keep_feasible=True
### Describe your issue.

I am creating a general app that can calculate the probability of failure of engineering systems based on a limit state function.  In its most simple form this requires the solution of a minimization problem with an equality constraint.  The app works on almost all tests but have found that one specifically (very academic) causes problems where the equality constraint is not satisfied. 

There is no error message but the result is not correct.  The result is close to x = [0,0] but it should be x = [0, 4].  The resulting equality constraint has a value equal to 1.0 and not 0.0 as required.  The correct result can be obtained if u0 is set to other values like [0,1] but not when set to [1,0].  My expectation was that the equality constraint would be always satisfied if the minimization was successful but it appears it is not, even when keep_feasible=True.  I can solve the issue by setting u0 to [0,1] but the goal of the app I am developing is to work for a general set of limit state functions (constraint functions) and so would like to set u0 = zeros for all cases.

I know this is a challenging problem because the gradient of both the objective and constraint functions is zero at u0 = [0,0].  But why does minimize also fail to satisfy the equality constraint when u0 = [1,0]?

I am using method="trust-constr" because my general app also requires bounds (linear constraints) on the variables and the objective and constraint functions can return gradient and hessian to speed up the solution in large problems.

### Reproducing Code Example

```python
import numpy as np
from scipy import optimize
from IPython.core.interactiveshell import InteractiveShell

def lsf(u):
    u1, u2 = u

    a, b, c, d = [0.0, 0.0, 3.0, 4.0]

    g = 1.0 + (((u1 - a) ** 2) / c**2 - (u2 - b) ** 2 / d**2)

    return g

def of(u):
    d = np.sum(u**2) # sum of squares
    return d

cons = optimize.NonlinearConstraint(lsf, 0, 0, keep_feasible=True)
u0 = [0.0, 0.0]

res = optimize.minimize(
            of,
            u0,
            method="trust-constr",
            constraints=cons,
        )

res
```


### Error message

```shell

```


### SciPy/NumPy/Python version and system information

```shell
1.11.1 1.25.1 sys.version_info(major=3, minor=11, micro=4, releaselevel='final', serial=0)
{
  "Compilers": {
    "c": {
      "name": "clang",
      "linker": "ld64",
      "version": "15.0.7",
      "commands": "arm64-apple-darwin20.0.0-clang"
    },
    "cython": {
      "name": "cython",
      "linker": "cython",
      "version": "0.29.35",
      "commands": "cython"
    },
    "c++": {
      "name": "clang",
      "linker": "ld64",
      "version": "15.0.7",
      "commands": "arm64-apple-darwin20.0.0-clang++"
    },
    "fortran": {
      "name": "gcc",
      "linker": "ld64",
      "version": "12.2.0",
      "commands": "/Users/runner/miniforge3/conda-bld/scipy-split_1687995684757/_build_env/bin/arm64-apple-darwin20.0.0-gfortran"
    },
    "pythran": {
      "version": "0.13.1",
      "include directory": "../../../_build_env/venv/lib/python3.11/site-packages/pythran"
    }
  },
  "Machine Information": {
    "host": {
      "cpu": "arm64",
      "family": "aarch64",
      "endian": "little",
      "system": "darwin"
    },
    "build": {
      "cpu": "x86_64",
      "family": "x86_64",
      "endian": "little",
      "system": "darwin"
    },
    "cross-compiled": true
  },
  "Build Dependencies": {
    "blas": {
      "name": "blas",
      "found": true,
      "version": "3.9.0",
      "detection method": "pkgconfig",
      "include directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/include",
      "lib directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/lib",
      "openblas configuration": "unknown",
      "pc file directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/lib/pkgconfig"
    },
    "lapack": {
      "name": "lapack",
      "found": true,
      "version": "3.9.0",
      "detection method": "pkgconfig",
      "include directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/include",
      "lib directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/lib",
      "openblas configuration": "unknown",
      "pc file directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/lib/pkgconfig"
    },
    "pybind11": {
      "name": "pybind11",
      "version": "2.10.4",
      "detection method": "pkgconfig",
      "include directory": "/Users/Carrasco/miniconda3/envs/carrasco_env311/include"
    }
  },
  "Python Information": {
    "path": "/Users/Carrasco/miniconda3/envs/carrasco_env311/bin/python",
    "version": "3.11"
  }
}
```

Hints from issue discussion:
Should clarify that the problem has two solutions at [0,4] and [0,-4] since equality constraint is a hyperbola.
Hi Cesar!

This is a great question. There are two things to note here.

First, keep_feasible is ignored when an equality constraint is used. This is described in the [documentation for Nonlinear Constraint](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.NonlinearConstraint.html). I believe the reason for this is that it is not practical to assure that every iteration is feasible for an equality constraint.

Second, your initial condition is on a saddle point. Because of this, the minimizer thinks your initial condition is a solution to the minimization problem. This is the expected behavior, and the solution for this is to change the y value in your initial condition. It also fails for u0 = [1,0] because $\frac{\partial of(x,y)}{\partial y} = 0$.

However, trust-constr incorrectly says "success: True" even though the constraint is not met. I am going to work on fixing this. 
Hello Ellie,
Thanks for the response.
Yes, this is a difficult problem and it is probably why it appears on some journal papers dealing with the topic of the app I am developing.  It is a completely academic problem but made it into my test unit and now I am stuck with it :)
Completely missed the documentation on keep_feasible for equality constraints and it of course makes sense.

Changing the result to "success=False" and maybe a "status=4" indicating that a minimum was found but  the constraint was not met within tolerance would be great.

Will probably modify my app to have the minimize function inside a while loop that first tries zeros as the initial condition but then adds a random disturbance when "success==False" and try a few times.

Thanks again for your response and look forward to the implementation of the "success==False" result.  I am new to python (about three months of experience) and do not feel proficient enough to attempt a truly pythonic modification of the code on my own.
Happy to help! This change is not going to show up until the next release of scipy, which will probably not be for a little while. In the mean time, you can use
```python 
success = res.success and (res.constr_violation < 1e-8)
```
as your success criterion.

Commit message: BUG/TST: optimize.minimize: fix trust-constr success=True when constraints are not satisfied
 
Reference issue
Fixes #18882.
 
What does this implement/fix?
Previously trust-constr returned success=True (i.e. status equal to 1 or 2) even when the constraints were not satisfied.  Now trust-constr returns success=False and status=4 with message="Constraints are not satisfied." A test is implemented in test_minimized_constrained.py to check success=False when constr_violation>gtol.
 
Additional information
This fix is based in part on the PR created by @ellieLitwack at https://github.com/scipy/scipy/pull/18890#issue-1806356775, but approaches the solution to the bug in a way that does not require a modification to the stop criteria.  The feedback from @andyfaff to this previous PR was also considered in the creation of this new fix.

Patch to review:
diff --git a/scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py b/scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py
index 0e1dd3746595..8b49140131c3 100644
--- a/scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py
+++ b/scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py
@@ -18,7 +18,8 @@
     0: "The maximum number of function evaluations is exceeded.",
     1: "`gtol` termination condition is satisfied.",
     2: "`xtol` termination condition is satisfied.",
-    3: "`callback` function requested termination."
+    3: "`callback` function requested termination.",
+    4: "Constraints are not satisfied."
 }
 
 
@@ -292,13 +293,14 @@ def _minimize_trustregion_constr(fun, x0, args, grad,
         Total execution time.
     message : str
         Termination message.
-    status : {0, 1, 2, 3}
+    status : {0, 1, 2, 3, 4}
         Termination status:
 
         * 0 : The maximum number of function evaluations is exceeded.
         * 1 : `gtol` termination condition is satisfied.
         * 2 : `xtol` termination condition is satisfied.
         * 3 : `callback` function requested termination.
+        * 4 : Constraints are not satisfied.
 
     cg_stop_cond : int
         Reason for CG subproblem termination at the last iteration:
@@ -540,6 +542,10 @@ def grad_and_jac(x):
             initial_constr_penalty, initial_tr_radius,
             factorization_method)
 
+    # Status 4 occurs when minimize is successful but constraints are not satisfied.
+    if result.status in (1, 2) and state.constr_violation > gtol:
+        result.status = 4
+
     # Status 3 occurs when the callback function requests termination,
     # this is assumed to not be a success.
     result.success = True if result.status in (1, 2) else False
diff --git a/scipy/optimize/tests/test_minimize_constrained.py b/scipy/optimize/tests/test_minimize_constrained.py
index ee700ec6e959..1f425b0442ce 100644
--- a/scipy/optimize/tests/test_minimize_constrained.py
+++ b/scipy/optimize/tests/test_minimize_constrained.py
@@ -763,6 +763,25 @@ def test_gh20665_too_many_constraints():
         minimize(rosen, x0, method='trust-constr', constraints=[g],
                  options={'factorization_method': 'SVDFactorization'})
 
+def test_issue_18882():
+    def lsf(u):
+        u1, u2 = u
+        a, b = [3.0, 4.0]
+        return 1.0 + u1**2 / a**2 - u2**2 / b**2
+
+    def of(u):
+        return np.sum(u**2)
+
+    with suppress_warnings() as sup:
+        sup.filter(UserWarning, "delta_grad == 0.0")
+        sup.filter(UserWarning, "Singular Jacobian matrix.")
+        res = minimize(
+            of,
+            [0.0, 0.0],
+            method="trust-constr",
+            constraints=NonlinearConstraint(lsf, 0, 0),
+        )
+    assert (not res.success) and (res.constr_violation > 1e-8)  
 
 class TestBoundedNelderMead:
 


Code review comments:
diff_hunk: @@ -18,7 +18,8 @@
     0: "The maximum number of function evaluations is exceeded.",
     1: "`gtol` termination condition is satisfied.",
     2: "`xtol` termination condition is satisfied.",
-    3: "`callback` function requested termination."
+    3: "`callback` function requested termination.",
+    4: "Constraints are not satisfied."
path: scipy/optimize/_trustregion_constr/minimize_trustregion_constr.py
text: @user1:
I think the status message can be improved here, currently it doesn't give enough information as to why the status was set to 4.

How about `Constraint violation exceeds 'gtol'`


@author:
I think you are right.  Your status message is more clear.  I'll go ahead and update it.

@author:
Quick question.  Should I rebase before attempting another commit? This is my first PR and am not familiar with the process to update an existing commit.

@user2:
Feel free to just commit the change on top of your current branch - a reviewer can always 'squash' commits down before merge.

@author:
Thanks, I just added a new commit with the updated status message.