#load plotting package
using PyPlot
#load optimization package
using Optim

#########Setting up the primitives of the model#############

#lower bound type distribution
const m = 0.5

#number of contestants
const k = 3

#cdf of type distribution; here uniform
function F(a::Float64)
    if m<= a<= 1.0
        return (a-m)/(1.0-m)
    elseif a<m
        return 0.0
    else 
        return 1.0
    end
end

#pdf of type distribution; here uniform
function f(a::Float64)
    if m<= a<= 1.0
        return 1.0/(1.0-m)
    else 
        return 0.0
    end
end    

#Bernoulli utility of prize money; here √x
function u(x::Float64)
    return sqrt(x)
end

#inverse of the Bernoulli utility function of prize money; here x^2
function uinv(v::Float64)
    return v^2
end

############Setting up the equilibrium strategies of the contestants###########
A(c) = (k-1)*quadgk(a-> (1-F(a))^(k-2) * f(a)/a,c,1)[1] #"quadgk" integrates the function (1st argument) using the 2nd and 3rd argument as boundaries
B(c) = (k-1)*quadgk(a-> (1-F(a))^(k-3) * ((k-1)*F(a)-1) * f(a)/a,c,1)[1]

b(W1,W2,c) = A(c)*W1 + B(c)*W2

####### Designer's problem ########
#this function returns the expected effort of a participant multiplied with -1 (as the optimization procedure below searches for the minimum)
function minusobjective(W1::Float64)
    function wrapped_objective(c::Float64)
        return b(W1,u(1.0-uinv(W1)),c)*f(c) #note that W2 = u(1.0-uinv(W1)) by the constraint V1+V2=1
    end
    return (-1.0)*quadgk(wrapped_objective,m,1)[1]
end

#result to optimization problem (minimizes "minusobjective" with W1 between u(0.5) and u(1.0); note that V1 cannot be less than 0.5)
r = optimize(minusobjective,u(0.5),u(1.0))

println(" ")
println("Optimal first prize is ", uinv(r.minimum))
println("Optimal second prize is ", 1.0-uinv(r.minimum))

#for plotting: returns average effort of a contestant if first prize is V1
function plotobjective(V1::Float64)
    return (-1.0)*minusobjective(u(V1))
end
function plotobjective(V1::FloatRange{Float64})
    return [(-1.0)*minusobjective(u(v)) for v in V1]
end

fig, ax = subplots()
ax[:set_xlabel]("first prize")
ax[:set_ylabel]("expected effort")
ax[:plot](0.0:0.01:1.0,plotobjective(0.0:0.01:1.0),"b-",linewidth=2,alpha=0.5)


((quadgk(c->A(c)*f(c),m,1)[1])^2) / ((quadgk(c->B(c)*f(c),m,1)[1])^2 )

(quadgk(c->A(c)*f(c),m,1)[1])^2 / ((quadgk(c->B(c)*f(c),m,1)[1])^2+(quadgk(c->A(c)*f(c),m,1)[1])^2 )