Pollen grains on the stamen of a flower.

Modeling Diffusion and Capture.

Collaborators: 

Andrew Bernoff (Harvey Mudd), 

Theo Kolokolnikov (Dalhousie University), 

Daniel Scmidt (Harvey Mudd) 

Ryan Spoonmore (Notre Dame), 

Justin Tzou, (U. British Columbia), 

Michael Ward (U. British Columbia) 

A.E. Lindsay

Applied and Computational Math & Stats

University of Notre Dame.

a.lindsay@nd.edu www.nd.edu/~alindsa1

NIST, 

June 13th, 2017.

Support:

NSF DMS - 1516753
Diffusion in heterogeneous and dynamic environments.

dXt = DdWt;

Main Questions: 

• How long till a diffusing particle finds a target?

• How does the distribution and target site mobility 

control possibility for, and time of capture?

square mean displacement

=

D x time

2D diffusion in bounded region. 3D diffusion in unbounded region.
• The cell nucleus - genetic material passes

 from the interior/exterior through small pores

celled Nuclear Pore Complexes (NPCs).

• Nucleus ~10% of cell volume.

• Roughly N = 2000 surface pores which 

occupy 2% of the surface area.

• Ref: Eilenberg et.al. Science 341(6146), 2013.

Maeshima et. al. Nature Struct. & Mol. Bio. 17 (2010)

• Pore density dynamic.

• Nuclear volume dynamic.

• Experimentalists measure pore

densities, not individual positions. 

Application I: Intracellular Transport.
• What is the probability of this binding occurring? 

• On average, how long does this take? 

• How does the distribution of the receptors affect this?

Antigens

T Cell

When an antigen (a toxin or a protein that promotes an immune response) 

binds to a receptor on a T-cell it can trigger the creation of antibodies.

Application II: Molecular Signaling
• General Aspects of the mathematical theory.

• Boundary Homogenization for spheres,

 disks and planes.

• A spectral boundary element method for 

capture problems.

• The roles of dynamics and cooperation

 factors in capture problems.

• Conclusions and Future Work.

TalkOutline
y = yn

y = yn+1

y = yn−1

xn−1 xn xn+1

p(xn, yn,t + ∆t) =

1

4 p(xn−1, yn,t) + p(xn+1, yn,t) + p(xn, yn+1,t) + p(xn, yn−1,t)

xn = nx

yn = ny

p(x, t)  probability that particle starting at x = (x, y) is free at time t

x = y ! 0

t ! 0

(x)2

4t

= D

The probability of absorption I 
The probability of absorption II 

In N dimensional free space:

Diffusion Kernel https://en.wikipedia.org/wiki/Brownian_motion
Absorbing (or Dirichlet).

Boundary Conditions.

⌦

Reflecting (or Neumann).

 ⌧ 1 Mostly Reflecting

  1 Mostly Absorbing

Partially Absorbing (or Robin)
The probability of Escape 
x 0 R

r=R

r=a

r=R

r=a

1D

2D

3D

The importance of dimension in 

diffusion and capture problems
x 0 R

Consider a delta-function release at x=R.

ut = Duxx x > 0, t> 0

u(x, 0) = (x  R) x > 0

u(0, t)=0 t > 0

u(x, t)

Capture in 1D

By method of images: u(x, t) = 1

2

p

⇡Dt[e(XR)2/4Dt  e(X+R)2/4Dt]

FPT: p(t) = Dux(0, t) = Probability of capture: ReR2/4Dt

2

p

⇡Dt3/2

(R)=1

t

p(t) Capture Time
Probability of capture: 

t

Capture Time

(R)=1 q(⌧, 3)

⇠

2

⌧ [ln ⌧ ]

2

r=R

Consider a d-function release at r=R and a circular trap of radius a.

r=a

ut = D



urr +

1

r

ur



r > 0, t> 0

u(r, 0) = 1

2⇡R(r  R) r>a

u(a, t)=0 t > 0

Capture in 2D

FPT: Obtained with a 

 Laplace Transform 

q(⌧, ) = 2

⇡

Z 1

s=0

ses2

⌧ J0(s)Y0(s)  J0(s)Y0(s)

[J0(s)]2 + [Y0(s)]2 ds
r=R

Consider a d-function release at r=R and a spherical trap of radius a.

r=a ut = D



urr +

2

r

ur



r > 0, t> 0

u(r, 0) = 1

4⇡R2 (r  R) r>a

u(a, t)=0 t > 0

Capture in 3D

q(⌧, 3)

FPT: 

Capture Time

q(⌧, ) = ✓1  1



◆ e

(1)2

4⌧

2

p⇡⌧ 3/2

r

(r)

(r) = a

r

p(t) = a2

d q(⌧, )  =

R

a

⌧ =

a2t

d

Probability of capture: 

τ
Problem: 3D Diffusion to a surface with small 

absorbing pores.

Problem 1: Pores on sphere. Problem 2: Pores on infinite plane.

Describe the rate of capture at the pores.
Diffusion to a structured 3D target with small 

absorbing pores.

Target

Sites
Berg & Purcell, Physics of Chemoreception (1977). 1500 Citations. Berg & Purcell - Physics of Chemoreception 1977 

(~1500 citations) 

Target

Sites

• Derived by interpolating between all absorbing and N independent pores.

 

• No information on influence of receptor spatial arrangements.

Limitations
Important Ingredient - Single Pore on a flat plane.

Ref: Sneddon, Mixed Boundary 

Value Problems, 1966.
Rain Drain Analogy - perimeter that matters in the 

biologically realistic limit!

Source: shutterstock.com
Literature since Berg & Purcell - How does 

clustering influence the capture rate?

• Exact analysis very difficult - classical potential theory utilizes

 separation of variables which is tricky for mixed boundary values problems

Boundary Homogenization or ‘effective medium theory’

Ref: Muratov, 

Shvartsmam, 

Berezhkovskii, 

SIAM MMS 2006.

Exact Planar capture time distribution
Literature since Berg & Purcell - How does 

clustering influence the capture rate?

• Exact analysis very difficult - classical potential theory utilizes

 separation of variables which is tricky for mixed boundary values problems

Boundary Homogenization or ‘effective medium theory’

Ref: Muratov, 

Shvartsmam, 

Berezhkovskii, 

SIAM MMS 2006.

• Depends on the absorbing area fraction.

• Depends on the particular receptor clustering.
Literature since Berg & Purcell - How does 

clustering influence the capture rate?

• Exact analysis very difficult - classical potential theory utilizes

 separation of variables which is tricky for mixed boundary values problems

Boundary Homogenization or ‘effective medium theory’

Ref: Muratov, 

Shvartsmam, 

Berezhkovskii, 

SIAM MMS 2006.

Using the Berg-Purcell 

Flux:

Ref: Shoup, 

Szabo 

BioPhys J. 

1982.
Literature since Berg & Purcell - How does 

clustering influence the capture rate?

• Exact analysis very difficult - classical potential theory utilizes

 separation of variables which is tricky for mixed boundary values problems

Boundary Homogenization or ‘effective medium theory’

Ref: Muratov, 

Shvartsmam, 

Berezhkovskii, 

SIAM MMS 2006.

Ref: 

Berezhkovskii 

2008, 2013, 

2016.

Particle simulations fit parameters. A = 1.62, 1.75, 1.37 and B = 1.36, 2.02, 2.59, 

for triangular, square and hexagonal lattices.
μ

0 0.5 1 1.5 2

g(μ)

0

2

4

6

8

10

Result requires an explicit surface Green’s Function. 

Lindsay, Ward, Bernoff - SIAM MMS, Vol. 15, No. 1 (2017) 

• From a detailed singular perturbation analysis with N pores of common radius.
Surface Green’s Function of Sphere.

Ref: Nemenman and Silbergleit, J. Appl. Phys., 86 (1999), pp. 614–624. 

Gravity Probe B 

• Launched in 2004 to verify 

predictions of General Relativity.

• Four Gyros measured the 

precessions over the period of a 

year.

• Gyroscopes the most perfect 

spheres ever manufactured at 

time.

A gyro for GP-B
Problem: Find global minimizers of

1

y

0

-1 x 1 0 -1

-1

1

0.5

0

-0.5

z

Smale’s 7th problem
Boundary Homogenization Result

0 500 1000 1500 2000

N

0

2

4

6

8

10

H

×105

Assuming a uniform distribution on points on sphere. 

Later: Numerical verification for 

realistic N (thousands).
The Equivalent 2D problem.

(`) = ⇡`

log ⇡`

2 Window/Pore size

Fixed ` = N, and N ! 1, gives e↵ective Robin problem

x

-1.5 -1 -0.5 0 0.5 1 1.5 2

y

-1.5

-1

-0.5

0

0.5

1

1.5

x

-1.5 -1 -0.5 0 0.5 1 1.5 2

y

-1.5

-1

-0.5

0

0.5

1

1.5

N =5& ` =

1

6 N = 10 & ` =

1

2

Lindsay, Tzou, Kolokolnikov

PRE, 91(3), 2015.
For N pores of radius a:

When the distance between pores, . dij  a

This tells us how to replace many pores with one.

Berg & Purcell 1977

Bernoff & Lindsay 2017

Jpore = 4aN  1

⇡

X

i6=j

8a2

dij

+

1

⇡2

X

j6=k

X

i6=j

16a3

dijdjk

Planar Homogenization Problem
Numerical Solution of the Capacitance Problem

Target

Sites

Popular Particle Based Solutions.
Step	1:	Project	Onto	Surface	

Step	2:	Check	if	par7cle	Has	

Been	Absorbed	

Step	3:	Diffuse	Par7cle	Off	

Surface	

Kinetic Monte Carlo for periodic domain.
Integral equation over 

the support of the pores.

Boundary Integral solution of the Capacitance 

Problem

p(y) =

x∈∂Ωa

G(x; y)q(x)dS. Neumann to

Dirichlet Map

Pores on a plane Pores on a sphere

Sphere Case

Planar Case
Integral equation over 

the support of the pores.

Boundary Integral solution of the Capacitance 

Problem

p(y) =

x∈∂Ωa

G(x; y)q(x)dS. Neumann to

Dirichlet Map

Pores on a plane Pores on a sphere

Main Difficulty: 

Flux has an (integrable) singularity 

on the edge of each each pore.

Solution: 

Expand surface potential and flux in 

basis which mimics the known 

singularity on boundary of pores.
The Zernike Polynomials.

Surface 

Potential

Surface

Flux

Z0 0 =

1

p2⇡

,

Z1 1 =

2

p⇡

r cos(✓), Z1 1 =

2

p⇡

r sin(✓),

Z2 2 =

r 6

⇡

r2 cos(2✓), Z2 0 =

r 3

⇡



2r2  1



, Z2 2 =

r 6

⇡

r2 sin(2✓).

At each pore, expand flux and potential:
• Dutch Physicist.

• 1953 Physics Nobel Prize for 

phase-contrast microscope. Imaging 

translucent samples by change in light 

phase, not intensity.

• Zernike Polynomial derived in 1934 and 

form a complete basis for square 

integrable functions with circular support. 

• Used extensively in beam optics for 

analyzing waveforms entering through 

circular apertures.

• Quantifying and correcting ocular 

aberrations (e.g. astigmatism) in 

optometry. 

Frits Zernike (1888-1966)
p(y) =

x∈∂Ωa

G(x; y)q(x)dS. Neumann to

Dirichlet Map

Boundary Integral solution Procedure

1. Expand flux and potential with M Zernike modes at the N pores. 

2. Gives N*(M+1)*(M+2)/2 unknowns for the system.

3. Form dense system for unknowns by projecting the flux onto Zernike modes.

4. Solve linear system and calculate the flux

Advantages

• Spectral accuracy. 

• Runs very quickly for low modes.

Disadvantages

• Relies on explicit Green’s Function. 

• Circular pore geometry assumed.
Numerical Results Sphere:

N=1:

N>1:

Relative Error - One Pore Agreement - Platonic Solids

-2 -1.5 -1 -0.5 0

-8

-6

-4

-2

0

log10 Erel

log10 σ

Relative error resolved to 1 part in 10^8!
d 1 1

Two Planar Pores

Ref: Strieder, 

J. Chem, Phys, 2012.

Accuracy as function of pore separation.

0 5 10 15

M

-10

-8

-6

-4

-2

log10

Erel

d = 2.1

d = 3

d = 4

d = 5

0 0.5 1 1.5 2

log10 d

-8

-7

-6

-5

-4

-3

-2

-1

log10

Erel

Slope -6

Convergence of numerics with exact solution Zernike Modes Taken
Validation of the homogenized sphere condition.

Fibonacci Spiral Points and 2% absorbing fraction.

Homogenized Formula 

 accurately predicts the 

 flux to the target in biologically 

realistic examples!!
Ongoing Work: The full First Passage Time (FPT) 

distribution.

∂Ω Ωε ε

Ω

x

The MFPT overestimates typical capture times.

Important to get the full capture time distribution

-6 -5 -4 -3 -2 -1 0 1

log t

0

0.5

1

1.5

2

2.5

3

F(t)

Capture Time Density

The Mean.
Ongoing Work II: Formulation of Capture 

problems in a dynamical setting.
Rotating Trap

Tzou & Kolokolnikov 

SIAM MMS 2015

Optimal radius for 

a given speed of 

rotation.
Ongoing Work II: Trap Cooperation Strategies I

Absorbing Boundaries

Random Motion

Random Motion

Absorbing Boundaries

x = x1(t)

x = x1(t)

x = x2(t)

x = x2(t)

Advection

x

t

Out-of Phase In Phase
z

0 50 100 150 200

χ(z)

-0.05

-0.04

-0.03

-0.02

-0.01

0

0.01

0.02

0.03

0.04

0.05

(z) =

cosh pz

2 sin pz2 + sinh pz2 cos pz2

cos p2z  cosh p2z . Trap interaction 

strength.

(!) > 0

In phase minimizes MFPT when

Out-of-phase minimizes MFPT when

(!) < 0

! - common trap frequency

Lindsay, Tzou, 

Kolokolnikov 2017

Ongoing Work: Trap Cooperation Strategies II
1. Showed explicitly how pore/receptor clustering influences diffusive sensing 

rates, resolving a long standing problem postulated by Berg & Purcell.

2. Derivation of macro scale capture laws from microscale clustering pattern. 

3. These averaged features can be used to give insight and reduce 

challenging computational aspects of these multi scale problems.

4. Developed a precision numerical tools for studying receptor clustering. 

Conclusions.

1. Formulate and solve problems for capture in more general dynamic 

environments. (e.g. moving domains, transient and growing pores).

2. On top of this, get the full distribution of capture times! Is the homogenized 

boundary condition verified here valid for the time dependent problem?

Future Work.
Thank you for your attention!!

References:

• Lindsay & Bernoff, Numerical approximation of diffusive capture rates by planar and spherical 

surfaces with absorbing pores, Submitted SIAM J. Applied Math. 2017.

• Lindsay, Tzou, Kolokolnikov, Optimization of first passage times by multiple cooperating

 mobile traps, To Appear, SIAM MMS (2017).

• Lindsay, Bernoff, Ward, First passage statistics for the capture of a Brownian particle by a 

structured spherical target with multiple surface traps, SIAM MMS, 15(1), pp. 74-109 (2017).

• Lindsay, Spoonmore, Tzou, Hybrid Asymptotic-Numerical Approach for estimating first passage 

densities of the two-dimensional narrow capture problem, PRE 94 042418, (2016)

• Lindsay, Tzou, Kolokolnikov, Narrow escape problem with a mixed trap and the effect of 

orientation. PRE, Vol. 91, No. 3, (2015).

Collaborators: 

Andrew Bernoff (Harvey Mudd) 

Theo Kolokolnikov (Dalhousie) 

Daniel Scmidt (Harvey Mudd) 

Ryan Spoonmore (Notre Dame) 

Justin Tzou, (UBC) 

Michael Ward (UBC) 

Support:

NSF DMS - 1516753 Preprints at: www.nd.edu/~alindsa1