1. Introduction
We consider the problem of describing the mean first passage time (MFPT) of two-dimensional Brownian motion in a bounded region to a small elliptical absorbing trap. The diffusive transport of molecules and individual agents from a source to a mobile or fixed target is a problem occurring in a variety of physical, biological and social systems [Reference Holcman and Schuss17, Reference Metzler, Oshanin and Redner36, Reference Redner42]. Ecological examples include the time required for an animal to find a mate or shelter [Reference Fauchald and Torkild12, Reference Kurella, Tzou, Coombs and Ward27, Reference Stepien, Zmurchok, Hengenius, Rivera, D’Orsogna and Lindsay43]. At the cellular scale, diffusion transports key cargoes within the cell [Reference Bressloff9, Reference Ghusinga, Dennehy and Singh13, Reference Gibbons, Chou and D’Orsogna14, Reference Grebenkov, Metzler and Oshanin16, Reference Lawley, Lindsay and Miles28], including fibroblasts to initiate wound healing [Reference Ambrosi, Gamba and Serini1], antigens for detection by T-cell receptors [Reference Morgan, Lindsay and Lavrik37, Reference Morgan, Pettmann, Dushek and Lindsay38], and material to and from the nucleus [Reference Leech, Hazel, Gatlin, Lindsay and Manhart30, Reference Windner, Manhart, Brown, Mogilner and Baylies45]. The scope of the target search problem has also been expanded to include features such as stochastic switching of target states, resetting [Reference Bressloff and Schumm7], extreme statistics, and homogenization [Reference Plunkett and Lawley41]. For an extensive review, we point the reader to the recent survey [Reference Grebenkov, Metzler and Oshanin15].
The MFPT
$u(\mathbf{x})$
describing the expected time to capture a diffusing particle initially at
$\mathbf{x}\in \Omega \setminus \Omega _{{\displaystyle \varepsilon }}$
solves the Poisson equation [Reference Bressloff and Newby6, Reference Holcman and Schuss17]
\begin{align} D \Delta u+1 &= 0, \qquad \mathbf{x}\in \Omega \setminus \Omega _{{\displaystyle \varepsilon }}; \\[-8pt] \nonumber \end{align}
\begin{align} D\nabla u \cdot {\hat {\textbf {n}}} &= 0, \qquad \mathbf{x} \in \partial \Omega ; \\[-8pt] \nonumber \end{align}
\begin{align} u &= 0, \qquad \mathbf{x} \in \partial \Omega _{{\displaystyle \varepsilon }}. \\[6pt] \nonumber \end{align}
The boundary conditions (1.1b) prescribe that the outer boundary
$\partial \Omega$
is reflecting and equation (1.1c) specifies that the trap
$\Omega _{{\displaystyle \varepsilon }}$
is absorbing. The aim of this paper is to construct a solution to (1.1) in the limit as
${\displaystyle \varepsilon }\to 0$
in the presence of an elliptical trap defined as
\begin{equation} \Omega _{{\displaystyle \varepsilon }} = \boldsymbol{\xi } + {\displaystyle \varepsilon } e^{i\phi } \mathcal{A}, \qquad \mathcal{A} = \Big \{ (y_1,y_2) \in \mathbb{R}^2 \ | \ \frac {y_1^2}{a^2} + \frac {y_2^2}{b^2} \lt 1\Big \}. \end{equation}
Schematic of the configuration of the domain
$\Omega$
with a single trap
$\Omega _{{\displaystyle \varepsilon }}$
as defined in (1.3). The trap is centred at a point
$\boldsymbol{\xi }\in \Omega$
located
$\mathcal{O}(1)$
from
$\partial \Omega$
and has semi-major and semi-minor axes
${\displaystyle \varepsilon } a$
and
${\displaystyle \varepsilon } b$
respectively. The semi-major axis of the trap is orientated at angle
$\phi$
with respect to the horizontal axis.
Here
$\boldsymbol{\xi }\in \Omega$
is the centring point of the trap such that
$\text{dist}(\boldsymbol{\xi },\partial \Omega ) = \mathcal{O}(1)$
as
${\displaystyle \varepsilon }\to 0$
,
${\displaystyle \varepsilon } a$
and
${\displaystyle \varepsilon } b$
are the semi-major and semi-minor axes, respectively, and
$\phi$
is the angle of orientation with respect to the horizontal axis (see Figure 1). The term
$e^{i\phi }$
corresponds to rotation by angle
$\phi$
in the counter clockwise direction. An important quantity, called the global MFPT (GMFPT), describes the overall capture rate based on a uniform distribution of start locations and is defined as
\begin{equation} \tau =\frac {1}{|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u(\mathbf{x})\, d\mathbf{x}. \end{equation}
Before outlining the rationale for this work and relationship to previous studies, we state our main result:
Principal Result: Consider equation (1.1) with a single elliptical trap centred at
$\boldsymbol{\xi }\in \Omega$
with semi-major and semi-minor axes
${\displaystyle \varepsilon } a$
and
${\displaystyle \varepsilon } b$
respectively (
$a\gt b$
) and the semi-major axis having elevation
$\phi$
from the horizontal. In the limit as
${\displaystyle \varepsilon }\to 0^{+}$
, a two term expansion of the solution to (1.1) and the GMFPT (1.3) is of form
\begin{equation} u(\mathbf{x}) = \frac {1}{D}\left [ u_0(\mathbf{x}) + {\displaystyle \varepsilon }^2 u_2(\mathbf{x}) + \mathcal{O}({\displaystyle \varepsilon }^4)\right ]; \qquad \tau = \frac {1}{D}\left [\tau _0 + {\displaystyle \varepsilon }^2 \tau _2 + \mathcal{O}({\displaystyle \varepsilon }^4)\right ]. \end{equation}
The terms in the above expansions are given explicitly as
\begin{align} u_0(\mathbf{x}) &= -|\Omega |\Big [ G(\mathbf{x};\boldsymbol{\xi }) - R(\boldsymbol{\xi };\boldsymbol{\xi }) \Big ] + \frac {|\Omega |}{2\pi \nu }; \\[-12pt] \nonumber \end{align}
\begin{align} u_2(\mathbf{x}) &= |\Omega | \Big [\frac 12 \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{\boldsymbol{\xi }} G(\mathbf{x};\boldsymbol{\xi })\big ) -2\pi \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi };\boldsymbol{\xi }) \cdot \mathcal{M} \nabla _{{\boldsymbol{\xi }}} G(\mathbf{x};\boldsymbol{\xi }) \Big ] + \chi _2. \\[8pt] \nonumber \end{align}
The logarithmic gauge function is
$\nu ({\displaystyle \varepsilon })= -1/\log ({\displaystyle \varepsilon } d_c)$
where
$d_c$
is the logarithmic capacitance which reflects the shape of the trap and is determined by (B.81). For an elliptical trap,
$d_c = \frac {1}{2}(a+b)$
. The constant
$\chi _2$
is given by
\begin{equation} \chi _2 = -|\Omega |\Big ( \mbox{Trace}\big (\mathcal{Q} \nabla _{{\boldsymbol{\xi }}}^2 R(\boldsymbol{\xi };\boldsymbol{\xi }) \big ) - 2\pi \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi };\boldsymbol{\xi })\cdot \mathcal{M} \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi };\boldsymbol{\xi }) \Big ) + \frac {a^2+b^2}{8}. \end{equation}
Here,
$\mathcal{Q}$
and
$\mathcal{M}$
are the quadrupole and moment polarization matrices and solve associated electrified disk problems (B.81) and (B.84), respectively. For the case of an elliptical trap, they are given explicitly by
\begin{equation} \mathcal{Q} = -\frac {a^2-b^2}{4} \begin{bmatrix} \cos 2\phi &\quad \phantom {-}\sin 2\phi \\[3pt] \sin 2\phi &\quad -\cos 2\phi \end{bmatrix}, \qquad \mathcal{M} = -\frac {(a+b)^2}{4}\mathcal{I} + \mathcal{Q}. \end{equation}
The terms of the GMFPT (1.4a) are given by
\begin{align} \tau _0 & = \frac {|\Omega |}{2\pi }\Big [ \frac {1}{\nu } + 2\pi R(\boldsymbol{\xi };\boldsymbol{\xi }) \Big ], \qquad \nu ({\displaystyle \varepsilon }) = \frac {-1}{\log ( {\displaystyle \varepsilon } d_c)};\\[-8pt]\nonumber \end{align}
\begin{align} \tau _2 & = \,\Big [ \frac {\pi ab }{|\Omega |}\tau _0 + \frac {a^2+b^2}{8} + \chi _2 \Big ]. \\[8pt] \nonumber \end{align}
We remark that the term
$\frac {\pi ab }{|\Omega |}\tau _0$
in (1.4g) arises due to fact that
$|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|= |\Omega | - \pi ab{\displaystyle \varepsilon }^2$
. In the above result,
$G(\mathbf{x};\boldsymbol{\xi })$
and
$R(\mathbf{x};\boldsymbol{\xi })$
are the Neumann Green’s function and its regular part, respectively, defined as the unique solution of
\begin{gather} \Delta G = \frac {1}{|\Omega |} - \delta (\mathbf{x}-\boldsymbol{\xi }), \quad \mathbf{x}\in \Omega ; \qquad \nabla G \cdot {\hat {\textbf {n}}} = 0,\quad \mathbf{x}\in \partial \Omega ; \\[-30pt] \nonumber \end{gather}
\begin{gather} \int _{\Omega } G(\mathbf{x};\boldsymbol{\xi }) d\mathbf{x} = 0; \qquad G(\mathbf{x};\boldsymbol{\xi }) = -\frac {1}{2\pi }\log | \mathbf{x}- \boldsymbol{\xi }| + R(\mathbf{x};\boldsymbol{\xi }). \\[-10pt] \nonumber \end{gather}
Before giving a detailed outline of the steps leading to the principal result, we review some recent works on related problems and motivations for this study. Over the past several decades, there has been extensive research in the asymptotic analysis of the two-dimensional Poisson problem (1.2) in the presence of small inhomogeneities [Reference Bressloff8, Reference Isaacson and Newby18, Reference Iyaniwura, Wong, Ward and Macdonald20, Reference Iyaniwura, Wong, Macdonald and Ward21, Reference Kolokolnikov, Titcombe and Ward24, Reference Kurella, Tzou, Coombs and Ward27, Reference Ward and Keller44], which serves as a canonical problem in the trafficking and delivery of small signalling molecules and cargoes [Reference Bressloff and Newby6]. More generally, there has been significant recent interest in the study of elliptic problems in punctured domains [Reference Campbell and Nazarov10, Reference Kolokolnikov and Lindsay23, Reference Kropinski, Lindsay and Ward26, Reference Lindsay, Quaife and Wendelberger33–Reference Lindsay, Ward and Kolokolnikov35], including higher-order corrections for eigenvalue problems in periodic domains [Reference Bernoff, Lindsay and Schmidt5, Reference Iyaniwura and Ward19, Reference Paquin-Lefebvre, Iyaniwura and Ward40].
Elliptical traps are of particular interest in cellular signalling problems due to the frequent observation of non-circularity in the cell itself [Reference Kaiyrbekov and Camley22] or other key organelles, such as the nucleus [Reference Leech, Dalwadi and Manhart29, Reference Windner, Manhart, Brown, Mogilner and Baylies45]. An oval or elliptical geometry accurately captures the aberrations to radial symmetry observed in these domains and hence it is natural to investigate how the capture rate of Brownian particles is modulated by non-circularity. This is one element in a broader mathematical effort to understand the contribution of geometry in cellular signalling [Reference Bernoff, Jilkine, Navarro Hernández and Lindsay3, Reference Bernoff and Lindsay4, Reference Kaiyrbekov and Camley22, Reference Lindsay, Bernoff and Hernández31, Reference Nakamura and Kobayashi39]
The leading order behaviour of the GMFPT (1.4f) as
${\displaystyle \varepsilon }\to 0$
(see [Reference Kurella, Tzou, Coombs and Ward27]) captures the effects of trap size and position. However, there is no information on how the orientation of the trap influences the solution. A correction to the leading order behaviour (1.4b) was derived in [Reference Lindsay, Kolokolnikov and Tzou32] that captures the effect of orientation,
\begin{equation} \tau = \frac {1}{D}\Big [\tau _0 - {\displaystyle \varepsilon }|\Omega |\,\big ( \, \textbf {d} \cdot \nabla _{\boldsymbol{\xi }} R(\boldsymbol{\xi };\boldsymbol{\xi })\big ) + \mathcal{O}({\displaystyle \varepsilon }^2) \Big ]. \end{equation}
In the result (1.6), the vector
$\textbf {d}$
is related to the dipole moment of the trap and is defined by associated problems that incorporate the shape and orientation of the trap. The contribution to the MFPT from the location of the trap in
$\Omega$
is captured by the quantity
$\nabla _{\boldsymbol{\xi }} R(\boldsymbol{\xi };\boldsymbol{\xi })$
, where the subscript reflects differentiation with respect to the source location.
For application of (1.6) to the case of an elliptical trap, which has two lines of symmetry, we find (see Appendix B.1) that the dipole term vanishes (d=0), thus (1.6) no longer describes the effect of orientation on the MFPT. Our refined result (1.4) describes the higher-order contribution to the MFPT due to the trap orientation. In particular, we can identify optimizing configurations by writing the GMFPT as
\begin{equation} \tau = \frac {\tau _0}{D} + \frac {{\displaystyle \varepsilon }^2}{D} \left [\frac {\pi ab}{|\Omega |}\tau _0 + \frac {a^2+b^2}{4} - \pi |\Omega | \frac {(a+b)^2}{2} (R_{\xi _1}^2 + R_{\xi _2}^2) + |\Omega | \frac {a^2-b^2}{4} \mathbf{p} \cdot \begin{bmatrix}\cos 2\phi \\ \sin 2\phi \end{bmatrix} \right ]. \end{equation}
The vector
$\mathbf{p}$
is found to be
\begin{equation} \mathbf{p} = \begin{bmatrix} R_{\xi _1\xi _1}-R_{\xi _2\xi _2} - 2\pi (R_{\xi _1}^2 -R_{\xi _2}^2) \\[4pt] 2R_{\xi _1\xi _2} - 4\pi R_{\xi _1} R_{\xi _2} \end{bmatrix}\, , \end{equation}
which gives the direction along which the trap should be orientated to optimize the correction term
$\tau _2$
of the GMFPT. We note from (1.4f), that when the trap centre
$\boldsymbol{\xi }$
is placed at a critical point of
$R(\boldsymbol{\xi };\boldsymbol{\xi })$
, so that
$\nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi };\boldsymbol{\xi }) = [R_{\xi _1},R_{\xi _2}]^T=[0,0]^T$
, the optimal orientation vector reduces to
$\mathbf{p} =[R_{\xi _1\xi _1}-R_{\xi _2\xi _2}, 2R_{\xi _1\xi _2}]^{T}$
.
The outline of the paper is as follows. In Section 2, we present a hierarchy of results, beginning with the solution of (1.1) in the reduced case of a circular trap located at the centre of a disk (Section 2.1). Following this, we derive the solution in the presence of an elliptical trap placed at the centre of a disk (Section 2.2). Finally, we present the corresponding result for the case of a general domain with an ellipse of arbitrary orientation (Section 2.3), which yields the principal result (1.4).
In Section 3, we first validate the asymptotic result on the unit disk domain where the regular part
$R(\mathbf{x};\boldsymbol{\xi })$
is known in closed form. In this unit disk case, we identify a bifurcation where for
$|\boldsymbol{\xi }|\gt r_c:=\sqrt {2-\sqrt {2}}$
, the GMFPT correction term
$\tau _2$
is minimized when the semi-major axis is orientated in the angular direction. Conversely, for
$|\boldsymbol{\xi }|\lt r_c$
, the GMFPT correction
$\tau _2$
is minimized when the semi-major axis is orientated in the radial direction. For certain regular domains, such as rectangles and ellipses, highly accurate series solutions for (1.5) are available, which allows us to determine
$\mathbf{p}$
. For these cases, we reveal similar bifurcations of the optimizing orientation depending on the centring point of the ellipse
$\boldsymbol{\xi }$
and proximity to
$\partial \Omega$
. Finally in Section 4, we discuss avenues for future research arising from this study.
2. Asymptotic analysis of the mean first passage time to a single elliptical trap
In this section, we perform the main asymptotic analysis on the MFPT problem (1.1). To guide the rationale for the higher order expansions, it is useful to first analyse two exactly solvable cases for a circular trap located at the centre of a disk and an elliptical trap located at the centre of a disk.
2.1. Unit disk with a circular trap at the origin
For a single circular trap at the origin, the MFPT (1.1) reduces to the ODE
\begin{equation} u_{rr} + \frac {1}{r} u_r = -\frac {1}{D}, \qquad {\displaystyle \varepsilon }\lt r\lt 1; \qquad u({\displaystyle \varepsilon }) = u'(1) = 0, \end{equation}
where
$r = |\mathbf{x}|$
. The exact solution of (2.9a) is
\begin{equation} u(r) = \frac {1}{2D} \left [ - \frac {r^2}{2} + \log \frac {r}{{\displaystyle \varepsilon }} + \frac {{\displaystyle \varepsilon }^2}{2} \right ]. \end{equation}
The corresponding GMFPT
\begin{align} \tau =& \frac {1}{|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u \, d\mathbf{x} = \frac {2\pi }{\pi (1-{\displaystyle \varepsilon }^2)} \int _{r={\displaystyle \varepsilon }}^1 u(r) r dr = \frac {1}{8D(1-{\displaystyle \varepsilon }^2)} \Big [ -3 + 4{\displaystyle \varepsilon }^2 - {\displaystyle \varepsilon }^4 - 4\log {\displaystyle \varepsilon } \Big ]\nonumber\\ & = \frac {1}{8D} \Big [ -3 - 4 \log {\displaystyle \varepsilon } + {\displaystyle \varepsilon }^2 (1 - 4 \log {\displaystyle \varepsilon }) + \mathcal{O}({\displaystyle \varepsilon }^4) \Big ]. \end{align}
The equation (2.9c) will be a useful case to validate solutions of (1.1) for more general configurations.
2.2. Unit disk with an elliptical trap at the origin
We now solve for the MFPT in the scenario of an elliptical trap at the origin, orientated along the horizontal axis (
$\phi =0$
) with semi-major and semi-minor axes
${\displaystyle \varepsilon } a$
and
${\displaystyle \varepsilon } b$
, respectively.
Outer Expansion We expand as
\begin{equation} u(r) = \frac {1}{D} \left [ u_0(r) + {\displaystyle \varepsilon }^2 u_2(r,\theta ) + \cdots \right ], \end{equation}
where
$\Delta u_0 +1 = 0$
with
$u_0'(1) = 0$
and
$\Delta u_2 = 0$
with
$u_2'(1) = 0$
. The general solutions for
$u_0$
and
$u_2$
in polar coordinates
$\mathbf{x} = re^{i\theta }$
are
\begin{align} u_0 &= \frac {1}{2}\log r - \frac {r^2}{4} + A_1, \\[-12pt] \nonumber \end{align}
\begin{align} u_2 &= B_2 \cos 2\theta ( r^2 + r^{-2} ) + B_1, \\[8pt] \nonumber \end{align}
for constants
$A_1,B_1, B_2$
to be determined. We remark that the choice of the
$\cos (2\theta )$
solution is even in
$\theta$
and generates the lowest order singularity. The general solution (2.11b) can accommodate a term of form
$\sin 2\theta \times \{r^2,r^{-2}\}$
in the case that the trap orientation moves off the perpendicular axis. Moving to a coordinate
$\mathbf{x} = {\displaystyle \varepsilon } \mathbf{y}$
, we find that
\begin{equation*} u \sim \frac {1}{2} \log |\mathbf{y}| + \frac {1}{2} \log {\displaystyle \varepsilon } + A_1 + B_2 \frac {\cos 2\theta }{|\mathbf{y}|^2} + {\displaystyle \varepsilon }^2 \Big [ B_1 -\frac {y^2}{4} \Big ] + \mathcal{O}({\displaystyle \varepsilon }^4). \end{equation*}
Inner Expansion: In the region
$\mathbf{x} = {\displaystyle \varepsilon }\mathbf{y}$
, the solution
$u(\mathbf{x}) = U(\mathbf{y})$
is expanded as
\begin{equation*} U(\mathbf{y}) = \frac {1}{D} \left [U_0(\mathbf{y}) + {\displaystyle \varepsilon }^2 U_2(\mathbf{y}) + \cdots \right ]. \end{equation*}
$\mathcal{O}({\displaystyle \varepsilon }^0)$
: The leading order problem satisfies
\begin{gather} \Delta _{\mathbf{y}} U_0 = 0, \quad \mbox{in} \quad \mathbb{R}^2\setminus \mathcal{A}; \qquad U_0 = 0 \quad \mbox{on} \quad \partial \mathcal{A}; \\[-30pt] \nonumber \end{gather}
\begin{gather} U_0 = \frac {1}{2} \log |\mathbf{y}| + \frac {1}{2} \log {\displaystyle \varepsilon } + A_1 + B_2 \frac {\cos 2\theta }{|\mathbf{y}|^2} + \cdots , \quad |\mathbf{y}|\to \infty , \\[4pt] \nonumber \end{gather}
where
$\mathcal{A}$
is the rescaled ellipse of semi-major and semi-minor axes
$a$
and
$b$
, respectively, with orientation
$\phi =0$
with respect to the origin. Using the solution
$v_{0c}(\mathbf{y})$
of the electrified disk problem derived in appendix B.1, we calculate that
\begin{equation} U_0 = \frac {1}{2}v_{0c}(\mathbf{y}) \sim \frac {1}{2} \Big [\log |\mathbf{y}| - \log \alpha - \frac {a^2-b^2}{4} \frac {\cos 2\theta }{|\mathbf{y}|^2} + \cdots \Big ], \qquad |\mathbf{y}| \to \infty , \end{equation}
where
$\alpha = (a+b)/2$
. Matching (2.12b) with (2.12c) yields that
\begin{equation} A_1 = \frac {1}{2\nu }, \qquad B_2 = -\frac {a^2-b^2}{8}, \qquad \nu = \frac {-1}{\log {\displaystyle \varepsilon }\alpha }. \end{equation}
$\mathcal{O}({\displaystyle \varepsilon }^2)$
: Proceeding to the next order, we have that
\begin{gather} \Delta _{\mathbf{y}} U_2 = -1, \quad \mbox{in} \quad \mathbb{R}^2\setminus \mathcal{A}; \qquad U_2 = 0 \quad \mbox{on} \quad \partial \mathcal{A}; \\[-25pt] \nonumber \end{gather}
\begin{gather} U_2 = B_1 - \frac {|\mathbf{y}|^2}{4} + \mathrm{\mathcal{O}}(1), \quad |\mathbf{y}|\to \infty . \\[8pt] \nonumber \end{gather}
We solve this problem in appendix B.3 by decomposing the solution as
$U_2 = -\frac {1}{4}|\mathbf{y}|^2 + U_{2h}$
where
$U_{2h}(\mathbf{y})$
solves a homogeneous problem. We find (see (B.91)) that the large argument behaviour of (2.14) is
\begin{equation} U_2 \sim - \frac {|\mathbf{y}|^2}{4} + \frac {a^2+b^2}{8} + \mathcal{O}(|\mathbf{y}|^{-2}), \qquad |\mathbf{y}|\to \infty . \end{equation}
Hence from comparison to (2.14b), we determine that
\begin{equation*} B_1 = \frac {a^2+b^2}{8}. \end{equation*}
We remark that the term
$\mathcal{B}_{11}$
in (B.91) is determined by the Hessian of the leading outer solution
$u_0$
, which vanishes in this case due to the trap being centred at the origin. This completes the derivation of the expansion (2.10) and yields the following expression for the MFPT
\begin{equation} u(r,\theta ) = \frac {1}{2D}\left [\log r - \frac {r^2}{2} + \frac {1}{\nu } + \frac {{\displaystyle \varepsilon }^2}{4} \Big ( (a^2+b^2) -(a^2-b^2)\big (r^2+r^{-2} \big )\cos 2\theta \Big ) \right ] + \cdots \end{equation}
as
${\displaystyle \varepsilon }\to 0$
. We remark that this expression reduces to the exact solution (2.9b) in the scenario
$a=b=1$
. Calculating the GMFPT, we must identify local and global contributions by introducing an intermediate scale
${\displaystyle \varepsilon } \ll \delta \ll 1$
.
\begin{equation*} \tau = \frac {1}{|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u \, d\mathbf{x} = \frac {1}{\pi (1- {\displaystyle \varepsilon }^2ab)} \Big [ \underbrace {\int _{|\mathbf{x}| \lt \delta } u d\mathbf{x}}_{I_1} + \underbrace {\int _{\delta \lt |\mathbf{x}|\lt 1} u d\mathbf{x}}_{I_2}\Big ]. \end{equation*}
We proceed to calculate both terms
$I_1$
and
$I_2$
to arrive at a final expression independent of
$\delta$
.
Calculation of I
1: Rescaling with
$\mathbf{x} = {\displaystyle \varepsilon }\mathbf{y}$
, we have that
\begin{equation} I_1 = {\displaystyle \varepsilon }^2 \int _{\substack {\mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}\\ |\mathbf{y}|\lt \delta /{\displaystyle \varepsilon }}} U\, d\mathbf{y} = \frac {{\displaystyle \varepsilon }^2}{D} \int _{|\mathbf{z}|=1}^{|\mathbf{z}| = {\frac {\delta }{{\displaystyle \varepsilon }\alpha }}} (\underbrace {U_0}_{I_{10}} + {\displaystyle \varepsilon }^2 \underbrace {U_1}_{I_{12}})|J|d\mathbf{z}. \end{equation}
In the above calculation of the integral
$I_1$
, the region exterior to the ellipse has been mapped to the exterior of the disk through the transformation
$\mathbf{y}\to \alpha \mathbf{z} + \beta /\mathbf{z}$
for
$\mathbf{z} = re^{i\theta }$
and Jacobian
$|J|$
given by
\begin{equation*} |J| = \begin{vmatrix} \alpha - \frac {\beta }{r^2}\cos 2\theta & - \frac {\beta }{r^2} \sin 2\theta \\[4pt] \frac {\beta }{r^2} \sin 2\theta & \alpha - \frac {\beta }{r^2}\cos 2\theta \end{vmatrix} = \alpha ^2- \frac {2\alpha \beta }{r^2}\cos 2\theta + \frac {\beta ^2}{r^4}. \end{equation*}
Using the fact that
$U_0(\mathbf{z}) = \frac {1}{2} \log |\mathbf{z}|$
, we have that
\begin{align} \nonumber I_{10} &\approx \frac {{\displaystyle \varepsilon }^2}{2} \int _{\theta =0}^{2\pi }\int _{r=1}^{r= \frac {\delta }{{\displaystyle \varepsilon }\alpha }} \log r \Big (\alpha ^2 - \frac {2\alpha \beta }{r^2} \cos 2\theta + \frac {\beta ^2}{r^4} \Big ) r dr d\theta = {\displaystyle \varepsilon }^2 \pi \int _{r=1}^{r= \frac {\delta }{{\displaystyle \varepsilon }\alpha }} \log r \Big (\alpha ^2 r + \frac {\beta ^2}{r^3} \Big ) dr \\[5pt] & \approx \frac {\pi }{4} \left [ {\displaystyle \varepsilon }^2(\alpha ^2+\beta ^2) + 2 \delta ^2 \log \frac {\delta }{{\displaystyle \varepsilon } \alpha } - \delta ^2 \right ]. \end{align}
The contribution from
$I_{12}$
is
$\mathcal{O}({\displaystyle \varepsilon }^4)$
. We now calculate the contribution from the outer solution.
Calculation of I 2:
\begin{align} \nonumber I_2 &= \frac {1}{D}\int _{|\mathbf{x}|\gt \delta } (u_0 + {\displaystyle \varepsilon }^2 u_2) d\mathbf{x} \\[3pt] \nonumber &=\frac {1}{2D} \int _{\theta =0}^{2\pi }\int _{r=\delta }^{1}\left [\log r - \frac {r^2}{2} + \frac {1}{\nu } + \frac {{\displaystyle \varepsilon }^2}{4} \Big ( (a^2+b^2) -(a^2-b^2)(r^2 +r^{-2})\cos 2\theta \Big ) \right ] r dr d\theta \\[3pt] \nonumber {} & = \frac {\pi }{D}\int _{r=\delta }^{1}\left [ \log r - \frac {r^2}{2} + \frac {1}{\nu } + \frac {{\displaystyle \varepsilon }^2}{4} \Big (a^2+b^2\Big ) \right ] rdr \\[3pt] {}& \approx \frac {\pi }{4D}\left [ -\frac {3}{2} + \frac {2}{\nu } + \frac {{\displaystyle \varepsilon }^2}{2}\Big (a^2+b^2\Big ) -2 \delta ^2 \log \delta + \delta ^2 - \frac {2\delta ^2}{\nu } \right ]. \end{align}
Combining the two terms (2.18) and (2.19), we have that
\begin{align} \nonumber \tau & = \frac {1}{\pi D(1-{\displaystyle \varepsilon }^2 ab)} \frac {\pi }{8} \left [ -3 + \frac {4}{\nu } + {\displaystyle \varepsilon }^2(a^2+b^2) + 2{\displaystyle \varepsilon }^2 (\alpha ^2+\beta ^2)\right ]\\[4pt] & = \frac {1}{8D} \left [ -3 + \frac {4}{\nu } + {\displaystyle \varepsilon }^2 \Big ( 2(a^2+ b^2) - 3 ab + \frac {4ab}{\nu } \Big ) \right ] + \mathcal{O}({\displaystyle \varepsilon }^4). \end{align}
As required, this expression is independent of
$\delta$
and substitution of
$a=b=1$
reduces (2.20) to (2.9c).
2.3. General case for a single trap
We will now determine the solution to the MFPT problem for a single elliptical trap, centred at
$\mathbf{x}=\boldsymbol{\xi }$
with semi-major and semi-minor axes
${\displaystyle \varepsilon } a,{\displaystyle \varepsilon } b$
, respectively, and orientation
$\phi$
with respect to the horizontal axis. The explicit form of the trap is given in (1.2).
In an outer region away from the elliptical trap, we expand the solution as
\begin{equation*} u(\mathbf{x}) = \frac {1}{D} \left [u_0(\mathbf{x}) + {\displaystyle \varepsilon } u_1(\mathbf{x}) + {\displaystyle \varepsilon }^2 u_2(\mathbf{x}) + \mathrm{\mathcal{O}}({\displaystyle \varepsilon }^2)\right ]. \end{equation*}
The outer problems
$u_j$
for
$j= 0,1,2,\ldots$
satisfy
\begin{align} \Delta u_j + \delta _{1j} &= 0, \quad \mathbf{x}\in \Omega \setminus \{\boldsymbol{\xi }\}; \\[-10pt] \nonumber \end{align}
\begin{align} \nabla u_j \cdot {\hat {\textbf {n}}} &= 0, \quad \mathbf{x}\in \partial \Omega . \\[8pt] \nonumber \end{align}
The local behaviour as
$\mathbf{x}\to \boldsymbol{\xi }$
is now established for each problem (2.21) through boundary layer analysis. In the vicinity of the trap, the solution is expanded in variables
\begin{equation} U = \frac {1}{D} \left [ U_0(\mathbf{y}) + {\displaystyle \varepsilon } U_1 (\mathbf{y}) + {\displaystyle \varepsilon }^2 U_2(\mathbf{y}) + \mathrm{\mathcal{O}}({\displaystyle \varepsilon }^2)\right ], \qquad \mathbf{y} = e^{-i\phi } \, \frac {\mathbf{x}-\boldsymbol{\xi }}{{\displaystyle \varepsilon }}. \end{equation}
Collecting terms at relevant orders gives a sequence of problems to be solved.
Inner Region
$\mathcal{O}({\displaystyle \varepsilon }^0)$
: The leading order problem for
$U_0$
satisfies
\begin{align} \Delta U_0 &= 0, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \\[-10pt] \nonumber \end{align}
\begin{align} U_0 &= 0, \quad \mathbf{y} \in \partial \mathcal{A}. \\[8pt] \nonumber \end{align}
In terms of the solution
$v_{0c}$
of the electrified disk problem obtained in Appendix B.1, we have that
\begin{equation} U_0(\mathbf{y}) = S\nu v_{0c}(\mathbf{y}), \qquad \nu = \frac {-1}{\log {\displaystyle \varepsilon }\alpha }. \end{equation}
Here,
$S$
is a constant to be determined in the matching process. The far field of equation (2.24) supplies the appropriate local behaviour for the outer solution. In (B.81), we establish that for an elliptical trap aligned in the horizontal direction (
$\phi =0$
), the far-field behaviour is
\begin{equation} v_{0c}(\mathbf{y}) = \log |\mathbf{y}| - \log \alpha + \frac {\mathbf{y}^T \tilde {\mathcal{Q}} \mathbf{y} }{|\mathbf{y}|^4}, \qquad \tilde {\mathcal{Q}} = - \alpha \beta \begin{bmatrix} 1 &\phantom {-}0 \\0& -1 \end{bmatrix} = - \frac {a^2-b^2}{4}\begin{bmatrix} 1 &\phantom {-}0 \\0& -1 \end{bmatrix}. \end{equation}
Hence, in terms of the outer coordinate
$\mathbf{y} = {\displaystyle \varepsilon }^{-1} e^{-i\phi } (\mathbf{x}-\boldsymbol{\xi })$
incorporating rotation by
$\phi$
with respect to the horizontal, and incorporating the far-field behaviour for
$|\mathbf{y}|\to \infty$
, we generate the local behaviour as
$\mathbf{x}\to \boldsymbol{\xi }$
\begin{equation} u\sim S\nu \left [ \log |\mathbf{x}-\boldsymbol{\xi }| + \frac {1}{\nu } + {\displaystyle \varepsilon }^2\frac {(\mathbf{x}-\boldsymbol{\xi })^T \mathcal{Q}(\mathbf{x}-\boldsymbol{\xi }) }{|\mathbf{x}-\boldsymbol{\xi }|^4}\right ] + \cdots . \end{equation}
The matrix
$\mathcal{Q}=e^{i\phi } \tilde {\mathcal{Q}}e^{-i\phi }$
is calculated as
\begin{align} \mathcal{Q} = -\frac {a^2-b^2}{4}\begin{bmatrix} \cos \phi &\quad -\sin \phi \\ \sin \phi &\quad \phantom {-}\cos \phi \end{bmatrix}\begin{bmatrix} 1&\quad \phantom {-}0 \\ 0 &\quad -1 \end{bmatrix}\begin{bmatrix} \phantom {-}\cos \phi &\quad \sin \phi \\ -\sin \phi &\quad \cos \phi \end{bmatrix} = -\frac {a^2-b^2}{4} \begin{bmatrix} \cos 2\phi &\quad \phantom {-}\sin 2\phi \\ \sin 2\phi &\quad -\cos 2\phi \end{bmatrix}.\nonumber\\[3pt] \end{align}
This behaviour is used to furnish terms in the outer expansion.
Outer Region
$\mathcal{O}({\displaystyle \varepsilon }^0)$
: The problem at leading order is
\begin{align} \Delta u_0 + 1 &= 0, \qquad \mathbf{x} \in \Omega \setminus \{\boldsymbol{\xi }\}; \\[-10pt] \nonumber \end{align}
\begin{align} \nabla u_0 \cdot {\hat {\textbf {n}}} &= 0, \qquad \mathbf{x}\in \partial \Omega ; \\[-10pt] \nonumber \end{align}
\begin{align} u_0 &\sim S\nu \log |\mathbf{x}-\boldsymbol{\xi }| + S + \cdots , \qquad \mathbf{x}\to \boldsymbol{\xi }. \\[8pt] \nonumber \end{align}
In terms of the Neumann’s Green’s function (1.5), we have that
\begin{equation} u_0 = -2\pi S \nu G(\mathbf{x};\boldsymbol{\xi }) + \tau _0, \end{equation}
where
$\tau _0$
is a constant. The expansion of
$u_0$
as
$\mathbf{x}\to \boldsymbol{\xi }$
gives the local behaviour
\begin{align} u_0 &= S\nu [\!\log |\mathbf{x}-\boldsymbol{\xi }| - 2\pi R(\mathbf{x}\,;\,\boldsymbol{\xi })] + \tau _0\\ \nonumber & \sim S\nu \log |\mathbf{x}-\boldsymbol{\xi }| + \tau _0 - 2\pi S \nu \Big [ R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) +\mathbf{a} \cdot (\mathbf{x} -\boldsymbol{\xi } ) + (\mathbf{x}-\boldsymbol{\xi })^T \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} (\mathbf{x}-\boldsymbol{\xi }) + \mathcal{O}(|\mathbf{x}-\boldsymbol{\xi }|^3) \Big ], \end{align}
where the coefficient terms for
$\mathbf{x}=(x_1,x_2)$
are given by
\begin{equation} \mathbf{a} : = \begin{bmatrix} \partial _{x_1} R(\mathbf{x}\,;\,\boldsymbol{\xi })\\[5pt] \partial _{x_2} R(\mathbf{x}\,;\,\boldsymbol{\xi }) \end{bmatrix}_{\mathbf{x} = \boldsymbol{\xi }}, \quad \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} = \frac {1}{2}\begin{bmatrix} \partial _{x_1x_1}R(\mathbf{x}\,;\,\boldsymbol{\xi }) &\quad \partial _{x_1x_2}R(\mathbf{x}\,;\,\boldsymbol{\xi }) \\[5pt] \partial _{x_1x_2}R(\mathbf{x}\,;\,\boldsymbol{\xi }) &\quad \partial _{x_2x_2}R(\mathbf{x}\,;\,\boldsymbol{\xi }) \end{bmatrix}_{\mathbf{x}=\boldsymbol{\xi }} = \frac 12 \begin{bmatrix} R_{11} &\quad R_{12}\\[5pt] R_{12} &\quad R_{22} \end{bmatrix}. \end{equation}
A system of two equations for unknowns
$(S,\tau _0)$
is found by both matching (2.29) with the local behaviour and integrating (2.28a). This yields that
\begin{equation} S = \frac {|\Omega | }{2\pi \nu }, \qquad \tau _0 = \frac {|\Omega |}{2\pi \nu } \Big [1 + 2\pi \nu R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \Big ]. \end{equation}
This result was established in [Reference Kolokolnikov, Titcombe and Ward24]. We now determine the correction to the inner expansion. In the local variable
$\mathbf{y} = {\displaystyle \varepsilon }^{-1}e^{-i\phi }(\mathbf{x}-\boldsymbol{\xi })$
, equation (2.30) yields the far-field behaviour
\begin{equation} U(\mathbf{y}) \sim S\nu \left [ \log |\mathbf{y}| - \log \alpha - 2\pi {\displaystyle \varepsilon }\, \mathbf{a} \cdot e^{i\phi } \mathbf{y} -2\pi {\displaystyle \varepsilon }^2 \, \mathbf{y}^{T} e^{-i\phi } \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} e^{i\phi } \mathbf{y} + \cdots \right ], \qquad |\mathbf{y}| \to \infty . \end{equation}
This reveals the leading order matching behaviour for the higher-order inner corrections. Specifically, we have that
\begin{align} U_1(\mathbf{y}) &\sim -2\pi S \nu \mathbf{a} \cdot e^{i\phi } \mathbf{y} + \cdots \\[-10pt] \nonumber \end{align}
\begin{align} U_2(\mathbf{y}) &\sim -2\pi S \nu \mathbf{y}^{T} e^{-i\phi } \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} e^{i\phi } \mathbf{y} + \cdots \\[8pt] \nonumber \end{align}
as
$|\mathbf{y}|\to \infty$
. This behaviour is now matched to corresponding inner problems.
Inner Region
$\mathcal{O}({\displaystyle \varepsilon }^1)$
: At this order, we must solve the exterior problem
\begin{gather} \Delta _{\mathbf{y}} {U}_{1} = 0\,, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad U_{1} = 0, \quad \mathbf{y} \in \partial {\mathcal{A}}; \\[-30pt] \nonumber \end{gather}
\begin{gather} U_{1} = -2\pi S \nu \mathbf{a} \cdot e^{i\phi } \mathbf{y} + \cdots , \qquad |\mathbf{y}|\to \infty ; \\[8pt] \nonumber \end{gather}
In appendix B.2, we introduce and solve the vector valued electrified disk problem
$\textbf {v}_{1c}$
\begin{gather} \Delta _{\mathbf{y}} \textbf {v}_{1c} = 0\,, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad \textbf {v}_{1c} = 0, \quad \mathbf{y} \in \partial \mathcal{A}; \\[-30pt] \nonumber \end{gather}
\begin{gather} \textbf {v}_{1c} = \mathbf{y} + \frac { \tilde {\mathcal{M}} \mathbf{y}}{|\mathbf{y}|^2} + \cdots \, \quad |\mathbf{y}|\to \infty ; \qquad \tilde {\mathcal{M}} = - \alpha \begin{pmatrix} a &\quad 0\\ 0&\quad b\end{pmatrix}. \\[8pt] \nonumber \end{gather}
In terms of this solution, we write that
\begin{equation} U_1 = -2\pi S \nu \, \mathbf{a} \cdot e^{i\phi } \textbf {v}_{1c}(\mathbf{y}). \end{equation}
Applying the far-field behaviour (2.36b) as
$|\mathbf{y}|\to \infty$
, while returning to outer coordinates with variable
$\mathbf{y} = {\displaystyle \varepsilon }^{-1} e^{-i\phi } (\mathbf{x}-\boldsymbol{\xi })$
, we generate the local behaviour
\begin{equation} U_1 \sim -\frac {2\pi S \nu }{{\displaystyle \varepsilon }} \, \mathbf{a} \cdot \left [ (\mathbf{x}-\boldsymbol{\xi }) + {\displaystyle \varepsilon }^2 \frac {\mathcal{M} (\mathbf{x}-\boldsymbol{\xi }) }{|\mathbf{x}-\boldsymbol{\xi }|^2} \right ] + \cdots , \qquad \mbox{as} \quad \mathbf{x}\to \boldsymbol{\xi }, \end{equation}
In the case of an elliptical trap, we calculate explicitly that
$\mathcal{M} = e^{i\phi } \tilde {\mathcal{M}} e^{-i\phi }$
has the form
\begin{align} \nonumber \mathcal{M} &= -\alpha \begin{bmatrix} \cos \phi &\quad -\sin \phi \\[5pt] \sin \phi &\quad \phantom {-}\cos \phi \end{bmatrix}\begin{bmatrix} a&\quad 0 \\[5pt] 0 &\quad b \end{bmatrix}\begin{bmatrix} \phantom {-}\cos \phi &\quad \sin \phi \\[5pt] -\sin \phi &\quad \cos \phi \end{bmatrix} \\[5pt] \nonumber &=-\alpha \begin{bmatrix} a\cos ^2\phi + b\sin ^2\phi &\quad (a-b)\sin \phi \cos \phi \\[5pt] (a-b)\sin \phi \cos \phi &\quad b\cos ^2\phi + a\sin ^2\phi \end{bmatrix}\\[5pt] \nonumber & = -\alpha ^2\begin{bmatrix} 1&\quad 0 \\[5pt] 0 &\quad 1 \end{bmatrix} -\alpha \beta \begin{bmatrix} \cos 2\phi &\quad \phantom {-}\sin 2\phi \\[5pt] \sin 2\phi &\quad -\cos 2\phi \end{bmatrix} \\[5pt] & = -\alpha ^2 \mathcal{I} + \mathcal{Q}, \end{align}
where
$\mathcal{Q}$
is the quadrupole matrix derived in (2.27).
Inner Region
$\mathcal{O}({\displaystyle \varepsilon }^2)$
: At
$\mathcal{O}({\displaystyle \varepsilon }^2)$
, we introduce the scalar valued problem
$U_2(\mathbf{y})$
where
\begin{gather} \Delta _{\mathbf{y}} U_2 = -1, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad U_2 = 0, \quad \mathbf{y} \in \partial \mathcal{A}; \\[-30pt] \nonumber \end{gather}
\begin{gather} U_2 = -2\pi S\nu \, \mathbf{y}^T e^{-i\phi } \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} e^{i\phi } \mathbf{y} + \cdots \, \quad |\mathbf{y}|\to \infty ; \\[8pt] \nonumber \end{gather}
To further decompose the far-field behaviour (2.39b), we write
\begin{align*} \mathcal{H}&=e^{-i\phi } \nabla _{\mathbf{x}}^2 R\mid _{\mathbf{x}=\boldsymbol{\xi }} e^{i\phi }\\[5pt] &= \frac {1}{2}\begin{bmatrix} \phantom {-}\cos \phi &\quad \sin \phi \\[5pt] -\sin \phi &\quad \cos \phi \end{bmatrix}\begin{bmatrix} R_{11}&\quad R_{12} \\[5pt] R_{12} &\quad R_{22} \end{bmatrix}\begin{bmatrix} \cos \phi &\quad -\sin \phi \\[5pt] \sin \phi &\quad \phantom {-}\cos \phi \end{bmatrix}\\[5pt] & = \frac {R_{11} + R_{22}}{4} \begin{bmatrix} 1 &\quad 0\\[5pt] 0 &\quad 1 \end{bmatrix} + \frac {1}{4}\begin{bmatrix} (R_{11}-R_{22})\cos 2\phi + 2 R_{12}\sin 2\phi &\quad \phantom {-}2R_{12}\cos 2\phi - (R_{11}-R_{22}) \sin 2\phi \\[5pt] 2R_{12}\cos 2\phi - (R_{11}-R_{22})\sin 2\phi &\quad -(R_{11}-R_{22})\cos 2\phi - 2 R_{12}\sin 2\phi \end{bmatrix}\\[5pt] &=\frac {R_{11} + R_{22}}{4} \begin{bmatrix} 1 &\quad 0\\[5pt] 0 &\quad 1 \end{bmatrix} + \begin{bmatrix} \mathcal{B}_{11} &\quad \phantom {-}\mathcal{B}_{12}\\[5pt] \mathcal{B}_{12} &\quad - \mathcal{B}_{11} \end{bmatrix}. \end{align*}
Here,
$\mathcal{B}$
is the matrix satisfying
$\mbox{Trace}(\mathcal{B}) = 0$
with components
\begin{align} \mathcal{B}_{11} &= \frac {1}{4} \Big [ (R_{11}-R_{22})\cos 2\phi + 2 R_{12}\sin 2\phi \Big ], \\[-10pt] \nonumber \end{align}
\begin{align} \mathcal{B}_{12} &= \frac {1}{4}\Big [ 2R_{12}\cos 2\phi - (R_{11}-R_{22})\sin 2\phi \Big ]. \\[8pt] \nonumber \end{align}
The far field behaviour is now in the form
\begin{equation} \mathbf{y}^T \mathcal{H} \mathbf{y} = \left [ (R_{11} + R_{22})\frac {y_1^2 + y_2^2}{4} + \mathcal{B}_{11} (y_1^2 - y_2^2) + 2 \mathcal{B}_{12} \, y_1 y_2 \right ]. \end{equation}
We remark that
$\Delta u_0 = \Delta ( S\nu \log |\mathbf{x}-\boldsymbol{\xi }| - 2\pi S\nu R(\mathbf{x}\,;\,\boldsymbol{\xi }) +\tau _0) = -1$
, which implies that
\begin{equation} \Delta R = R_{11} + R_{22} = \frac {1}{|\Omega |}. \end{equation}
After applying these reductions, together with
$2\pi S\nu = |\Omega |$
from (2.32), we can restate equation (2.39) as
\begin{gather} \Delta _{\mathbf{y}} U_2 = -1, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad U_2 = 0, \quad \mathbf{y} \in \partial \mathcal{A}; \\[-30pt] \nonumber \end{gather}
\begin{gather} U_2 \sim -\frac {|\mathbf{y}|^2}{4}- 2\pi S\nu \, \mathbf{y}^{T} \mathcal{B}\mathbf{y} + \cdots , \qquad |\mathbf{y}|\to \infty . \\[6pt] \nonumber \end{gather}
In Appendix B.3, we state and solve the canonical problem (2.43) and obtain the refined behaviour
\begin{equation} U_2 \sim -\frac {|\mathbf{y}|^2}{4} -2\pi S\nu \mathbf{y}^{T} \mathcal{B} \mathbf{y} + d_{2c} + \mathcal{O}(|\mathbf{y}|^{-2}), \quad \mbox{as} \quad |\mathbf{y}|\to \infty , \end{equation}
where the constant term is
\begin{equation} d_{2c} = \frac {\alpha ^2 + \beta ^2}{4} + 4\pi S\nu \alpha \beta \, \mathcal{B}_{11}. \end{equation}
If we incorporate the value of
$\mathcal{B}_{11}$
shown in (2.40a), this term reduces to
\begin{equation} d_{2c} = \frac {\alpha ^2 + \beta ^2}{4} - \pi S\nu \mbox{Trace}\big (\mathcal{Q} \nabla _{\mathbf{x}}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ), \end{equation}
where we have used the identity
\begin{equation*} \mbox{Trace}\big (\mathcal{Q} \nabla _{\mathbf{x}}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) = \mathcal{Q}_{11} (R_{11} - R_{22}) + 2\mathcal{Q}_{12} R_{12}; \quad \mathcal{Q} = \begin{bmatrix} \mathcal{Q}_{11} &\quad \phantom {-}\mathcal{Q}_{12}\\[5pt] \mathcal{Q}_{12} &\quad -\mathcal{Q}_{11} \end{bmatrix} = -\alpha \beta \begin{bmatrix} \cos 2\phi &\quad \phantom {-}\sin 2\phi \\[5pt] \sin 2\phi &\quad -\cos 2\phi \end{bmatrix}. \end{equation*}
This completes the solution of the inner expansion (2.22) the inner problem up to
$\mathcal{O}({\displaystyle \varepsilon }^2)$
. A combination of equations (2.26), (2.38) and (2.44) yields the local behaviour
\begin{align} u &\sim S\nu \log |\mathbf{x}-\boldsymbol{\xi }| + S \\[5pt] \nonumber &+ {\displaystyle \varepsilon }^2 \left (S\nu \left [ \frac {(\mathbf{x}-\boldsymbol{\xi })^T \mathcal{Q}(\mathbf{x}-\boldsymbol{\xi }) }{|\mathbf{x}-\boldsymbol{\xi }|^4} - 2\pi \,\mathbf{a} \cdot \frac {\mathcal{M} (\mathbf{x}-\boldsymbol{\xi }) }{|\mathbf{x}-\boldsymbol{\xi }|^2}\right ] + d_{2c} \right ) \quad \mbox{as} \quad \mathbf{x}\to \boldsymbol{\xi }. \end{align}
We now return to the outer expansion.
Outer region
$\mathcal{O}({\displaystyle \varepsilon }^1)$
: At this order, the problem is given by
\begin{align} \Delta u_1 &= 0, \qquad \mathbf{x} \in \Omega \setminus \{\boldsymbol{\xi }\}; \\[-10pt] \nonumber \end{align}
\begin{align} \nabla u_1 \cdot {\hat {\textbf {n}}} &= 0, \qquad \mathbf{x}\in \partial \Omega ; \\[-10pt] \nonumber \end{align}
\begin{align} u_1 &\to 0, \qquad \mathbf{x}\to \boldsymbol{\xi }. \\[8pt] \nonumber \end{align}
The unique solution of (2.47) is
$u_1 \equiv 0$
.
Outer region
$\mathcal{O}({\displaystyle \varepsilon }^2)$
: At this order, we have that
\begin{align} \Delta u_2 &= 0, \qquad \mathbf{x} \in \Omega \setminus \{\boldsymbol{\xi }\}; \\[-10pt] \nonumber \end{align}
\begin{align} \nabla u_2 \cdot {\hat {\textbf {n}}} &= 0, \qquad \mathbf{x}\in \partial \Omega ; \\[-10pt] \nonumber \end{align}
\begin{align} u_2 & \sim S \nu \left [ \frac {(\mathbf{x}-\boldsymbol{\xi })^T \mathcal{Q}(\mathbf{x}-\boldsymbol{\xi })}{|\mathbf{x}-\boldsymbol{\xi }|^4} -2\pi \mathbf{a} \cdot \frac {\mathcal{M} (\mathbf{x}-\boldsymbol{\xi })}{|\mathbf{x}-\boldsymbol{\xi }|^2}\right ] + d_{2c} + \cdots , \qquad \mathbf{x}\to \boldsymbol{\xi }. \\[8pt] \nonumber\end{align}
To express the solution of (2.48) in terms of the Green’s function, we first notice by direct calculation that
\begin{align} \nabla _{\boldsymbol{\xi }} \log {|\mathbf{x}-\boldsymbol{\xi }|} &= -\frac {\mathbf{x}-\boldsymbol{\xi }}{|\mathbf{x}-\boldsymbol{\xi }|^2}, \\[-10pt] \nonumber \end{align}
\begin{align} -\frac {1}{2}\mbox{Trace} \big (\mathcal{Q}\nabla _{{\boldsymbol{\xi }}}^2 \log |\mathbf{x}-\boldsymbol{\xi }| \big ) &= \frac {(\mathbf{x}-\boldsymbol{\xi })^T \mathcal{Q}(\mathbf{x}-\boldsymbol{\xi }) }{|\mathbf{x}-\boldsymbol{\xi }|^4}, \qquad \mathcal{Q} = \begin{bmatrix} Q_{11} & \phantom {-}Q_{12} \\ Q_{12} & -Q_{11} \end{bmatrix}. \\[8pt] \nonumber \end{align}
The solution of (2.48) can then be written as
\begin{equation} u_2(\mathbf{x}) = S\nu \Big [\pi \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi })\big ) -4\pi ^2\mathbf{a} \cdot \mathcal{M} \nabla _{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi }) \Big ] + \chi _2, \end{equation}
where
$\chi _2$
is a constant to be determined. In the formulation (2.50), the derivatives with respect to the source location
$\mathbf{y}=(y_1,y_2)$
are
\begin{equation*} \nabla _{\boldsymbol{\xi }} = \begin{bmatrix} \partial _{\xi _1}\\[5pt] \partial _{\xi _2}\end{bmatrix}, \qquad \nabla ^2_{\boldsymbol{\xi }} = \begin{bmatrix} \partial ^2_{\xi _1 \xi _1} &\quad \partial ^2_{\xi _1\xi _2}\\[5pt] \partial ^2_{\xi _2 \xi _1} &\quad \partial ^2_{\xi _2 \xi _2} \end{bmatrix}. \end{equation*}
This in particular leads to the identity
\begin{equation} \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{\boldsymbol{\xi }} G(\mathbf{x}\,;\,\mathbf{y})\big ) = Q_{11} (G_{y_1y_1} -G_{y_2y_2}) + 2Q_{12} G_{y_1y_2}. \end{equation}
To complete the matching, we evaluate (2.50) as
$\mathbf{x}\to \boldsymbol{\xi }$
. Since
$\mathbf{a} = \nabla _{\mathbf{x}} R (\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) = \nabla _{\boldsymbol{\xi }}R (\boldsymbol{\xi }\,;\,\boldsymbol{\xi })$
, we calculate
\begin{equation*} u_2 \sim \chi _2 + S\nu \left [ \frac {(\mathbf{x}-\boldsymbol{\xi })^T \mathcal{Q} (\mathbf{x}-\boldsymbol{\xi })}{|\mathbf{x}-\boldsymbol{\xi }|^4} -2\pi \mathbf{a} \cdot \frac {\mathcal{M} (\mathbf{x}-\boldsymbol{\xi })}{|\mathbf{x}-\boldsymbol{\xi }|^2} + \pi \mbox{Trace}\big (\mathcal{Q} \nabla _{{\mathbf{y}}}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) - 4\pi ^2 \mathbf{a} \cdot \mathcal{M} \mathbf{a} \right ], \quad \mathbf{x}\to \boldsymbol{\xi }. \end{equation*}
Matching with (2.48c), we have that
\begin{equation*} \chi _2 = -S\nu \Big ( \pi \mbox{Trace}\big (\mathcal{Q} \nabla _{{\mathbf{y}}}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) - 4\pi ^2 \mathbf{a} \cdot \mathcal{M} \mathbf{a} \Big ) + d_{2c}. \end{equation*}
In the case particular to the elliptical trap, we apply (2.38b) so that
$\mathcal{M} = - \alpha ^2 \mathcal{I} + \mathcal{Q}$
. In addition, we apply the value of
$d_{2c}$
given in (2.45), together with the symmetry relationship
$\nabla _{\mathbf{x}}^2R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) = \nabla _{\boldsymbol{\xi }}^2R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi })$
, and
$S\nu = |\Omega |/2\pi$
from (2.32) to further reduce
$\chi _2$
to
\begin{align} \nonumber \chi _2 &= -S\nu \Big ( 2\pi \mbox{Trace}\big (\mathcal{Q} \nabla _{\boldsymbol{\xi }}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) + 4\pi ^2 \alpha ^2 |\mathbf{a}|^2 - 4\pi ^2 \mathbf{a} \cdot \mathcal{Q} \mathbf{a} \Big ) + \frac {a^2+b^2}{8}\\[5pt] &= -|\Omega |\Big ( \mbox{Trace}\big (\mathcal{Q} \nabla _{\boldsymbol{\xi }}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) + 2\pi \alpha ^2 |\mathbf{a}|^2 - 2\pi \mathbf{a} \cdot \mathcal{Q} \mathbf{a} \Big ) + \frac {a^2+b^2}{8}. \end{align}
In summary, the solution of (1.1) admits the expansion
$u = \frac {1}{D}[u_0 + {\displaystyle \varepsilon }^2 u_2+\cdots ]$
where the terms are
\begin{align} u_0(\mathbf{x}\,;\,\boldsymbol{\xi }) &= -|\Omega | \Big [ G(\mathbf{x}\,;\,\boldsymbol{\xi }) - R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \Big ] + \frac {|\Omega |}{2\pi \nu }, \\[-10pt] \nonumber \end{align}
\begin{align} u_2(\mathbf{x}\,;\,\boldsymbol{\xi }) &= |\Omega |\Big [ \frac 12 \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{\boldsymbol{\xi }} G(\mathbf{x}\,;\,\boldsymbol{\xi })\big ) -2\pi \mathbf{a} \cdot \mathcal{M} \nabla _{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi }) \Big ] + \chi _2, \\[8pt] \nonumber \end{align}
where
$\chi _2$
is the constant given in (2.52).
2.4. Calculation of the global MFPT
In this section, we calculate the GMFPT defined as
\begin{equation*} \tau = \frac {1}{|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|}\int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}}u \, d\mathbf{x} = \frac {1}{D|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|}\int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 \, d\mathbf{x}. \end{equation*}
The challenge as before is to determine the correct expansion accurate to
$\mathcal{O}({\displaystyle \varepsilon }^2)$
by accounting for contributions from the inner expansion. We first decompose the region of integration
$\Omega \setminus \Omega _{{\displaystyle \varepsilon }} := (\Omega \setminus B_{\delta })\cup (B_{\delta }\setminus \Omega _{{\displaystyle \varepsilon }})$
for the disk
$B_{\delta } = \{ \mathbf{x}\in \mathbb{R}^2 \ | \ |\mathbf{x}-\boldsymbol{\xi }|\lt \delta \}$
and then apply the limit
$\delta \to 0$
. This analysis is valid in the regime
${\displaystyle \varepsilon }\ll \delta \ll 1$
and results in a final answer independent of
$\delta$
. The integral becomes
\begin{equation*} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u \, d\mathbf{x} = \underbrace {\int _{B_{\delta }\setminus \Omega _{\displaystyle \varepsilon }} u \, d\mathbf{x}}_{I_1} + \underbrace {\int _{\Omega \setminus B_{\delta }} u \, d\mathbf{x}}_{I_2}. \end{equation*}
The integral
$I_1$
is evaluated by first transforming to the coordinate
$\mathbf{y} = e^{-i\phi }(\mathbf{x}-\boldsymbol{\xi })/{\displaystyle \varepsilon }$
followed by
$\mathbf{y} = \alpha \mathbf{z} + \beta /\mathbf{z}$
to find that
\begin{equation} I_1 = {\displaystyle \varepsilon }^2 \int _{\substack {\mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}\\ |\mathbf{y}|\lt \delta /{\displaystyle \varepsilon }}} U\, d\mathbf{y} = \frac {{\displaystyle \varepsilon }^2}{D} \int _{|\mathbf{z}|=1}^{|\mathbf{z}| = {\frac {\delta }{{\displaystyle \varepsilon }\alpha }}} (\underbrace {U_0}_{I_{10}} + \underbrace { {\displaystyle \varepsilon } U_1}_{{I_{11}}} + \underbrace { {\displaystyle \varepsilon }^2 U_2}_{I_{12}} + \cdots )|J|d\mathbf{z}. \end{equation}
From equation (2.24), we have that
$U_0(\mathbf{z})= S\nu \log |\mathbf{z}|$
so that for
$\mathbf{z} = r e^{i \theta }$
,
\begin{align} \nonumber I_{10} &= S\nu {\displaystyle \varepsilon }^2\int _{\theta = 0}^{2\pi }\int _{r=1}^{\frac {\delta }{{\displaystyle \varepsilon }\alpha } } \log r \left [ \alpha ^2- \frac {2\alpha \beta }{r^2}\cos 2\theta + \frac {\beta ^2}{r^4} \right ] r dr d\theta = 2\pi S\nu {\displaystyle \varepsilon }^2 \int _{r=1}^{r=\frac {\delta }{{\displaystyle \varepsilon }\alpha }}\log r \left [ \alpha ^2 r + \frac {\beta ^2}{r^3} \right ] dr \\[5pt] & = 2\pi S\nu {\displaystyle \varepsilon }^2 \Big [\alpha ^2 \Big (\frac {r^2}{2}\log r -\frac {r^2}{4}\Big ) - \beta ^2 \Big ( \frac {1+ 2\log r}{4 r^2} \Big ) \Big ]_{r=1}^{r=\frac {\delta }{{\displaystyle \varepsilon }\alpha }} \approx \frac {S\nu \pi }{2} \left [ 2\delta ^2\log \frac {\delta }{{\displaystyle \varepsilon }\alpha }-\delta ^2 + {\displaystyle \varepsilon }^2 (\alpha ^2 + \beta ^2)\right ].\nonumber\\[3pt] \end{align}
In the final calculation of (2.55), the
$\mathcal{O}({\displaystyle \varepsilon }^4)$
terms have been ignored. Following on to the term
$I_{11}$
, we apply from (2.37) that
$U_1 = -2\pi S \nu \mathbf{a} \cdot e^{-i\phi } v_{1c}$
where
$v_{1c}(\mathbf{z}) = (\mathbf{z} - \mathbf{z}/|\mathbf{z}|^2)$
. This leads to
$I_{11} = {\displaystyle \varepsilon }^3 \int _{|\mathbf{z}|=1}^{|\mathbf{z}|= {\delta /({\displaystyle \varepsilon }\alpha )}} U_1 |J| d\mathbf{z} = 0$
. The contribution from
$I_{22} = \mathcal{O}({\displaystyle \varepsilon }^4)$
.
The contribution from the outer region is now
\begin{equation} I_2 = \int _{\Omega \setminus B_{\delta }} u \, d\mathbf{x} = \frac {1}{D} \underbrace {\int _{\Omega \setminus B_{\delta }} u_0 \, d\mathbf{x}}_{I_{20}} + \frac {{\displaystyle \varepsilon }^2}{D} \underbrace {\int _{\Omega \setminus B_{\delta }} u_2 \, d\mathbf{x}}_{I_{22}}. \end{equation}
The first term is calculated with
$u_0 = -2\pi S\nu G(\mathbf{x}\,;\,\boldsymbol{\xi }) + \tau _0$
with
$\tau _0 = S(1+ 2\pi \nu R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }))$
to determine
\begin{align} \nonumber I_{20} &= \int _{\Omega \setminus B_{\delta } } u_0\, d \mathbf{x} = \int _{\Omega } u_0\,d\mathbf{x} - \int _{B_{\delta }} u_0\, d\mathbf{x} \\[4pt] \nonumber {} &= |\Omega | \tau _0 - \int _{B_{\delta } } [S\nu \log |\mathbf{x}-\boldsymbol{\xi }| +S]\, d\mathbf{x} -2\pi S\nu \int _{B_{\delta } } (R(\mathbf{x}\,;\,\boldsymbol{\xi })- R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) )\, d\mathbf{x}\\[4pt] {} &= |\Omega | \tau _0 - S\nu \pi \delta ^2\Big (\log \delta - \frac 12\Big ) - S\pi \delta ^2 + \mathcal{O}(\delta ^4). \end{align}
Following on to the second term in (2.56), we calculate that
\begin{align*} I_{22} = \int _{\Omega \setminus B_{\delta } } u_2 \, d\mathbf{x} = \int _{\Omega } u_2\,d\mathbf{x} - \int _{B_{\delta }} u_2\, d\mathbf{x}. \end{align*}
From (2.48b), the contributions of the second term have vanishing average, hence
\begin{equation*} \int _{B_{\delta }} u_2\, d\mathbf{x} = {\pi \delta ^2 d_{2c}}, \end{equation*}
where
$d_{2c}$
is defined in equation (2.43). Combining (2.55) and (2.57), we have that
\begin{align} \nonumber \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 \, d\mathbf{x} &= |\Omega |\tau _0 + {\displaystyle \varepsilon }^2 \Big ( \frac {S\nu \pi }{2} (\alpha ^2 + \beta ^2) + \int _{\Omega } u_2 \, d\mathbf{x} \Big ) + \mathcal{O}({\displaystyle \varepsilon }^3),\\ &= |\Omega | \tau _0 + {\displaystyle \varepsilon }^2 |\Omega | \Big ( \frac {a^2+b^2}{8} + \chi _2 + \int _{\Omega } \frac {1}{2} \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi })\big ) \\ \nonumber &\quad-2\pi \mathbf{a} \cdot \mathcal{M} \nabla _{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi }) d\mathbf{x} \Big )+ \mathcal{O}({\displaystyle \varepsilon }^3). \end{align}
Applying the identity (2.51), we note that the integral terms in (2.58) vanish since
\begin{equation} \int _{\Omega } (G_{\xi _1\xi _1} - G_{\xi _2\xi _2}) d \mathbf{x} = 0, \quad \int _{\Omega } G_{\xi _1\xi _2} d \mathbf{x} = 0, \quad \int _{\Omega } G_{\xi _1} d \mathbf{x} = 0, \quad \int _{\Omega } G_{\xi _2} d \mathbf{x} =0. \end{equation}
Hence, we obtain the final expression that incorporates the local and global contributions to the GMFPT,
\begin{equation} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 \, d\mathbf{x} = |\Omega | \tau _0 + {\displaystyle \varepsilon }^2 |\Omega | \Big ( \frac {a^2+b^2}{8} + \chi _2 \Big ). \end{equation}
Recalling the value of
$\chi _2$
given in (2.52), equation (2.58) then reduces to
\begin{align} \int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 \, d\mathbf{x} = |\Omega | \tau _0 + {\displaystyle \varepsilon }^2 |\Omega | \left ( \frac {a^2+b^2}{4} - |\Omega |\Big ( \mbox{Trace}\big (\mathcal{Q} \nabla _{\boldsymbol{\xi }}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) + \pi \frac {(a+b)^2}{2} |\mathbf{a}|^2 - 2\pi \mathbf{a} \cdot \mathcal{Q} \mathbf{a} \Big)\right)\!.\nonumber\\[-1pt] \end{align}
Finally, using
$|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}| = |\Omega | -\pi {\displaystyle \varepsilon }^2 ab$
, we have that
\begin{align} \nonumber \tau &= \frac {1}{D|\Omega \setminus \Omega _{{\displaystyle \varepsilon }}|}\int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 d\mathbf{x} = \frac {1}{D|\Omega |}\left [ 1 + {\displaystyle \varepsilon }^2 \frac {\pi a b}{|\Omega |} \right ]\int _{\Omega \setminus \Omega _{{\displaystyle \varepsilon }}} u_0 + {\displaystyle \varepsilon }^2 u_2 \, d\mathbf{x}\\[5pt] &= \frac {\tau _0}{D} + \frac {{\displaystyle \varepsilon }^2}{D} \left ( \frac {\pi ab}{|\Omega |}\tau _0 + \frac {a^2+b^2}{4} - |\Omega |\Big ( \mbox{Trace}\big (\mathcal{Q} \nabla _{\boldsymbol{\xi }}^2 R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) + \pi \frac {(a+b)^2}{2} |\mathbf{a}|^2 - 2\pi \mathbf{a} \cdot \mathcal{Q} \mathbf{a} \Big ) \right ). \end{align}
3. Results
In this section, we demonstrate the validity of the expansion and investigate how trap orientation modulates the GMFPT. In terms of minimizing the GMFPT, we remark that trap location is felt in the leading order through the term
$R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi })$
. Hence, when considering the role of orientation, we primarily focus on the correction terms
$u_2(\mathbf{x})$
and
$\tau _2$
. A globally optimizing configuration of the GMFPT would first locate the trap at the critical points of
$\tau _0$
followed by choosing the orientation that optimizes
$\tau _2$
.
Example 3.1. Convergence in the disk geometry with a single elliptical trap.
In this case where the enclosing geometry
$\Omega$
is the unit disk, we have exact formula for the Green’s function, and we can obtain the necessary derivatives. From the equations (A.73), we calculate that
\begin{align} \mathbf{a} = \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) &= \frac {1}{2\pi }\frac { 2 - |\boldsymbol{\xi }|^2}{1-|\boldsymbol{\xi }|^2} \boldsymbol{\xi }; \quad \mbox{Trace}\big (\mathcal{Q} \nabla ^2_{\boldsymbol{\xi }} R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \big ) = \frac {1}{\pi } \frac {\boldsymbol{\xi } \cdot \mathcal{Q} \boldsymbol{\xi }}{(1-|\boldsymbol{\xi }|^2)^2}\, ; \quad \mathbf{a} \cdot \mathcal{Q} \mathbf{a} =\frac {1}{4\pi ^2} \frac {(2-|\boldsymbol{\xi }|^2)^2}{(1-|\boldsymbol{\xi }|^2)^2}\, \boldsymbol{\xi } \cdot \mathcal{Q} \boldsymbol{\xi }; \nonumber \\[5pt] S\nu &= \frac {|\Omega | }{2\pi } = \frac {1 }{2}; \qquad \tau _0 = \frac {|\Omega |}{2\pi \nu } \Big [1 + 2\pi \nu R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \Big ] = \frac {1}{2} \Big [ \frac {1}{\nu } - \log (1-|\boldsymbol{\xi }|^2) + |\boldsymbol{\xi }|^2 -\frac 34 \Big ]. \end{align}
Applying (3.63) to (2.62) together with
$|\Omega | = \pi$
and
$D=1$
yields the GMFPT
\begin{align} \nonumber \tau & = \tau _0 + {\displaystyle \varepsilon }^2 \left [ ab\tau _0 + \frac {a^2+b^2}{4} - \frac {(a+b)^2}{8} \left (\frac {2-|\boldsymbol{\xi }|^2}{1-|\boldsymbol{\xi }|^2}\right )^2 |\boldsymbol{\xi }|^2 - \frac {2- (2- |\boldsymbol{\xi }|^2)^2}{2(1-|\boldsymbol{\xi }|^2)^2} \boldsymbol{\xi } \cdot \mathcal{Q} \boldsymbol{\xi } \right ]\\[5pt] & = \tau _0 + {\displaystyle \varepsilon }^2 \left [ab\tau _0 + \frac {a^2+b^2}{4} - \frac {(a+b)^2}{8} \left (\frac {2-|\boldsymbol{\xi }|^2}{1-|\boldsymbol{\xi }|^2}\right )^2 |\boldsymbol{\xi }|^2 + \frac {a^2-b^2}{4} \frac {2 - (2- |\boldsymbol{\xi }|^2)^2}{2(1-|\boldsymbol{\xi }|^2)^2} \begin{bmatrix} \xi _1^2 - \xi _2^2\\[5pt] 2\xi _1 \xi _2 \end{bmatrix} \cdot \begin{bmatrix} \cos 2\phi \\[5pt] \sin 2\phi \end{bmatrix}\right ].\nonumber\\[3pt] \end{align}
In Figure 2, we show agreement between finite element solutions of (1.1) and the asymptotic result (3.64). The common variables in this validation are the disk geometry, the trap centre
$\boldsymbol{\xi } = (0.3,0.4)$
and the semi-axes dimensions
$(a,b)=(3,1)$
. In Figures 2a–Figures 2b, we show agreement of the correction terms
$u_2 = (u(\mathbf{x})- u_0(\mathbf{x}))/{\displaystyle \varepsilon }^2$
at the point
$\mathbf{x} = (\!-0.2,-0.4)$
and
$\tau _2 = (\tau - \tau _0)/{\displaystyle \varepsilon }^2$
as the trap orientation
$\phi$
varies with
$ {\displaystyle \varepsilon } = 0.03$
fixed. A schematic of the domain is shown in Figure 3a.
Convergence of the asymptotic approximation (1.4) in the disk case with a trap centred at
$\boldsymbol{\xi } = (0.3,0.4)$
. Panel (a): Agreement between the solution correction
$u_2(\mathbf{x}) = {\displaystyle \varepsilon }^{-2}(u(\mathbf{x})-u_0(\mathbf{x}))$
for
$\mathbf{x} = (\!-0.2,-0.4)$
. Panel (b): The GMFPT correction
$\tau _2= {\displaystyle \varepsilon }^{-2}(\tau -\tau _0)$
from numerical and asymptotic approximations for
${\displaystyle \varepsilon } = 0.03$
as orientation
$\phi$
varies. Panels (c – d): Convergence as
${\displaystyle \varepsilon } \to 0$
of the relative errors between numerical and asymptotic approximations (leading and correction) of
$u(\mathbf{x})$
for
$\mathbf{x} = (\!-0.2,-0.4)$
(c), and
$\tau$
with fixed
$\phi = \pi /6$
(d). Straight lines are of slope
$2$
(blue) and
$4$
(red) indicating convergence rates. Domain schematic shown in Figure 3a.
To directly confirm the convergence of the asymptotic expansion, we calculate a sequence of relative errors
\begin{equation} \mathcal{E}_{\textrm {rel}}[z] = \left |\frac {z_{\text{approx}}({\displaystyle \varepsilon }) - z_{\text{true}}}{z_{\text{true}}}\right |, \end{equation}
at a range of
$\displaystyle \varepsilon$
values. In (3.65) the true value is calculated from finite element simulations of (1.1) and the approximate values come from taking one or two terms in the expansion (1.4a). In Figures 2c–Figures 2d, we show that the rate of convergence in
$\mathcal{E}_{\textrm {rel}}$
as
${\displaystyle \varepsilon }\to 0^{+}$
is in agreement with the principal result (1.4), and in particular, the two term approximation (3.64) of the GMFPT has error
$\mathcal{O}({\displaystyle \varepsilon }^4)$
.
Example 3.2. Optimization of trap orientation for the GMFPT in the unit disk domain.
To study the optimizing trap orientations, we set
$\boldsymbol{\xi } = r e^{i\theta }$
in (3.64) which reduces the GMFPT to
\begin{align} &\tau \sim \tau _0 + {\displaystyle \varepsilon }^2 \left [ ab\tau _0 +\frac {a^2+b^2}{4} - \frac {(a+b)^2}{8} r^2\left (\frac {2-r^2}{1-r^2}\right )^2 + \frac {a^2-b^2}{4} g(r) \cos 2(\theta -\phi ) \right ]\!,\nonumber\\[5pt] &\quad g(r)= \frac {2 - (2- r^2)^2}{2(1-r^2)^2}r^2.\nonumber\\[1pt] \end{align}
Minimization of the GMFPT correction in the disk with a single elliptical trap. Panel (a): Domain with a single elliptical trap at
$\boldsymbol{\xi } = (0.3,0.4)$
, axes
${\displaystyle \varepsilon }(a,b) = {\displaystyle \varepsilon }(3,1)$
and orientation
$\phi = \pi /6$
. The highlighted point (black dot) is
$\mathbf{x}=(\!-0.2,-0.4)$
. Panel (b): The function
$g(r)$
and the critical radius
$r=r_c$
. For
$r\lt r_c$
, the GMFPT correction is minimized when the ellipse has major axis pointed towards the centre of the disk.
The function
$g(r)$
has the following (see Figure 3b) simple properties,
\begin{equation*} g(r) \to 0^{-},\quad \text{as}\quad r \to 0^{+}; \qquad g(r)\to +\infty , \quad \text{as}\quad r \to 1^{-}; \qquad g(r_c) = 0, \quad r_c = \sqrt {2-\sqrt {2}}, \end{equation*}
and hence we can conclude that the correction
$\tau _2$
of the GMFPT is minimized when the semi-major axis of the trap is aligned in the radial direction, i.e.
$\phi = \theta$
, provided
$g(r)\lt 0$
or
$(2-r^2)^2\gt 2$
. The GMFPT minimizing configuration is therefore
\begin{equation} \phi = \left \{ \begin{array}{rl} \theta , & r \in (0,r_c),\\[5pt] \theta +\frac {\pi }{2}, & r \in (r_c,1); \end{array} \right . \qquad r_c = \sqrt {2-\sqrt {2}} \approx 0.7654. \end{equation}
In Figure 3b, we plot the function
$g(r)$
together with corresponding trap orientations that minimize the GMFPT. We also remark that the trap centre, which minimizes the leading order term
$\tau _0$
in GMFPT is
$\boldsymbol{\xi } = (0,0)$
. For a trap centred at this location, there is no contribution (
$g(0)=0$
) to the GMFPT from the trap orientation as expected by symmetry considerations.
Example 3.3. The MFPT from the origin to an elliptical trap in a disk.
In this example, we calculate the MFPT for a particle in the disk domain initially at
$\mathbf{x}=0$
to arrive at an elliptical trap centred at
$\boldsymbol{\xi } = re^{i\theta }$
. From the main result (1.4c), we have the correction term
\begin{equation*} u_2(\mathbf{x}\,;\,\boldsymbol{\xi }) = |\Omega |\Big [ \frac 12 \mbox{Trace} \big ( \mathcal{Q} \nabla ^2_{\boldsymbol{\xi }} G(\mathbf{x}\,;\,\boldsymbol{\xi })\big ) -2\pi \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi }) \cdot \mathcal{M} \nabla _{{\boldsymbol{\xi }}} G(\mathbf{x}\,;\,\boldsymbol{\xi }) \Big ] + \chi _2, \end{equation*}
which describes the role that trap orientation plays in the capture rate. In Appendix A.1, we calculate the relevant terms to be
\begin{align} \text{Trace}(\mathcal{Q} \nabla _{\boldsymbol{\xi }}^2G(\textbf {0}\,;\,\boldsymbol{\xi })) &= \frac {1}{\pi }\frac {\boldsymbol{\xi } \cdot \mathcal{Q} \boldsymbol{\xi }}{|\boldsymbol{\xi }|^4}; \\[-10pt] \nonumber \end{align}
\begin{align} \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi }\,;\,\boldsymbol{\xi })\cdot \mathcal{M}\nabla _{\boldsymbol{\xi }}G(\textbf {0}\,;\,\boldsymbol{\xi }) &= -\frac {1}{4\pi ^2}\frac {2-|\boldsymbol{\xi }|^2}{|\boldsymbol{\xi }|^2} \left [ \boldsymbol{\xi }\cdot \mathcal{Q} \boldsymbol{\xi } - \frac {(a+b)^2}{4} |\boldsymbol{\xi }|^2 \right ]. \\[8pt] \nonumber \end{align}
We then have that
\begin{equation} u_2(\textbf {0}) = \frac {a^2+b^2}{8} - \frac {(a+b)^2}{8} \frac {(2-r^2)}{(1-r^2)^2} +\frac {a^2-b^2}{8} \frac {2r^4-1}{r^2(1-r^2)^2} \cos \big ( 2(\phi -\theta )\big ). \end{equation}
Favourable agreement between this result and numerical simulations is shown in Figure 4. We conclude that
$u_2(\textbf {0})$
is minimised by orientating the trap in the radial direction (
$\theta = \phi$
) for
$2r^4-1\lt 0$
or
$0\lt r\lt 2^{-\frac 14}\approx 0.8409$
. For
$2^{-\frac 14}\lt r\lt 1$
,
$u_2(\textbf {0})$
is minimised by orientating the semi-major axis of the trap parallel to the boundary.
The effects of trap orientation on the MFPT staring at the centre of a disk. The correction
$u_2 = {\displaystyle \varepsilon }^{-2}(u-u_0)$
from (3.68c) with a single trap of extent
${\displaystyle \varepsilon } = 0.01$
, semi-major axes
$(a,b)=(3,1)$
centred at
$\boldsymbol{\xi }=(r,0)$
. Curves shown for orientations
$\phi =\pi /2$
and
$\phi = 0$
.
Example 3.4. Here we consider a general domain
$\Omega$
with a single elliptical trap centred at
$\mathbf{x}=\boldsymbol{\xi }$
of semi-major and semi-minor axes
${\displaystyle \varepsilon } a$
and
${\displaystyle \varepsilon } b$
, respectively, and orientation
$\phi$
with respect to the horizontal.
The GMFPT (2.61) reduces to the form
\begin{equation} \tau = \tau _0 + {\displaystyle \varepsilon }^2 \left [\frac {\pi ab}{|\Omega |}\tau _0 + \frac {a^2+b^2}{4} - \pi |\Omega | \frac {(a+b)^2}{2} (R_{\xi _1}^2 + R_{\xi _2}^2) + |\Omega | \frac {a^2-b^2}{4} \mathbf{p} \cdot \begin{bmatrix}\cos 2\phi \\[5pt] \sin 2\phi \end{bmatrix} \right ]. \end{equation}
The vector
$\mathbf{p}$
is given as
\begin{equation} \mathbf{p} = \begin{bmatrix} R_{\xi _1\xi _1}-R_{\xi _2\xi _2} - 2\pi (R_{\xi _1}^2 -R_{\xi _2}^2) \\[5pt] 2R_{\xi _1\xi _2} - 4\pi R_{\xi _1} R_{\xi _2} \end{bmatrix} = \begin{bmatrix} R_{\xi _1\xi _1}-R_{\xi _2\xi _2} \\[5pt] 2R_{\xi _1\xi _2} \end{bmatrix} - 2\pi \begin{bmatrix} R_{\xi _1}^2 -R_{\xi _2}^2 \\[5pt] 2 R_{\xi _1} R_{\xi _2} \end{bmatrix}. \end{equation}
For some regular domains such as ellipses and rectangles, the function
$R(\mathbf{x}\,;\,\boldsymbol{\xi })$
is available in the form of rapidly convergent series (see Appendices A.2 and A.3). The necessary derivatives
$R_{\xi _1\xi _1},R_{\xi _2\xi _2},R_{\xi _1\xi _2},R_{\xi _1}$
and
$R_{\xi _2}$
in (3.70) can then be calculated by finite differences. We remark that if the trap centre is chosen to optimize
$\tau _0$
, then we have
$R_{\xi _1} = R_{\xi _2} = 0$
and so
\begin{equation} \tau = \tau _0 + {\displaystyle \varepsilon }^2 \left [\frac {\pi ab}{|\Omega |}\tau _0 + \frac {a^2+b^2}{4} + |\Omega | \frac {a^2-b^2}{4} \begin{bmatrix} R_{\xi _1\xi _1}-R_{\xi _2\xi _2} \\[4pt] 2R_{\xi _1\xi _2}. \end{bmatrix} \cdot \begin{bmatrix}\cos 2\phi \\[5pt] \sin 2\phi \end{bmatrix} \right ]. \end{equation}
Minimisation of
$\tau _2$
for a single elliptical trap placed in the rectangular domain
$\Omega = [0,L]\times [0,d]$
for
$d = 1$
and
$L = 1$
(a),
$L=1.01$
, (b)
$L=1.02$
, (c)
$L=1.04$
, (d)
$L=1.1$
, (e)
$L=1.5$
. The directional arrow indicates the direction on which the semi-major axis should be aligned to minimize
$\tau _2$
, the higher-order GMFPT correction term.
In Figure 5 and Figure 6, we plot vector fields of the orientation vector
$\textbf {p}$
for the example of rectangular and elliptical domains, respectively. In each case, we form the necessary derivatives of
$R(\mathbf{x}\,;\,\boldsymbol{\xi })$
by combining the rapidly convergent series stated in Appendix A together with centred finite difference approximations for the first and second derivatives. As the disk domain deforms into an ellipse, the bifurcation of the minimizing orientation noted in Figure 3 is smoothed out. However, a generic observation remains that for traps centred close to a smooth boundary,
$\tau _2$
is minimized by orienting the semi-major axis of the trap parallel to
$\partial \Omega$
.
In the case of the rectangular domain
$[0,L]\times [0,d]$
, we observe similar discontinuous structures that deform from the square case (
$L=d$
) as the rectangle elongates
$(L\gt d)$
. Curiously, the minimizing orientation of
$\tau _2$
at the corners is observed to be when the semi-major axis is aligned into the corners.
Minimization of
$\tau _2$
for circular (a) and elliptical domains (b – c) at various locations. The directional arrow indicates the direction along which the semi-major axis should be aligned so that the correction term to the GMFPT is minimized. In panel (a), the dashed blue line is the disk of radius
$r_c = \sqrt {2-\sqrt {2}} \approx 0.7654$
where the optimal orientation flips.
Example 3.5. The limit of an infinitely thin ellipse to a slit.
In this example, we consider the GMFPT in the limit as
$b\to 0$
as the elliptical trap tends towards a thin slit. The formula for the GMFPT is uniformly valid in this limit and we find from (3.69) that
\begin{equation} \lim _{b\to 0^+}\tau = \tau _0 + {\displaystyle \varepsilon }^2 \left [ \frac {a^2}{4} - \frac {\pi a^2|\Omega |}{2} (R_{\xi _1}^2 + R_{\xi _2}^2) + \frac {a^2|\Omega | }{4} \mathbf{p} \cdot \begin{bmatrix}\cos 2\phi \\[5pt] \sin 2\phi \end{bmatrix} \right ]. \end{equation}
In the case of a rectangular domain, we use this formula to explore the effect of trap orientation on the GMFPT. In Figure 7, we display
$\tau _2$
for a trap with extent
${\displaystyle \varepsilon }=0.2$
, centred at
$\boldsymbol{\xi } = [0.3,0.4]$
inside the rectangular domain
$\Omega = [0,1]\times [0,0.8]$
. The two curves plotted are for trap orientations
$\phi = \{\pi /2,\pi /6\}$
and for
$a=1$
and varying
$0\lt b\lt 1$
. As expected, we see no effects of orientation when
$a=b=1$
and smooth behaviour as
$b\to 0$
.
The effects of trap orientation and ellipticity on the high-order correction to the GMFPT in the limit as
$b\to 0$
. The correction
$\tau _2 = {\displaystyle \varepsilon }^{-2}(\tau -\tau _0)$
to the GMFPT for a rectangular domain
$\Omega = [0,1]\times [0,0.8]$
with a single trap of extent
${\displaystyle \varepsilon } = 0.2$
, semi-major axis
$a=1$
centred at
$\boldsymbol{\xi }=[0.3,0.4]$
and varying semi-minor axis
$b$
. Curves shown for orientations
$\phi =\pi /2$
and
$\phi = \pi /6$
, which coincide for circular traps (
$b=1$
).
4. Discussion
In this work, we have developed a high-order matched asymptotic expansion for the MFPT of a Brownian walker to a small trap, enclosed in a two dimensional domain
$\Omega$
. The high-order correction term describes the effect that the orientation of the trap has on the capture rate. We investigated the role that trap orientation has on the MFPT and the GMFPT, observing a sensitive dependence on the centring point of the trap in the domain. In the specific case where the enclosing domain is a disk, we found a bifurcation where the correction to the GMFPT is minimized by orientating the semi-major axis of the trap in the radial direction when the centring point satisfies
$0\lt |\boldsymbol{\xi }|\lt \sqrt {2-\sqrt {2}}$
and is minimized by orientating in the angular direction when
$\sqrt {2-\sqrt {2}}\lt |\boldsymbol{\xi }|\lt 1$
(Figure 3). The discontinuous nature of this transition appears to be related to the symmetries of the domain, and similar effects are noted in rectangular domains (Figure 5). In domains with smooth boundaries, such as ellipses, we observe that the discontinuity in the vector field of optimal directions is smoothed out (Figure 6). However, we observe that the GMFPT correction is generally minimized when the semi-major axis of the trap is aligned parallel to the boundary.
Our exposition of the role of trap orientation has focused on the case of an elliptical trap. This simple geometry allowed for explicit calculation of several key quantities, in particular the logarithmic capacitance
$d_c$
, the quadrupole matrix
$\mathcal{Q}$
and the moment polarization tensor
$\mathcal{M}$
. The determination of these key quantities was facilitated by a complex variable approach that utilized the known mapping between the disk and ellipses to solve three variations of Laplace problems (Appendix B). However, our Principal Result (1.4) is valid for any trap geometry that features two lines of symmetry, for example those with rectangular or dumbbell shape. We hypothesize that the correction term derived in (1.4) would vanish for trap geometries with additional lines of symmetries, such as equilateral triangles or square. In such cases, further corrections would be necessary to describe the effects of orientation. To apply the present results to these geometries, it would be necessary to calculate the three previously mentioned key quantities, either by a suitable complex transformation [Reference Kolokolnikov, Titcombe and Ward24] or numerical method [Reference Baddoo and Trefethen2]. A similar set of obstacles would arise in extending this analysis to the case of a Robin boundary condition
$\partial _n u+ \kappa u =0$
on each trap. It is likely the relevant inner problems can be more easily solved by an application of the elliptical coordinate system.
In future work, we aim to derive the MFPT in the presence of multiple absorbing elliptical bodies. This will involve generalizing the outer solution to a superposition of Green’s functions as performed in [Reference Kolokolnikov, Titcombe and Ward24] followed by matching of gradient and Hessian terms to local problems [Reference Lindsay, Kolokolnikov and Tzou32]. This situation is relevant in the consideration of fly muscle cells, which contain multiple well spaced nuclei of elliptical shape and correlated orientations [Reference Windner, Manhart, Brown, Mogilner and Baylies45]. This work shows the steps required to derive higher-order corrections in a variety of narrow capture problems.
Financial support
AEL acknowledges support from the NSF under award DMS 2052636.
Competing interests
The authors’ declare no conflicts of interest.
Appendix A. Neumann Green’s functions for disks, rectangles and ellipses
Here we state some known expressions for the Neumann Green’s function (1.5) for the disk, ellipse and rectangle domains.
A.1. Neumann Green’s function for a disk
In the case of the disk domain
$\Omega = \{ \mathbf{x} = (x_1,x_2) \ | \ x_1^2 + x_2^2 \leq 1\}$
and source
$\boldsymbol{\xi } = (\xi _1,\xi _2)$
, we have [Reference Kolokolnikov, Titcombe and Ward24] that
\begin{equation*} G = \frac {-1}{2\pi } \log |\mathbf{x}-\boldsymbol{\xi }| + R(\mathbf{x};\boldsymbol{\xi }); \qquad R(\mathbf{x};\boldsymbol{\xi }) = -\frac {1}{2\pi } \left [ \frac 12 \log (1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi }) - \frac {1}{2}(|\boldsymbol{\xi }|^2 + |\mathbf{x}|^2) + \frac {3}{4} \right ]. \end{equation*}
We calculate gradients
$\nabla _{\mathbf{x}} = (\partial _{x_1},\partial _{x_2})$
,
$\nabla _{\boldsymbol{\xi }} = (\partial _{\xi _1},\partial _{\xi _2})$
as
\begin{gather} \nabla _{\mathbf{x}}G = \frac {-1}{2\pi } \frac {\mathbf{x} - \boldsymbol{\xi }}{|\mathbf{x}-\boldsymbol{\xi }|^2} + \nabla _{\mathbf{x}}R ; \qquad \nabla _{\mathbf{x}}R = \frac {-1}{2\pi }\left [ \frac {|\boldsymbol{\xi }|^2 \mathbf{x} - \boldsymbol{\xi }}{1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi }} - \mathbf{x} \right ], \end{gather}
\begin{gather} \nabla _{\boldsymbol{\xi }}G = \frac {1}{2\pi } \frac {\mathbf{x} - \boldsymbol{\xi }}{|\mathbf{x}-\boldsymbol{\xi }|^2} + \nabla _{\boldsymbol{\xi }}R; \qquad \nabla _{\boldsymbol{\xi }}R = \frac {-1}{2\pi }\left [ \frac {|\mathbf{x}|^2 \boldsymbol{\xi } - \mathbf{x}}{1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi }} - \boldsymbol{\xi } \right ]. \end{gather}
As
$\mathbf{x}\to \boldsymbol{\xi }$
, we have
\begin{equation*} R(\boldsymbol{\xi };\boldsymbol{\xi }) =\frac {-1}{2\pi } \left [ \log (1- |\boldsymbol{\xi }|^2) - |\boldsymbol{\xi }|^2 + \frac 34 \right ], \qquad \nabla _{\boldsymbol{\xi }}R(\boldsymbol{\xi };\boldsymbol{\xi }) = \frac {1}{2\pi }\left [ \frac {2-|\boldsymbol{\xi }|^2}{1-|\boldsymbol{\xi }|^2} \right ] \boldsymbol{\xi } . \end{equation*}
The second derivatives for the Hessian are
\begin{align} \frac {\partial ^2 R}{\partial {\xi _j}^2} = \frac {-1}{2\pi } \left [ \frac {|\mathbf{x}|^2(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2\mathbf{x}\cdot \boldsymbol{\xi }) - 2(|\mathbf{x}|^2 \xi _j - x_j)^2 }{(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi })^2} -1\right ],\quad \frac {\partial ^2 R}{\partial {\xi _1\partial \xi _2}} = \frac {1}{\pi } \frac {(|\mathbf{x}|^2 \xi _1-x_1)(|\mathbf{x}|^2 \xi _2 - x_2)}{(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi })^2},\nonumber\\[3pt] \end{align}
\begin{align} \frac {\partial ^2 R}{\partial {x_j}^2} = \frac {-1}{2\pi } \left [ \frac {|\boldsymbol{\xi }|^2(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2\mathbf{x}\cdot \boldsymbol{\xi }) - 2(|\boldsymbol{\xi }|^2 x_j - \xi _j)^2 }{(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi })^2} -1\right ],\quad \frac {\partial ^2 R}{\partial x_1 \partial x_2} = \frac {1}{\pi } \frac {(|\boldsymbol{\xi }|^2 x_1-\xi _1)(|\boldsymbol{\xi }|^2 x_2 - \xi _2)}{(1+ |\mathbf{x}|^2 |\boldsymbol{\xi }|^2 - 2 \mathbf{x}\cdot \boldsymbol{\xi })^2}.\nonumber\\[3pt] \end{align}
The terms
$(R_{\xi _1\xi _1} - R_{\xi _2 \xi _2})$
and
$R_{\xi _1\xi _2}$
as
$\mathbf{x}\to \boldsymbol{\xi }$
are then
\begin{equation} \lim _{\mathbf{x}\to \boldsymbol{\xi }} \left ( \frac {\partial ^2 R}{\partial {\xi _1}^2}-\frac {\partial ^2 R}{\partial {\xi _2}^2} \right ) = \frac {1}{\pi } \frac {\xi _1^2 - \xi _2^2}{(1-|\boldsymbol{\xi }|^2)^2}, \qquad \lim _{\mathbf{x}\to \boldsymbol{\xi }} \frac {\partial ^2 R}{\partial \xi _1 \partial \xi _2} = \frac {1}{\pi } \frac {\xi _1\xi _2}{(1-|\boldsymbol{\xi }|^2)^2}. \end{equation}
A.2. Neumann Green’s function for a rectangle
The Green’s function for a rectangle
$\Omega = [0,L]\times [0,d]$
is known [Reference Chen and Ward11, Reference Kolokolnikov, Ward and Wei25] in the form of a rapidly convergent series. For
$\mathbf{x} = (x_1,x_2)$
and
$\boldsymbol{\xi } = (\xi _1,\xi _2)$
, we have
\begin{align} \nonumber \kern-8pt R(\mathbf{x};\boldsymbol{\xi }) &= \frac {-1}{2\pi } \sum _{n=0}^{\infty } \log ( |1-q^n z_{+,+}||1-q^n z_{+,-}||1-q^n z_{-,+}||1-q^n \zeta _{+,+}||1-q^n \zeta _{+,-}||1-q^n \zeta _{-,+}||1-q^n \zeta _{-,-}|)\\ &\quad -\frac {1}{2\pi } \log \frac {|1- z_{-,-}|}{|r_{-,-}|} + \frac {L}{d}\left [ \frac 13 - \frac {\max (x_1,\xi _1)}{L} + \frac {x_1^2+\xi _1^2}{2L^2}\right ] - \frac {1}{2\pi } \sum _{n=1}^{\infty } \log | 1- q^n z_{-,-}|. \end{align}
where
\begin{align} z_{\pm ,\pm } &\equiv e^{\mu r_{\pm ,\pm }/2}, \qquad \zeta _{\pm ,\pm } \equiv e^{\mu \rho _{\pm ,\pm }/2}, \qquad \mu \equiv \frac {2\pi }{d}, \quad q \equiv e^{-\mu L},\\[-10pt]\nonumber \end{align}
\begin{align} r_{+,\pm } &\equiv -|x_1 + \xi _1| + i (x_2 \pm \xi _2), \qquad r_{-,\pm } \equiv -|x_1 - \xi _1| + i (x_2 \pm \xi _2),\\[-6pt]\nonumber \end{align}
\begin{align} \rho _{+,\pm } &\equiv |x_1 + \xi _1| + i (x_2 \pm \xi _2) -2L, \qquad \rho _{-,\pm } \equiv |x_1 - \xi _1| + i (x_2 \pm \xi _2) - 2L. \end{align}
The self-interaction term
$R(\mathbf{x};\mathbf{x})$
is given by
\begin{align} \nonumber R(\mathbf{x};\mathbf{x}) &= \frac {-1}{2\pi } \sum _{n=0}^{\infty } \log ( |1-q^n z^0_{+,+}||1-q^n z^0_{+,-}||1-q^n z^0_{-,+}||1-q^n \zeta ^0_{+,+}||1-q^n \zeta ^0_{+,-}||1-q^n \zeta ^0_{-,+}||1-q^n \zeta ^0_{-,-}|)\\[5pt] {}&\quad + \frac {L}{d} \left ( \frac 13 - \frac {x_1}{L} + \frac {x_1^2}{L^2} \right ) -\frac {1}{2\pi }\log \left ( \frac {\pi }{d}\right ) - \frac {1}{2\pi }\sum _{n=1}^{\infty } \log (1- q^n ). \end{align}
where
\begin{align} z^0_{+,+} &\equiv e^{\mu (-x_1 + i x_2)}, \qquad z^0_{+,-} \equiv e^{-\mu x_1 }, \qquad z^0_{-,+} \equiv e^{\mu i x_2 },\\[-6pt]\nonumber \end{align}
\begin{align} \zeta ^0_{+,+} & \equiv e^{\mu (x_1 -L + i x_2)}, \qquad \zeta ^0_{-,+} \equiv e^{\mu (-L + i x_2)}, \qquad \zeta ^0_{+,-} \equiv e^{\mu (x_1 -L)}, \qquad \zeta ^0_{-,-} \equiv e^{- \mu L}. \end{align}
A.3. Neumann Green’s function for an ellipse
A rapidly convergent series for the solution of (1.5) in the elliptical domain
$\Omega = \{ \mathbf{x} = (x_1,x_2) \; | \; (x_1/a)^2 + (x_2/b)^2 \leq 1\}$
was derived in [Reference Iyaniwura, Wong, Macdonald and Ward21]. For completeness, we restate the final result here. The first step is to introduce the transformation,
\begin{equation} x_1 = f \cosh \xi \cos \eta , \qquad x_2 = f \sinh \xi \sin \eta , \qquad f = \sqrt {a^2-b^2}, \end{equation}
which maps
$\mathbf{x} = (x_1,x_2)\in \Omega$
to the rectangle
$0\leq \xi \leq \xi _b$
and
$0\leq \eta \leq 2\pi$
where
$a = f \cosh \xi _b$
and
$b = f\sinh \xi _b$
so that
\begin{equation} f = \sqrt {a^2-b^2}, \qquad \xi _b = \tanh ^{-1} \frac {b}{a} = -\frac {1}{2} \log \gamma , \qquad \gamma = \left (\frac {a-b}{a+b}\right ). \end{equation}
For a pair
$(x_1,x_2)$
, the corresponding
$(\xi ,\eta )$
satisfy
\begin{equation*} \xi = \frac 12 \log \Big (1-2s + s\sqrt {s^2-s} \Big ), \qquad s = \frac {-\mu - \sqrt {\mu ^2 + 4 f^2 y^2}}{2f^2}, \qquad \mu = x_1^2 + x_2^2 - f^2. \end{equation*}
For
$\eta _{\ast } = \sin ^{-1} (\sqrt {p})$
, the value of
$\eta$
is given by
\begin{align} \eta = \left \{ \begin{array}{lr} \eta _{*}, &\quad \text{for } x_1 \geq 0,\ x_2 \geq 0\\[5pt] \pi -\eta _{*}, &\quad \text{for } x_1 \lt 0,\ x_2 \geq 0\\[5pt] \pi +\eta _{*}, &\quad \text{for } x_1 \leq 0,\ x_2 \lt 0\\[5pt] 2\pi -\eta _{*}, &\quad \text{for } x_1 \gt 0,\ x_2 \lt 0\\ \end{array} \right . \, , \qquad \text{where} \qquad p = \frac {-\mu + \sqrt {\mu ^2 + 4f^2 y^2}}{2f^2}. \end{align}
For points
$\mathbf{x}= (x_1,x_2)$
and
$\mathbf{y} = (y_1,y_2)$
, the Green’s function
$G(\mathbf{x};\mathbf{y})$
for
$\mathbf{x}\neq \mathbf{y}$
is given by
\begin{align} G(\mathbf{x};\mathbf{y}) &= \frac {1}{4|\Omega |}(|\mathbf{x}|^2 + |\mathbf{y}|^2) - \frac {3}{16|\Omega |}(a^2 + b^2) - \frac {1}{4\pi }\log \gamma - \frac {1}{2\pi }\max (\xi ,\xi _0)\nonumber\\[5pt] &\quad- \frac {1}{2\pi }\sum _{n=0}^{\infty } \log \left ( \prod _{j=1}^8 |1- \gamma ^{2n} z_j| \right )\!, \end{align}
where
$|\Omega | = \pi ab$
. The complex constants
$z_1,\ldots , z_8$
are defined in terms of
$(\xi ,\eta )$
,
$(\xi _0,\eta _0)$
and
$\xi _b$
by
\begin{gather*} z_1 = e^{-|\xi -\xi _0| + i (\eta -\eta _0)}, \quad z_2 = e^{|\xi -\xi _0| - 4\xi _b + i (\eta -\eta _0)}, \quad z_3 = e^{(\xi +\xi _0) - 2\xi _b + i (\eta -\eta _0)},\\[5pt] z_4 = e^{(\xi +\xi _0) - 2\xi _b + i (\eta -\eta _0)}, \quad z_5 = e^{(\xi +\xi _0) - 4\xi _b + i (\eta +\eta _0)}, \quad z_6 = e^{-(\xi +\xi _0) + i (\eta +\eta _0)},\\[5pt] z_7 = e^{|\xi +\xi _0| - 2\xi _b + i (\eta +\eta _0)}, \quad z_8 = e^{-|\xi +\xi _0| - 2\xi _b + i (\eta +\eta _0)}. \end{gather*}
The point
$(x_1,x_2)$
is mapped to
$(\xi ,\eta )$
while the source point
$(y_1,y_2)$
is mapped to
$(\xi _0,\eta _0)$
by the transformation (A.76). The quantity
$R(\mathbf{y};\mathbf{y})$
is given by
\begin{align*} R(\mathbf{y};\mathbf{y}) &= \frac {|\mathbf{y}|^2}{2|\Omega |} - \frac {3}{16|\Omega |}(a^2 + b^2) + \frac {1}{2\pi } \log (a+b) - \frac {\xi _0}{2\pi } + \frac {1}{4\pi } \log \big ( \cosh ^2 \xi _0 - \cos ^2 \eta _0 \big )\\[5pt] &- \frac {1}{2\pi } \sum _{n=1}^{\infty } \log (1- \gamma ^{2n}) - \frac {1}{2\pi } \sum _{n=0}^{\infty } \log \left ( \prod _{j=2}^8 |1-\gamma ^{2n} z_j^0| \right ). \end{align*}
Here the constants
$z_j^0$
for
$j=2,\ldots ,8$
are
\begin{gather*} z_2^0 = \gamma ^2, \quad z_3^0 = \gamma e^{-2\xi _0}, \quad z_4^0 = \gamma e^{2\xi _0}, \quad z_4^0 = \gamma ^2 e^{2\xi _0+2i\eta _0},\\[4pt] z_6^0 = e^{-2\xi _0 + 2i \eta _0}, \quad z_7^0 = \gamma e^{2i\eta _0}, \quad z_8^0 = \gamma e^{2i \eta _0}, \qquad \gamma = \frac {a-b}{a+b}. \end{gather*}
Appendix B. Inner problems for the exterior of the ellipse
We solve a variety of Laplace equations posed in the exterior of the elliptical domain
$\mathcal{A} = \{ \mathbf{y} = (y_1,y_2) \ | \ y_1^2/a^2 + y_2^2/b^2 \lt 1 \}$
. In each case, we make use of the complex transformation
\begin{equation} \mathbf{y} = \alpha \, \mathbf{z} + \frac {\beta }{\mathbf{z}}, \end{equation}
which maps the unit disk to the ellipse
$\mathcal{A}$
with semi-major and semi-minor axes
$a$
and
$b$
, respectively, with the semi-major axis aligned on the horizontal axis. On the unit disk
$\mathbf{z} = e^{i\theta }$
, we have that
\begin{equation*} \mathbf{y} = \alpha e^{i\theta } + \beta e^{-i\theta } = (\alpha + \beta ) \cos \theta + i(\alpha - \beta )\sin \theta = a \cos \theta + i b \sin \theta . \end{equation*}
The mapping parameters
$\alpha$
and
$\beta$
are then
\begin{equation} \alpha = \frac {a+b}{2}, \qquad \beta = \frac {a-b}{2}. \end{equation}
The Laplace equations to be considered will be solved exterior to the unit disk then mapped to the ellipse by the inverse transformation of (B.78). The large argument behaviour of the inverse transform of (B.78) is calculated as follows for
$|\mathbf{y}|\gg 1$
\begin{equation} \mathbf{z} \sim \frac {\mathbf{y}}{\alpha } - \frac {\beta }{\mathbf{y}} = \frac {\mathbf{y}}{\alpha } - \frac {\beta \bar {\mathbf{y}}}{|\mathbf{y}|^2}, \qquad |\mathbf{z}| = \frac {|\mathbf{y}|}{\alpha }\left [1 - \frac {\alpha \beta }{|\mathbf{y}|^4} \mathbf{y}^T \begin{pmatrix} 1 & \phantom {-}0\\ 0& -1\end{pmatrix} \mathbf{y} + \mathcal{O}(|\mathbf{y}|^{-4}) \right ] , \qquad |\mathbf{y}| \to \infty . \end{equation}
B.1. The order
$\mathcal{O}({\displaystyle \varepsilon }^0)$
problem
The leading order inner problem is given by
\begin{align} \Delta _{\mathbf{y}} v_{0c} = 0\,, \qquad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad v_{0c} = 0, \qquad \mathbf{y} \in \partial \mathcal{A};\nonumber\\[-24pt] \end{align}
\begin{align} v_{0c} = \log |\mathbf{y}| - \log d_c + \frac {\mathbf{d} \cdot \mathbf{y}}{|\mathbf{y}|^2} + \frac {\mathbf{y}^{T} \tilde {\mathcal{Q}} \mathbf{y}}{|\mathbf{y}|^4} + \cdots \, \quad |\mathbf{y}|\to \infty . \end{align}
In classic potential theory, the logarithmic term is the monopole,
$\mathbf{d}$
is the dipole vector and
$\tilde {\mathcal{Q}}$
is the quadrupole matrix. Our goal is to obtain the solution of (B.81) and identify the logarithmic capacitance
$d_c$
and quadrupole matrix
$\tilde {\mathcal{Q}}$
. In the scenario of the unit disk (
$a=b=1$
), the unique solution is
$\log |\mathbf{z}|$
. For the problem (B.81), we apply the mapping (B.80) to obtain the far-field behaviour
\begin{equation} v_{0c}(\mathbf{y}) = \log |\mathbf{z}| \sim \log |\mathbf{y}| - \log \alpha - \frac {\alpha \beta }{|\mathbf{y}|^4} \mathbf{y}^T \begin{pmatrix} 1 & \phantom {-}0\\ 0& -1\end{pmatrix} \mathbf{y} + \mathcal{O}(|\mathbf{y}|^{-4}). \end{equation}
Hence in comparing (B.81b) with (B.82) and applying (B.79), we establish that
\begin{equation} d_c = \log \frac {a+b}{2}, \qquad \mathbf{d} = \mathbf{0}, \qquad \tilde {\mathcal{Q}} = -\alpha \beta \begin{pmatrix} 1 & \phantom {-}0\\ 0& -1\end{pmatrix}= - \frac {a^2 -b^2}{4} \begin{pmatrix} 1 & \phantom {-}0\\ 0& -1\end{pmatrix}. \end{equation}
The vanishing dipole vector
$\mathbf{d} = 0$
is consistent with the two lines of symmetry of the elliptical domain. Extending this analogy, domains with four lines of symmetries (e.g. square) would have a vanishing quadrupole matrix requiring an even higher order of expansion.
B.2. The order
$\mathcal{O}({\displaystyle \varepsilon }^1)$
problem
The first-order inner problem is vector valued and given by
\begin{align} \Delta _{\mathbf{y}} \textbf {v}_{1c} = 0\,, \quad \mathbf{y}\in \mathbb{R}^2\setminus \mathcal{A}; \qquad \textbf {v}_{1c} = 0, \quad \mathbf{y} \in \partial \mathcal{A}; \end{align}
\begin{align} \textbf {v}_{1c} = \mathbf{y} + \frac { \tilde {\mathcal{M}} \mathbf{y}}{|\mathbf{y}|^2} + \cdots \, \quad |\mathbf{y}|\to \infty . \end{align}
For the case of the unit disk, we have that
$\textbf {v}_{1d}= \mathbf{z} - \mathbf{z}/|\mathbf{z}|^2$
. Hence, for the elliptical domain under the transformation (B.80), we have that
\begin{equation} \textbf {v}_{1c}(\mathbf{y}) = \alpha \textbf {v}_{1d}(\mathbf{y}) = \mathbf{y} + \frac {\tilde {\mathcal{M}}\mathbf{y}}{|\mathbf{y}|^2} +\cdots \qquad |\mathbf{y}|\to \infty ; \qquad \tilde {\mathcal{M}} = - \alpha \begin{pmatrix} a &\quad 0\\[5pt] 0&\quad b\end{pmatrix}. \end{equation}
B.3. The order
$\mathcal{O}({\displaystyle \varepsilon }^2)$
problem
In this subsection, we derive the solution of the higher-order correction problem
\begin{align} \Delta _{\mathbf{y}} v_{2c} = -1, \qquad \mathbf{y} \in \mathbb{R}^2\setminus \mathcal{A}; \qquad v_{2c} = 0, \quad \mathbf{y}\in \partial \mathcal{A}; \end{align}
\begin{align} v_{2c} = -\frac {|\mathbf{y}|^2}{4} + \mathbf{y}^{T} \mathcal{B} \mathbf{y} + \cdots \, \quad |\mathbf{y}|\to \infty , \qquad \mathcal{B} = \begin{bmatrix} \mathcal{B}_{11} &\quad\phantom {-}\mathcal{B}_{12} \\[5pt] \mathcal{B}_{12} &\quad -\mathcal{B}_{11}\end{bmatrix}. \end{align}
where
$\mbox{Trace}(\mathcal{B}) = 0$
. The general solution takes the form
$v_{2c} = -\frac {1}{4}|\mathbf{y}|^2 + v_{2h}$
where the homogeneous solution satisfies
\begin{align} \Delta _{\mathbf{y}} v_{2h} = 0, \qquad \mathbf{y} \in \mathbb{R}^2\setminus \mathcal{A}; \qquad v_{2h} = \frac {|\mathbf{y}|^2}{4}, \qquad \mathbf{y}\in \partial \mathcal{A}; \end{align}
\begin{align} v_{2h} \sim \mathbf{y}^{T} \mathcal{B} \mathbf{y} + \cdots \quad |\mathbf{y}|\to \infty . \end{align}
As with previous solutions of inner problems, we solve the corresponding problem on the disk and use the complex transformation (B.78) to map the solution to
$\mathcal{A}$
. Recalling that when
$\mathbf{z}=e^{i\theta }$
, we have
$|\mathbf{y}|^2 = (\alpha ^2 + \beta ^2) + 2\alpha \beta \cos 2\theta$
. Hence, the homogeneous solution is expressed in terms of the disk solution in complex form as
\begin{equation} v_{2h} = a_1 + \mbox{Re}[ b_1\mathbf{z}^2 +b_2 \mathbf{z}^{-2}] + \mbox{Im} [c_1\mathbf{z}^2 +c_2 \mathbf{z}^{-2}]. \end{equation}
On the boundary
$\mathbf{z}=e^{i\theta }$
, we have the conditions
\begin{equation*} \frac 14(\alpha ^2 + \beta ^2) + \frac 12 \alpha \beta \cos 2\theta = a_1 + (b_1 + b_2)\cos 2\theta + (c_1 + c_2)\sin 2\theta , \end{equation*}
which yields the conditions
\begin{equation*} a_1 = \frac 14 (\alpha ^2 + \beta ^2), \qquad b_1 + b_2 = \frac 12\alpha \beta , \qquad c_1+c_2 = 0. \end{equation*}
To establish the behaviour as
$|\mathbf{y}| \to \infty$
, we consider from (B.80) that
\begin{equation*} \mathbf{z}^2 \sim \left (\frac {\mathbf{y}}{\alpha } - \frac {\beta }{\mathbf{y}}\right )^2 = \frac {\mathbf{y}^2}{\alpha ^2} - 2\frac {\beta }{\alpha } + \frac {\beta ^2}{\mathbf{y}^2}, \quad \mbox{as} \quad |\mathbf{y}|\to \infty . \end{equation*}
The large argument behaviour of
$v_{2h}$
given in (B.88) is
\begin{align*} v_{2h} &\sim \frac {b_1}{\alpha ^2} \mbox{Re}[\mathbf{y}^2] + \frac {c_1}{\alpha ^2} \mbox{Im}[\mathbf{y}^2] + \frac 14 (\alpha ^2 + \beta ^2) - 2\frac {\beta }{\alpha } b_1 + \mathcal{O}(|\mathbf{y}|^{-2}). \end{align*}
Now, comparing with (B.87b), we see that
\begin{equation} b_1 = \alpha ^2 \mathcal{B}_{11}, \qquad c_1 = \alpha ^2 \mathcal{B}_{12}. \end{equation}
Hence, we obtain the large argument behaviour of (B.87) to be
\begin{align} v_{2h} &\sim \mathbf{y}^T \mathcal{B}\mathbf{y} + \frac {1}{4}(\alpha ^2 + \beta ^2) - 2\alpha \beta \, \mathcal{B}_{11} + + \mathcal{O}(|\mathbf{y}|^{-2}), \end{align}
and finally establish the large argument behaviour of (B.86) to be
\begin{equation} v_{2c} = -\frac {|\mathbf{y}|^2}{4} + \mathbf{y}^{T} \mathcal{B} \mathbf{y} + d_{2c} + \mathcal{O}(|\mathbf{y}|^{-2}), \quad \mbox{as} \quad |\mathbf{y}|\to \infty ; \qquad d_{2c} = \frac {1}{4}( \alpha ^2 + \beta ^2) - 2\alpha \beta \, \mathcal{B}_{11} . \end{equation}
Open access
