reset;

option randseed'';

set A := 1..100;                                # set of potential prices

param T default 10;								# time horizon

param a1 := 20; #round(Uniform(1,100));			# lower bound for two-bound strategy
param a2 := 30; #round(Uniform(a1,100));		# upper bound for two-bound strategy

param S1 {p in A}:= max(if p<a1 or p>a2 then a2 else p-1,1); # strategy player 1
param S2 {a in A}:= max(a-1,1); 							 # strategy player 2

param P{a in A,p in A} := 1-min(a,p)/100;       # probabilities for a sale (for period of length 0.5)


#####  Simulation/Evaluation of Price Responses 

param pt {t in 0..T by 0.5,k in 1..2} :=        # evaluation of mutual price reactions
         if t=0 then (if k=1 then 20 else 20)	
	else if k=2 then (if t<>floor(t) then S2[pt[t-0.5,1]] else pt[t-0.5,2])
	else if k=1 then (if t= floor(t) then S1[pt[t-0.5,2]] else pt[t-0.5,1]);

param Rt {t in 0..T by 0.5,k in 1..2} default	# for profit accumulation of both players
         if t=0 then 0 else Rt[t-0.5,k];

for {t in 0..T by 0.5} {						# simulate profit realizations for both players
         if Uniform(0,1)<P[pt[t,1],pt[t,2]]
   then {
        if pt[t,1]<pt[t,2]  then let Rt[t,1]:=Rt[t,1]+pt[t,1];
   else if pt[t,1]>pt[t,2]  then let Rt[t,2]:=Rt[t,2]+pt[t,2];
   else if Uniform(0,1)<0.5 then let Rt[t,1]:=Rt[t,1]+pt[t,1];
   else                          let Rt[t,2]:=Rt[t,2]+pt[t,2];
         }};
	
	
display pt;                     # display price trajectories
display Rt;                     # display accumulated profits

display Rt[T,1],Rt[T,2];        # final total profits

display a1,a2;                  # used bounds for strategy of player 1


#end;

#####   Optimal Strategy (Finite Horizon)  #####

param delta default 1;   	# discount factor

let T:=200;              	# number of iterations / length of time horizon

param P1{a in A,p in A} := P[a,p] *        # probability that player 1 sells an item in (t,t+0.5)
                          (if a<p then 1 else if a=p then 1/2);     # (no sale in (t+0.5,t+1))

param R{a in A,p in A} := a;    # profit for player 1 in case of a sale

					 
param V{t in 0..T,p in A}:= if t=T then 0 else max {a in A}         # value function
	(P1[a,p] * R[a,p] + delta * V[t+1,S2[a]] );

param at{t in 0..T-1,p in A}:= min{a in A: V[t,p] = if t<T then     # arg max strategy
	(P1[a,p] * R[a,p] + delta * V[t+1,S2[a]] )} a;

	
param Vt{t in 0..min(50,T)} := V[t,50]/T;	# final 50 iterations for state p=50 (check convergence)

param V0{p in A} := V[0,p];                 # expected profit in t=0
param a0{p in A} := at[0,p];	            # price response in t=0

param VS1{t in 0..T,p in A}:= if t=T then 0 else max {a in A:a=S1[p]} # evaluate S1 (no maximization)
	(P1[a,p] * R[a,p] + delta * VS1[t+1,S2[a]] );
param VS10{p in A} := VS1[0,p];             # expected profit in t=0

display Vt;                                 # check convergence

display V0;                                 # compare total exp. profits S1 vs. optimal (T<00)
display VS10;

display a0;                                 # compare S1 vs. optimal strategy (T<00)
display S1;

#end;


##########  Optimal Strategy (Infinite Horizon)  ##########

printf"Infinite Horizon Case\n\n";

let delta := 0.999;                        # discount factor (for 1 period)

var value_opt {p in A} := 2000;            # variables for value function
var value_heu {p in A} := 2000;            # variables for heuristic expected profits


#####   expected results of the heuristic strategy S1 for player 1  #####
								
subject to NB_heu {p in A}:  value_heu[p] = max {a in A:a=S1[p]}  (  P1[a,p]
                                * R[a,p] + delta * value_heu[S2[a]]) ;

								
#####   optimal strategy for player 1 against player 2  #####

subject to NB_opt {p in A}:   value_opt[p] = max {a in A}         (  P1[a,p]
                                * R[a,p] + delta * value_opt[S2[a]]) ;

option solver minos;

solve; 

#####   get optimal strategy for player 1 (arg max)  #####

param V00{p in A} default 0;

for {p in A} let V00[p] := value_opt[p]; 	# save value function as parameter

param p00 {p in A} := 		# optimal price response strategy player 1 (arg max)

              min {a in A :     0.000001 > abs( V00[p]  - (P1[a,p]
                                * R[a,p] + delta * V00[S2[a]]) ) } a;

display p00;
#display S1;
display value_opt[20],value_heu[20];    # comparison of results for state p=20 (T=00)


end;

printf"" > duopoly_out.txt;
for{p in A diff{0}} {
	printf"%1.3f %1.6f %1.6f\n", p, p00[p],value_opt[p] >> duopoly_out.txt;
					};