NMinimize violates constraints

Question

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?

Answer 1

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 *)

Answer 2

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}}*)