Solving partial differential equations using pdesolver

Hello, I am a college student learning MATLAB. I have a problem while studying, so I ask a question.

I want to solve this partial differential equation.

𝜕𝐶/𝜕𝑡=𝐷 (𝜕^2 C)/(𝜕z^2 ) − 𝑢/((1+𝐹𝐻)) 𝜕𝐶/𝜕𝑧

%% fuction script
function [j,f,s] = pde1(z,t,c,dcdz)
D =0.01;
e=0.6;
F=(1-e)/e;
H=1.88;
u=0.05; 
j = 1;
f = (D/(1+(F*H)))*dcdz;
s = -(u/(1+(F*H)))*dcdz;
%p= F*dQdt;
end

%% boundary condition script
function [pl,ql,pr,qr]=pde1BC(zl, cl,zr,cr,t)
pl=zl;
ql=0.2;
pr=zr;
qr=0;
end

%% initial condition script
function c0=pde1ic(z,t)
c0=0.2;
end

%% main script (for 'run')
z=linspace(0,10,100)
t=linspace(0,200,2000)
m=0;
%eqn= @ (z,t,c,q,dcdz,dqdt)
sol=pdepe(m,@pde1, @pde1ic, @pde1BC,t,z)

When I 'run', the following error appears.

An error occurred while using: Failed to discretize pdepe space. Discretization supports only parabolic and elliptic equations with spatial derivative related flux terms.

An error occurred: pde1_main (line 6)

sol=pdepe(m,@pde1, @pde1ic, @pde1BC,t,z)

how can i solve this error

Thanks for reading this long post.

Answer 1

[Bill Greene]

First, in your call to pdepe, you have the t and z arguments switched. It should be

sol=pdepe(m,@pde1, @pde1ic, @pde1BC,z,t);

Second, your boundary condition at the right end of the region(z=10) is incorrect. What is the boundary condition you want to impose at this point? If, for example, it is the Dirichlet condition c(10)=0, you should replace

pr=zr;

with

pr=cr;

[Torsten]

@인범 황

And - if your equation from the beginning is correct -

f = D*dcdz;

instead of

f = (D/(1+(F*H)))*dcdz;

Note that you have to modify your boundary condition setting

pl=zl;
ql=0.2;

accordingly.

[인범 황]

Thanks for the two replies.

I have one more question.

If the boundary condition is as follows, how should I write the boundary condition script?

I've tried very hard, but I can't get it to work, so I'm asking.

Thank you.

[Bill Greene]

I modified my answer to reflect your actual outlet boundary conditon. The code for the inlet boundary condition is simply

ql=0;

pl=cl-ci0;

I am assuming you have determined that your domain length of ten is sufficiently large that it approximates infiniity for this particular problem.

Apparently you found the pdepe documentation for defining boundary conditions confusing. Can you say more about what confused you? I have written a document describing pdepe boundary condition definitions in a more tutorial fashion. Can you please look at this document and let me know if you find it more or less confusing? Thanks.

[Torsten]

@인범 황

There is no mass sink in your equation. So the boundary condition c=0 at x=Inf seems inadequate to me. Better use dc/dx = 0.