function demoCoef
% compute the transport coefficients from the Jenkins-Richman theory

% vary volume fraction phi, fix coefficient of restitution
e_R = [1 0.9 0.8 0.7];
phi = [0:0.01:0.5]';

Ne = length(e_R);
N = length(phi);

% Jenkins & Richman (Arch. Ration. Mech. Anal., 1985)
etaJR_e = fetaJR(phi, e_R);     % normalized shear viscosity
kJR_e = fKJR(phi, e_R);       % normalized thermal conductivity
dissipJR_e = fDissipJR(phi, e_R);      % normalized volumetric dissipation rate

linesJR = {'r-', 'g-', 'b-', 'm-', 'c', 'y-', 'k-'};

% shear viscosity
figure;
legendtxt = [];
for i = 1:Ne
    plot(phi, etaJR_e(:,i), char(linesJR(i)));
    hold on;
    legendtxt = strvcat(legendtxt, strcat('e = ', num2str(e_R(i))));
end
legend(legendtxt);
xlabel('volume fraction \phi');
ylabel('normalized shear viscosity \eta*(e, \phi)');

% thermal conductivity
figure;
for i = 1:Ne
    plot(phi, kJR_e(:,i), char(linesJR(i)));
    hold on;
end
legend(legendtxt);
xlabel('volume fraction \phi');
ylabel('normalized heat conductivity \kappa*(e, \phi)');

% collisional dissipation rate
figure;
for i = 1:Ne
    plot(phi, dissipJR_e(:,i), char(linesJR(i)));
    hold on;
end
legend(legendtxt);
xlabel('volume fraction \phi');
ylabel('normalized dissipation rate \gamma*(e, \phi)');


end
        

