1D зависящий от времени массоперенос и реакция в очень маленьком газовом фильтре

Question

Я пытаюсь смоделировать PDE, включая адвекцию-диффузию и условия реакции. Я хочу видеть профиль концентрации после долгого времени (2*60*60 сек) в очень крошечной пластинке (5 микрометров). Я использовал метод линии (ode15s), чтобы решить этот PDE, но мой результат не совпадает с экспериментальным результатом. Потому что концентрация на выходе получает максимальную концентрацию сразу (как видно на моем графике). У кого-нибудь есть предложение?

clc
clear 
close all
%%%%%%%%%%%parameters%%%%%%%%%%%%%%%%
global a b  g  Ki dx cin

Lf=5e-6;     %microM
dm=8.6*10^-6;   %m2/min

area=0.65;       % m^-1
cin=20; % 
ki=1;           %mg/m3 min
Ki=0.24;        %m3/mg
ux=0.5;       %m/min
%%%%%%%%%%%%%%%%%%constant%%%%%%%%%%%%%%%%%%%
a=dm;
b=ux;
g=ki*Ki;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
tf=2*60*60;
nx=50;
xf=Lf;
dx=xf/nx;
nt=1000;
dt=tf/nt;

x = 0:dx:xf;
t = 0:dt:tf;
N = numel(x);
u0 = zeros(N,1);
y0 = u0 ;

[t,y] = ode15s(@(t,y)fun(t,y,x,N),t,y0);
u = y(:,1:2:end);
x = (1:N)/(N+1);
plot(t/60,u(:,end))

function dydt=fun(~,y,~,N)
  global a b  g dx cin
  bc=cin;
  dydt = zeros(N,1);
  i = 1;
  dydt(i,:)=a*(y(i+1,:)-2*y(i,:)+bc)/(dx^2)-b*(y(i+1,:)-bc)./(2*dx)-g.*y(i+1,:);

  i = 2:N-1;
  dydt(i,:)=a*(y(i+1,:)-2*y(i,:)+y(i-1,:))/(dx^2)-b*(y(i+1,:)-y(i-1,:))./(2*dx)-g.*y(i,:);

  i = N;
  dydt(i,:)=a*(y(N,:)-2*y(i,:)+y(i-1,:))/(dx^2)-b*(y(N,:)-y(i-1,:))./(2*dx)-g.*y(i,:);

end