next up previous
Next: About this document ...

7. Absorbing boundaries


LAST TIME
Tim Duty described branching processes e.g. birth and death processes:
System of n objects of type A.
Elementary processes: Discussed methods to calculate survival probability and cascade size distribution.


State with no A- particles absorbing state. Once in such a state system trapped and there is no escape!
Today want to discuss how to treat absorbing states in continuum limit.

LIMBO STATE and FICTITIOUS STATE


Master equation for one step processes:

\begin{displaymath}\frac{\partial P(n,t)}{\partial t}=(E-1)d(n)P(n,t)-(E^{-1}-1)b(n)P(n,t)\end{displaymath}

If b(0)=0, call n=0 limbo state.
P*(t)=probability of limbo state.
If $b(n)\Rightarrow 0$ as $n \Rightarrow 0$ natural boundary.
Otherwise n=0 artificial boundary.


In continuum limit must avoid treating n=0 as separate case.
If boundary artificial introduce fictitious state n=0. Boundary condition P(0,t)=0.
Normalization:

\begin{displaymath}P^*(t)+\sum_{n=0}^\infty P(n,t)=1\end{displaymath}

If $b\propto n$ and $P(n,t)\Rightarrow const.$ as $n \Rightarrow 0$ there is probability current into limbo state!

Example: DIFFUSION NEAR A CLIFF

\begin{figure}
\epsfxsize=400pt
\epsffile{cliff.eps}
\end{figure}
Jumps at rate $\lambda$ distance a in $\pm x$ and $\pm y$ directions.
S absorbing boundary. Master equation:

\begin{displaymath}\frac{\partial P(n,m,t)}{\partial t}=\lambda [E_x+E_x^{-1}+E_y+E_y^{-1}-4]P(n,m,t)\end{displaymath}


\begin{displaymath}x=na;\;\;y=ma\end{displaymath}

a "smallness parameter"
L= macroscopic length scale
$\Omega=L/a$ system size parameter
S= artificial boundary!


\begin{displaymath}\frac{1}{a^2}\Pi(x,y,t)dx dy=P(n,m,t)\end{displaymath}


\begin{displaymath}E_x-1=a\frac{\partial}{\partial x}+a^2\frac{\partial^2}{\partial x^2}
+\cdots\end{displaymath}


\begin{displaymath}E_y-1=a\frac{\partial}{\partial y}+a^2\frac{\partial^2}{\partial y^2}
+\cdots\end{displaymath}

Diffusion equation:

\begin{displaymath}\frac{\partial\Pi(x,y,t)}{\partial t}=D\left[\frac{\partial^2}{\partial x^2}+
\frac{\partial^2}{\partial y^2}\right]\Pi\end{displaymath}

where $D=\lambda a^2$
Boundary condition $\Pi=0$ on $S
\Rightarrow$ eigenvalue problem:

\fbox{\parbox{9cm}{
\begin{displaymath}\Pi(x,y,t)=\sum_{i=0}^\infty \alpha_i e^...
...displaymath}\Pi(x,y,0)=\sum_{i=0}^\infty\alpha_i\psi_i(x,y)\end{displaymath}
}}


Example:
Square with corners at (0,0),(0,L),(L,L),(L,0) and absorbing sides

\begin{displaymath}\Pi(x,y,t)=\sum_{n,m=1}^\infty \alpha_{nm}\psi_{nm}\exp
(-\lambda_{nm}t)\end{displaymath}


\begin{displaymath}\psi_{nm}=\frac{2}{L}\sin\frac{n\pi}{L}\sin\frac{m\pi}{L}\end{displaymath}


\begin{displaymath}\lambda_{nm}=\frac{D(n^2+m^2)\pi^2}{L^2}\end{displaymath}



\begin{displaymath}\alpha_{nm}=\int_0^Ldx\int_0^Ldy\psi_{nm}(x,y)\Pi(x,y,0)\end{displaymath}


For long times smallest eigenvalue (ground state) dominates!

\begin{displaymath}\lambda_{11}=\frac{2D\pi^2}{L^2}\end{displaymath}

Asymptotic formula for survival probability

\begin{displaymath}\Pi_s(t)\approx\frac{8L\alpha_{11}}{\pi^2}\exp(-\frac{2D\pi^2 t}{L^2})\end{displaymath}



Example of problem with natural boundaries:
GENETIC DRIFT
A gene is expressed as allele a or A.
The allele is neutral: there is no competititve advantage to being either A or a.
Stable population size: Population pairs off.
Each pair produces two offspring.

\begin{displaymath}\begin{array}{cccll}
A+A&\Rightarrow&A+A&{\rm probability}&1...
...&\Rightarrow&a+a&&1/4\\
&\Rightarrow&A+A&&1/4\\
\end{array}\end{displaymath}

Only last two processes change composition of population. Rate of
contributing processes:

\begin{displaymath}\beta n(\Omega-n)/\Omega\end{displaymath}

$\Omega =$ size of population.
n = number of individuals of type a.


Master equation:

\begin{displaymath}\frac{\partial P(n,t)}{\partial t}=(E_n+E_n^{-1}-2)\frac{\beta n(\Omega-n)}{\Omega}P(n,t)\end{displaymath}

Macroscopic variables $x=\frac{n}{\Omega};\;\tau=\frac{t}{\Omega}.$

\begin{displaymath}\frac{1}{\Omega}\Pi(x,\tau)=P(n,t)\end{displaymath}

Fokker-Planck equation

\fbox{\parbox{9cm}{
\begin{displaymath}\frac{\partial\Pi(x,\tau)}{\partial\tau}=\beta\frac{\partial^2}{\partial x^2}x(1-x)
\Pi(x,\tau)\end{displaymath}
}}

Initial condition $\Pi(x,0)=\Pi_0(x)$.
Try solution on form

\begin{displaymath}\Pi(x,\tau)=\sum_{m=0}^\infty \alpha_me^{-\lambda_m\tau}\psi_m(x)\end{displaymath}


\begin{displaymath}-\lambda_m \psi_m(x)=\beta\frac{d^2}{d x^2}x(1-x)\psi_m(x)\end{displaymath}

Differential equation singular at x=0,x=1! Cannot impose arbitrary boundary condition.

Require that eigenfunctions well behaved at
x=0,x=1. Try

\begin{displaymath}\psi=\sum_{l=0}a_kx^k\end{displaymath}

Get

\begin{displaymath}a_{k+1}=\left(1-\frac{\lambda}{\beta(k+2)(k+1)}\right)a_k\end{displaymath}

Radius of convergence of power series x=1. For solution to be well behaved at x=1 require that series breaks off! Eigenvalues:

\begin{displaymath}\lambda_m=\beta (m+2)(m+1);\;\;m=0,1,2\cdots\end{displaymath}

Lowest eigenvalue $\lambda_0=2\beta$. Asymptotically, probability that both alleles can be found in population decays as

\begin{displaymath}P_s\propto \exp(-2\beta\tau)=\exp(-\frac{2\beta t}{\Omega})\end{displaymath}

Note that lifetime $\propto \Omega$, population size!
Equation for $\psi$ special case of hypergeometric equation. General solution can be expressed in terms of Gegenbauer polynomials! For details see Kimura [1955]

KRAMERS ESCAPE RATE
Consider an overdamped partice moving in a potential well near a "cliff".
\begin{figure}
\epsfysize=200pt
\epsffile{pot.eps}
\end{figure}
Force

\begin{displaymath}f(x)=-\frac{dU(x)}{dx}\end{displaymath}



Fokker-Planck equation

\begin{displaymath}\frac{\partial P(x,t)}{\partial t}= -\frac{D}{kT}\frac{\partial}{\partial x}f(x) P(x,t)+D\frac{\partial^2P(x,t)}{\partial x^2}\end{displaymath}


Probability current

\begin{displaymath}J=\frac{D}{kT}f(x) P(x,t) -D\frac{\partial P(x,t)}{\partial x}\end{displaymath}


\begin{displaymath}=-D\exp(-\frac{U(x)}{kT})\frac{\partial}{\partial x}\exp(\frac{U(x)}{kT})P(x,t)\end{displaymath}


Integrate from b to e with P(e,t)=0

\begin{displaymath}D\exp(\frac{U(b)}{kT})P(b,t)=J\int_{b}^{e}dx\exp(\frac{U(x)}{kT})\end{displaymath}


Assume that near b

\begin{displaymath}P(x,t)\approx P(b,t)\exp(-\frac{U(x)-U(b)}{kT})\end{displaymath}



Probability that particle is "near" b

\begin{displaymath}p=P(b,t)\exp(\frac{U(b)}{kT})\int_a^c\exp(-\frac{U(x)}{kT})\end{displaymath}


Escape rate

\fbox{\parbox{10cm}{
\begin{displaymath}r=\frac{J}{p}=\frac{D}{\int_a^c dx\exp(\frac{-U(x)}{kT})
\int_b^edy\exp(\frac{-U(y)}{kT})}\end{displaymath}
}}

The integrals in the boxed expression above can be evaluated numerically to arbitrary accuracy. In order to get a simple and easy to remember expression we note that the dominant contribution from first integral comes from near minimum. Expand near b

\begin{displaymath}U(x)\approx U(b)+\frac{1}{2}kT(x-b)^2\alpha\end{displaymath}

Second integral dominated by contribution from near maximum

\begin{displaymath}U(y)\approx U(b)-\frac{1}{2}kT (y-d)^2\beta\end{displaymath}

Get Kramers' escape rate

\fbox{\parbox{7cm}{
\begin{displaymath}r=\frac{D\sqrt{\alpha\beta}}{2\pi}\exp(-\frac{\Delta U}{kT})\end{displaymath}
}}

This is an asymptotic formula valid in the limit $\Delta U/kT\Rightarrow
\infty$. See the book of Risken for low order correction terms!

SUMMARY

Click here for Return to title page
Click here for 8. Transformations and blow torches

 
next up previous
Next: About this document ...
Birger Bergersen
1998-10-03