PDE propagating from point source

2015 7 11
1 个回答
回答已采纳
更新时间：2018 2 2
5 次查看（30 天）

0 个投票

I want to solve numerically a nonlinear diffusion equation from an instantaneous point source. Thus, I have initial conditions, but not boundary conditions. How should I go about writing a code to solve circular propagation from a point?
Thanks!!

2 个评论
显示 无 隐藏 无

Torsten
Torsten 2015-7-13
What is the equation you try to solve (because you are talking about a nonlinear diffusion equation) ?
Best wishes
Torsten
María Jesús
María Jesús 2015-7-14
$\frac{\partial C}{\partial t}=r^{1-s}\frac{\partial}{\partial r}[r^{s-1}D\frac{\partial C}{\partial r})]$ where $s$ is constant and $D=D_0(\frac{\partial C}{\partial C_0})^n$ and $n>0$

请先登录，再进行评论。

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

 采纳的回答

Torsten
Torsten 2015-7-14

0 个投票

I assume you want the point source appear at r=0.
Choose the interval of integration as [0:R] where R is big enough to ensure that C=C(t=0) throughout the period of integration.
As initial condition, choose an approximation to the delta function.
As boundary conditions, choose dC/dr = 0 at both ends.
Best wishes
Torsten.

6 个评论
显示 4更早的评论 隐藏 4更早的评论

Torsten
Torsten 2015-7-14
... and try MATLAB's pdepe to solve.
Best wishes
Torsten.
María Jesús
María Jesús 2015-7-23
编辑：María Jesús 2015-7-23
thanks! just one question: how do I specify the delta function? Also, do you mean Kronecker or dirac delta?
Torsten
Torsten 2015-7-23
Take a look at the paragraph "Approximations to the identity" under
https://en.wikipedia.org/wiki/Dirac_delta_function
Best wishes
Torsten.
María Jesús
María Jesús 2015-7-25
I've written the code, but it gives errors... can anyone help with this? The code is
clear all
m = 1;
r = linspace(0, 1, 20);
t = linspace(0, 2, 10);
u = pdepe(m, @pat2, @ic1, @bc1, r, t);
surf(r, t, u);
title('Surface plot of solution');
xlabel('Distance r');
ylabel('Time t');
figure
plot(r, u(end,:))
title('Solution at t = 2')
xlabel('Distance r')
ylabel('u(r,2)')
function [c, b, s, D] = pat2(r, t, u, DuDr)
n = 1;
D_0 = 1;
D = D_0/(KronD(r, 0))^n;
c = 1;
b = D*u^n*DuDr;
s = 0;
function [d] = KronD(j,k)
if nargin < 2
error('Too few inputs. See help KronD')
elseif nargin > 2
error('Too many inputs. See help KronD')
end
if j == k
d = Inf;
else
d = 0;
end
function [pl, ql, pr, qr] = bc1(rl, ul, rr, ur, t)
pl = ul;
ql = 0;
pr = ur;
qr = 0;
function value = ic1(r)
value = KronD(r, 0);
And the error is
Warning: Failure at t=0.000000e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (7.905050e-323) at time t. > In ode15s (line 668) In pdepe (line 289) In pattle2_0 (line 7) Warning: Time integration has failed. Solution is available at requested time points up to t=0.000000e+00. > In pdepe (line 303) In pattle2_0 (line 7) Error using surf (line 57) Z must be a matrix, not a scalar or vector.
Error in pattle2_0 (line 8) surf(r, t, u);
Thanks again!!!
Torsten
Torsten 2015-7-27
1. You will have to work with a numerical approximation to the delta function. I gave you a suitable link.
2. Your boundary conditions are incorrect. You will have to set
pl=0, ql=1, pr=0, qr=1
3. I don't understand your definition of D. The setting
D = D_0/(KronD(r, 0))^n;
doesn't make sense.
Best wishes
Torsten.
Torsten, regarding the earlier answer you provided, whats the reasoning behind using the dirac delta approximation for the point source rather than just setting the initial condition to the source concentration/density as u0(x==0) = initial condition?

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Eigenvalue Problems 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by