kill(all);

/* Model */
u:a*(q1+q2)-b/2*(q1^2+q2^2)-d*(q1*q2)-p1*q1-p2*q2$

solve( [ diff(u,q1)=0, diff(u,q2)=0 ],[q1,q2] )[1]$
Demand:%$

piU1:w1*q1+F1$
piU2:w2*q2+F2$

piD1:(p1-w1)*q1-F1$
piD2:(p2-w2)*q2-F2$

/* Equilibrium with passive beliefs */

piD1$
ev(%,Demand)$
diff(%,p1)$
FOCp1PB:%$

piD2$
ev(%,Demand)$
diff(%,p2)$
FOCp2PB:%$

solve( [FOCp1PB=0,FOCp2PB=0],[p1,p2] )[1]$
factor(%)$
Outcomep1p2PB:%;

piD1$
ev(%,Demand)$
solve(%=0,F1)$
factor(%)$
OutcomeF1PB:%;

piD2$
ev(%,Demand)$
solve(%=0,F2)$
factor(%)$
OutcomeF2PB:%;

piU1$
ev(%,Demand)$
ev(%,OutcomeF1PB)$
ev(%,OutcomeF2PB)$
ev(%,Outcomep1p2PB[1])$
diff(%,w1)$
ev(%,Outcomep1p2PB)$
factor(%)$
FOCw1PB:%$

piU2$
ev(%,Demand)$
ev(%,OutcomeF1PB)$
ev(%,OutcomeF2PB)$
ev(%,Outcomep1p2PB[2])$
diff(%,w2)$
ev(%,Outcomep1p2PB)$
factor(%)$
FOCw2PB:%$

solve( [FOCw1PB=0,FOCw2PB=0],[w1,w2] )[1]$
factor(%)$
Outcomew1w2PB:%;

[p1,p2,piU1,piU2]$
ev(%,Demand)$
ev(%,OutcomeF1PB)$
ev(%,OutcomeF2PB)$
ev(%,Outcomep1p2PB)$
ev(%,Outcomew1w2PB)$
factor(%);

/* Equilibrium with symmetric beliefs */

piD1$
ev(%,Demand)$
diff(%,p1)$
ev(%, [p2=p1])$
solve( %=0, p1 )$
Outcomep1SB:%;

piD2$
ev(%,Demand)$
diff(%,p2)$
ev(%, [p1=p2])$
solve( %=0, p2 )$
Outcomep2SB:%;

piD1$
ev(%,Demand)$
ev(%,[p2=p1])$
ev(%,Outcomep1SB)$
solve(%=0,F1)$
factor(%)$
OutcomeF1SB:%;

piD2$
ev(%,Demand)$
ev(%,[p1=p2])$
ev(%,Outcomep2SB)$
solve(%=0,F2)$
factor(%)$
OutcomeF2SB:%;

piU1$
ev(%,Demand)$
ev(%,OutcomeF1SB)$
ev(%,OutcomeF2SB)$
ev(%,Outcomep1SB)$
ev(%,Outcomep2SB)$
diff(%,w1)$
factor(%)$
FOCw1SB:%$

piU2$
ev(%,Demand)$
ev(%,OutcomeF1SB)$
ev(%,OutcomeF2SB)$
ev(%,Outcomep1SB)$
ev(%,Outcomep2SB)$
diff(%,w2)$
factor(%)$
FOCw2SB:%$

solve( [FOCw1SB=0,FOCw2SB=0],[w1,w2] )[1]$
factor(%)$
Outcomew1w2SB:%;

[p1,p2,piU1,piU2]$
ev(%,Demand)$
ev(%,OutcomeF1SB)$
ev(%,OutcomeF2SB)$
ev(%,Outcomep1SB)$
ev(%,Outcomep2SB)$
ev(%,Outcomew1w2SB)$
factor(%);