我在看对流扩散方程的解析解
初始条件为
和狄利克雷边界条件
纽曼边界条件
但是,当对上述解析解进行编码时,我会得到 C 的负值,这是不现实的。
function AnalyticalSoln2()
format long
R = 1;
D = 900;
v = 200;
L = 60;
co = 4100;
ci = 4000;
X = linspace(0,60,10);
t = 0:0.001:2;
sol=[];
for pos = 1:length(X)
x = X(pos);
A1 = 0.5*erfc((x.*R-t.*v)./(2*(t.*D*R).^0.5));
A2 = 0.5*exp(x.*v/D).*erfc((x.*R+t.*v)./(2*(D*t.*R).^0.5));
A31 = 0.5*(2+v*(2*L-x)/D + (t.*v^2)/(D*R))*exp(v*L/D);
A32 = erfc((R*(2*L-x)+t.*v)./(2*(D*t.*R).^0.5));
A41 = -((t.*v^2)./(pi*D*R)).^0.5;
A42 = exp((v*L/D)-(R./(4*t.*D)).*(2*L -x + t.*v/R).^2);
A = A1 + A2 + A31.*A32 + A41.*A42;
if t==0
C = ci + (co - ci)*A'
else
C = ci + (co - ci)*A' - co*A';
end
sol = horzcat(sol,C);
end
sol(1:100:end,:)
end
而从 pdepe 求解器获得的数值解是非负的。
这是使用 pdepe 获得的解决方案
function DiffusionConvectionMATLAB
format short
global D
m = 0;
x = linspace(0,60,10);
t = 0.1:0.1:2;
D = 900;
sol = pdepe(m,@pdefun,@icfun,@bcfun,x,t)
function [g,f,s] = pdefun(x,t,c,DcDx)
v = 200;
g = 1;
f = D*DcDx;
s = -v*DcDx;
end
function c0 = icfun(x)
c0 = 4000;
end
function [pl,ql,pr,qr] = bcfun(xl,cl,xr,cr,t)
pl = cl-4100;
ql = 0;
pr = 0;
qr = 1;
end
end
有人可以建议我从分析表达式计算解决方案的方式是否有任何问题?