backtraces for nonlinear expressions

[mlubin]

They don't currently work.

[mlubin]

Related to #320

[tobiassalz]

Below is a code example. It is supposed to solve an asymmetric auction with collocation methods. In this case the inverses of the optimal bidding strategies are approximated with simple polynomials. The code terminates with

WARNING: Ipopt finished with status Invalid_Number_Detected
WARNING: Not solved to optimality, status: Error

One possible source of this is the fact that a truncated weibull distribution (the distribution of cost draws) is evaluated. But the expression for the pdf/cdf might not be defined for x outside of [ll,uu]. That is why I am imposing some constraints on the argument (in this case constraints on non-linear expressions) but I am not sure whether these constraints are allowed to be violated, in which case this might be the source of the error.

Not sure whether this is a useful example.

[tobiassalz]

using JuMP

Hypergeometric Distribution

function hyperprob(NL,NH,m,k)

v1 = binomial(NL,k)
v2 = binomial(NH,m-k)
v3 = binomial(NL+NH,m)

return (v1*v2)/v3

end;

Define parameters of weibull distributions + number of firms

npoints = 300; # number of points for which we measure deviation at FOC
ll = 0.1 # weibull lower truncation point
uu = 7.0 # weibull higher truncation point
aparamH = 1.5
lparamH = 1.
aparamL = 1.5
lparamL = 2.
npoly = 5

M = 5 # M = number of firms + 1
NH = 2 # number of firms of type H
NL = 2 # number of firms of type L

Define probabilities of firm composition

cL = zeros(NL+NH-1,NL+NH)
cH = zeros(NL+NH-1,NL+NH)
for m = 1:NL+NH-1
for k = 0:m
cL[m,k+1] = hyperprob(NL-1,NH,m,k)
cH[m,k+1] = hyperprob(NH-1,NL,m,k)
end;
end;

define probability weights on number of competitors

h = rand(M-1)
h = h./sum(h)
h = [0.0 1/2 1/2 0.0]

define jump model

m = Model()

@defvar(m, 0.0000 <= alphaL[1:npoly] <= 5.0);
@defvar(m, 0.0000 <= alphaH[1:npoly] <= 5.0);
@defvar(m, ll <= plowerL <= 10.2);
@defvar(m, ll <= plowerH <= 10.2);

for i = 1:npoly
setValue(alphaL[i],0.0)
setValue(alphaH[i],0.0)
end
setValue(plowerL,ll+0.1)
setValue(plowerH,ll+0.1)

define endogenous grid of player L and H

@defNLExpr(psetL[j=1:npoints],plowerL + ((j-1)(uu-plowerL))/(npoints-1));
@defNLExpr(psetH[j=1:npoints],plowerH + ((j-1)(uu-plowerH))/(npoints-1));

define the inverse bidding function as well as it's derivative of player L

@defNLExpr(betainvL[j=1:npoints],plowerL + alphaL[1]+sum{alphaL[k](psetL[j]-plowerL)^(k-1),k=2:npoly});
@defNLExpr(dbetainvL[j=1:npoints],sum{alphaL[k](k-1)*(psetL[j]-plowerL)^(k-2),k=2:npoly});

define the inverse bidding function as well as it's derivative of player H

@defNLExpr(betainvH[j=1:npoints],plowerH + alphaH[1]+sum{alphaH[k](psetH[j]-plowerH)^(k-1),k=2:npoly});
@defNLExpr(dbetainvH[j=1:npoints],sum{alphaL[k](k-1)*(psetH[j]-plowerH)^(k-2),k=2:npoly});

@addNLConstraint(m, mbetainvLC[i=1:npoints], ll <= betainvL[i] <= uu);
@addNLConstraint(m, mbetainvHC[i=1:npoints], ll <= betainvH[i] <= uu);

truncated weibull pdf evaluated at inverse bidding function

@defNLExpr(wL[i=1:npoints],((aparamL/lparamL)_((betainvL[i]/lparamL)^(aparamL-1))exp(-((betainvL[i]/lparamL)^aparamL)))/(exp(-((ll/lparamL)^aparamL))-exp(-((uu/lparamL)^aparamL))));
@defNLExpr(wH[i=1:npoints],((aparamH/lparamH)((betainvH[i]/lparamH)^(aparamH-1))_exp(-((betainvH[i]/lparamH)^aparamH)))/(exp(-((ll/lparamH)^aparamH))-exp(-((uu/lparamH)^aparamH))));

truncated weibull cdf evaluated at inverse bidding function

@defNLExpr(WL[i=1:npoints],(exp(-((ll/lparamL)^aparamL))-exp(-((betainvL[i]/lparamL)^aparamL)))/(exp(-((ll/lparamL)^aparamL))-exp(-((uu/lparamL)^aparamL))));
@defNLExpr(WH[i=1:npoints],(exp(-((ll/lparamH)^aparamH))-exp(-((betainvH[i]/lparamH)^aparamH)))/(exp(-((ll/lparamH)^aparamH))-exp(-((uu/lparamH)^aparamH))));

expressions that make up FOC of player L and H

@defNLExpr(psi1L[i=1:npoints],sum{h[m]sum{cL[m,k+1](WL[i]^k)(WH[i]^(m-k)),k=0:m},m=1:NL+NH-1});
@defNLExpr(psi1H[i=1:npoints],sum{h[m]sum{cH[m,k+1](WH[i]^k)(WL[i]^(m-k)),k=0:m},m=1:NL+NH-1});

@defNLExpr(psi2L[i=1:npoints],(psetL[i]-betainvL[i])_sum{h[m]sum{(k_wL[i]cL[m,k+1](WL[i]^(k-1))(WH[i]^(m-k)))/(0.0001+dbetainvL[i]),k=1:m},m=1:NL+NH-1});
@defNLExpr(psi2H[i=1:npoints],(psetH[i]-betainvH[i])_sum{h[m]sum{(k_wH[i]cH[m,k+1](WH[i]^(k-1))(WL[i]^(m-k)))/(0.0001+dbetainvH[i]),k=1:m},m=1:NL+NH-1});

@defNLExpr(psi3L[i=1:npoints],(psetL[i]-betainvL[i])_sum{h[m]_sum{((m-k)wH[i]cL[m,k+1](WL[i]^(k))(WH[i]^(m-k-1)))/(0.0001+dbetainvH[i]),k=0:m-1},m=1:NL+NH-1});
@defNLExpr(psi3H[i=1:npoints],(psetH[i]-betainvH[i])_sum{h[m]_sum{((m-k)wL[i]cH[m,k+1](WH[i]^(k))(WL[i]^(m-k-1)))/(0.0001+dbetainvL[i]),k=0:m-1},m=1:NL+NH-1});

@defNLExpr(betainvHighL,plowerL + alphaL[1]+sum{alphaL[k](uu-plowerL)^(k-1),k=2:npoly});
@defNLExpr(betainvHighH,plowerH + alphaH[1]+sum{alphaH[k](uu-plowerH)^(k-1),k=2:npoly});

@defNLExpr(betainvLowL,plowerL + alphaL[1]);
@defNLExpr(betainvLowH,plowerH + alphaH[1]);

expressions that make up FOC of player H

@setNLObjective(m,Min, (sum{(psi1L[i]-psi2L[i]-psi3L[i])^2,i=1:npoints}+sum{(psi1H[i]-psi2H[i]-psi3H[i])^2,i=1:npoints}+ npoints_(betainvHighL-uu)^2 + npoints_(betainvHighH-uu)^2 + npoints_(betainvLowL-ll)^2 + npoints_(betainvLowH-ll)^2)/10000000)
solve(m)

[mlubin]

I think this issue went away with the NLP rewrite?