my nonlinear constraints HAVE to be satisfied, what to do next?

I am working on a sophisticated anisotropic tomographic modeling project. I recently learned that the fmincon function in MATLAB does not always honor nonlinear constraints. I am aware that there are codes on GitHub that claim to always honor nonlinear constraints, but I do not have much confidence in GitHub code. My nonlinear constraints are quite complicated: they are 4th order polynomials with 4 variables for my current work. In the future, I may have 15th order polynomials with 20 variables.

I have two questions:

1. Is there any MATLAB program that always honors nonlinear constraints?
2. If there is no MATLAB program that always honors nonlinear constraints, I may have to write my own code that will always honor nonlinear constraints. I am thinking of boxing or bounding the mth order polynomials with n variables first (of course, this is a complicated task), so that I turn the nonlinear constraints into bound constraints next. Then, I could apply fmincon. Do you have any suggestions?

Thank you very much for your help!

[Walter Roberson]

In the future, I may have 15th order polynomials with 20 variables.

Yikes! Unless the range of inputs is -1 to +1, 15th order polynomials are likely to involve a lot of numeric noise!

4th order polynomials are not so bad, but establishing bounds for them can be sensitive to numeric noise, especially at the place where the branches "cross".

[Jon]

When you say "always honors nonlinear constraints", do you mean are the constraints always enforced as it iterates to a solution (that is the path to the solution is always feasible) or do you mean does the algorithm always converge to a feasible (non-linear constraints honored) solution?

I think either way there are options within MATLAB that will do this.

[Torsten]

You wrote in another post about the consequences if the constraints are not satisfied:

The problem is not physical. For example, (1) energy becomes negative; (2) Snell's law is violated.

But does it really matter if during the iteration process, energy becomes negative and/or Snell's law is violated ? My guess is it would be adjusted a few iterations later.

[Frank]

To Jon: Thank you very much for your reply. I meant both. I was told in another post that MATLAB can't do this in fmincon.

[Frank]

To Torsten: Thank you very much! It matters. If energy becomes negative or Snell's law is violated, that means I can't finish computing the objective function (In other words, I have made the program returns if energy becomes negative or Snell's law is violated).

[Walter Roberson]

@Jon

fmincon does not check nonlinear constraints before running the objective function for every set of model parameters. (For one thing, it does not check them during the initial estimate of the jacobian, if I recall correctly.)

The poster is looking for a way to ensure that every trial model parameters passed to the objective function already satisfy all constraints.

[Frank]

To Walter: Thanks! For now, I don't think that is a big issue. I can always loose the bounds if I am stuck. I am not seeking a completely rigorous solution.

[Torsten]

@Frank

To Torsten: Thank you very much! It matters. If energy becomes negative or Snell's law is violated, that means I can't finish computing the objective function.

Then first project the infeasible points onto the feasible region defined by your constraints and evaluate your objective function with the result obtained after this projection:

https://de.mathworks.com/matlabcentral/answers/1921445-fmincon-any-way-to-enforce-linear-inequality-constraints-at-intermediate-iterations?s_tid=srchtitle

[Frank]

@Torsten Thank you very much! Could you please tell me what do you mean "project"?

[Torsten]

Project means:

Solve the subproblem of finding the nearest feasible point to the one suggested by fmincon:

min: sum((x-x_suggested_by_fmincon).^2)

s.t.

your constraints on x.

in each call of your objective function.

But maybe evaluating your constraints gives similar problems as evaluating your objective for infeasible points - I don't know.

I think Matt J has more experience than I with optimization (see his answer below and the link I included to a similar discussion).

Answer 1

[Matt J]

(1) It is not possible to enforce equality constraints (e.g. Snell's Law) exactly , whether linear or non-linear. This is true both during the iterative process or at convergence. This is true for any optimization software app you might pursue. It is a fundamental limitation of finite precision digital computing.

(2) Inequality constraints, whether linear or nonlinear can be enforced at convergence by introducing slack, e.g., instead of c(x)<=0 impose c(x)<=-Delta<0. Set the ConstraintTolerance parameter to be less than Delta. Note: This will not work if your inequality constrained region has no interior, e.g., the constraint x^2<=0, but in this case, it means your inequality constraint should really be written as an equality constraint anyway.

(3) If for some reason, the computation of the objective or nonlinear constraint function cannot be completed for a certain input x, the recommended approach is to just return NaN or Inf. Certain fmincon algorithms like sqp or interior-point can recover from NaN or Inf values encountered during the iterative search.

[Frank]

@Matt J Thank you very much! I Googled but didn't find the answer. I am quite familiar with classical optimization algorithms. Unfortunately, fmincon source code is too complicated for me to understand. Could you please tell me further how fmincon interior-point algorithm recover from NaN or Inf values?

[Matt J]

I don't know how the recovery mechanism works, or if there is any documentation describing it. However, it isn't something you're supposed to need to worry about as the user. It's handled under the hood.

[Walter Roberson]

(Note: some optimization algorithms can recover from inf but not nan)

Answer 2

[Walter Roberson]

If you use surrogateopt then there is no explicit non-linear equality of inequality function parameter. However, you can have the function return a struct in which the field Fval is to be minimized, and the function will try to ensure that the field Ineq <= 0 .

The idea is that the function will be passed a set of model parameters that satisfy the bounds and linear equalities and inequalities, and the function will be responsible for testing nonlinear inequalities itself. If they pass, then compute a valid Fval; if they do not pass then return a high value such as Inf. So you would never compute with the locations that do not pass, because you tested for yourself.

[Frank]

@Walter Roberson Thanks! I believe surrogateopt is a derivative-free optimization algorithm. My algorithm requires derivatives.