1
$\begingroup$ 
g1[a0_, b0_, a1_, b1_, a2_, b2_] :=Piecewise[{{(a0*a1*a2)/(a1*a2*(-b0 + b1)+ a0*(a1*(1 + a2 + a2*b0 - b2) + a2*(-b1 + b2))), 0 <= (a0*a1*a2)/(a0*a1 + a0*a1*a2 - a1*a2*b0 + a0*a1*a2*b0 - a0*a2*b1 + a1*a2*b1 - a0*a1*b2 + a0*a2*b2) <= a2}}, 1]

g2[a0_, b0_, a1_, b1_, a2_, b2_] :=Piecewise[{{((a0*b1 - a1*((-1 + a0)*b0 + b1))^2 - 2*a0*(a1*(b0 - b1) + a0*(a1^2 - a1*b0 + b1))*b2 + a0^2*(1 + a1^2)*b2^2)/(2*(a1*(-b0 + b1) + a0*(a1*(-1 + b0) - b1 + b2))*(a1*(-b0 + b1) + a0*(a1 + a1*b0 - b1 + b2))), a2 <= (a0*a1*(-1 + b2)*(a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2)))/((a1*(-b0 + b1) + a0*(a1*(-1 + b0) - b1 + b2))*
  (a1*(-b0 + b1) + a0*(a1 + a1*b0 - b1 + b2))) <= a1}}, 1]

g3[a0_, b0_, a1_, b1_, a2_, b2_] :=Piecewise[{{(a0*a1*(2*a1*(-b0 + b1)*b2 - a0*(a1 + 2*(-(a1*b0) + b1)*b2 + (-2 + a1)*b2^2)))/(2*(a1*(-b0 + b1) + a0*(a1*(-1 + b0) - b1 + b2))*(a1*(-b0 + b1) + a0*(a1 + a1*b0 - b1 + b2))), a2 <= (a0*a1*(-(a0*a1) + b2*(a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2))))/
 (-(a0^2*a1^2) + (a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2))^2) <= a1}}, 1]

g4[a0_, b0_, a1_, b1_] :=Piecewise[{{((-1 + a0)^2*b0^2 + 2*(-1 + a0)*b0*b1 + b1*(a0^2*(-2 + b1) + b1))/(2*(a0*(-1 + b0) - b0 + b1)*(a0 - b0 + a0*b0 + b1)), a1 <= (a0*(-b0 + a0*b0 + b1 + b0*b1 - a0*b0*b1 - b1^2))/(a0^2 - b0^2 + 2*a0*b0^2 - a0^2*b0^2 + 2*b0*b1 - 2*a0*b0*b1 - b1^2) <= a0}}, 1]

g5[a0_, b0_, a1_, b1_] :=Piecewise[{{(a0*(-a0 + 2*(-1 + a0)*b0*b1 - (-2 + a0)*b1^2))/(2*(a0*(-1 + b0) - b0 + b1)*(a0 - b0 + a0*b0 + b1)), a1 <= (a0*(-a0 - b0*b1 + a0*b0*b1 + b1^2))/(-a0^2 + b0^2 - 2*a0*b0^2 + a0^2*b0^2 - 2*b0*b1 + 2*a0*b0*b1 + b1^2) <= a0}}, 1]

g6[a0_, b0_] :=Piecewise[{{b0/(1 + b0), a0 <= b0/(1 + b0) <= 1}}, 1]

fmymin[a0_, b0_, a1_, b1_, a2_, b2_] :=Min[1/2, (a0^2 - 2*(-1 + a0)*b0 - (-1 + a0)^2*b0^2)/2, ((2 - a0)*a0 + (-1 + a0)^2*b0^2)/2, (a1^2 + 2*b1 - (2*a1*((-1 + a0)*b0 + b1))/a0 - (b1 - (a1*((-1 + a0)*b0 + b1))/a0)^2)/2, (2*a1 - a1^2 + (b1 - (a1*((-1 + a0)*b0 + b1))/a0)^2)/2, (a2^2 + 2*b2 - (2*a2*(a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2)))/(a0*a1) - (b2 - (a2*(a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2)))/(a0*a1))^2)/2, (2*a2 - a2^2 + (b2 - (a2*(a1*((-1 + a0)*b0 + b1) + a0*(-b1 + b2)))/(a0*a1))^2)/2, g1[a0, b0, a1, b1, a2, b2], g2[a0, b0, a1, b1, a2, b2], g3[a0, b0, a1, b1, a2, b2], g4[a0, b0, a1, b1], g5[a0, b0, a1, b1], g6[a0, b0]]

fmyminorg[a0_, b0_, a1_, b1_, a2_, b2_] :=Min[(a0^2 + 2*(1 - a0)*b0 - (1 - a0)^2*b0^2)/2, (2*a0 - a0^2 + (1 - a0)^2*b0^2)/2, (a1*(a0^2*a1 - 2*(-1 + a0)*a0*b0 - (-1 + a0)^2*a1*b0^2) + 2*(a0 - a1)*(a0 + (-1 + a0)*a1*b0)*b1 - (a0 - a1)^2*b1^2)/(2*a0^2), (a1*(-(a0^2*(-2 + a1)) + (-1 + a0)^2*a1*b0^2) + 2*(-1 + a0)*a1*(-a0 + a1)*b0*b1 + (a0 - a1)^2*b1^2)/(2*a0^2), (a0 + a1 - a0*a1 + ((-1 + a0)*b0*(-(a0*b1) + a1*((-1 + a0)*b0 + b1)))/a0)/2, (a0^2*(a1 - a1*b0^2 + b0*(-1 + b1)) - a1*(-1 + b0)*(b0 - b1) + a0*(2*a1*b0^2 + b1 - b0*(-1 + a1 + b1 + a1*b1)))/(2*a0), (a0 + a2 - a0*a2 + ((-1 + a0)*b0*(a1*a2*((-1 + a0)*b0 + b1) - a0*a1*b2 + a0*a2*(-b1 + b2)))/(a0*a1))/2, (a1 + a2 - a1*a2 + ((-(a0*b1) + a1*((-1 + a0)*b0 + b1))*(a1*a2*((-1 + a0)*b0 + b1) - a0*a1*b2 + a0*a2*(-b1 + b2)))/(a0^2*a1))/2, (2 + a2^2 + 2*a2*((b0 - a0*b0 - b1)/a0 + (b1 - b2)/a1 + (-1 + b2)/a2) - (a1*a2*((-1 + a0)*b0 + b1) - a0*a1*b2 + a0*a2*(-b1 + b2))^2/(a0^2*a1^2))/2, (2*a2 - a2^2 + (a1*a2*((-1 + a0)*b0 + b1) - a0*a1*b2 + a0*a2*(-b1 + b2))^2/(a0^2*a1^2))/2, (-(a1^2*a2*(b0 - b1)^2) + a0*a1*(b0 - b1)*(a1 + a2 + 2*a1*a2*b0 - 2*a2*b1 - a1*b2 + a2*b2) + a0^2*(a1^2*(a2 - a2*b0^2 + b0*(-1 + b2)) - a2*(-1 + b1)*(b1 - b2) + a1*(b1 + a2*b0*(-1 + 2*b1 - b2) + b2 - b1*b2)))/(2*a0^2*a1), (-(a1*a2*(-1 + b0)*(b0 - b1)) + a0^2*(a1*(a2 - a2*b0^2 + b0*(-1 + b2)) + a2*b0*(b1 - b2)) + a0*(-(a2*(-1 + b0)*(b1 - b2)) + a1*(2*a2*b0^2 + b2 - b0*(-1 + a2 + a2*b1 + b2))))/(2*a0*a1)]

max = NMaximize[{fmymin[a0, b0, a1, b1, a2, b2] - fmyminorg[a0, b0, a1, b1, a2, b2], {a1 <= a0 <= 1, 0 <= b0 <= 1, a2 <= a1 <= a0, b0 <= b1 <= b2, 0 <= a2 <= 1, b1 <= b2 <= 1}}, {a0, b0, a1, b1, a2, b2}]
(* {0.0218682, {a0 -> 0.612102, b0 -> 0.0022207, a1 -> 0.144964, 
                b1 -> 0.211932, a2 -> 0.144964, b2 -> 0.545072}} *)

max = NMaximize[{fmymin[a0, b0, a1, b1, a2, b2] - fmyminorg[a0, b0, a1, b1, a2, b2], {0 <= a0 <= 1, 0 <= b0 <= 1, a2 <= a1 <= a0, b0 <= b1 <= b2, 0 <= a2 <= a1, b1 <= b2 <= 1}}, {a0, b0, a1, b1, a2, b2}]
(* {0.0231773, {a0 -> 0.705064, b0 -> 0., a1 -> 0.20014, b1 -> 0.306156, 
                a2 -> 0.00654739, b2 -> 0.842462}} *)

max = NMaximize[{fmymin[a0, b0, a1, b1, a2, b2] - fmyminorg[a0, b0, a1, b1, a2, b2], {0 <= a0 <= 1, 0 <= b0 <= b1, a2 <= a1 <= a0, b0 <= b1 <= b2, 0 <= a2 <= a1, 0 <= b2 <= 1}}, {a0, b0, a1, b1, a2, b2}]
(* {0.0232671, {a0 -> 0.691563, b0 -> 0., a1 -> 0.192949, b1 -> 0.298011,
                a2 -> 0.144842, b2 -> 0.95349}} *)

max = NMaximize[{fmymin[a0, b0, a1, b1, a2, b2] - fmyminorg[a0, b0, a1, b1, a2, b2], {a1 <= a0 <= 1, 0 <= b0 <= b1, a2 <= a1 <= a0, b0 <= b1 <= b2, 0 <= a2 <= 1, 0 <= b2 <= 1}}, {a0,b0, a1, b1, a2, b2}]
(* {0.0216001, {a0 -> 0.620419, b0 -> 0.0905576, a1 -> 0.212008, 
                b1 -> 0.290074, a2 -> 5.25103*10^-9, b2 -> 0.679186}} *)

Just as a side note, all $4$ (same) optimization problems above have exactly the same constraints but written slightly differently. Namely we have:

$$1\geq a_0\geq a_1\geq a_2\geq 0$$ $$1\geq b_2\geq b_1\geq b_0\geq 0$$

Improve this question
$\endgroup$

2 Answers 2

Reset to default
4
$\begingroup$

This is not a satisfactory answer but maybe a start using a brute force approach.

First the online help states:

  • NMaximize always attempts to find a global maximum of f subject to the constraints given.

  • If f and cons are linear, NMaximize can always find global maxima,
    over both real and integer values.

  • Otherwise, NMaximize may sometimes find only a local maximum.

Note the word "attempts" which applies in your case. So the brute force approach would be to use FindMaximum which only purports to find a local maximum. One could use random starting values for the parameters (subject to the constraints) or select starting values from a fine enough grid (again subject to the constraints).

Here's the random starting value approach:

SeedRandom[12345];
Do[
 startingValues = 
  Flatten[{Sort[RandomReal[1, 3]], Sort[RandomReal[1, 3]]}];
 parms = {a2, a1, a0, b0, b1, b2};
 init = Transpose[{parms, startingValues}];
 max = FindMaximum[{fmymin[a0, b0, a1, b1, a2, b2] - 
     fmyminorg[a0, b0, a1, b1, a2, b2], {0 < a2 <= a1 <= a0 < 1, 
     0 < b0 <= b1 <= b2 < 1}}, init];
 If[i == 1, best = max, If[max[[1]] > best[[1]], best = max]],
 {i, 1000}]
best
(* {0.0226371, {a2 -> 1.58994*10^-17, a1 -> 0.266554, a0 -> 0.941873,
                b0 -> 0.260884, b1 -> 0.333375, b2 -> 0.767125}} *)

Here's a grid search approach:

n = 6;
max = ConstantArray[0, (1/6 (-1 + n) n (1 + n))^2];
k = 0;
Do[Do[Do[Do[Do[Do[
      k = k + 1;
      startingValues = {i2/n, i1/n, i0/n, j0/n, j1/n, j2/n};
      parms = {a2, a1, a0, b0, b1, b2};
      init = Transpose[{parms, startingValues}];
      max[[k]] = 
       FindMaximum[{fmymin[a0, b0, a1, b1, a2, b2] - 
          fmyminorg[a0, b0, a1, b1, a2, b2], {0 < a2 <= a1 <= a0 < 1, 
          0 < b0 <= b1 <= b2 < 1}}, init];
      If[max[[k, 1]] > best[[1]], best = max[[k]]],
      {i2, 1, i1}], {i1, 1, i0}], {i0, 1, n - 1}], {j0, 1, j1}], {j1, 1, j2}], {j2, 1, n - 1}]
best
(* {0.0231124, {a2 -> 3.33598*10^-9, a1 -> 0.194261, a0 -> 0.827781,
                b0 -> 0.213226, b1 -> 0.218863, b2 -> 0.672657}} *)
Improve this answer
$\endgroup$
4
  • $\begingroup$ Thank you very much for the answer. Brute force search is just okay. The point is that I am still not able to get the global optimum. I know that the answer is at least $0.0233$ and I obtained with Nmaximize $0.0232$ with {a0 -> 0.699024, b0 -> 0.00182182, a1 -> 0.194679, b1 -> 0.296958, a2 -> 0.17766, b2 -> 0.92299}} $\endgroup$
    – Seyhmus Güngören
    Commented Aug 25, 2018 at 7:10
  • $\begingroup$ btw. I am running the second grid search method on my computer and it has been some hours already. How much more should I wait? $\endgroup$
    – Seyhmus Güngören
    Commented Aug 25, 2018 at 11:36
  • $\begingroup$ It took a little over 5 minutes with n = 6 on my laptop. $\endgroup$
    – JimB
    Commented Aug 25, 2018 at 15:56
  • $\begingroup$ still running on my computer. I think there should be something wrong. $\endgroup$
    – Seyhmus Güngören
    Commented Aug 25, 2018 at 16:10
2
$\begingroup$

Here is a more direct approach: playing with Method.

NMaximize[{fmymin[a0, b0, a1, b1, a2, b2] - fmyminorg[a0, b0, a1, b1, a2, b2],
  {0 < a2 <= a1 <= a0 < 1, 0 <= b0 <= b1 <= b2 < 1}}, {a0, b0, a1, b1, a2, b2},
 Method -> {"RandomSearch", "SearchPoints" -> 500}]
(* {0.0233256, {a0 -> 0.90198, b0 -> 0., a1 -> 0.711278,
                b1 -> 3.28891*10^-6, a2 -> 0.200922, b2 -> 0.307387}} *)
Improve this answer
$\endgroup$
1
  • $\begingroup$ This is also currently what I am running on my machine with some different seeds. I am still unsure whether I have the global optimum. I would prefer to use your other answer but It doesnt work. Could you please add a pure brute force search method so that I can accept the answer? I mean all parameters can take values from $\{0.01,0.02,...,1\}$. It is like 6 times for loop as in your answer with DODODO.. but no need for Findmaximum, It must just calculate the difference of those two functions and keep the best result at each iteration. Thanks for the answers once again. $\endgroup$
    – Seyhmus Güngören
    Commented Aug 26, 2018 at 15:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.