T = 1273; 
E = 8.617e-5 * T;

% Neutral vacancies
D0_0 = 0.066;
E0_0 = 3.44;

D0 = D0_0 * exp(-E0_0/E);

% Negative vacancies
Dn_0 = 22.9;
En_0 = 4.1;

Dn = Dn_0 * exp(-En_0/E);

% Temperature dependence
n_D = 10^20;
n_i = ni(T);

D = D0 + (n_D/n_i) * Dn;

