function [Dgs_moldrup] = lijun_co2_code (CO2, Moist, T)

CO2=CO2;
Moist=moist;
T=temp;

%parameters
R=8.314;
D0=1.47*10^(-5);
P=1.013*10^5;
z=[0.02 0.08 0.16];
depth=[0 0.05 0.12];
porosity_t=0.735  ;
Fcp=0.88;
%
L=length(CO2);
Coeff=zeros(L,3);
CO2_mol=zeros(L,3);
elta=zeros(L,3);
Dgs_penman=zeros(L,3);
Dgs_millington=zeros(L,3);
Dgs_moldrup=zeros(L,3);
F_penman=zeros(L,3);
F_millington=zeros(L,3);
F_moldrup=zeros(L,3);

%
for i=1:L
  T_soil=mean(T(i,:));
  D0=D0*((273.15+T_soil)/(273.15+20))^1.5;
  CO2_mol(i,:)=CO2(i,:)*P/(R*(273.15+T_soil));
  CO2_gradient=CO2_mol(i,:);
  Coeff(i,:)=polyfit(z,CO2_gradient,2);
  gradient(i,:)=2*Coeff(i,1)*depth+Coeff(i,2);
  elta(i,:)=porosity_t-Moist(i,:);
  Dgs_penman(i,:)=D0*0.66.*elta(i,:);
  Dgs_millington(i,:)=D0*elta(i,:).^(10/3).*porosity_t^(-2);
  Dgs_moldrup(i,:)=D0*porosity_t^2.*(elta(i,:).*porosity_t^(-1)).^(2.9*Fcp);
  F_penman(i,:)=Dgs_penman(i,:).*gradient(i,:);
  F_millington(i,:)=Dgs_millington(i,:).*gradient(i,:);
  F_moldrup(i,:)=Dgs_moldrup(i,:).*gradient(i,:);
  D0=1.47*10^(-5);
end;

end

%{
fid = fopen('F_penmman.txt','w');
fprintf(fid,'6.2%f %3f\n',F_penman);
fclose(fid);
%}