PDEPE boundary condition outputting odd results.
1 回表示 (過去 30 日間)
古いコメントを表示
For the boundary condition I've set the left (pl) to a constant value of 280 with the idea being that this value will always be 280, however when i output the results the value is instead 4e10 and appears to be changing. What am I doing incorrectly so that my boundary condition is a) not constant and b) way above the value I've stated.
I've gone over the pdepe page a few times trying to figure it out but unfortauntely haven't made much progress.
Thank you for any help provided.
m = 0;
tspan = linspace(0,365*20*86400,21);
xmesh = linspace(0,20,201);
sol = pdepe(m,@pdefun,@icfun,@bcfun,xmesh,tspan);
u = sol(:,:,1);
figure, plot(u(2,:),xmesh)
ylabel('Depth (m)')
xlabel('Oxygen Concentration (g/m3)');
set(gca, 'XAxisLocation', 'top')
set(gca, 'YDir','reverse')
function[c,f,s]=pdefun(x,t,u,DuDx)
c = 0.18;
f = 1e-8*DuDx;
s = -1e-9*u;
end
function u=icfun(x)
u = 0;
end
function[pl,ql,pr,qr]=bcfun(xl,ul,xr,ur,t)
pl = 280;
ql = 1;
pr = ur;
qr = 0;
end
4 件のコメント
採用された回答
Torsten
2019 年 3 月 22 日
編集済み: Torsten
2019 年 3 月 22 日
function main
m = 0;
tspan = linspace(0,365*20*86400,21);
xmesh = linspace(0,20,201);
sol = pdepe(m,@pdefun,@icfun,@bcfun,xmesh,tspan);
u = sol(:,:,1);
figure, plot(u(2,:),xmesh)
ylabel('Depth (m)')
xlabel('Oxygen Concentration (g/m3)');
set(gca, 'XAxisLocation', 'top')
set(gca, 'YDir','reverse')
end
function[c,f,s]=pdefun(x,t,u,DuDx)
c = 0.18;
f = 1e-8*DuDx;
s = -1e-9*u;
end
function u=icfun(x)
u = 0;
end
function[pl,ql,pr,qr]=bcfun(xl,ul,xr,ur,t)
u0 = 280;
pl = ul - u0;
ql = 0;
pr = ur;
qr = 0;
end
2 件のコメント
Torsten
2019 年 3 月 22 日
Yes, call pdepe twice: first from t=0 to t=t_switch, then from t=t_switch to t=tend. And use the solution obtained at the end of the first time interval as initial profile for the second time interval.
その他の回答 (2 件)
Scott Lines
2019 年 3 月 22 日
1 件のコメント
Torsten
2019 年 3 月 22 日
編集済み: Torsten
2019 年 3 月 22 日
function main
m = 0;
xmesh = linspace(0,20,101);
tspan = linspace(0,365*5*86400,201);
u0 = 280.0;
icarg = @(x) 0.01*ones(size(x));
sol = pdepe(m,@pdefun,@(x)icfun(x,icarg),@(xl,ul,xr,ur,t)bcfun(xl,ul,xr,ur,t,u0),xmesh,tspan);
w = sol(end,:);
plot(xmesh,w)
tspan2 = linspace(tspan(end),365*20*86400,201);
u0 = 0.0;
icarg = @(x)interp1(xmesh,w,x);
sol2 = pdepe(m,@pdefun,@(x)icfun(x,icarg),@(xl,ul,xr,ur,t)bcfun(xl,ul,xr,ur,t,u0),xmesh,tspan2);
w2 = sol2(1,:);
hold on
plot(xmesh,w2)
end
function [c,f,s] = pdefun(xmesh,tspan,u,DuDx)
c = 2;
f = 1e-8*DuDx;
s = -1e-7*u;
end
function u = icfun(x,icarg)
u = icarg(x);
end
function [pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t,u0)
pl = ul - u0;
ql = 0;
pr = ur;
qr = 0;
end
参考
カテゴリ
Help Center および File Exchange で Environment and Settings についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!