# NSolve does not return a value

I'm interested in obtaining the KKT multipliers in a constrained optimization problem. I've checked the NMaximize documentation, but I've been unable to find how to recover the multipliers. As an alternative, I've tried NSolve on the FOCs and Mathematica stalls (even though NMaximize does return a solution).

The code is as follows

ClearAll["Global*"]
T = {{2, 3, 1}, {3, 2, 1}, {1, 4, 2}, {4, 1, 3}};
techniques = Dimensions[T][[1]];
factors = Dimensions[T][[2]];
αVec = ConstantArray[4/6, techniques];
λVec = Array[λ, techniques];
γVec = Array[γ, factors];
μVec = Array[μ, techniques];
onesVec = ConstantArray[1, techniques];
zerosVec = ConstantArray[0, techniques];
needs = Transpose[T].λVec^(1/αVec);
FactorNeeds[mix_] := Transpose[T].mix^(1/αVec);

λcritical = Table[λVec /.
Flatten[
Minimize[{needs[[j]], onesVec.λVec == 1,
onesVec >= λVec >= zerosVec}, λVec, Reals]][[
2 ;; techniques + 1]], {j, 1, factors}];

NMaximize[{λVec.onesVec,needs <= FactorNeeds[Transpose[λcritical].{0, 1, 0}],onesVec >= λVec >= zerosVec}, λVec, Reals]

MatrixForm[Flatten[{D[λVec.onesVec + γVec.(needs -
FactorNeeds[
Transpose[λcritical].{0, 1,
0}]) + μVec.λVec, {λVec}],Table[γVec[[s]] (needs -
FactorNeeds[Transpose[λcritical].{1, 0, 0}])[[s]], {s,
1, factors}],Table[μVec[[s]] λVec[[s]], {s, 1, techniques}]}]]

NSolve[Flatten[{D[λVec.onesVec + γVec.(needs -
FactorNeeds[
Transpose[λcritical].{0, 1,
0}]) + μVec.λVec, {λVec}],
Table[γVec[[
s]] (needs -
FactorNeeds[Transpose[λcritical].{1, 0, 0}])[[
s]], {s, 1, factors}],
Table[μVec[[s]] λVec[[s]], {s, 1, techniques}]}] ==ConstantArray[0, 2*techniques + factors],Flatten[{λVec, γVec, μVec}], Reals]


I'm sure my coding can be improved upon, so feel free to comment on it as well. In particular, I think I should be able to do without the definition of either needs or FactorNeeds[].

Try this workaround:

eq = Flatten[{D[λVec.onesVec + γVec.(needs -FactorNeeds[Transpose[\
[Lambda]critical].{0, 1,0}]) + μVec.λVec, {λVec}],Table[\
[Gamma]Vec[[s]] (needs -FactorNeeds[Transpose[λcritical].{1, 0, 0}])[[s]],
{s, 1,factors}],
Table[μVec[[s]] λVec[[s]], {s, 1, techniques}]}] ==
ConstantArray[0, 2*techniques + factors]


NMinimize has a very robust equation solver

NMinimize[{1, eq}, Flatten[{λVec, γVec, μVec}]]
(*{1., {λ[1] -> 0.504511, λ[2] ->0.707354, λ[3] -> 0.176838,
λ[4] ->0.163564, γ[1] -> -4.7945*10^-9, γ[2] -> -0.145919,
γ[3] -> -0.500829, μ[1] -> 5.87935*10^-8, μ[2] ->6.94144*10^-8,
μ[3] ->-5.16111*10^-7, μ[4] -> -6.38508*10^-7}}*)

• Thanks. I was expecting to obtain the same results (for the $\lambda$s) than NMaximize does. Maybe I've made a mistake in my FOCs, I'll check. Also, I don't understand your use of NMinimizehere. Why do you use it and what is the purpose of the 1? Commented Jan 9, 2019 at 8:54
• It is a kind of misuse of NMinimize, which try to minimize the constant 1 thereby fullfilling the constraint eq . NMinimize searchs a global real solution. Commented Jan 9, 2019 at 11:21
• Ok. Thank you so much Commented Jan 9, 2019 at 11:25

Try this (look at the command Thread[] and to the change $$\{0.,1.,0.\},\{1.,0.,0.\}$$ instead of $$\{0,1,0\},\{1,0,0\}$$)

NMaximize[{λVec.onesVec, Thread[needs <= FactorNeeds[Transpose[λcritical].{0., 1., 0.}]],Thread[onesVec >= λVec >= zerosVec]}, λVec, Reals]

FactorNeeds[
Transpose[λcritical].{0., 1., 0.}]) + μVec.λVec, {\
[Lambda]Vec}],
Table[γVec[[s]] (needs - FactorNeeds[Transpose[λcritical].{1., 0., 0.}])[[s]], {s, 1, factors}],
Table[μVec[[s]] λVec[[s]], {s, 1, techniques}]}] == ConstantArray[0, 2*techniques + factors]],
Flatten[{λVec, γVec, μVec}], Reals]


giving the result

{{λ[1] -> 0.167225, λ[2] -> 0.376255, λ[3] ->
0.0940638, λ[4] -> 1.50502, μ[1] -> 0, μ[2] ->
0, μ[3] -> 0, μ[4] -> 0, γ[1] ->
0, γ[2] -> -0.543422, γ[3] -> 0}, {λ[1] ->
0.671919, λ[2] -> 0.671919, λ[3] ->
0.16798, λ[4] -> 0.0746576, μ[1] -> 0, μ[2] ->
0, μ[3] -> 0, μ[4] -> 0, γ[1] -> 0, γ[2] ->
0, γ[3] -> -0.813299}, {λ[1] ->
0.504511, λ[2] -> 0.707353, λ[3] ->
0.176838, λ[4] -> 0.163564, μ[1] -> 0, μ[2] ->
0, μ[3] -> 0, μ[4] -> 0, γ[1] ->
0, γ[2] -> -0.145918, γ[3] -> -0.50083}, {λ[
1] -> 0, λ[2] -> 0, λ[3] -> 0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] ->
0}, {λ[1] -> 0, λ[2] -> 0, λ[3] ->
0, λ[4] ->
0, μ[1] -> -1., μ[2] -> -1., μ[3] -> -1., μ[
4] -> -1., γ[1] -> 0, γ[2] -> 0, γ[3] -> 0}}

• Why do we need to Thread` the inequality instead of using vector inequalities? Commented Jan 9, 2019 at 8:12
• @Patricio Is NSolve returning now the due result? Commented Jan 9, 2019 at 8:56
• No, it still stalls Commented Jan 9, 2019 at 8:59
• @Patricio Did you noticed that in my script I represent $\{0.,1.,0.\}$ instead $\{0,1,0\}$ and $\{1.0,0.0,0.0\}$ instead $\{1,0,0\}$ ? Try that. Commented Jan 9, 2019 at 9:04