$\displaystyle+\widetilde{\bm{Q}}_{-n_{0_{1}}+1+2\kappa_{1}N_{1},...,-n_{0_{d}}+1+2\kappa_{d}N_{d}}(n_{1},...,n_{d},z)\bigg{\\}},$

we obtain, from Eq. (3), the reflective propagator

$\displaystyle\widetilde{\bm{P}}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},z)=$
(S28)
$\displaystyle\frac{1}{2\pi}\int_{-\pi}^{\pi}...\int_{-\pi}^{\pi}\sum_{\kappa_{1}=-\infty}^{\infty}...\sum_{\kappa_{d}=-\infty}^{\infty}\left(\left[e^{-i(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\bm{\xi}}+e^{-i(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\bm{\xi}}\right]e^{i\bm{\kappa}^{\intercal}\cdot
2\bm{N}\bm{\xi}}\right)\left[\mathbb{I}-z\bm{\lambda}\left(\xi_{1},...,\xi_{d}\right)\right]^{-1}\cdot\bm{U_{m_{0}}}\text{d}^{d}\xi.$

Using Eq. (S20) with $\mu_{d}=\frac{\pi}{N_{d}}$, one obtains

$\displaystyle\widetilde{\bm{P}}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},z)=$
(S29)
$\displaystyle\frac{1}{2\pi}\int_{-\pi}^{\pi}...\int_{-\pi}^{\pi}\sum_{\kappa_{1}=-\infty}^{\infty}...\sum_{\kappa_{1}=-\infty}^{\infty}\bigg{(}\delta\left(\bm{\xi}-\frac{\bm{\kappa}}{\bm{N}}\pi\right)\left[e^{-i(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\bm{\xi}}+e^{-(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\bm{\xi}}\right]\bigg{)}\left[\mathbb{I}-z\bm{\lambda}\left(\xi_{1},...,\xi_{d}\right)\right]^{-1}\cdot\bm{U_{m_{0}}}\text{d}\xi,$

where
$\frac{\bm{k}}{\bm{N}}=\left[\frac{\kappa_{1}}{N_{1}},...,\frac{\kappa_{d}}{N_{d}}\right]$
which is directly solvable via the same procedure as Eq. (S22) leading to

$\displaystyle\widetilde{\bm{P}}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},z)=$
(S30)
$\displaystyle\frac{1}{(2N)^{d}}\sum_{\kappa_{1}=0}^{2N_{1}-1}...\sum_{\kappa_{d}=0}^{2N_{d}-1}\left(\big{[}e^{-i\pi(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}+e^{-i\pi(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}\big{]}\right)\left[\mathbb{I}-z\bm{\lambda}\left(\frac{\pi\kappa_{1}}{N_{1}},...,\frac{\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\cdot\bm{U_{m_{0}}},$

which, after some $z$-inverting, utilising that all the functions are periodic
to shift the summation and splitting the summation into parts and isolating
the $(\kappa_{1},...,\kappa_{d})=(0,...,0)$-th term, allows passage to

$\displaystyle\bm{P}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},t)=\frac{1}{(2N)^{d}}\underset{\kappa_{1}+\kappa_{2}+...\kappa_{d}>0}{\sum_{\kappa_{1}=0}^{N_{1}-1}...\sum_{\kappa_{d}=0}^{N_{d}-1}}\bigg{\\{}\left(\big{[}e^{-i\pi(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}+e^{-i\pi(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}\big{]}\right)\bm{\lambda}\left(\frac{\pi\kappa_{1}}{N_{1}},...,\frac{\pi\kappa_{d}}{N_{d}}\right)^{t}\cdot\bm{U_{m_{0}}}$
(S31)
$\displaystyle+\left(\big{[}e^{i\pi(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}+e^{i\pi(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}\big{]}\right)\bm{\lambda}\left(\frac{-\pi\kappa_{1}}{N_{1}},...,\frac{-\pi\kappa_{d}}{N_{d}}\right)^{t}\cdot\bm{U_{m_{0}}}\bigg{\\}}+\frac{2}{(2N)^{d}}\bm{\lambda}(0,...,0)^{t}\cdot\bm{U_{m_{0}}}.$

To proceed we note that while
$\bm{\lambda}\left(\frac{\pi\kappa_{1}}{N_{1}},...,\frac{\pi\kappa_{d}}{N_{d}}\right)$
and
$\bm{\lambda}\left(\frac{-\pi\kappa_{1}}{N_{1}},...,\frac{-\pi\kappa_{d}}{N_{d}}\right)$
govern different state level dynamics, the site level dynamics are equivalent
as the change of sign in the matrix exponent merely alters the state that is
entered on each movement direction. Since the propagation is symmetric along
all axes (ensured by the uniform initial weighting), then the above reasoning
allows one to realise that
$\sum_{j=1}^{2d}\left[\bm{\lambda}\left(\frac{\pi\kappa_{1}}{N_{1}},...,\frac{\pi\kappa_{d}}{N_{d}}\right)\cdot\bm{U_{m_{0}}}\right]_{j}=\sum_{j=1}^{2d}\left[\bm{\lambda}\left(\frac{-\pi\kappa_{1}}{N_{1}},...,\frac{-\pi\kappa_{d}}{N_{d}}\right)\cdot\bm{U_{m_{0}}}\right]_{j}$.
Hence, upon summation over the states to find the site level propagator

$P^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},t)=\sum_{m=1}^{M}\bm{P}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},t)_{m},$
(S32)

one can factor out the matricial terms in the first two lines of Eq. (S31)
leading directly to

$\displaystyle\sum_{m=1}^{M}\bm{P}^{(r)}_{n_{0_{1}},...,n_{0_{d}}}(n_{1},...,n_{d},t)_{m}=\sum_{m=1}^{M}\bigg{[}$
$\displaystyle\frac{1}{(2N)^{d}}\underset{\kappa_{1}+\kappa_{2}+...\kappa_{d}>0}{\sum_{\kappa_{1}=0}^{N_{1}-1}...\sum_{\kappa_{d}=0}^{N_{d}-1}}\bigg{\\{}\bigg{[}\big{[}e^{-i\pi(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}+e^{-i\pi(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}\big{]}$
(S33) $\displaystyle+$
$\displaystyle\big{[}e^{i\pi(\bm{n}-\bm{n_{0}})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}+e^{i\pi(\bm{n}+\bm{n_{0}}-\bm{1})^{\intercal}\cdot\frac{\bm{\kappa}}{\bm{N}}}]\bigg{]}\bm{\lambda}\left(\frac{\pi\kappa_{1}}{N_{1}},...,\frac{\pi\kappa_{d}}{N_{d}}\right)^{t}\cdot\bm{U_{m_{0}}}\bigg{\\}}+$
$\displaystyle\frac{2}{(2N)^{d}}\bm{\lambda}(0,...,0)^{t}\cdot\bm{U_{m_{0}}}\bigg{]}_{m},$

which via some trigonometric manipulation of the scalar terms in Eq. (S33)
allows for direct passage to Eq. (14) with $\gamma=r$ of the main text.

### The Squeezing Method

We detail here, via a schematic, the method of transforming a periodic domain
into one with reflecting boundaries. The method, first introduced in ref. [45]
for the standard LRW case can be applied also to the CRW. The essence of the
method can be evinced from Fig. S1. As shown, if the walker occupies site 1 of
the periodic lattice, to move off the site its only options are sites 2 and 6.
Thus, if we combine sites 2 and 6 together, that is we squeeze the lattice, to
create a new site (labelled $2+6$ in the far right panel), we have precisely
the reflecting boundary condition considered by Chandrasekhar [50].

Figure S1: (Colour Online). In the style of Fig. 1b in ref. [45], a schematic
depicting the method to turn a periodic domain into a reflecting domain by
squeezing the internal sites, while the boundary points (sites 1 and 4) remain
unaffected by the transformation. For visual clarity we omit arrows depicting
the state level dynamics, that is the directions from which a walker may enter
a site, depicted in Fig. 1 of the main text. As shown in the right-hand panel,
it is clear that after the transformation the walker may enter both internal
states from either direction.

## The propagator in a two-dimensional Square lattice

The generating function of Eq. (14) when $d=2$ is given by

$\widetilde{\bm{P}}^{(\gamma)}_{n_{01},n_{02}}(n_{1},n_{2},t)=\frac{1}{N^{2}}\sum_{\kappa_{1}=0}^{N-1}\sum_{\kappa_{2}=0}^{N-1}\left[\prod_{i=1}^{2}g_{\kappa_{i}}^{(\gamma)}(n_{i},n_{0_{i}})\right]\left[\mathbb{I}-z\bm{\lambda}\left(\pi\mathcal{N}^{(\gamma)}_{\kappa_{1}},\pi\mathcal{N}^{(\gamma)}_{\kappa_{2}}\right)\right]^{-1}\cdot\bm{U}_{\bm{m_{0}}}.$
(S34)

In general, one can invert the matrix in Eq. (S34), however with arbitrary
parameters it is a fourth order polynomial, and writing the roots explicitly
may not be convenient. There are, however, two non-trivial cases for which the
matrix inversion is greatly simplified. They occur when $c^{(2)}=0$ and with a
uniform localised initial condition, that is
$\bm{U}_{\bm{m}_{0}}=[1/4,1/4,1/4,1/4]^{\intercal}$, and we present them
below.

### The case with $b=\ell$

By inverting the matrix in Eq. (S34), one finds

$\displaystyle\widetilde{P}^{(\gamma)}_{n_{0_{1}},n_{0_{2}}}(n_{1},n_{2},z)=\frac{1}{2N^{2}}\sum_{\kappa_{1}=0}^{N-1}\sum_{\kappa_{2}=0}^{N-1}\left[\prod_{i=1}^{2}g_{\kappa_{i}}^{(\gamma)}(n_{i},n_{0_{i}})\right]$
(S35)
$\displaystyle\times\frac{2+3z(4\ell-1)\psi(\kappa_{1},\kappa_{2})+2z^{2}(1-4\ell)^{2}\sigma(\kappa_{1},\kappa_{2})-z^{3}(1-4\ell)^{3}\psi(\kappa_{1},\kappa_{2})}{1-2z(1-3\ell)\psi(\kappa_{1},\kappa_{2})+2z^{2}(1-6\ell+8\ell^{2})\sigma(\kappa_{1},\kappa_{2})+2z^{3}(1-4\ell)^{2}(\ell-1)\psi(\kappa_{1},\kappa_{2})-z^{4}(4\ell-1)^{3}},$

where
$\psi(\kappa_{1},\kappa_{2})=\cos\left(\pi\mathcal{N}_{\kappa_{1}}^{(\gamma)}\right)+\cos\left(\pi\mathcal{N}_{\kappa_{2}}^{(\gamma)}\right)$,
and
$\sigma(\kappa_{1},\kappa_{2})=1+2\cos\left(\pi\mathcal{N}_{\kappa_{1}}^{(\gamma)}\right)\cos\left(\pi\mathcal{N}_{\kappa_{2}}^{(\gamma)}\right)$.

### The case with $f=\ell$

To study this case, we make use of the theory of matrix resolvent [104], for
which

$\det\left[z\mathbb{I}-\bm{\lambda}\left(\pi\mathcal{N}^{(\gamma)}_{\kappa_{1}},\pi\mathcal{N}^{(\gamma)}_{\kappa_{2}}\right)\right]=\big{[}(1-4\ell)^{2}-z^{2}\big{]}\Big{\\{}1-4\ell-z^{2}+2\ell
z\left[\cos\left(\pi\mathcal{N}_{\kappa_{1}}^{(\gamma)}\right)+\cos\left(\pi\mathcal{N}_{\kappa_{2}}^{(\gamma)}\right)\right]\Big{\\}},$
(S36)

to find the eigenvalues of
$\bm{\lambda}\left(\pi\mathcal{N}^{(\gamma)}_{\kappa_{1}},\pi\mathcal{N}^{(\gamma)}_{\kappa_{2}}\right)$
as

$\lambda_{1}\\!=\\!1\\!-\\!4\ell,\quad\lambda_{2}\\!=\\!-1\\!+\\!4\ell,\quad\lambda_{\\{3,4\\}}\\!=\\!\ell\left[\cos\left(\pi\mathcal{N}_{\kappa_{1}}^{(\gamma)}\right)\\!+\\!\cos\left(\pi\mathcal{N}_{\kappa_{2}}^{(\gamma)}\right)\right]\pm\sqrt{1\\!-\\!4\ell\\!+\\!\ell^{2}\left[\cos\left(\pi\mathcal{N}_{\kappa_{1}}^{(\gamma)}\right)\\!+\\!\cos\left(\pi\mathcal{N}_{\kappa_{2}}^{(\gamma)}\right)\right]^{2}}.$
(S37)

By finding the corresponding eigenvectors, the matrix power in Eq. (S34) can
be computed exactly and after some laborious algebra, one arrives at Eq. (15)
of the main text.

## Derivation of Correlated Hexagonal Propagator

With the sojourn probability, $c^{(4)}=1-(f+b+\sum_{i=1}^{4}\ell_{i})$, we
write the Master equation as

$\displaystyle Q($ $\displaystyle
n_{1},n_{2},1,t+1)=\bigg{[}fQ(n_{1}+1,n_{2},1,t)+bQ(n_{1}+1,n_{2},2,t)+\ell_{4}Q(n_{1}+1,n_{2},3,t)+$
(S38)
$\displaystyle\ell_{2}Q(n_{1}+1,n_{2},4,t)+\ell_{3}Q(n_{1}+1,n_{2},5,t)+\ell_{1}Q(n_{1}+1,n_{2},6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},1,t),$
$\displaystyle Q($ $\displaystyle
n_{1},n_{2},2,t+1)=\bigg{[}bQ(n_{1}-1,n_{2},1,t)+fQ(n_{1}-1,n_{2},2,t)+\ell_{2}Q(n_{1}-1,n_{2},3,t)+$
(S39)
$\displaystyle\ell_{4}Q(n_{1}-1,n_{2},4,t)+\ell_{1}Q(n_{1}-1,n_{2},5,t)+\ell_{3}Q(n_{1}-1,n_{2},6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},2,t),$
$\displaystyle Q($ $\displaystyle
n_{1},n_{2},3,t+1)=\bigg{[}\ell_{1}Q(n_{1},n_{2}+1,1,t)+\ell_{3}Q(n_{1},n_{2}+1,2,t)+fQ(n_{1},n_{2}+1,3,t)+$
(S40) $\displaystyle
bQ(n_{1},n_{2}+1,4,t)+\ell_{4}Q(n_{1},n_{2}+1,5,t)+\ell_{2}Q(n_{1},n_{2}+1,6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},3,t),$
$\displaystyle Q($ $\displaystyle
n_{1},n_{2},4,t+1)=\bigg{[}\ell_{3}Q(n_{1},n_{2}-1,1,t)+\ell_{1}Q(n_{1},n_{2}-1,2,t)+bQ(n_{1},n_{2}-1,3,t)+$
(S41) $\displaystyle
fQ(n_{1},n_{2}-1,4,t)+\ell_{2}Q(n_{1},n_{2}-1,5,t)+\ell_{4}Q(n_{1},n_{2}-1,6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},4,t),$
$\displaystyle Q($ $\displaystyle
n_{1},n_{2},5,t+1)=\bigg{[}\ell_{2}Q(n_{1}-1,n_{2}+1,1,t)+\ell_{4}Q(n_{1}-1,n_{2}+1,2,t)+\ell_{1}Q(n_{1}-1,n_{2}+1,3,t)+$
(S42)
$\displaystyle\ell_{3}Q(n_{1}-1,n_{2}+1,4,t)+fQ(n_{1}-1,n_{2}+1,5,t)+bQ(n_{1}-1,n_{2}+1,6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},5,t),$
$\displaystyle Q($ $\displaystyle
n_{1},n_{2},6,t+1)=\bigg{[}\ell_{4}Q(n_{1}+1,n_{2}-1,1,t)+\ell_{2}Q(n_{1}+1,n_{2}-1,2,t)+\ell_{3}Q(n_{1}+1,n_{2}-1,3,t)+$
(S43)
$\displaystyle\ell_{1}Q(n_{1}+1,n_{2}-1,4,t)+bQ(n_{1}+1,n_{2}-1,5,t)+fQ(n_{1}+1,n_{2}-1,6,t)\bigg{]}+c^{(4)}Q(n_{1},n_{2},6,t).$

With the initial condition
$\bm{U}_{m_{0}}=\delta_{n_{1},n_{0_{1}}}\delta_{n_{2},n_{0_{2}}}[\alpha_{1},...,\alpha_{6}]^{\intercal}$,
and following the techniques in Sec. Derivation of the Unbounded Propagator, ,
it is possible to write the lattice Green’s function as

$\widetilde{\bm{Q}}^{(\mathcal{H})}_{n_{0}}(\bm{n},z)=\frac{1}{4\pi^{2}}\int_{-\pi}^{\pi}\int_{-\pi}^{\pi}e^{-i\bm{\xi}(\bm{n}-\bm{n}_{0})}\bm{\lambda}^{(\mathcal{H})}(\xi_{1},\xi_{2})\cdot\bm{U_{m_{0}}}\text{d}\xi_{1}\text{d}\xi_{2}$
(S44)

where

$\displaystyle\bm{\lambda}^{(\mathcal{H})}(\xi_{1},\xi_{2})=$ (S45)
$\displaystyle\hbox{}\;\vbox{\kern 63.55554pt\hbox{$\kern
319.68983pt\kern-4.75pt\left[\kern-319.68983pt\vbox{\kern-63.55554pt\vbox{\halign{$#$\hfil\kern
2\p<EMAIL_ADDRESS>2.0pt\kern
4.75pt&\nwarrow&\searrow&\nearrow&\swarrow&\longrightarrow&\longleftarrow\crcr\kern
2.0pt\cr\nwarrow$\hfil\kern 2.0pt\kern
4.75pt&fe^{-i\xi_{1}}+c^{(4)}&be^{-i\xi_{1}}&\ell_{4}e^{-i\xi_{1}}&\ell_{2}e^{-i\xi_{1}}&\ell_{3}e^{-i\xi_{1}}&\ell_{1}e^{-i\xi_{1}}\cr\searrow$\hfil\kern
2.0pt\kern
4.75pt&be^{i\xi_{1}}&fe^{i\xi_{1}}+c^{(4)}&\ell_{2}e^{i\xi_{1}}&\ell_{4}e^{i\xi_{1}}&\ell_{1}e^{i\xi_{1}}&\ell_{3}e^{i\xi_{1}}\cr\nearrow$\hfil\kern
2.0pt\kern
4.75pt&\ell_{1}e^{-i\xi_{2}}&\ell_{3}e^{-i\xi_{2}}&fe^{-i\xi_{2}}+c^{(4)}&be^{-i\xi_{2}}&\ell_{4}e^{-i\xi_{2}}&\ell_{2}e^{-i\xi_{2}}\cr\swarrow$\hfil\kern
2.0pt\kern
4.75pt&\ell_{3}e^{i\xi_{2}}&\ell_{1}e^{i\xi_{2}}&be^{i\xi_{2}}&fe^{i\xi_{2}}+c^{(4)}&\ell_{2}e^{i\xi_{2}}&\ell_{4}e^{i\xi_{2}}\cr\longrightarrow$\hfil\kern
2.0pt\kern
4.75pt&\ell_{2}e^{i(\xi_{1}-\xi_{2})}&\ell_{4}e^{i(\xi_{1}-\xi_{2})}&\ell_{1}e^{i(\xi_{1}-\xi_{2})}&\ell_{3}e^{i(\xi_{1}-\xi_{2})}&fe^{i(\xi_{1}-\xi_{2})}+c^{(4)}&be^{i(\xi_{1}-\xi_{2})}\cr\longleftarrow$\hfil\kern
2.0pt\kern
4.75pt&\ell_{4}e^{-i(\xi_{1}-\xi_{2})}&\ell_{2}e^{-i(\xi_{1}-\xi_{2})}&\ell_{3}e^{-i(\xi_{1}-\xi_{2})}&\ell_{1}e^{-i(\xi_{1}-\xi_{2})}&be^{-i(\xi_{1}-\xi_{2})}&fe^{-i(\xi_{1}-\xi_{2})}+c^{(4)}\crcr\cr}}\kern-12.0pt}\,\right]$}},$

with the elements representing the different movement options are labeled, in
Eq. (S45), via the arrows.

Figure S2: (Colour Online). A schematic representing two of the six
permissible movement steps in the correlated hexagonal random walk where the
solid arrows represent the movement probabilities if the walker had entered
the site via the dotted arrow. The movements shown on the left side case are
shown in the $3^{\text{rd}}$ column of Eq. (S45) and the movements shows in
the right panel are encoded in the $2^{\text{nd}}$ column of Eq. (S45).

### The Hexagonal Periodic Boundary Condition

A valid image set for the standard (diffusive) periodic hexagonal propagator
is known as [55]

$\widetilde{\bm{P}}^{(\mathcal{H})}_{\bm{n}_{0}}(n_{1},n_{2},z)=\sum_{\kappa_{1}=-\infty}^{\infty}\sum_{\kappa_{2}=-\infty}^{\infty}\widetilde{\bm{Q}}_{\bm{n}_{0}}^{(\mathcal{H})}(n_{1}+\hat{n}_{1},n_{2}+\hat{n}_{2},z),$
(S46)

where

$\begin{bmatrix}\widehat{n}_{1}\\\
\widehat{n}_{2}\end{bmatrix}=\begin{bmatrix}-R\kappa_{1}+(2R+1)\kappa_{2}\\\
-(R+1)\kappa_{1}-R\kappa_{2}\end{bmatrix},$ (S47)

and since again, the correlated periodic boundary conditions are identical to
their diffusive counterpart, Eq. (S44) can be solved explicitly following the
same procedure as the solution to the LRW with periodic boundary conditions
derived in [55]. This leads to Eq. (20) of the main text.

As mentioned in Sec. IV.2 of the main text, the non-orthogonality of the
coordinates leads to a ‘shift’ in the position of the walker when traversing
the periodic boundary. For a deeper discussion we again refer the interested
reader to [55], and provide here a simple explanation of the dynamics seen in
Fig. 3 of the main text. We describe here the path taken by the probability
pulse that is near North-East of the domain in the ballistic case Fig. 3(a)
(other pulses may then be understood similarly due to the rotational symmetry
of the system).

With initial condition at the origin, the pulse in the North-East first
travelled due West and reached the West most corner of the hexagon, at $t=8$.
Crossing the boundary at $t=9$ the hexagonal periodicity brings the walker
back into the domain at the North-East corner. Since the pulse is
ballistically travelling due West, it will remain travelling on this
trajectory along the top boundary until $t=17$, passing the point
$\bm{n}=(-6,-2,8)$ at time $t=11$ as shown in Fig. 3(a).

As mentioned in Sec IV.2, this shift leads allows the persistent walker to
visit every site in the domain. We show this formally by considering the
following [42]. First taking the $z$-transfrom of Eq. (20) such that

$\widetilde{\bm{P}}^{(\mathcal{H})}_{\bm{n}_{0}}(\bm{n},z)=\frac{1}{\Omega}\left[\mathbb{I}-z\bm{\lambda}^{(\mathcal{H})}(0,0)\right]^{-1}+\widetilde{\bm{P}}^{(\kappa\neq
0)}_{\bm{n}_{0}}(\bm{n},z),$ (S48)

where $\widetilde{P}^{(\kappa\neq 0)}_{\bm{n}_{0}}(\bm{n},z)$ is the
generating function of the second term in Eq. (20). Since
$\bm{\lambda}^{(\mathcal{H})}(0,0)=\mathbb{I}$ for
$(\kappa_{1},\kappa_{2})=(0,0)$ only, it is clear that at $z=1$, each element
$\widetilde{\bm{P}}^{(\kappa\neq 0)}_{\bm{n}_{0}}(\bm{n},z)$ is holomorphic.
Therefore,
$\widetilde{\bm{P}}^{(\mathcal{H})}_{\bm{n}_{0}}(\bm{n},z)\to\infty$ as $z\to
1$ for any $\bm{n}$ meaning it is certain all sites will be hit, irrespective
of the movement parameters.

## Derivations of Specific First Passage Cases

### One State Case

#### Survival Arguments

With a single target we have a set $\mathcal{S}$ of targets with cardinality
one, namely $\mathdutchcal{S}=\left\\{(\bm{s},m_{\bm{s}})\right\\}$. Starting
from Eq. (33) of Appendix C in the main text, we find

$\displaystyle\widetilde{P}^{(a)}_{\bm{n}_{0},m_{0}}(\bm{s},m,z)=\alpha_{m_{0}}\widetilde{P}_{\bm{n}_{0},m_{0}}(\bm{s},m,z)+\frac{\rho}{\rho-1}\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m,z)\widetilde{P}^{(a)}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z).$
(S49)

To proceed we set here $(\bm{s},m)=(\bm{s},m_{\bm{s}})$, solve for the unknown
function $\widetilde{P}^{(a)}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z)$ in
terms of the known functions,

$\widetilde{P}^{(a)}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z)=\frac{\alpha_{m_{0}}\widetilde{P}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z)}{1+\frac{\rho}{\rho-1}\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)}$
(S50)

and then substitute Eq. (S50) into Eq. (S49), which leads to

$\widetilde{P}^{(a)}_{\bm{n}_{0},m_{0}}(\bm{s},m,z)=\alpha_{m_{0}}\widetilde{P}_{\bm{n}_{0},m_{0}}(\bm{s},m,z)-\rho\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m,z)\frac{\alpha_{m_{0}}\widetilde{P}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z)}{1-\rho+\rho\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)}.$
(S51)

Performing a summation over $m$ and $m_{0}$, we find the occupation
probability with one partially absorbing trap as

$\widetilde{P}^{(a)}_{\bm{n}_{0}}(\bm{s},z)=\widetilde{P}_{\bm{n}_{0}}(\bm{s},z)-\rho\alpha_{m_{\bm{s}}}\sum_{m=1}^{M}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m,z)\frac{\sum_{j=1}^{M}\alpha_{m_{0_{j}}}\widetilde{P}_{\bm{n}_{0},m_{0_{j}}(\bm{s},m_{\bm{s}},z)}}{1-\rho+\rho\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)}.$
(S52)

Following the same procedure outlined in Sec. V to pass from Eq. (21) to Eq.
(22) of the main text we find the first-absorption probability with one
absorbing state as

$\widetilde{F}_{\bm{n}_{0}}(\bm{s},m_{\bm{s}},z)=\rho\alpha_{m_{\bm{s}}}\frac{\sum_{j=1}^{M}\alpha_{m_{0_{j}}}\widetilde{P}_{\bm{n}_{0},m_{0_{j}}}(\bm{s},m_{\bm{s}},z)}{1-\rho+\rho\alpha_{m_{\bm{s}}}\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)},$
(S53)

which, in the $\rho\to 1$ limit, gives Eq. (24) of the main text.

#### Renewal Arguments

We start with a set of two coupled equations, which relate the probability of
being found in the specific target state and the corresponding first passages

$\displaystyle\alpha_{m_{0_{1}}}P_{\bm{n}_{0},1}(\bm{s},m_{\bm{s}},t)$
$\displaystyle=\alpha_{m_{0_{1}}}\delta_{m_{0},1}\delta_{t,0}\delta_{\bm{s},\bm{n}_{0}}+\sum_{t^{\prime}=0}^{t}F_{\bm{n}_{0},1}(\bm{s},m_{\bm{s}},t^{\prime})P_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},t-t^{\prime}),$
(S54) $\displaystyle\alpha_{m_{0_{2}}}P_{\bm{n}_{0},2}(\bm{s},m_{\bm{s}},t)$
$\displaystyle=\alpha_{m_{0_{2}}}\delta_{m_{0},2}\delta_{t,0}\delta_{\bm{s},\bm{n}_{0}}+\sum_{t^{\prime}=0}^{t}F_{\bm{n}_{0},2}(\bm{s},m_{\bm{s}},t^{\prime})P_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},t-t^{\prime}),$
$\displaystyle\;\;\vdots$
$\displaystyle\alpha_{m_{0_{M}}}P_{\bm{n}_{0},M}(\bm{s},m_{\bm{s}},t)$
$\displaystyle=\alpha_{m_{0_{M}}}\delta_{m_{0},M}\delta_{t,0}\delta_{\bm{s},\bm{n}_{0}}+\sum_{t^{\prime}=0}^{t}F_{\bm{n}_{0},M}(\bm{s},m_{\bm{s}},t^{\prime})P_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},t-t^{\prime}).$

Assuming $\bm{n}\neq\bm{n}_{0}$ and taking the generating function we find

$\displaystyle\alpha_{m_{0_{1}}}\widetilde{P}_{\bm{n}_{0},1}(\bm{s},m_{\bm{s}},z)$
$\displaystyle=\widetilde{F}_{\bm{n}_{0},1}(\bm{s},m_{\bm{s}},z)\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z),$
(S55)
$\displaystyle\alpha_{m_{0_{2}}}\widetilde{P}_{\bm{n}_{0},2}(\bm{s},m_{\bm{s}},z)$
$\displaystyle=\widetilde{F}_{\bm{n}_{0},2}(\bm{s},m_{\bm{s}},z)\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z),$
$\displaystyle\;\;\vdots$
$\displaystyle\alpha_{m_{0_{M}}}\widetilde{P}_{\bm{n}_{0},M}(\bm{s},m_{\bm{s}},z)$
$\displaystyle=\widetilde{F}_{\bm{n}_{0},M}(\bm{s},m_{\bm{s}},z)\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z).$

For each equation in Eq. (S55) one may solve for

$\widetilde{F}_{\bm{n}_{0},i}(\bm{s},m_{\bm{s}},z)=\frac{\alpha_{m_{0_{i}}}\widetilde{P}_{\bm{n}_{0},i}(\bm{s},m_{\bm{s}},z)}{\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)},$
(S56)

independently of the other equations and upon summation
$\widetilde{F}_{\bm{n}_{0}}(\bm{s},m_{\bm{s}},z)=\sum_{i=1}^{M}\widetilde{F}_{\bm{n}_{0},i}(\bm{s},m_{\bm{s}},z)$,
we find Eq. (24) of the main text. Evidently one also sees here that if
$\alpha_{m_{0_{i}}}=1$ then one recovers the known cases of ref. [35] and ref.
[34].

### Two State Case

#### Survival Arguments

As explained in Sec. V of the main text, the first-passage probability can be
found from (22) as

$\widetilde{F}_{n_{0}}(\mathdutchcal{S},z)=\widetilde{F}_{\bm{n}_{0}}(\mathdutchcal{S},z)=\sum_{j=1}^{M}\sum_{i=1}^{S}\alpha_{m_{\bm{s}_{i}}}\frac{\det[\mathbb{H}^{(i)}(\bm{n}_{0},m_{0_{j}},z)]}{\det[\mathbb{H}(z)]}.$
(S57)

For the cases with low numbers of targets, Eq. (S57) can be expanded
explicitly. An important case is when there are two targets, which is
pertinent to the one-dimensional case with a target in both states of $s$.
With two targets, one finds

$\widetilde{F}_{n_{0}}(\mathdutchcal{S},z)\\!=\\!\frac{\alpha_{m_{\bm{s}_{1}}}\det[\mathbb{H}^{(1)}(\bm{n}_{0},1,z)]\\!+\\!\alpha_{m_{\bm{s}_{2}}}\det[\mathbb{H}^{(2)}(\bm{n}_{0},1,z)]}{\det[\mathbb{H}(z)]}\\!+\\!\frac{\alpha_{m_{\bm{s}_{1}}}\det[\mathbb{H}^{(1)}(\bm{n}_{0},2,z)]\\!+\\!\alpha_{m_{\bm{s}_{2}}}\det[\mathbb{H}^{(2)}(\bm{n}_{0},2,z)]}{\det[\mathbb{H}(z)]}.$
(S58)

where
$\mathbb{H}(z)_{l,k}=\alpha_{m_{\bm{s}_{k}}}\widetilde{P}_{\bm{s}_{k},m_{\bm{s}_{k}}}(\bm{s}_{l},m_{\bm{s}_{l}},z)$,
$\mathbb{H}(z)_{k,k}=\alpha_{m_{\bm{s}_{k}}}\widetilde{P}_{\bm{s}_{k},m_{\bm{s}_{k}}}(\bm{s}_{k},m_{\bm{s}_{k}},z)$
and $\mathbb{H}^{(i)}(\bm{n}_{0},j,z)$ is the same, but with the
$i^{\text{th}}$ column replaced with
$\alpha_{m_{0_{j}}}\left[\widetilde{P}_{\bm{n}_{0},m_{0_{j}}}(\bm{s}_{1},m_{\bm{s}_{1}},z),...,\widetilde{P}_{\bm{n}_{0},m_{0_{j}}}(\bm{s}_{S},m_{\bm{s}_{S}},z)\right]^{\intercal}$.
One may write this more explicitly as

$\displaystyle\widetilde{F}_{n_{0}}(\mathdutchcal{S},z)=$ (S59)
$\displaystyle\alpha_{m_{\bm{s}_{1}}}\frac{\begin{vmatrix}\alpha_{m_{0_{1}}}\widetilde{P}_{\bm{n}_{0},m_{0_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{0_{1}}}\widetilde{P}_{\bm{n}_{0},m_{0_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}+\begin{vmatrix}\alpha_{m_{0_{2}}}\widetilde{P}_{\bm{n}_{0},m_{0_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{0_{2}}}\widetilde{P}_{\bm{n}_{0},m_{0_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}}{\begin{vmatrix}\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}}$
$\displaystyle+\alpha_{m_{\bm{s}_{2}}}\frac{\begin{vmatrix}\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{0_{1}}}\widetilde{P}_{\bm{n}_{0},m_{0_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{0_{1}}}\widetilde{P}_{\bm{n}_{0},m_{0_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}\\!+\\!\begin{vmatrix}\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{0_{2}}}\widetilde{P}_{\bm{n}_{0},m_{0_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{0_{2}}}\widetilde{P}_{\bm{n}_{0},m_{0_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}}{\begin{vmatrix}\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{1},m_{\bm{s}_{1}},z)\\\
\alpha_{m_{\bm{s}_{1}}}\widetilde{P}_{\bm{s}_{1},m_{\bm{s}_{1}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)&\alpha_{m_{\bm{s}_{2}}}\widetilde{P}_{\bm{s}_{2},m_{\bm{s}_{2}}}(\bm{s}_{2},m_{\bm{s}_{2}},z)\\\
\end{vmatrix}}.$

Performing the determinants and some trivial simplification leads to Eq. (25)
of the main text.

#### Renewal Arguments

We generalise here the approach taken in ref. [34] to the case of two target
sites. We follow the same approach as for Eq. (S54) but now include two target
sites. We note that for simplicity in the algebra we assume without loss of
generality that the system is one-dimensional with two internal states
$m=\\{1,2\\}$. The equations relating the state specific occupation
probability and first passage are now coupled, namely:

$\displaystyle\alpha_{m_{0_{1}}}P_{n_{0},1}(s,1,t)$
$\displaystyle=\alpha_{m_{0_{1}}}\delta_{m_{0},1}\delta_{t,0}\delta_{n,n_{0}}+\sum_{t^{\prime}=0}^{t}F_{n_{0},1}(s,1,t^{\prime})P_{s,1}(s,1,t-t^{\prime})+F_{n_{0},1}(s,2,t^{\prime})P_{s,2}(s,1,t-t^{\prime}),$
(S60) $\displaystyle\alpha_{m_{0_{2}}}P_{n_{0},2}(s,1,t)$
$\displaystyle=\alpha_{m_{0_{2}}}\delta_{m_{0},1}\delta_{t,0}\delta_{n,n_{0}}+\sum_{t^{\prime}=0}^{t}F_{n_{0},2}(s,1,t^{\prime})P_{s,1}(s,1,t-t^{\prime})+F_{n_{0},2}(s,2,t^{\prime})P_{s,2}(s,1,t-t^{\prime}),$
$\displaystyle\alpha_{m_{0_{1}}}P_{n_{0},1}(s,2,t)$
$\displaystyle=\alpha_{m_{0_{1}}}\delta_{m_{0},1}\delta_{t,0}\delta_{n,n_{0}}+\sum_{t^{\prime}=0}^{t}F_{n_{0},1}(s,1,t^{\prime})P_{s,1}(s,2,t-t^{\prime})+F_{n_{0},1}(s,2,t^{\prime})P_{s,2}(s,2,t-t^{\prime}),$
$\displaystyle\alpha_{m_{0_{2}}}P_{n_{0},2}(s,2,t)$
$\displaystyle=\alpha_{m_{0_{2}}}\delta_{m_{0},1}\delta_{t,0}\delta_{n,n_{0}}+\sum_{t^{\prime}=0}^{t}F_{n_{0},2}(s,1,t^{\prime})P_{s,1}(s,2,t-t^{\prime})+F_{n_{0},2}(s,2,t^{\prime})P_{s,2}(s,2,t-t^{\prime}).$

Proceeding in the $z$-domain, and assuming $n\neq n_{0}$, we find

$\displaystyle\alpha_{m_{0_{1}}}\widetilde{P}_{n_{0},1}(s,1,z)$
$\displaystyle=\widetilde{F}_{n_{0},1}(s,1,z)\widetilde{P}_{s,1}(s,1,z)+\widetilde{F}_{n_{0},1}(s,2,z)\widetilde{P}_{s,2}(s,1,z),$
(S61) $\displaystyle\alpha_{m_{0_{2}}}\widetilde{P}_{n_{0},2}(s,1,z)$
$\displaystyle=\widetilde{F}_{n_{0},2}(s,1,z)\widetilde{P}_{s,1}(s,1,z)+\widetilde{F}_{n_{0},2}(s,2,z)\widetilde{P}_{s,2}(s,1,z),$
$\displaystyle\alpha_{m_{0_{1}}}\widetilde{P}_{n_{0},1}(s,2,z)$
$\displaystyle=\widetilde{F}_{n_{0},1}(s,1,z)\widetilde{P}_{s,1}(s,2,z)+\widetilde{F}_{n_{0},1}(s,2,z)\widetilde{P}_{s,2}(s,2,z),$
$\displaystyle\alpha_{m_{0_{2}}}\widetilde{P}_{n_{0},2}(s,2,z)$
$\displaystyle=\widetilde{F}_{n_{0},2}(s,1,z)\widetilde{P}_{s,1}(s,2,z)+\widetilde{F}_{n_{0},2}(s,2,z)\widetilde{P}_{s,2}(s,2,z),$

which one may formulate as a matrix equation

$\begin{bmatrix}\widetilde{F}_{n_{0},1}(s,1,z)\\\
\widetilde{F}_{n_{0},1}(s,2,z)\\\ \widetilde{F}_{n_{0},2}(s,1,z)\\\
\widetilde{F}_{n_{0},2}(s,2,z)\end{bmatrix}=\begin{bmatrix}\widetilde{P}_{s,1}(s,1,z)&\widetilde{P}_{s,2}(s,1,z)&0&0\\\
0&0&\widetilde{P}_{s,1}(s,1,z)&\widetilde{P}_{s,2}(s,1,z)\\\
\widetilde{P}_{s,1}(s,2,z)&\widetilde{P}_{s,2}(s,2,z)&0&0\\\
0&0&\widetilde{P}_{s,1}(s,2,z)&\widetilde{P}_{s,2}(s,2,z)\end{bmatrix}^{-1}\cdot\begin{bmatrix}\alpha_{m_{0_{1}}}\widetilde{P}_{n_{0},1}(s,1,z)\\\
\alpha_{m_{0_{2}}}\widetilde{P}_{n_{0},2}(s,1,z)\\\
\alpha_{m_{0_{1}}}\widetilde{P}_{n_{0},1}(s,2,z)\\\
\alpha_{m_{0_{2}}}\widetilde{P}_{n_{0},2}(s,2,z)\end{bmatrix}.$ (S62)

Upon evaluating Eq. (S62), one finds Eq. (25) of the main text by performing
$\sum_{m_{0}=1}^{2}\sum_{m=1}^{2}\widetilde{F}_{n_{0},m_{0}}(s,m,z)$, namely
performing a summation over the column vector of state level first passages.

The generalisation of the above procedure to an arbitrary number of states
allows passage to a solution equivalent to Eq. (22) of the main text, while to
study the directional first passage, the first-passage time to one of the
states becomes zero and Eq. (S60) reduces to the form of Eq. (S54).

## The Mean First-Passage Time

### Derivation of the Multi Target Mean First-Passage Time

As mentioned in the main text, the general procedure of finding the multi-
target mean first-passage time (MFPT) may be found in ref. [63] and we refer
the interested reader there. Here, we present the extension from that work
that is required to find the MFPT with correlated random variables.

We start from Eq. (22) of the main text and since we assume here that the
$\bm{U_{m_{0}}}$ is uniformly distributed, one may factor out $\alpha=1/2d$ as
it is a constant factor across all elements of each determinant and rewrite
the first passage as

$\widetilde{F}_{\bm{n}_{0}}(\mathdutchcal{S},z)=\frac{1}{2d}\sum_{j=1}^{M}\sum_{i=1}^{S}\frac{\det[\mathbb{H}^{(i)}(\bm{n}_{0},m_{0_{j}},z)]}{\det[\mathbb{H}(z)]},$
(S63)

where
$\mathbb{H}(z)_{l,k}=\widetilde{P}_{\bm{s}_{k},m_{\bm{s}_{k}}}(\bm{s}_{l},m_{\bm{s}_{l}},z)$,
$\mathbb{H}(z)_{k,k}=\widetilde{P}_{\bm{s}_{k},m_{\bm{s}_{k}}}(\bm{s}_{k},m_{\bm{s}_{k}},z)$
and $\mathbb{H}^{(i)}(\bm{n}_{0},m_{0},z)$ is the same, but with the
$i^{\text{th}}$ column replaced with
$\left[\widetilde{P}_{\bm{n}_{0},m_{0_{j}}}(\bm{s}_{1},m_{\bm{s}_{1}},z),...,\widetilde{P}_{\bm{n}_{0},m_{0_{j}}}(\bm{s}_{S},m_{\bm{s}_{S}},z)\right]^{\intercal}$.
Using the linearity of the outer summation, one may perform
$\mathdutchcal{F}_{\bm{n_{0}}\to\\{\mathdutchcal{S}\\}}=\left.\frac{\partial\widetilde{F}_{\bm{n_{0}}}(\mathdutchcal{S},z)}{\partial
z}\right|_{z=1}$ on each element and following ref. [63], we find

$\mathdutchcal{F}_{\bm{n_{0}}\to\\{\mathdutchcal{S}\\}}=\frac{1}{2d}\sum_{j=1}^{2d}\frac{\det[\mathbb{T}_{m_{0_{j}}}]}{\det[\mathbb{T}_{1}]-\det[\mathbb{T}]},$
(S64)

where $\mathbb{T}_{i,i}=0$,
$\mathbb{T}_{i,k}=\mathdutchcal{F}_{\bm{s}_{k},m_{\bm{s}_{k}}\to\bm{s}_{i},m_{\bm{s}_{i}}}$
($i\neq k$), $\mathbb{T}_{{1}_{i,k}}=\mathbb{T}_{i,k}-1$ and the
$(i,k)^{\text{th}}$ element of $\mathbb{T}_{m_{0_{j}}}$ is
$\mathbb{T}_{i,k}-\mathdutchcal{F}_{\bm{n}_{0},m_{0_{j}}\to\bm{s}_{i},m_{\bm{s}_{i}}}$.
Repeated use of the multi-linearity property of the determinant allows one to
bring inside the summation over the initial states allowing us to arrive at
Eq. (26) of the main text.

### Derivation of the State Level Mean First-Passage Time

We derive here the MFPT from one state to another state, as needed to use Eq.
(26) of the main text. Namely, we are interested in calculating the quantity
$\mathdutchcal{F}_{\bm{n}_{0},m_{0}\to\bm{s}_{i},m_{\bm{s}_{i}}}=\frac{\partial}{\partial
z}\left.\frac{\widetilde{P}_{\bm{n}_{0},m_{0}}(\bm{s},m_{\bm{s}},z)}{\widetilde{P}_{\bm{s},m_{\bm{s}}}(\bm{s},m_{\bm{s}},z)}\right|_{z=1}$.
We give the procedure here for the periodic propagator but note that the
derivation is general as long as three assumptions are held. Namely, (i) one
may use the propagator as the defect-free dynamics on which to place traps
(see Sec. V of the main text for discussion), (ii) the walk is recurrent and
(iii) is not confined to a genuine sublattice.

We begin by writing the first-passage as

$\widetilde{F}^{(p)}_{\bm{n}_{0},m_{0}}(\bm{s}_{i},m_{\bm{s}_{i}},z)=\frac{\bm{e}_{m}^{\intercal}\cdot\left\\{\left[\mathbb{I}-z\bm{\lambda}\left(0,...,0\right)\right]^{-1}+\sum_{\bm{\kappa}\neq\bm{0}}\left(\prod_{i=1}^{d}g^{(p)}_{\kappa_{i}}(s_{i},n_{0_{i}})\right)\left[\mathbb{I}-z\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\right\\}\cdot\bm{e}_{m_{0}}}{\bm{e}_{m}^{\intercal}\cdot\left\\{\left[\mathbb{I}-z\bm{\lambda}\left(0,...,0\right)\right]^{-1}+\sum_{\bm{\kappa}\neq\bm{0}}\left[\mathbb{I}-z\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\right\\}\cdot\bm{e}_{m}^{\intercal}},$
(S65)

and multiplying the top and bottom of Eq. (S66) by
$\det[\mathbb{I}-z\bm{\lambda}\left(0,...,0\right)]$ leads to

$\displaystyle\widetilde{F}^{{}^{(p)}}_{\bm{n}_{0},m_{0}}(\bm{s}_{i},m_{\bm{s}_{i}},z)=$
(S66)
$\displaystyle\frac{\bm{e}_{m}^{\intercal}\cdot\left\\{\text{Adj}\left[\mathbb{I}-z\bm{\lambda}\left(0,...,0\right)\right]+\det[\mathbb{I}-z\bm{\lambda}\left(0,...,0\right)]\sum_{\bm{\kappa}\neq\bm{0}}\left(\prod_{i=1}g^{(p)}_{\kappa_{i}}(s_{i},n_{0_{i}})\right)\left[\mathbb{I}-z\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\right\\}\cdot\bm{e}_{m_{0}}}{\bm{e}_{m}^{\intercal}\cdot\left\\{\text{Adj}\left[\mathbb{I}-z\bm{\lambda}\left(0\right)\right]+\det[\mathbb{I}-z\bm{\lambda}\left(0\right)]\sum_{\bm{\kappa}\neq\bm{0}}\left[\mathbb{I}-z\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\right\\}\cdot\bm{e}_{m}^{\intercal}},$

where $\text{Adj}\left[\mathbb{A}\right]$ denotes the adjoint of $\mathbb{A}$.
The recurrence of the walk ensures
$\det[\mathbb{I}-\bm{\lambda}\left(\bm{0}\right)]=0$, and we utilise this
property upon differentiating to find

$\displaystyle\mathdutchcal{F}^{{}^{(p)}}_{\bm{n}_{0},m_{0}\to\bm{s}_{i},m_{\bm{s}_{i}}}\\!=\\!\frac{1}{\left[\bm{e}_{m}^{\intercal}\cdot\text{Adj}\left[\mathbb{I}-\bm{\lambda}\left(0,...,0\right)\right]\cdot\bm{e}_{m}\right]^{2}}\bigg{\\{}\bm{e}_{m}^{\intercal}\cdot\text{Adj}\left[\mathbb{I}\\!-\\!\bm{\lambda}\left(0,...,0\right)\right]\cdot\bm{e}_{m}\bigg{[}\bm{e}_{m}^{\intercal}\cdot\bigg{(}\left.\frac{\partial\text{Adj}[\mathbb{I}\\!-\\!z\bm{\lambda}(0,...,0)]}{\partial
z}\right|_{z=1}$ (S67)
$\displaystyle\\!+\\!\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,...,0)]}{\partial
z}\right|_{z=1}\sum_{\bm{\kappa}\neq\bm{0}}\left(\prod_{i=1}^{d}g^{(p)}_{\kappa_{i}}(s_{i},n_{0_{i}})\right)\\!\left[\mathbb{I}\\!-\\!\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\\!\bigg{)}\\!\cdot\\!\bm{e}_{m_{0}}\bigg{]}\\!-\\!\bm{e}_{m}^{\intercal}\\!\cdot\\!\text{Adj}\left[\mathbb{I}-\bm{\lambda}\left(0,...,0\right)\right]\\!\cdot\\!\bm{e}_{m_{0}}$
$\displaystyle\times\bigg{[}\bm{e}_{m}^{\intercal}\cdot\bigg{(}\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,...,0)]}{\partial
z}\right|_{z=1}+\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,...,0)]}{\partial
z}\right|_{z=1}\sum_{\bm{\kappa}\neq\bm{0}}\left[\mathbb{I}-\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N_{1}},...,\frac{2\pi\kappa_{d}}{N_{d}}\right)\right]^{-1}\bigg{)}\cdot\bm{e}_{m}^{\intercal}\bigg{]}\bigg{\\}}.$

It is the knowledge of the structure function $\bm{\lambda}(0,....,0)$ that
allows simplification of Eq. (S67). In one dimension we have

$\text{Adj}[\mathbb{I}-\bm{\lambda}(0)]=b\mathbb{J},\quad\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0)]}{\partial
z}\right|_{z=1}=\begin{bmatrix}b-1&b\\\
b&b-1\end{bmatrix}\;\;\text{and}\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0))}{\partial
z}\right|_{z=1}=-2b,$ (S68)

where $\mathbb{J}$ is the all ones matrix, and substitution into Eq. (S67)
leads to, with $g^{(p)}_{\kappa_{i}}(s_{i},n_{0_{i}})$ taken from the main
text,

$\displaystyle\mathdutchcal{F}^{(p)}_{n_{0},m_{0}\to s,m_{s}}=$
$\displaystyle\bm{e}_{m_{s}}^{\intercal}\cdot\begin{bmatrix}0&\frac{1}{b}\\\
\frac{1}{b}&0\end{bmatrix}\cdot\bm{e}_{m_{0}}+$ (S69) $\displaystyle
2\sum_{\kappa=1}^{N-1}\left[\bm{e}_{m_{s}}^{\intercal}\cdot\left[\mathbb{I}-\bm{\lambda}\left(\frac{2\pi\kappa}{N}\right)\right]^{-1}\cdot\bm{e}_{m_{s}}-\exp\left(\frac{-2\pi
i\kappa(s-n_{0})}{N}\right)\bm{e}_{m_{s}}^{\intercal}\cdot\left[\mathbb{I}-\bm{\lambda}\left(\frac{2\pi\kappa}{N}\right)\right]^{-1}\cdot\bm{e}_{m_{0}}\right],$

where

$\displaystyle\big{[}\mathbb{I}-\bm{\lambda}(x)\big{]}^{-1}=\frac{1}{2f(b+f)(1-\cos(x))}\begin{bmatrix}b+f-fe^{-ix}&be^{ix}\\\
be^{-ix}&b+f-fe^{ix}\end{bmatrix}.$ (S70)

The first term in Eq. (S69) corresponds to the dynamics within the site and
leads to, in the degenerate $N=1$ case, a MFPT of zero if $m=m_{0}$ and $b$ if
$m\neq m_{0}$, while the summation corresponds to the dynamics between the
sites.

To simplify Eq. (S69) further, one can exploit the structure of the matrices
and the property that when $m_{0}=m_{s}$ ($m_{0}\neq m_{s}$) a diagonal (anti-
diagonal) element will be chosen. To illustrate, consider the four options,
namely (I) $m_{0}=m_{s}=1$, (II) $m_{0}=m_{s}=2$, (III) $m_{0}=1,m_{s}=2$ and
(IV) $m_{0}=2,m_{s}=1$. In turn, these options give

(I)

$\mathdutchcal{F}^{(p)}_{n_{0},1\to
s,1}=\sum_{\kappa=1}^{N-1}\frac{\left(b+f-fe^{\frac{-2\pi
i\kappa}{N}}\right)\left(1-e^{\frac{-2\pi
i\kappa(s-n_{0})}{N}}\right)}{f(b+f)\left(1-\cos\left(\frac{2\pi\kappa}{N}\right)\right)},$
(S71)

(II)

$\mathdutchcal{F}^{(p)}_{n_{0},2\to
s,2}=\sum_{\kappa=1}^{N-1}\frac{\left(b+f-fe^{\frac{2\pi
i\kappa}{N}}\right)\left(1-e^{\frac{-2\pi
i\kappa(s-n_{0})}{N}}\right)}{f(b+f)\left(1-\cos\left(\frac{2\pi\kappa}{N}\right)\right)},$
(S72)

(III)

$\mathdutchcal{F}^{(p)}_{n_{0},2\to
s,1}=\frac{1}{b}+\sum_{\kappa=1}^{N-1}\frac{\left(b+f-fe^{\frac{2\pi
i\kappa}{N}}\right)-be^{\frac{-2\pi i\kappa(s-n_{0})}{N}}e^{\frac{-2\pi
i\kappa}{N}}}{f(b+f)\left(1-\cos\left(\frac{2\pi\kappa}{N}\right)\right)},$
(S73)

and (IV)

$\mathdutchcal{F}^{(p)}_{n_{0},1\to
s,2}=\frac{1}{b}+\sum_{\kappa=1}^{N-1}\frac{\left(b+f-fe^{\frac{-2\pi
i\kappa}{N}}\right)-be^{\frac{-2\pi i\kappa(s-n_{0})}{N}}e^{\frac{2\pi
i\kappa}{N}}}{f(b+f)\left(1-\cos\left(\frac{2\pi\kappa}{N}\right)\right)}.$
(S74)

Upon combining into one equation, one finds

$\displaystyle\mathdutchcal{F}_{n_{0},m_{0}\rightarrow
s,m_{s}}^{(p)}=\frac{1\\!-\\!\delta_{m_{s},m_{0}}}{b}$ (S75)
$\displaystyle+\sum_{\kappa=1}^{N-1}\frac{\left(b\\!+\\!f\\!-\\!fe^{\frac{2(\delta_{m_{s},2}-\delta_{m_{s},1})\pi
i\kappa}{N}}\right)\\!-\\!e^{\frac{-2\pi
i\kappa(s-n_{0})}{N}}\left[\delta_{m_{s},m_{0}}\left(b\\!+\\!f\\!-\\!fe^{\frac{2(\delta_{m_{s},2}\\!-\\!\delta_{m_{s},1})\pi
i\kappa}{N}}\right)\\!+\\!(1\\!-\\!\delta_{m_{s},m_{0}})be^{\frac{2\pi
i\kappa(m_{s}-m_{0})}{N}}\right]}{f(b+f)(1-\cos\left(\frac{2\pi\kappa}{N}\right))}.$

The degenerate $N=1$ case, with $f=b=\frac{q}{2}$ here, corresponds exactly to
the symmetric walk in an $N=2$ reflective domain, where the MFPT is precisely
zero if $n=n_{0}$ and $\frac{q}{2}$ if $n\neq n_{0}$ [6].

Since the reflective propagator, namely Eq. (11) in the main text, is simply a
sum of periodic propagators, one follows the same procedure used to obtain Eq.
(S69) making it straightforward to find the MFPT in the reflective case as

$\displaystyle\mathdutchcal{F}^{(r_{s})}_{n_{0},m_{0}\to s,m_{s}}$
$\displaystyle=\frac{1-\delta_{m_{s},m_{0}}}{b}+2\sum_{\kappa=1}^{2N-3}\left\\{\bm{e}_{m}^{\intercal}\cdot\left[1+\mu(s)\exp\left(\frac{2\pi
i\kappa(N-s)}{N-1}\right)\right]\left[\mathbb{I}-\bm{\lambda}\left(\frac{\pi\kappa}{N-1}\right)\right]^{-1}\cdot\bm{e}_{m}\right.$
(S76)
$\displaystyle\left.-\bm{e}_{m}^{\intercal}\cdot\left[\exp\left(\frac{\pi
i\kappa(s-n_{0})}{N-1}\right)+\mu(s)\exp\left(\frac{\pi\kappa(2N-s-
n_{0})}{N-1}\right)\right]\left[\mathbb{I}-\bm{\lambda}\left(\frac{\pi\kappa}{N-1}\right)\right]^{-1}\cdot\bm{e}_{m_{0}}\right\\}.$

In the square lattice,

$\displaystyle\mathdutchcal{F}^{{}^{(p)}}_{\bm{n}_{0},m_{0}\to\bm{s}_{i},m_{\bm{s}_{i}}}=\frac{1}{\left[\bm{e}_{m}^{\intercal}\cdot\text{Adj}\left[\mathbb{I}-\bm{\lambda}\left(0,0\right)\right]\cdot\bm{e}_{m}\right]^{2}}\bigg{\\{}\bm{e}_{m}^{\intercal}\cdot\text{Adj}\left[\mathbb{I}-\bm{\lambda}\left(0,0\right)\right]\cdot\bm{e}_{m}\bigg{[}\bm{e}_{m}^{\intercal}\cdot\bigg{(}\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,0)]}{\partial
z}\right|_{z=1}$ (S77)
$\displaystyle+\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,0))}{\partial
z}\right|_{z=1}\sum_{\bm{\kappa}\neq\bm{0}}\exp\left(-\frac{2\pi
i[\kappa_{1}(n_{1}-n_{0_{1}})+\kappa_{2}(n_{2}-n_{0_{2}})]}{N}\right)\left[\mathbb{I}-\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N},\frac{2\pi\kappa_{2}}{N}\right)\right]^{-1}\bigg{)}\cdot\bm{e}_{m_{0}}\bigg{]}$
$\displaystyle-\bm{e}_{m}^{\intercal}\cdot\text{Adj}\left[\mathbb{I}-\bm{\lambda}\left(0,0\right)\right]\cdot\bm{e}_{m_{0}}\bigg{[}\bm{e}_{m}^{\intercal}\cdot\bigg{(}\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,0)]}{\partial
z}\right|_{z=1}+$
$\displaystyle\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,0))}{\partial
z}\right|_{z=1}\sum_{\bm{\kappa}\neq\bm{0}}\left[\mathbb{I}-\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N},\frac{2\pi\kappa_{2}}{N}\right)\right]^{-1}\bigg{)}\cdot\bm{e}_{m}^{\intercal}\bigg{]}\bigg{\\}},$

and we find

$\displaystyle\text{Adj}[\mathbb{I}-\bm{\lambda}(0,0)]=4\ell(b+\ell)^{2}\mathbb{J},\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,0))}{\partial
z}\right|_{z=1}=-16\ell(b+\ell)^{2}\;\;\text{and}\quad\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,0)]}{\partial
z}\right|_{z=1}$ (S78)
$\displaystyle=\left[\begin{smallmatrix}-2(b+\ell)[b-\ell+6(sf+\ell^{2}-c\ell)]&2(b+\ell)^{2}(-1+6\ell)&-4\ell(b+\ell)(-2+3c+3f+3\ell)&-4\ell(b+\ell)(-2+3c+3f+3\ell)\\\
2(b+\ell)^{2}(-1+6\ell)&-2(b+\ell)[b-\ell+6(sf+\ell^{2}-c\ell)]&-4\ell(b+\ell)(-2+3c+3f+3\ell)&-4\ell(b+\ell)(-2+3c+3f+3\ell)\\\
-4\ell(b+\ell)(-2+3c+3f+3\ell)&-4\ell(b+\ell)(-2+3c+3f+3\ell)&-2(b+\ell)[b-\ell+6(sf+\ell^{2}-c\ell)]&2(b+\ell)^{2}(-1+6\ell)\\\
-4\ell(b+\ell)(-2+3c+3f+3\ell)&-4\ell(b+\ell)(-2+3c+3f+3\ell)&2(b+\ell)^{2}(-1+6\ell)&-2(b+\ell)[b-\ell+6(sf+\ell^{2}-c\ell)]\end{smallmatrix}\right],$

which upon substitution and simplification, gives

$\displaystyle\mathdutchcal{F}^{(p)}_{\bm{n}_{0},m_{0}\to\bm{s}_{i},m_{\bm{s}_{i}}}=\delta_{|m-m_{0}|>1}\frac{3\ell+b}{2\ell(b+\ell)}+\delta_{|m-m_{0}|,1}\left((1-\delta_{m_{0}\in\\{2,3\\}}\delta_{m\in\\{2,3\\}})\frac{2}{b+\ell}\\!+\\!\delta_{m_{0}\in\\{2,3\\}}\delta_{m\in\\{2,3\\}}\frac{3\ell+b}{2\ell(b+\ell)}\right)$
(S79)
$\displaystyle+4\underset{\kappa_{1}+\kappa_{2}>0}{\sum_{\kappa_{1}=0}^{N-1}\sum_{\kappa_{2}=0}^{N-1}}\left[\\!\bm{e}_{m}^{\intercal}\cdot\left[\mathbb{I}\\!-\\!\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N},\frac{2\pi\kappa_{2}}{N}\right)\right]^{-1}\cdot\bm{e}_{m}-\exp\left(\frac{-i2\pi[\bm{\kappa}\\!\cdot\\!(\bm{s}-\bm{n}_{0})]}{N}\right)\bm{e}_{m}^{\intercal}\cdot\left[\mathbb{I}\\!-\\!\bm{\lambda}\left(\frac{2\pi\kappa_{1}}{N},\frac{2\pi\kappa_{2}}{N}\right)\\!\right]^{-1}\\!\cdot\\!\bm{e}_{m_{0}}\\!\right].$

In the hexagonal case, we find

$\displaystyle\text{Adj}[\mathbb{I}-\bm{\lambda}(0,0)]=48\ell^{2}(b+2\ell)^{3}\mathbb{J},\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,0))}{\partial
z}\right|_{z=1}=-288\ell^{2}(b+2\ell)^{3},$ (S80)

and we omit explicitly giving
$\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,0)]}{\partial
z}\right|_{z=1}$, owing to its large size. However, upon simplification, many
terms cancel out and one obtains

$\displaystyle\mathdutchcal{F}^{(\mathcal{H})}_{\bm{n}_{0},m_{0}\to\bm{s}_{i},m_{\bm{s}_{i}}}=\delta_{|m-m_{0}|>1}\frac{5\ell+b}{2\ell(b+2\ell)}+\delta_{|m-m_{0}|,1}\left(\delta^{(m,m_{0})}\frac{5\ell+b}{2\ell(b+2\ell)}+(1-\delta^{(m,m_{0})})\frac{3}{b+2\ell}\right)$
(S81)
$\displaystyle+6\sum_{r=0}^{R-1}\sum_{s=0}^{3r+2}\Bigg{\\{}\bm{e}_{m}^{\intercal}\cdot\left(\left[\mathbb{I}-\bm{\lambda}^{(\mathcal{H})}\left(\frac{2\pi\kappa_{1}}{\Omega},\frac{2\pi\kappa_{2}}{\Omega}\right)\right]^{-1}+\left[\mathbb{I}-\bm{\lambda}^{(\mathcal{H})}\left(\frac{-2\pi\kappa_{1}}{\Omega},\frac{-2\pi\kappa_{2}}{\Omega}\right)\right]^{-1}\right)\cdot\bm{e}_{m}$
$\displaystyle-\bm{e}_{m}^{\intercal}\cdot\left(e^{\frac{2\pi
i\bm{k}\cdot(\bm{s}-\bm{n}_{0})}{\Omega}}\left[\mathbb{I}-\bm{\lambda}^{(\mathcal{H})}\left(\frac{2\pi\kappa_{1}}{\Omega},\frac{2\pi\kappa_{2}}{\Omega}\right)\right]^{-1}+e^{\frac{-2\pi
i\bm{\kappa}\cdot(\bm{s}-\bm{n}_{0})}{\Omega}}\left[\mathbb{I}-\bm{\lambda}^{(\mathcal{H})}\left(\frac{-2\pi\kappa_{1}}{\Omega},\frac{-2\pi\kappa_{2}}{\Omega}\right)\right]^{-1}\right)\cdot\bm{e}_{m_{0}}\Bigg{\\}},$

where
$\delta^{(m,m_{0})}=\delta_{m\in\\{2,3\\}}\delta_{m_{0}\in\\{2,3\\}}+\delta_{m\in\\{5,6\\}}\delta_{m_{0}\in\\{5,6\\}}$,
where we have lightened the notation from $\kappa_{i}(r,s)$ (defined just
below Eq. (20) in the main text) to $\kappa_{i}$.

We note here that the procedure can be carried out in arbitrary dimensions,
however numerical procedures may be required to find
$\text{Adj}[\mathbb{I}-\bm{\lambda}(0,...,0)],\quad\left.\frac{\partial\text{Adj}[\mathbb{I}-z\bm{\lambda}(0,...,0)]}{\partial
z}\right|_{z=1}\;\;\text{and}\left.\frac{\partial\det[\mathbb{I}-z\bm{\lambda}(0,...,0))}{\partial
z}\right|_{z=1}$ in higher dimensions.

### The Global Mean First-Passage Time

Here we make connection with ref. [36], where they consider the global mean
first-passage time (GMFPT), that is the MFPT averaged over all initial
conditions. To do so, we first consider the persistence length $l_{p}$,
defined in ref. [36] as the mean number of steps in a ballistic excursion.
Since our formalisms allow for the inclusion of a sojourn probability, we
slightly redefine this quantity as the mean excursion time before a turn takes
place. More precisely, this is the mean of a Geometric distribution i.e.,
$l_{p}=\sum_{t=0}^{\infty}t\left[1-\left(f+c^{(\mathfrak{d})}\right)\right]\left(f+c^{(\mathfrak{d})}\right)^{(t-1)}$,
$(t\geq 1)$ such that
$p_{l}=\left[1-\left(f+c^{(\mathfrak{d})}\right)\right]^{-1}$. Using the
knowledge that the movement parameters sum up to one (see Sec. IV), one may
rewrite this as

$p_{l}=\frac{1}{b+(Z-2)\ell}.$ (S82)

In the diffusive case ($b=\ell=q/Z$) this becomes
$p_{l}=Z\left[q(Z-1)\right]^{-1}$, where $0<q\leq 1$ is the probability of
changing lattice site, with $q=1$ a walker which always moves [6]. In the
square lattice the movement probabilities are $q/4$ to move in any direction
and $1-q$ to remain on the current lattice site, implying $p_{l}=4(3q)^{-1}$.
Eq. (S82) allows for the inclusion of a sojourn probability and does not
demand that $b=\ell$ or $f=\ell$. The definition via the coordination number
of the lattice, $Z$, rather than the spatial dimension means that it is valid
on all lattices rather than just hypercubic.

By performing a summation over $\bm{n}_{0}$ of Eq. (26) of the main text, and
repeatedly using the multi-linearity property of the determinant, we find the
GMFPT as

$\mathdutchcal{G}_{\\{\mathdutchcal{S}\\}}=\frac{\det[\mathbb{G}]}{\det[\mathbb{T}_{1}]-\det[\mathbb{T}]},$
(S83)

where
$\mathbb{G}_{i,j}=\mathbb{T}_{i,j}-(2dN)^{-1}\sum_{\bm{n}_{0}}\sum_{j=1}^{2d}\mathdutchcal{F}_{\bm{n}_{0},m_{0_{j}}\to\bm{s}_{i},m_{\bm{s}_{i}}}$,
while $\mathbb{T}_{1}$ and $\mathbb{T}$ are defined above in Sec. Derivation
of the Multi Target Mean First-Passage Time. In certain situations such as the
one target cases in infinite space, treated in [35], and in periodic domains,
treated in [36], the GMFPT is independent of the target location(s), however
in general this is not the case. In fact,
$\mathdutchcal{G}_{\\{\mathdutchcal{S}\\}}$ is generally dependent on the
location of the target(s), the size of the domain and the movement parameters.

In Fig. S3, we plot, in panel (a), the GMFPT in a reflecting domain with two
targets, with one fixed target at $s_{1}=30$ and varying the position of the
second $s_{2}$. In panel (b), we plot a 2D periodic domain with one target and
obtain minima as a function of $p_{l}$ as observed also in ref. [36]. To do
so, we fix a very low $b=0.01$ and vary $p_{l}$ by changing both $f$ and
$\ell$.

Figure S3: Global mean first-passage time, rescaled by $N$, plotted using Eq.
(S83). (a) By employing Eq. (S76) in (S83), we show the case of two targets in
a one-dimensional reflective domain with one fixed target at $s_{1}=30$ as a
function of the placement of the second target $s_{2}$ for varying persistence
lengths. (b) Using Eq. (S79) we treat the case of an $N\times N$ periodic
domain and plot the GMFPT as a function of the rescaled persistence length for
different domain sizes. In panel (a) when $s_{1}=s_{2}$ the calculation
reduces to a one target case, which explains the slight kink in the plots for
this case. In (b) because of the periodicity of the lattice there is no
dependence on the target location and we drop the subscript
$\\{\mathdutchcal{S}\\}$ on the $y$-label. To vary $p_{l}$ in this case, we
fix $b=0.01$ and vary $f\in[0.02,0.96]$ from low to high values, while
simultaneously varying $\ell\in[0.01,0.58]$ from high to low values to ensure
normalisation with $c^{(2)}=0.01$. The vertical dashed lines show $p_{l}/N$
values for the diffusive cases.

In the one-dimensional reflection case, panel (a), an increase in the
persistence length lowers the MFPT for all target placement, which is expected
due an increase in persistence always aiding coverage in one-dimension. We
also see strong dependence on the placement of the second target with the
higher GMFPTs occurring when $s_{1}$ and $s_{2}$ are close together. This is
expected as any trajectories, which head away from $s_{1}$ do not find a
target until a reflection has occurred. We note the higher GMFPT when $s_{2}$
is to the right of $s_{1}$. To explain this feature, consider the situation
when $n_{0}<s_{1}<s_{2}$, $s_{2}$ is inaccessible, and likewise when
$s_{2}<s_{1}<n_{0}$. With the static target at $s_{1}>N/2$, the case where
$n_{0}<s_{1}<s_{2}$ leads to, on average, longer MFPTs and also occupies a
larger proportion of the contributions to the GMFPT than the
$s_{2}<s_{1}<n_{0}$ case. If $s_{1}=N/2$, one would see a symmetric plot
around $s_{2}$.

In panel (b), we recover similar behaviour to the results reported in ref.
[36], whereby the average search time may be minimised as a function of the
persistence length. With larger domains, the minima occurs at a smaller value
of $p_{l}$, a feature which is expected as the average number of sideways
steps required to reach the target grows with the domain size. In all cases
the minima falls within the persistent regime.

## Stochastic Simulations

Discrete space-time stochastic simulations are fairly straightforward, and the
inclusion of an internal variable, which keeps track of the short term memory
of the system, requires only a marginal extension of standard computational
routines. To simulate the movement process, one must keep track of the
internal variable, which is updated based on the movement direction of the
previous step. Depending on this internal variable, and the global movement
variables ($f,b,\ell,c$, etc), a new movement direction is chosen. The
position of the walker is then updated and, if necessary, one imposes a
boundary condition. We note that when we have considered the ‘squeezed’
boundary condition, the simulations are independent of any periodic dynamics,
and we simply impose the Chandrasekhar boundary condition described in the
main text.

Then, disregarding the internal variable, we check whether the walker is
located at a target site $\bm{s}$, and if so, the first passage event is
recorded with some probability $\rho_{\bm{s}}$, where for the fully absorbing
targets $\rho_{\bm{s}}=1$. As an exception, when one simulates the directional
first-passage, the internal variable is required to check whether a first-
passage event has occurred, since one must know the direction the walker
entered the site from.