# NMinimize violates constraints

I have a rather convoluted objective function to minimize, subject to some constraints. The results violates the constraints. The following example highlights the issue:

ObjectFunc = 1 - 0.000911933 s1 - 0.000911933 s2 - 0.000911933 s3;
Mat = {{1. - 0.0027358 s1, 0. + 0.0192484 s1, 0. + 0.022066 s1}, {0. + 0.0192484 s1, 1. - 0.135427 s1 - 0.138162 s2, 0. - 0.15525 s1 - 0.111996 s2}, {0. + 0.022066 s1, 0. - 0.15525 s1 - 0.111996 s2, 1 - 0.177976 s1 - 0.090786 s2 - 0.0198013 s3}};

NMinimize[{ObjectFunc, Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0}, {s1, s2, s3}]


Specifically, if I evaluate Det[Mat] with the optimised parameters, I get a negative number. Admittedly, it is a very small negative number, but negative nonetheless. I have looked at changing the WorkingPrecision and PrecisionGoal.

Any ideas why NMinimize fails to appropriately constrain the objective function?

You should not expect high precision results using machine precision numbers.

Start by converting all values to exact numbers.

ObjectFunc =
1 - 0.000911933 s1 - 0.000911933 s2 - 0.000911933 s3 //
Rationalize[#, 0] & // Simplify;

Mat = {{1. - 0.0027358 s1, 0. + 0.0192484 s1,
0. + 0.022066 s1}, {0. + 0.0192484 s1, 1. - 0.135427 s1 - 0.138162 s2,
0. - 0.15525 s1 - 0.111996 s2}, {0. + 0.022066 s1,
0. - 0.15525 s1 - 0.111996 s2,
1 - 0.177976 s1 - 0.090786 s2 - 0.0198013 s3}} // Rationalize[#, 0] & //
Simplify;


Unfortunately, there is no solution that satisfies the constraints.

({min, arg} =
Minimize[{ObjectFunc, Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0}, {s1, s2,
s3}])

(* Minimize::natt: The minimum is not attained at any point satisfying the given constraints.

{-∞, {s1 -> Indeterminate, s2 -> Indeterminate,
s3 -> Indeterminate}} *)


Consequently, the closest you can get is an approximate numeric solution on the boundary of the constrained region.

({min, arg} =
NMinimize[{ObjectFunc, Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0}, {s1,
s2, s3}, WorkingPrecision -> 50]) // N

(* {0.953946, {s1 -> 0., s2 -> 0., s3 -> 50.5017}} *)

Det[Mat] /. arg // N

(* -5.91946*10^-25 *)


Increasing the WorkingPrecision will bring the Det closer to zero.

({min, arg} =
NMinimize[{ObjectFunc, Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0}, {s1,
s2, s3}, WorkingPrecision -> 100]) // N

(* {0.953946, {s1 -> 0., s2 -> 0., s3 -> 50.5017}} *)

(Det[Mat] /. arg) // N

(* -1.26523*10^-49 *)

• Sorry @Bob Hanlon , this is not a problem of solution at boundaries, but the fact, that NMinimize simply finds a local minimum. As Minimize says, minimum is at - Infinity with si -> Infinity. Constraints give two distinct 3 dim regions, as log-Plot of it shows. Nminimize gets catched in the small region at si near zero. See rp = RegionPlot3D[ red /. {s1 -> Exp[ls1] - 1, s2 -> Exp[ls2] - 1, s3 -> Exp[ls3] - 1}, {ls1, 0, 15}, {ls2, 0, 15}, {ls3, 0, 15}, PlotPoints -> 50, PlotStyle -> Directive[Opacity[0.7], Yellow], Mesh -> False];  and see next comment. Commented Nov 23, 2020 at 8:46
• of[ls1_, ls2_, ls3_] = ObjectFunc /. {s1 -> Exp[ls1] - 1, s2 -> Exp[ls2] - 1, s3 -> Exp[ls3] - 1}; cp = ContourPlot3D[ of[ls1, ls2, ls3], {ls1, 0, 20}, {ls2, 0, 20}, {ls3, 0, 20}, Contours -> {0.953946, -20, -2000}, ContourStyle -> Blue];  Show[rp, cp]  Commented Nov 23, 2020 at 8:46
• red = Reduce[ Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0, {s1, s2, s3}];  Commented Nov 23, 2020 at 8:55

First, let us convert decimals to rationals.

ObjectFunc = Rationalize[1 - 0.000911933 s1 - 0.000911933 s2 - 0.000911933 s3,10^-20];
Mat = Rationalize[{{1. - 0.0027358 s1, 0. + 0.0192484 s1,
0. + 0.022066 s1}, {0. + 0.0192484 s1,
1. - 0.135427 s1 - 0.138162 s2,
0. - 0.15525 s1 - 0.111996 s2}, {0. + 0.022066 s1,
0. - 0.15525 s1 - 0.111996 s2,
1 - 0.177976 s1 - 0.090786 s2 - 0.0198013 s3}}, 10^-20]


Second, let us find the exact solution of the problem under consideration.

Minimize[{ObjectFunc,Det[Mat] > 0 && s1 > 0 && s2 > 0 && s3 > 0}, {s1, s2, s3}]


Minimize::natt: The minimum is not attained at any point satisfying the given constraints.

(*{-\[Infinity], {s1 -> Indeterminate, s2 -> Indeterminate, s3 -> Indeterminate}}*)


Indeed,

FindInstance[ObjectFunc == -1000 && Det[Mat] > 0 && s1 > 0 && s2 > 0 &&
s3 > 0, {s1, s2, s3}, Reals]
(*{{s1 -> 90996849686/82903, s2 -> 1, s3 -> 37}}*)