kill(all);

/* (Inverse) Demand and Profit */

solve(
[v12*u1-p1=v12*u2-p2, v2*u2-p2=0]
,[v12,v2]
)[1]$
Defv12v2:%;

[q1=v0-v12, q2=v12-v2]$
ev(%,Defv12v2)$
Demand:%;

solve(
[q1=rhs(Demand[1]), q2=rhs(Demand[2])]
,[p1,p2]
)[1]$
InverseDemand:%;

pi1:p1*q1-u1^2/2$
pi2:p2*q2-u2^2/2$

/* II(i) Bertrand competition */

/* p1 and p2 */

pi1$
ev(%,Demand)$
diff(%,p1)$
FOCp1B:%$

pi2$
ev(%,Demand)$
diff(%,p2)$
FOCp2B:%$

solve(
[FOCp1B=0, FOCp2B=0]
,[p1,p2]
)[1]$
Outcomep1p2B:%;

[pi1,pi2]$
ev(%,Demand)$
ev(%,Outcomep1p2B)$
[%[1]+u1^2/2, %[2]+u2^2/2]$
factor(%)$
[%[1]-u1^2/2, %[2]-u2^2/2];

/* u1 and u2 */

pi1$
ev(%,Demand)$
ev(%,Outcomep1p2B)$
diff(%,u1)$
%+u1$
factor(%)$
%-u1$
FOCu1B:%;

pi2$
ev(%,Demand)$
ev(%,Outcomep1p2B)$
diff(%,u2)$
%+u2$
factor(%)$
%-u2$
FOCu2B:%;

FOCu1B$
factor(%)$
ev(%,v0=V^(1/2))$
solve([%=0],V)$
factor(%)$
DefVB:%;

FOCu2B$
ev(%,v0=V^(1/2))$
factor(%)$
ev(%,DefVB)$
factor(%);
num(%);
ev(%, u1=m*u2)$
factor(%)/(u2^3);
allroots(%=0);
[%[3]]$
OutcomemB:%;

FOCu1B$
ev(%,u2=u1/m)$
ev(%,OutcomemB)$
factor(%)$
solve( [%=0],u1 )$
float(%)$
Outcomeu1B:%;

FOCu2B$
ev(%,u1=m*u2)$
ev(%,OutcomemB)$
factor(%)$
solve( [%=0],u2 )$
float(%)$
Outcomeu2B:%;

/* SOC */

pi1$
ev(%,Demand)$
ev(%,Outcomep1p2B)$
diff(%,u1)$
diff(%,u1)$
%+1$
factor(%)$
%-1;

pi2$
ev(%,Demand)$
ev(%,Outcomep1p2B)$
diff(%,u2)$
diff(%,u2)$
%+1$
factor(%)$
%-1;

/* II(ii) Cournot competition */

/* q1 and q2 */

pi1$
ev(%,InverseDemand)$
diff(%,q1)$
FOCq1C:%$

pi2$
ev(%,InverseDemand)$
diff(%,q2)$
FOCq2C:%$

solve(
[FOCq1C=0, FOCq2C=0]
,[q1,q2]
)[1]$
Outcomeq1q2C:%;

[p1,p2]$
ev(%,InverseDemand)$
ev(%,Outcomeq1q2C)$
factor(%);

[pi1,pi2]$
ev(%,InverseDemand)$
ev(%,Outcomeq1q2C)$
[%[1]+u1^2/2, %[2]+u2^2/2]$
factor(%)$
[%[1]-u1^2/2, %[2]-u2^2/2];

/* u1 and u2 */

pi1$
ev(%,InverseDemand)$
ev(%,Outcomeq1q2C)$
diff(%,u1)$
%+u1$
factor(%)$
%-u1$
FOCu1C:%;


pi2$
ev(%,InverseDemand)$
ev(%,Outcomeq1q2C)$
diff(%,u2)$
%+u2$
factor(%)$
%-u2$
FOCu2C:%;

FOCu1C$
factor(%)$
ev(%, v0=V^(1/2))$
solve( [%=0],V)$
factor(%)$
DefVC:%;

FOCu2C$
factor(%)$
ev(%, v0=V^(1/2))$
ev(%, DefVC)$
factor(%);
num(%);
ev(%,u1=k*u2)$
%/(u2^4)$
factor(%);
allroots(%=0);
[%[4]]$
OutcomekC:%;

FOCu1C$
ev(%,u2=u1/k)$
ev(%,OutcomekC)$
factor(%)$
solve( [%=0], u1)$
float(%);

FOCu2C$
ev(%,u1=k*u2)$
ev(%,OutcomekC)$
solve( [%=0], u2)$
float(%);