# Optimization Problem: NMaximize goes out of the constraints

I am having issues with the code:

CC = p0^2/(p0 - L)
Dbuy = (1 - CC)*(L/p0) + (CC - p0)*Log[CC/p0]
Dwait = (1 - p1)*(L/(p1 - p0))
Dlock = (1 - p1)*(1 - L/p0 - L/(p1 - p0)) + (p1 - CC)*(1 - L/p0) - L*Log[(p1 - p0)/(CC - p0)]
RevLocking1 = Table[NMaximize[{x*Dbuy*p0 + x*Dwait*0.4*p1 + x*Dlock*(p0 + 0.4*L), x*(Dbuy + Dwait + Dlock) <= 1,  0 <= L <= p0 <= p1 <= 1}, {p0, p1, L}, MaxIterations -> 10000], {x, 1, 6, 0.25}]


The code works, but it ignores the constraint:

0 <= L <= p0 <= p1 <= 1


Since it reaches values like:

L -> 0.303712, p0 -> -0.0156877, p1 -> 0.900987


Any ideas?

• Possible duplicate: mathematica.stackexchange.com/questions/59706/… – Michael E2 Oct 19 '18 at 21:21
• I used the "Re" trick, but didn't really work. Do you have a better suggestion? – Ovunc Oct 19 '18 at 23:21
• @ovunc: You should use && in between constraints. You are using ,, which is not correct. NMaximize[{obj, const1&&const2&&...},{vars}] – Tugrul Temel Oct 19 '18 at 23:49

It might be that the starting values are too far away from those that maximize the function. Rather than using Table one could put everything into a Do loop and use the previous estimates as the starting values using FindMaximum rather than NMinimize. Using (or stealing) some of the code from @TugrulTemel:

ClearAll[LL, CC, Dbuy, Dwait, Dlock, xx, xmin, xmax, n];
CC = p0^2/(p0 - LL);
Dbuy = (1 - CC)*(LL/p0) + (CC - p0)*Log[CC/p0];
Dwait = (1 - p1)*(LL/(p1 - p0));
Dlock = (1 - p1)*(1 - LL/p0 - LL/(p1 - p0)) + (p1 - CC)*(1 - LL/p0) - LL*Log[(p1 - p0)/(CC - p0)];

max[xx_, LLInit_, p0Init_, p1Init_] :=
FindMaximum[{xx*(Dbuy*p0 + Dwait*0.4*p1 + Dlock*(p0 + 0.4*LL)),
xx*(Dbuy + Dwait + Dlock) <= 1 && 0 <= LL <= p0 && p0 <= p1 <= 1},
{{LL, LLInit}, {p0, p0Init}, {p1, p1Init}}]

xmin = 1.;
xmax = 6;
n = 20;
RevLocking1 = ConstantArray[{0, 0, 0, 0, 0}, n + 1];
Do[xx = xmin + (xmax - xmin) i/n;
sol = If[i == 0, max[xx, 0.28, 0.95, 0.99],
max[xx, RevLocking1[[i, 3]], RevLocking1[[i, 4]], RevLocking1[[i, 5]]]];
RevLocking1[[i + 1]] = {xx, sol[[1]], LL, p0, p1} /. sol[[2]],
{i, 0, n}]
RevLocking1 // TableForm


$$\left( \begin{array}{ccccc} 1. & 1.12066 & 0.373178 & 0.981281 & 1. \\ 1.25 & 1.10575 & 0.334785 & 0.979521 & 1. \\ 1.5 & 1.09425 & 0.304713 & 0.978577 & 1. \\ 1.75 & 1.0851 & 0.280353 & 0.978128 & 1. \\ 2. & 1.07763 & 0.260118 & 0.977986 & 1. \\ 2.25 & 1.07141 & 0.242977 & 0.97804 & 1. \\ 2.5 & 1.06615 & 0.228229 & 0.978219 & 1. \\ 2.75 & 1.06164 & 0.215377 & 0.978478 & 1. \\ 3. & 1.05773 & 0.204055 & 0.978788 & 1. \\ 3.25 & 1.0543 & 0.193992 & 0.979129 & 1. \\ 3.5 & 1.05127 & 0.184976 & 0.979488 & 1. \\ 3.75 & 1.04857 & 0.176844 & 0.979855 & 1. \\ 4. & 1.04614 & 0.169466 & 0.980226 & 1. \\ 4.25 & 1.04396 & 0.162735 & 0.980594 & 1. \\ 4.5 & 1.04198 & 0.156567 & 0.980958 & 1. \\ 4.75 & 1.04017 & 0.15089 & 0.981316 & 1. \\ 5. & 1.03852 & 0.145644 & 0.981665 & 1. \\ 5.25 & 1.037 & 0.140781 & 0.982006 & 1. \\ 5.5 & 1.0356 & 0.136258 & 0.982337 & 1. \\ 5.75 & 1.03431 & 0.132039 & 0.982659 & 1. \\ 6. & 1.03311 & 0.128093 & 0.982972 & 1. \\ \end{array} \right)$$

• This is a promising start. I made a few changes in the maximization (there were some issues in the model) and it works just fine. Can we use NMaximize instead of FindMaximum though? – Ovunc Oct 20 '18 at 5:44
• How are the models @TugrulTemel and I used different from your changed model? – JimB Oct 20 '18 at 13:25
• I made a change in a constraint: Dbuy = (1 - CC)*(LL/p0) + (CC - p0) - p0*Log[CC/p0]; I also added a new constraint into the FindMaximum: p1*(p0 - LL) >= p0^2 Things work very good right now, and numbers make much more sense. I will go use my original method with these new constraints and let you know about it. Thanks so much! – Ovunc Oct 20 '18 at 16:53
• Another change I made is "p0 + 0.4*LL" became "0.4*p0+LL", which is quite important:P – Ovunc Oct 20 '18 at 17:17
• can you tell me how the following works? Do[xx = xmin + (xmax - xmin) i/n; sol = If[i == 0, max[xx, 0.28, 0.95, 0.99], max[xx, RevLocking1[[i, 3]], RevLocking1[[i, 4]], RevLocking1[[i, 5]]]]; RevLocking1[[i + 1]] = {xx, sol[[1]], LL, p0, p1} /. sol[[2]], {i, 0, n}] – Ovunc Oct 20 '18 at 17:53

You forgot to add x to the variable list. Also avoid using single capital letters as variable names since they may interfere with Mathematica built-in functions. I made some changes in your variable names, such x replaced with xx.

Try this:

ClearAll[LL, CC, Dbuy, Dwait, Dlock];
CC = p0^2/(p0 - LL)
Dbuy = (1 - CC)*(LL/p0) + (CC - p0)*Log[CC/p0]
Dwait = (1 - p1)*(LL/(p1 - p0))
Dlock = (1 - p1)*(1 - LL/p0 - LL/(p1 - p0)) + (p1 - CC)*(1 - LL/p0) -
LL*Log[(p1 - p0)/(CC - p0)]
RevLocking1 =
NMaximize[{xx*(Dbuy*p0 + Dwait*0.4*p1 + Dlock*(p0 + 0.4*LL)),
xx*(Dbuy + Dwait + Dlock) <= 1 && 0 <= LL <= p0 &&
p0 <= p1 <= 1}, {p0, p1, LL, xx}]


giving the following result:

{0.147246, {LL -> 0.3322, p0 -> 0.571466, p1 -> 0.960381,
xx -> 0.457388}}


and it gives the following warning:

NMaximize::nrnum: "The function value -10.5461+1.01463\ I is not a real number at {LL,p0,p1,xx} = {0.21037,0.534614,0.529012,0.991907}. \!$$\*ButtonBox[\">>\", Appearance->{Automatic, None}, BaseStyle->\"Link\", ButtonData:>\"paclet:ref/NMaximize\", ButtonNote->\"NMaximize::nrnum\"]$$"

• x is not a variable, but a parameter. So I basically solve the same optimization problem for different values of x (between 1 and 6, with increments of 0.25). Thanks for the suggestions though, they are certainly helpful! – Ovunc Oct 20 '18 at 4:31
• @Ovunc: You may want to try Manipulate for visualizing the solution with varying x. – Tugrul Temel Oct 20 '18 at 10:08