Collision time of Brownian particles

Question

Let's assume two spherical particles $p_1$ and $p_2$ of finite radius $r_1$ and $r_2$, which are at locations $(\pm\frac{d}{2},0,0)$ a distance $d$ apart at initial time $t$. These particles diffuse with coefficients $D_1$ and $D_2$, respectively. How can I obtain the probability distribution of the collision time (that is, if the time at which they collide is $t + \Delta t$, I would like to know the probability density function of the random variable $\Delta t$), or at least some of its moments, in these two cases:

1. when the domain is infinite, e.g., $\mathbb{R}^3$, and
2. when it is bounded (inside a sphere with surface area $A$ and the origin $(0,0,0)$ as its center)?

I tried to look up for an answer in Van Kampen's book, but I couldn't really come up with an answer.

Answer 1

The previous answers do not exhaust this question, so I will add a solution to the problem as posted. The question is: what is the probability distribution for the meeting time of two diffusing spheres of radii $r_1$ and $r_2$ starting at initial separation R. By transforming to relative coordinates (as Peter Shor said), you reduce the problem to one diffusing point particle which is absorbed on a sphere of radius $A=r_1+r_2$.

The following is true in 3d

• The particles will not certainly collide. The probability of eventual collision is A/R. The probability of reaching infinity is $1-A/R$.
• If you ask the question--- do the particles collide first or reach a separation of B first, the probability of colliding first is ${ {1\over R} - {1\over B} \over {1\over A} - {1\over B} }$
• There is the usual 3d miracle of Laplacians, that there is a change of variables which normalizes the radial part to a free 1d Laplacian, which allows a closed form expression for the probability density of colliding at any finite time. This gives a complete answer to the question: the probability of not having collided at time t is given by $ P(t) = (1-{A\over R}) + {A\over R} S(\Delta,t)$, where $\Delta= R-A$, and $S(\Delta,t)$ is the probability that a 1 dimensional Brownian motion starting at position $\Delta$ has not hit the origin by time t.

$S(\Delta,t)$ is a simple function which is essentially just the error function:

$$ S(\Delta,t) = ( {2\over \sqrt{\pi}} \mathrm{erf}({\Delta \over \sqrt{2t}})) $$

It's only a function of ${\Delta \over\sqrt{t}}$ by scaling. The probability density of a collision at time t is found by differentiating the above expression, and it only involves elementary functions:

$$\rho(t) = {A\over R} {1\over 4\sqrt{\pi}}{\Delta e^{-{\Delta^2\over 2t}} \over t^{3\over 2}}$$

The answer is strange, because up to the factor of A/R which accounts for the fact that not all walks collide, it is exactly the same as the probability density for a 1d Brownian motion to collide with the origin at time t.

Image charge magnitude is collision probability

Diffusion is by the Laplacian. If you introduce a steady stream of particles at position R, they will come to a steady equilibrium where some are running off to infinity, and some are absorbed by the sphere of radius A. If you introduce these particles at a random angular position on the sphere of radius R, you get a spherically symmetric solution of Laplace's equation at steady state $\phi(r)$, which represents the steady-state density of particles.

The gradient of $\phi$ is the particle steady-state flow, it is the electric field in an electrostatic analog. The source on the sphere at R is a charge Q uniformly distributed on the sphere, which produces an outgoing electric field, an outgoing flow of particles. The particle sink on the sphere of radius A is a potential 0 surface, a metal grounded to infinity.

The solution to this places an image sphere at $A^2/R$ which has charge as required to cancel the potential. The charge is QA/R. This means that the flux into the metal is a fraction A/R of the flux introduced at R, which means that of all the outgoing flux, a fraction A/R is captured by the sphere. This gives the collision probability.

In one and two dimensions, the image charge magnitude in a sphere is always the same as the charge magnitude, and this says that 2d walks are recurrent.

Time-independent Problems are Easy

The solution to the problem of collision with an inner sphere vs. an outer sphere can be solved in the same way. The sphere at A is now is grounded on the outer sphere at a radius B, so that, matching the potentials:

$${Q \over R} - {Q_i \over A} = {Q-Q_i \over B}$$

Where $Q$ is the charge and $Q_i$ is the image charge. The image charge magnitude is

$$ Q_i = Q { {1\over R} - {1\over B} \over {1\over A} - {1\over B} }$$

This is the probability of reaching A before reaching B. It asymptotes to the correct answer for $B=\infty$, of course. This is a time-independent problem, so it has a simple answer using grounded conductors.

This gives you a qualitative indication of what the distribution for first colllsion is, By using B as a surrogate for $\sqrt{t}$. Expanding in B to leading order, and substituting $\sqrt{t}$ for B gives that the probability of surviving is

$$P(t) = {A\over R} (1 - {R-A\over \sqrt{t}})$$

Which is to be compared to the exact answer below.

Time dependent problem

The time dependent problem, the original question, also has a closed form solution. The first thing to note is that the distance between the two particles is undergoing a diffusion too. This distance prefers to go out rather than in, according to the extra phase space volume of being at a larger radius.

The form of the effective one dimensional diffusion can be worked out using a simple trick--- the uniform 3d distribution must be the Boltzmann distribution for diffusion in the effective potential on r. The resulting equation for the radial diffusion is:

$${\partial\rho\over\partial t} = {\partial \over\partial r} ( {\partial\rho\over \partial r} - {(d-1)\over r} \rho )$$

Where $d=3$ is the case of interest, and the boundary conditions you enforce are $\rho(A)=0$ and $\rho(\infty)=0$. This has a simple solution in the usual 3d way: you write $\rho(r,t) = rf(r,t)$ and then

$${\partial f\over \partial t} = {\partial^2 f\over \partial^2 r}$$

Which is the usual miracle of three dimensions, radial Laplacians become free 1d Laplacians once you extract a 1/r factor. In this case, since the distribution $\rho$ is already mutliplied by $r^2$ relative to the radially symmetric solution of Laplace's equation, you extract an r instead.

I will digress to clear up a few confusions:

• How is it possible for a 3d Laplacian, which conserves the 1d integral $r^2\rho(r)$, to be equivalent to a 1d Laplacian, which conserves the integral $f(r)$ without powers of r?

The reason is that diffusion conserves two different moments: the total probability, and the center of mass. The center of mass for f is the total probability for $\rho$. This is important, because f is not to be interpreted as a probability directly, it's first moment is the total probability.

• The location of images in 1d and in 3d are totally different! In 1d, if you have a Dirichlet condition, you always reflect at that point to find the image position. In 3d, you have to reflect in a sphere. How could the two problems be related?

There is an important note regarding this: if you just map a free diffusion in 3d to 1d, the 1d problem has absorbing boundary at the origin. The reason is that the 1d diffusion must be 0 at 0, because it comes from a 3d diffusion which has a finite density at 0 (or from a 1d radial diffusion whose density vanishes as $r^2$ near 0--- same thing). The correct images to use when there is a Dirichlet boundary condition on a sphere are the 1d images. These give the correct 1d solution. Although this is a trivial point, considering the previous section, it is confusing when you take the limits to extract the behavior at long times.

The solution for a delta function at r=R is

$$ f_D(r,t) = {1\over \sqrt{2\pi t}} (e^{-{(r-R)^2\over 2t}} - e^{-{(r-(2A-R))^2\over 2t}}) $$

Which is the image-charge method for producing a solution which is zero at $r=A$. This solution will decay, because of the negative probability coming from the image, or in real life, because of the absorption of the probability into the boundary at r=A, i.e. because of collisions. The leftover probability density is the answer to the question.

So the answer is the first moment of $f_D(r,t)$:

$$ {1\over N} \int_A^\infty r f_D(r,t) = {1\over N} \int_A (r-A) f_D(r,t) + A f_D(r,t) $$

Where the normalization constant N is so that the answer will be 1 at $t=0$. $$ N = \int_A^\infty r f_D(r,0) = \int_A^\infty r \delta(r-R) = R$$ Using the variable $u=r-A$, the first integral is symmetric between u and -u, and becomes half the integral from minus infinity to infinity of something simple to do.

$${1\over 2R} \int_{-\infty}^\infty u f_d(u,t) = {A\over R}$$

The second integral $\int_0^\infty A f_D(r,t)$ is simply the probability of a 1d walk to avoid the origin up to time t, and it is multiplied by ${A\over R}$.

Answer 2

Marek suggested I post my comment (which doesn't completely answer the question) as an answer. Here it is:

Suppose you have two Brownian motions with diffusion coefficients $D_1$ and $D_2$, which start at the same point at $t=0$. Let $x_i$ be the average displacement vector for particle $i$ after time $t$. Then, $\langle x_i^2 \rangle = 2D_it$, where $x_i^2$ means the square of the length of vector $x_i$. Now, consider the relative position of these two random walks after time $t$ (assuming they start out at the same place). The expected value of the relative displacement will be $$\langle (x_1 - x_2)^2 \rangle = \langle x_1^2 \rangle + \langle x_2^2 \rangle = 2 (D_1 + D_2)t,$$ where you can ignore the cross terms since the expected inner product of the two displacement vectors is 0.

Hopefully, this argument convinces you that the relative position of two points undergoing Brownian motion is still Brownian. You can make it completely rigorous if you want.

So your question boils down to: what is the probability distribution for the time it takes Brownian motion starting at distance $d$ from the origin to first reach a disc of radius $r=r_1+r_2$ around the origin. It's straightforward to find the probability that Brownian motion is in that disc after time $t$, but computing the probability distribution for the first time it's reached is harder. I am sure this has been done, but right now I don't have time to do a literature search for this.

Answer 3

Here is an answer in the first case (unbounded domain).

As already noted, the difference between the positions of the two particles performs a Brownian motion starting from point $(d,0,0)$ with diffusion coefficient $D=D_1+D_2$ and one is interested in the time $\tau$ of the first hitting of the ball $B(r)$ around zero with radius $r=r_1+r_2$ when $d>r$.

It happens that in three or more dimensions, the Brownian motion is transient. This means that the particle goes to infinity almost surely, in the sense that for any radius $R$, there will exist an almost surely finite time $L_R$ such that after time $L_R$, the particle is not in $B(R)$. The same holds if one replaces $B(R)$ by any bounded domain: after a random time, finite with full probability, the particle wanders out of the domain.

Another consequence is that, starting from a point not in $B(r)$, for example from the point $(d,0,0)$, the particle has a nonzero probability to never hit $B(r)$. Thus, it happens with positive probability that the two particles the OP is considering never collide in $\mathbb{R}^3$ (and of course this probability depends on $d$, $r_1$ and $r_2$, but for any $d>r_1+r_2$ it is positive). In particular:

The collision time is infinite with nonzero probability. All its (positive) moments are infinite.

Two more remarks.

First, in one dimension the collision time $\tau$ would be almost surely finite, with the well known density $$ P(\tau\in\mathrm{d}t)=\frac{a}{\sqrt{2\pi}}\mathrm{e}^{-a^2/2t}\frac{\mathrm{d}t}{t^{3/2}},\qquad a=d-r_1-r_2. $$ Note that, since the density above is equivalent to a multiple of $t^{-3/2}$ when $t\to+\infty$, the mean of $\tau$ is infinite. More precisely, for positive $\gamma$, $\tau^\gamma$ is integrable if $\gamma<1/2$ but not if $\gamma\ge1/2$. Beware that we followed the mathematicians' convention for Brownian motion, that is, the formula above applies to a one dimensional Brownian motion $(B_t)$ with transition probabilities $$ P(B_t\in\mathrm{d}x\vert B_0=y)=\frac1{\sqrt{2\pi t}}\mathrm{e}^{-(x-y)^2/2t}\mathrm{d}x. $$ By scaling, it should be straightforward to deduce the distribution of $\tau$ for a one dimensional Brownian motion with diffusion constant $D$ since $\tau$ for a diffusion constant $D$ is distributed like $\tau^{(1)}/D$ where $\tau^{(1)}$ corresponds to the first hitting time for a diffusion constant $1$.

Second, in case 2 (bounded domain), one can expect the collision time to be almost surely finite and integrable (and even exponentially integrable) but to make it possible to answer the question, one should explain what happens to the particles when one of them hits the boundary of the domain before their collision time.