% Thermo-electric current

% Physical constants
hbar = 1.054e-34;
q = 1.602e-19;

% Parameters (eV)
kBT_1 = 0.025;
kBT_2 = 0.026;
mu = 0;
gamma_1 = 0.005;
gamma_2 = gamma_1;
gamma_sum = gamma_1 + gamma_2;

% Channel energy levels, varying between -0.25eV and 0.25eV
epsilon = linspace(-0.25, 0.25, 101);

% Energy grid
E = linspace(-1, 1, 501);

% Contact fermi functions
f_1 = 1 ./ (1 + exp((E - mu)./kBT1));
f_2 = 1 ./ (1 + exp((E - mu)./kBT2));

% Iterate through channel energy levels
for n = 1:length(epsilon)
    % Compute energy level density functions - integral normalized to unity
    D = (gamma./(2*pi))./((E-cal_E(count)).ˆ2+((gamma./2).ˆ2));
    D = D./(dE*sum(D));


end