Clear[R,P,Q,f,f0,f1,A,B,k,m,s,β,ζ,γ0,γ1,ϵ,t];

$Assumptions=(m∈Integers)&&s∈Integers&&s≥2&&k≥2&&k>s&&(0<ϵ<0.5)&&(0<β<1/k)&&(γ0>0)&&(γ1>0);

----------------Define the parameters----------------

\mu_{s-1} = \beta, \mu_s = \beta * e^{\gamma_0}.

---------------------------------------------------------------

β;

ζ=1-β*(E^(γ0));

A=(s-1)*t+(s-1)*γ0;

B=k*t+(s-1)*γ0;

----------------Define the measure μ----------------

Clear[j,μ];
μ[j_]:=Which[j==0,1-(s-1)*β-(k-s+1)*β*(E^(γ0)),
1≤j≤s-1,β,
j≥s,β*(E^(γ0))];

----------------Define the functions μf----------------
---------------- j is a notation for e_j, j=0,1,...,k, 
---------------- m is the index for player, m = 1,...,k, 
---------------- t is time
----------------------------------------------------------------

μf[j_,m_,t_]:=Which[t≤0&&1≤m≤s-1&&j==0,μ[0]*(E^(-A)),
t≤0&&1≤m≤s-1&&1≤j≤s-1,β*(E^(γ0))*(E^t)*(E^(-A))*Boole[j≠m],
t≤0&&1≤m≤s-1&&s≤j≤k,β*(E^(γ0))*(E^(-A)),
t≤0&&m≥s,0,
t>0&&1≤m≤k&&j==0,μ[0]*(E^(-B)),
t>0&&1≤m≤k&&1≤j≤k,β*(E^(γ0))*(E^t)*(E^(-B))*Boole[j≠m]];

μf0[j_,m_,t_]:=Which[t>0,(1-ϵ*ζ)*μf[j,m,t],
t≤0&&m≠s&&j≠s,(1-ϵ*ζ)*(E^(-t))*μf[j,m,t],
t≤0&&m≠s&&j==s,(1-ϵ*ζ)*μf[j,m,t],
t≤0&&m==s&&j==0,(1-ϵ*ζ)*μ[0]*(E^(-t))*(E^(-A)),
t≤0&&m==s&&1≤j≤s-1,(1-ϵ*ζ)*β*(E^(γ0))*(E^(-A)),
t≤0&&m==s&&j==s,0,
t≤0&&m==s&&s+1≤j≤k,(1-ϵ*ζ)*β*(E^(γ0))*(E^(-t))*(E^(-A))];

μf1[j_,m_,t_]:=Which[m==s,0,
t≤0&&m≠s&&j≠s,(1-ϵ*β*(E^(γ0)))*μf[j,m,t],
t≤0&&m≠s&&j==s,(1+ϵ*ζ)*μf[j,m,t],
t>0&&m≠s&&j≠s,(1-ϵ*β*(E^(γ0)))*(E^t)*μf[j,m,t],
t>0&&m≠s&&j==s,(1+ϵ*ζ)*μf[j,m,t]];

----------------Define the functions f----------------

f[m_,t_]:=Which[t≤0&&s≤m≤k,0,
t≤0&&1≤m≤s-1,(1-(s-1)*β+(s-2)*β*(E^(γ0))*(E^t))*(E^(-A)),
t>0,(1-(s-1)*β-(k-s+1)*β*(E^(γ0))+(k-1)*β*(E^(γ0))*(E^t))*(E^(-B))];

f0[m_,t_]:=Which[t≤0&&s+1≤m≤k,0,
t≤0&&1≤m≤s,(1-ϵ*(1-β*(E^(γ0))))*((1-β*(E^(γ0))-(s-1)*β)*(E^(-t))+(s-1)*β*(E^(γ0)))*(E^(-A)),
t>0,(1-ϵ*(1-β*(E^(γ0))))*(1-(s-1)*β-(k-s+1)*β*(E^(γ0))+(k-1)*β*(E^(γ0))*(E^t))*(E^(-B))];

f1[m_,t_]:=Which[t≤0&&s≤m≤k,0,
t≤0&&1≤m≤s-1,(1+β*(E^(γ0))*(1-ϵ*β*(E^(γ0)))*((s-2)*(E^t)-(s-1)*(E^(-γ0))))*(E^(-A)),
t>0&&m==s,0,
t>0&&m≠s,(1+β*(E^(γ0))*(1-ϵ*β*(E^(γ0)))*((k-2)*(E^t)-(s-1)*(E^(-γ0))-k+s))*(E^t)*(E^(-B))];

----------------To compute f^2 term, Result is recorded in S2----------------

Clear[R,P,Q,S2];

R[m_,t_]:=(((f[m,t]))^2)-((((f0[m,t]))^2)+(((f1[m,t]))^2))/2;

P=Simplify[Sum[Integrate[R[m,t],{t,0,γ1},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

Q=Simplify[Sum[Integrate[R[m,t],{t,-γ0,0},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

S2=Simplify[P+Q,(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)];

Print[S2];

----------------To compute μf[0]^2 term, Result is recorded in S02----------------

Clear[R,P,Q,S02];

R[m_,t_]:=(((μf[0,m,t]))^2)-((((μf0[0,m,t]))^2)+(((μf1[0,m,t]))^2))/2;

P=Simplify[Sum[Integrate[R[m,t],{t,0,γ1},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

Q=Simplify[Sum[Integrate[R[m,t],{t,-γ0,0},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

S02=Simplify[P+Q,(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)];

Print[S02];

----------------To compute sum of μf[0]*f term, Result is recorded in S0f----------------

Clear[R,P,Q,S0f];
R[m_,t_]:=μf[0,m,t]*f[m,t]-(μf0[0,m,t]*f0[m,t]+μf1[0,m,t]*f1[m,t])/2;
P=Simplify[Sum[Integrate[R[m,t],{t,0,γ1},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];
Q=Simplify[Sum[Integrate[R[m,t],{t,-γ0,0},
Assumptions->(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];
S0f=Simplify[P+Q,(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)];

Print[S0f];

----------------To compute sum of μf[j]^2 term, Result is recorded in SS2----------------

Clear[R,P,Q,SS2];

R[j_,m_,t_]:=(((μf[j,m,t]))^2)-((((μf0[j,m,t]))^2)+(((μf1[j,m,t]))^2))/2;

P=Simplify[Sum[Sum[Integrate[R[j,m,t]*Boole[j≠m],{t,0,γ1},
Assumptions->(j∈Integers)&&(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->(j∈Integers)&&(m∈Integers)&&s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{j,1,k},
Assumptions->(j∈Integers)&&(m∈Integers)&&s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

Q=Simplify[Sum[Sum[Integrate[R[j,m,t]*Boole[j≠m],{t,-γ0,0},
Assumptions->(j∈Integers)&&(m∈Integers)&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{m,1,k},
Assumptions->(j∈Integers)&&(m∈Integers)&&s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)],
{j,1,k},
Assumptions->(j∈Integers)&&(m∈Integers)&&s∈Integers&&s>2&&k>s&&(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)]];

SS2=Simplify[P+Q,(0<ϵ<0.5)&&(0<β<1)&&(γ0>0)&&(γ1>0)];

Print[SS2];

----------------To compute the final result for Internal, Result is recorded in IC----------------

Clear[R,IC];

R = ((1/(2*(1-β)))+((k-3)/(2*(1-2*β))))*S2-(k/(2*(1-k*β)))*S02+((1/(2*(1-2*β)))-((k-1)/(2*β)))*SS2+(1/(1-2*β))*S0f;
IC =Simplify[R,k≠s+1];

Simplify[IC/.{ γ0->0,γ1->0,ϵ->0}]

Clear[g];
g[ϵ_]:=IC;
Simplify[D[g[ϵ],ϵ,NonConstants->{γ0,γ1}]/.{D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]

Simplify[D[g[ϵ],{ϵ,2},NonConstants->{γ0,γ1}]/.{D[γ0,{ϵ,2},NonConstants->{γ0,γ1}]->1-2*β,D[γ1,{ϵ,2},NonConstants->{γ0,γ1}]->2*β-1,D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]

Simplify[D[g[ϵ],{ϵ,3},NonConstants->{γ0,γ1}]/.{D[γ0,{ϵ,3},NonConstants->{γ0,γ1}]->2-10*β+8*(β^2),D[γ1,{ϵ,3},NonConstants->{γ0,γ1}]->2+6*(β^2),D[γ1,{ϵ,2},NonConstants->{γ0,γ1}]->2*β-1,D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]

----------------To compute the final result for External, Result is recorded in ICExt----------------

Clear[R,ICExt];

R = (1/(2*(1-β)))*S2-(1/(2*(1-k*β)))*S02-(1/(2*β))*SS2;
ICExt =Simplify[R,k≠s+1];

Simplify[ICExt/.{ γ0->0,γ1->0,ϵ->0}]

Clear[g];
g[ϵ_]:=ICExt;
Simplify[D[g[ϵ],ϵ,NonConstants->{γ0,γ1}]/.{D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]

Simplify[D[g[ϵ],{ϵ,2},NonConstants->{γ0,γ1}]/.{D[γ0,{ϵ,2},NonConstants->{γ0,γ1}]->1-2*β,D[γ1,{ϵ,2},NonConstants->{γ0,γ1}]->2*β-1,D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]

Simplify[D[g[ϵ],{ϵ,3},NonConstants->{γ0,γ1}]/.{D[γ0,{ϵ,3},NonConstants->{γ0,γ1}]->2-10*β+8*(β^2),D[γ1,{ϵ,3},NonConstants->{γ0,γ1}]->2+6*(β^2),D[γ1,{ϵ,2},NonConstants->{γ0,γ1}]->2*β-1,D[γ0,ϵ,NonConstants->{γ0,γ1}]->1,D[γ1,ϵ,NonConstants->{γ0,γ1}]->1, 
γ0->0,γ1->0,ϵ->0}]