load dos_si_per.out;

% trapzoidal rule
ss = dos_si_per(1,2)*dos_si_per(1,1)/2.;
for i=2:size(dos_si_per,1)-1,
 ss = ss+(dos_si_per(i,2)+dos_si_per(i-1,2))/2.*(dos_si_per(i,1)-dos_si_per(i-1,1));
end;

fprintf (1,'the sum for Si in perfect crystal is %e \n',ss);

load dos_c_per.out;

% trapzoidal rule
ss = dos_c_per(1,2)*dos_c_per(1,1)/2.;
for i=2:size(dos_c_per,1)-1,
 ss = ss+(dos_c_per(i,2)+dos_c_per(i-1,2))/2.*(dos_c_per(i,1)-dos_c_per(i-1,1));
end;

fprintf (1,'the sum for C in perfect crystal is %e \n',ss);

load dos_si_def.out;

% trapzoidal rule
ss = dos_si_def(1,2)*dos_si_def(1,1)/2.;
for i=2:size(dos_si_def,1)-1,
 ss = ss+(dos_si_def(i,2)+dos_si_def(i-1,2))/2.*(dos_si_def(i,1)-dos_si_def(i-1,1));
end;

fprintf (1,'the sum for Si in defect crystal is %e \n',ss);

load dos_c_def.out;

% trapzoidal rule
ss = dos_c_def(1,2)*dos_c_def(1,1)/2.;
for i=2:size(dos_c_def,1)-1,
 ss = ss+(dos_c_def(i,2)+dos_c_def(i-1,2))/2.*(dos_c_def(i,1)-dos_c_def(i-1,1));
end;

fprintf (1,'the sum for C in defect crystal is %e \n',ss);

