0
$\begingroup$

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[].

$\endgroup$

2 Answers 2

0
$\begingroup$

Perhaps NSolve needs more information about the solution range.

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}}*)
$\endgroup$
3
  • $\begingroup$ 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? $\endgroup$
    – Patricio
    Commented Jan 9, 2019 at 8:54
  • $\begingroup$ 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. $\endgroup$ Commented Jan 9, 2019 at 11:21
  • $\begingroup$ Ok. Thank you so much $\endgroup$
    – Patricio
    Commented Jan 9, 2019 at 11:25
0
$\begingroup$

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]


NSolve[Thread[Flatten[{D[λVec.onesVec + γVec.(needs - 
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}}
$\endgroup$
4
  • $\begingroup$ Why do we need to Thread the inequality instead of using vector inequalities? $\endgroup$
    – Patricio
    Commented Jan 9, 2019 at 8:12
  • $\begingroup$ @Patricio Is NSolve returning now the due result? $\endgroup$
    – Cesareo
    Commented Jan 9, 2019 at 8:56
  • $\begingroup$ No, it still stalls $\endgroup$
    – Patricio
    Commented Jan 9, 2019 at 8:59
  • $\begingroup$ @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. $\endgroup$
    – Cesareo
    Commented Jan 9, 2019 at 9:04

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.