function [sX,vx] = fContestEquilibriumTullock(vN)
  % This function computes the equilibria in sequential contests as
  % described in "Optimal sequential contests" by Toomas Hinnosaar, 
  % avalable at http://toomas.hinnosaar.net/contests.pdf
  % Matlab code written by Toomas Hinnosaar (toomas@hinnosaar.net)
  % Inputs: vN     = [n_1,...,n_T]
  [sXapprox,vxapprox] = fContestEquilibriumLinear(vN,1);
  bUseapprox = 0;
  sT     = max(size(vN));
  mS     = zeros(sT,sT);   % Matrix to hold all information measures
  mG     = zeros(sT,sT+2); % Matrix to hold g_k functions, which are 
                           % polynomials of degrees 2,...,T+1
  mF     = zeros(sT+1,sT+2); % Matrix to hold f_T;...;f_0                         
  mS(1,1)    = vN(sT); % Information before the last period is just S=(n_T) 
  mG(1,sT)   = -1;  % g_1(X) = X(1-X) = -1 X^2 + 1 X + 0
  mG(1,sT+1) = 1;
  mF(1,sT+1) = 1; % f_T(X) = X
  for k=2:sT
    % S_k(n^{t-1}) = S_k(n^t) + n_t S_{k-1}(n^t)  
    mS(k,:) = mS(k-1,:)+[1 mS(k-1,1:sT-1)]*vN(sT+1-k);
    % g_k(X) = -g_k'(X) X (1-X) = X^2 g_k'(X) - X g_k'(X)  
    vPD = polyder(mG(k-1,:));
    mG(k,sT+1-k:sT)   = vPD;
    mG(k,sT+2-k:sT+1) = mG(k,sT+2-k:sT+1)-vPD;
    % f_t = X - \sum_{k=1}^T S_k(n^t) g_k(X)
    mF(k,:) = mF(1,:)-mS(k-1,:)*mG;
    sX = max(roots(mF(k,:)));
    if max(abs(1-sX),abs(sXapprox-sX))<1e-6 
      bUseapprox = 1;
      disp('Approximation seems to be good enough ');
      break
    end
  end 
  if bUseapprox>0
    vx = vxapprox;
    sX = sXapprox;
  else 
    mF(k+1,:) = mF(1,:)-mS(k,:)*mG; % f_0(X)
    sX = max(roots(mF(sT+1,:)));
    vx = zeros(sum(vN),1);
    sfkX = sX;
    for k=1:sT
      t = sT+1-k;  
      sfkX1 = polyval(mF(k+1,:),sX);
      vx(sum(vN(1:t-1))+[1:vN(t)]) = (sfkX-sfkX1)/vN(t);
      sfkX  = sfkX1;
    end
  end
end
