function fK = fKJR(phi, e_R)
% This function calculates the normalized heat conductivity K from the 
% Jenkins & Richman (1985) theory.
% The normalized conductivity is defined as
% 
% fK = K / [(25*pi^(1/2)/128)* rho_s * d * T^(1/2)]
%
% syntax:
%   fK = fKJR(phi, e_R)
%
% inputs: volume fraction phi and the coefficient of restitution e_R, 
%
% output: normalized conductivity fK, its size is [length(phi), length(e_R)]
%
Nphi = length(phi);
Ne = length(e_R);
feta = zeros(Nphi, Ne);

chi = chiCS(phi);

for i=1:Ne
    ei = e_R(i);
    E1 = 32 / ((1+ei)*(49-33*ei));
    E2 = (96/5)*(1+ei+2*ei*ei) / (49-33*ei);
    E3 = 64*(1+ei)/(25*(49-33*ei))*( 9*(1+ei)*(2*ei-1)+(4/pi)*(49-33*ei));
    
    fK(:,i) = E1./chi + E2*phi + E3*(phi.^2).*chi;
end
