Using PDEPE for Ficks second law ?

16 次查看（过去 30 天）
Stefan
Stefan 2021-6-10
回答： Stefan 2021-6-14
Hello together,
i tried to use the PDEPE-Solver for Fick's second law of diffusion:
But I realy got confused with the initial and boundary conditions.
Is it possible to set the initial condition as c(x>0,0)=0
and the two boundaries as c(x=0,t)=10 und dc(x=l,t)/dt=0 ?
clc
clear
global dcoeff flux
L = 0.01; % m length of slab
tend = 36000; % s total of 10 hours
dcoeff = 10^-10; % m^2/s
flux = 10; % surface concentration wt%
m = 0;
x = linspace(0,L,20);
t = linspace(0,tend,20);
sol = pdepe(m,@pdefun,@icfun,@bcfun,x,t);
u = sol(:,:,1);
surf(x,t,u)
xlabel('distance')
ylabel('time')
figure
ptot(x,u(:,:))
xlabel('distance')
ylabel('wt%')
figure
plot(t,u(:,:))
xlabel('time')
ylabel('wt%')
function [c,f,s] = pdefun(x,t,u,DuDx)
global dcoeff
c = 1/dcoeff;
f = DuDx;
s = 0;
end
function u0 = icfun(x)
if x > 0
u0 = 0;
else
u0 = 10,
end
end
function [pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)
pl = ul-10;
ql = 0;
pr = 0;
qr = 1;
end
Thank you for every hint.
  2 个评论
Bill Greene
Bill Greene 2021-6-11
You say you want a BC . Is this a typo? This is a strange BC that would have no other effect than to set the value of c there equal to the initial condition.
Stefan
Stefan 2021-6-14
Hello Bill,
thank you for your comment! To be honest, I tried to reproduce following example (Fick's 2. Law) as an easy exercise for the beginning. I thought the BC: would set a constant c if x >= l ?
The biggest issue at the moment is that the If-Loop as IC isn't working right, it always sets u0=10 for all x at t=0 but i would like to set it as described in the example described above.
(To avoid confusio u in the code corresponds to c above)
function u0 = icfun(x)
if x > 0
u0 = 0;
else
u0 = 10,
end
Thank you very much in advance

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Stefan
Stefan 2021-6-14
Ok I could solve the problem by myself!
I just mixed the IC and BC up! The If-Loop wasn't necessary.
I added the Code if someone is face the same topic!
Regards, Stefan
clc
clear
global dcoeff flux
L = 10; % m length of slab
tend = 50; % s total of 10 hours
dcoeff = 1; % m^2/s
flux = 10; % surface concentration wt%
m = 0;
x = linspace(0,L,100);
t = linspace(0,tend,100);
sol = pdepe(m,@pdefun,@icfun,@bcfun,x,t);
u = sol(:,:,1);
surf(x,t,u)
xlabel('distance')
ylabel('time')
figure
contourf(x,t,u,[1:9],'ShowText','on')
ax = gca;
ax.YDir = 'reverse';
xlabel('x')
ylabel('time')
function [c,f,s] = pdefun(x,t,u,DuDx)
global dcoeff
c = 1/dcoeff;
f = DuDx;
s = 0;
end
function u0 = icfun(x)
u0 = 0;
end
function [pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)
pl = ul-10;
ql = 0;
pr = 0;
qr = 1;
end

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 1-D Partial Differential Equations 的更多信息

标签

产品


版本

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by