This is a follow up to my previous post

I've been trying to compare the diffusion time obtained from theoretical derivation(answered in my previous post) and what is obtained computationally, for a 1D diffusion equation.

Theoretically, the diffusion time is given by

$$t_D = \frac{l^2}{2D}$$ Calculating $t_D$ , l = 5 nm and D = 500 $nm^2/min$ gives $t_D$ = 0.025 min.

The following is the MATLAB code that simulates 1D diffusion system using pdepe solver.

function sol=so()
format short
global D nnode init_co find_index
m = 0;
xend = 5; 
D = 500; 
x = 0:1:xend;
find_index  = 0:1:xend;
t = 0:0.00001:0.5;
init_co = 1*ones(length(x),1);
nnode = length(x);
sol = pdepe(m,@pdefun,@icfun,@bcfun,x,t)
function [g,f,s] = pdefun(x,t,c,DcDx)
g = 1;
f = D*DcDx;
s = 0;

function c0 = icfun(x)
c0 = init_co(find(find_index==x));

function [pl,ql,pr,qr] = bcfun(xl,cl,xr,cr,t)
% Dirichlet at left (concentration at left boundary = 2, Neumann at right(dC/dx = 0)
    pl = cl - 3;
    ql = 0;
    pr = 0;
    qr = 1;

Result: enter image description here From the plot of C vs time, all the curves reach a steady-state at 0.1 min.

Could someone explain why there is a difference in what is computed theoretically and computationally? Is it appropriate to make this comparison?

  • 1
    $\begingroup$ “All curves reach a steady state at 0.1”. That’s simply not true. $\endgroup$ – David Ketcheson Jan 10 '20 at 15:18
  • 1
    $\begingroup$ My comment had nothing to do with the labels. It still stands. The curves exponentially approach the steady state. They certainly do not reach it at 0.1. $\endgroup$ – David Ketcheson Jan 11 '20 at 15:25
  • 1
    $\begingroup$ I think what you really need is to understand what "diffusion time" means. This is not a computational question at all; it's a mathematical definition. You seem to interepret it as "the time at which a plot of the solution will look like it is very close to the steady state", but that's too imprecise to be a mathematical definition. $\endgroup$ – David Ketcheson Jan 12 '20 at 13:15
  • 1
    $\begingroup$ "I intend to check for the time that a substance takes to diffuse through a length L" -- Again, you need a precise definition first. Checking when the graph of the value looks close to the steady state is too imprecise for any kind of mathematical or scientific work. A useful definition would be something like "the time at which the concentration at a distance L reaches 50% of the concentration at the point of origin". The 50% is not important, but you need something precisely measurable. $\endgroup$ – David Ketcheson Jan 12 '20 at 14:14
  • 1
    $\begingroup$ I think @DavidKetcheson is right here. For that 50%, I suggest to replace it with $100 \times (1-e^{-1})$%, which is approximately 63%. It comes from the fact that most popular kinetics theory to describe a diffusion is exponential kinetics or change of concentration with respect to time must be: $C(t) = C_{\infty} (1-\exp(-\frac{t}{\tau}))$, where $C_{\infty}$ is your concentration at $t = \infty$ and $\tau$ is your relaxation time or kinetics parameter that is directly related to diffusion coefficient. Now, when you want to know when system will reach steady state: $t = \tau$. $\endgroup$ – Alone Programmer Jan 13 '20 at 3:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.