MATLAB: Solve 1-D interfacial mass transfer using pdepe

mass transferpdepe

Hi Everyone,
I am trying to use pdepe to solve a diffusion problem and Im having issues trying to set my left side boundary condition.
dP=0.04; %thickness polymer layer [cm]
d=100; %times in days [d]
tt=d*86400; %time in seconds [s]
x = linspace(0,dP,5);
t = linspace(0,tt,100);
m = 0;
sol = pdepe(m,@equa,@IC,@BC,x,t);
u =sol(:,:,1)
function [c,f,s,algo] = equa(x,t,u,DuDx)
c = 1;
f = D*DuDx;
s = 0;
function u0 = IC(x)
u0 = Cpo; %Initial concentration[microg/cm^3]
function [pl,ql,pr,qr] = BC(xl,ul,xr,ur,t)
pl =0;
ql =1;
pr =h*((ur/K)-Cinf);
qr =D;
the first boundary condition (0,t)
the second boundary condition (x,t)
Rigth now looks like it is working using Cinf as a constant but actually Cinf should change and increase with time (accumulation).
Cinf=A/V*integral(Cp/x (t) dt, 0,t)
I really dont know how to solve it this way
Could someone please guide me?

Best Answer

First of all one error in your actual code:
qr=1 instead of qr=D.
If you set qr=D, your boundary condition reads
Concerning your question about Cinf there are two ways to solve this problem:
1. Discretize your PDE equation in space, add the equation
dCinf/dt = A/V*dCp/dx(@x=dP)
to the system and solve all the equations using ODE15S. Look up "method-of-lines" for more details.
2. Iterative procedure:
Fix a list of constant output times t_i.
a) Solve your equation for a constant Cinf.
b) From the solution, evaluate dCp/dx(@x=dP) at the output times t_i and use "cumtrapz" to build an approximation to Cinf(t_i)=A/V*integral_{t=0}^{t=t_i} dCp/dx(@x=dP) dt.
c) Now hand these two vectors (t_i and C_inf(t_i)) to BC and use "interp1" to interpolate Cinf at the time t requested by the solver.
d) Compare with the first solution. If error<eps, stop, else go to step b)
Best wishes