A geometric criterion for the optimal spreading of active polymers in porous media

Efficient navigation through disordered, porous environments poses a major challenge for swimming microorganisms and future synthetic cargo-carriers. We perform Brownian dynamics simulations of active stiff polymers undergoing run-reverse dynamics, and so mimic bacterial swimming, in porous media. In accord with experiments of Escherichia coli, the polymer dynamics are characterized by trapping phases interrupted by directed hopping motion through the pores. Our findings show that the spreading of active agents in porous media can be optimized by tuning their run lengths, which we rationalize using a coarse-grained model. More significantly, we discover a geometric criterion for the optimal spreading, which emerges when their run lengths are comparable to the longest straight path available in the porous medium. Our criterion unifies results for porous media with disparate pore sizes and shapes and for run-and-tumble polymers. It thus provides a fundamental principle for optimal transport of active agents in densely-packed biological and environmental settings.

passive counterparts (e.g., ordinary colloids). Yet, their ability to self-propel might not suffice to make them generally suitable for performing complex tasks in biomedical and environmental settings, where, for example, they may be expected to deliver drugs to a specific target [19], penetrate the porous structure of tumors [20,21], or find and induce degradation of contaminants [16,22]. In fact, randomizing the swimming direction of these autonomous agents could be an efficient strategy for reaching a target. To date, however, experimental realizations of controlled reorientation of self-propelled synthetic agents are sparse [23,24].
Experiments of biological microswimmers [25][26][27][28] and synthetic, active agents [29,30] in confined, disordered environments are often concerned with near-surface motility. These studies display a range of unusual phenomena, ranging from the circular swimming motion of bacteria near walls [25,26], hydrodynamic trapping [29,30], and enhancement of bacterial transport near surfaces due to the presence of obstacles [27]. Similarly, theoretical studies on active transport in crowded environments mainly focus on 2D models for active, point-like particles moving in periodic structures [31], on a lattice with obstacles [32], or disordered environments [33]. Accounting for elongated shapes, Brownian dynamics simulations of self-propelled flexible polymers have revealed subdiffusive motion in 2D porous media [34]. Quantitative studies of active transport in 3D porous media, however, are sparse and it was only recently that the hop and trap mechanism of individual E. coli cells moving in a 3D porous structure was identified [4,5].
While such studies shed light on how pore-scale confinement influences bacterial motility, it is still unclear what motility patterns are optimal for spreading in porous media. A clue comes from the seminal work of Wolfe and Berg [1], who studied the spreading of engineered bacteria, which lacked the ability to sense chemical gradients and whose tumbling rate could be controlled chemically. Their experiments indicated that smooth swimming strains of E. coli get stuck in the porous structure of the semi-solid agar, similarly to incessantly tumbling cells, while at an intermediate tumbling rate bacterial transport appeared more efficient. However, the underlying optimal transport mechanism, dictated by the interplay of swimming characteristics and geometric features of the 3D porous medium, remains an open question.
In this work, we elucidate the spreading of selfpropelled stiff polymers, as model systems for elongated microorganisms, in porous media by performing Brownian dynamics simulations and developing a coarsegrained theory. We demonstrate that reorientation mechanisms are indispensable for efficient dispersion through porous media, as intrinsic reversals enhance the overall diffusivities by up to two orders in magnitude. We identify a competition between the pore length, a direct measurement of straight pathways available in the pore space, and the run length of self-propelled polymers. In particular, the hopping lengths of rarely and frequently reversing polymers are constrained by the pore length and their intrinsic run length, respectively. Most importantly, maximal spreading occurs when the intrinsic run length of the polymers is comparable to the longest pore length of the porous medium, which allows us to introduce a simple but robust geometric criterion for optimal transport. Subsequently, we rationalize the non-monotonic transport behavior in terms of a renewal theory.
While such a non-monotonic behavior is predicted in Refs. [32,35], our study unravels the underlying mechanism and demonstrates that this large-scale nonmonotonic behavior persists irrespective of the pore shapes and is dictated by the maximal pore length only. These findings together with our geometric criterion provide fundamental, physical insight into earlier experimental observations [1] and thereby should guide the future design of synthetic cargo-carriers, applicable in biomedical and environmental settings.

I. RESULTS
A. Model: run-reverse polymer in a porous environment We model the elongated shape of a bacterial cell by a stiff polymer with aspect ratio L/σ = 5, where L denotes the polymer length and σ the diameter of the individual monomers. We employ Brownian dynamics simulations of the discretized polymer, where the monomers are connected via stiff springs according to the well-established bead-spring model (see Fig. 1a and Methods). The stiffness of the polymer is characterized by a persistence length p , which exceeds its contour length L, p /L 1. The self-propulsion of the polymer is modeled by active forces of magnitude F p acting on each monomer in the direction tangential to its backbone [36] [Fig. 1a]. The velocity along the contour of the polymer is then determined by the friction coefficient ζ of a single monomer, v c = F p /ζ, and each monomer is subject to translational diffusion with diffusivity D 0 = k B T /ζ, where k B is the Boltzmann constant and T the temperature. The diffusive time scale of a single monomer is τ 0 = σ 2 /D 0 .
In addition, the active polymer reverses its swimming direction at exponentially distributed times with reversal rate λ [ Fig. 1b]. At these 'reversal events' the polymer instantaneously changes its swimming direction and moves along the opposite direction. The run length of the polymer is then defined as the length the polymer moves before it reverses: run ≡ v c /λ. A run-reverse mechanism is employed by several microorganisms [9,12], yet for studying the large-scale spreading of active agents in a porous environment we anticipate that it can be used to model run-and-tumble bacteria [3,8] since often the only way to escape narrow pores is by reversing the swimming direction [4,5]. As such the reversal rate for run-reverse motion should be chosen smaller than the real tumbling rate of bacteria. However, if one is interested in the details of how bacteria explore the pore space and perform foraging, it is important to take into account details of the reorientation mechanism.
The swimming characteristics of run-reverse polymers can be described by two dimensionless parameters: (1) the Péclet number Pe ≡ v c L/D 0 , which measures the self-propulsion strength relative to diffusion, and (2) the reversal rate λ with respect to the characteristic diffusive time scale of a single monomer τ 0 , i.e., λτ 0 .
To model the porous medium, we generate a disordered, monodisperse, porous structure composed of N overlapping spheres of size 4σ within a cubic box of length 30σ [Figs. 1(c) and (d) for N = 10 3 ]. The microarchitecture of the porous medium can be characterized by the distribution of the straight paths, referred to as chord lengths, available in the medium, ϕ Lc ( ) [37]. It is shown for different N in Fig. 1e. We further introduce the maximal chord length L c,max via ϕ Lc (L c,max ) = 10 −5 /σ. Throughout the manuscript we mainly focus on very densely packed environments and hence N = 10 3 , unless stated otherwise. For this case, the maximal chord length is L c,max 20σ and the medium is characterized by narrow channels of average pore diameter σ p comparable in size to the individual monomer σ p /σ 3.5 (or the polymer σ p /L 0.7), as is typical of many natural environments [38] (see Methods).

B. Dynamics: mean-square displacement
To investigate the dynamics of a self-propelled polymer in a porous environment, we measure the meansquare displacement (MSD) of the center monomer FIG. 1. Model set-up of a run-reverse polymer in a porous environment. a Sketch of a self-propelled stiff polymer composed of monomers (gray spheres) with active forces acting on each monomer in the direction tangential to the backbone of the polymer. (Zoom) The polymer chain is modeled using a bead-spring model, where ri is the position of bead i with diameter σ, and ti,i+1 is the tangent vector between bead i and i + 1. The beads are connected by elastic springs. Furthermore, F (i) p = Fp(ti−1,i + ti,i+1) denotes the active force acting on bead i. b Schematic of the run-reverse mechanism of the active polymer. At the run-reverse event the polymer reverses its swimming direction, in particular, the active forces change sign, F c Simulation snap-shots of several active stiff polymers immersed in a porous environment composed of overlapping spheres. |∆r c (t)| 2 with ∆r c (t) = r c (t) − r c (0) and ∆r c (t) = [∆x c (t), ∆y c (t), ∆z c (t)] T . We keep the Péclet number fixed, Pe = 50, and tune the rate of reversal events, λ [ Fig. 2]. At short times t τ diff ≡ D 0 /v 2 c , the MSDs of active agents with different reversal rates λ collapse and display a linear increase reflecting their diffusive motion at short times, which remains independent of the porous structure of the environment. At intermediate times, the directed motion of the polymers dominates and the MSDs exhibit a superdiffusive increase (for λτ 0 5 · 10 −4 ), which varies for different λ. The MSDs eventually cross over to a linear regime, which characterizes the effective diffusive behavior of the run-reverse polymers in a porous medium.
In contrast, non-reversing polymers can only reorient due to the interplay of thermal fluctuations, activity, and conformational chain dynamics [36]. The corresponding rotational relaxation time of the polymer is denoted by τ 0 rot and can be extracted by measuring the fluctuation of the end-to-end direction of the polymer, e(t) = [r 5 (t) − r 1 (t)]/|r 5 (t) − r 1 (t)|, which evolves as |∆e(t)| 2 = 2 − 2 exp(−t/τ 0 rot ) [39]. The MSD for nonreversing polymers, λ = 0 (and polymers with λτ 0 = 5 · 10 −4 ), displays a subdiffusive behavior at intermediate times, t ∼ τ 0 rot , where the effect of the porous environment becomes apparent as the polymer slows down. At long times, the polymer has merely moved a distance   comparable to its own length over the whole simulation time, which is a fingerprint of confined transport and indicates that motion in tight narrow spaces becomes hindered if active agents cannot reverse efficiently.
Moreover, we find that the cross-over time to long-time diffusion is determined by the faster of the two times: the reorientation time due to reversing 1/λ and the orientational relaxation time τ 0 rot (indicated in Fig. 2). To quantify the long-time behavior, we extract the effective diffusivities via Our findings demonstrate that the long-time effective diffusivities can be enhanced by more than two orders of magnitude (over the whole simulation time of 10 3 τ 0 ) upon introducing a reversal mechanism, D sim eff (λτ 0 = 0.5)/D sim eff (λτ 0 = 0) 3 · 10 2 . This is in stark contrast to the motion of run-reverse polymers in a dilute environment, whose long time transport becomes suppressed upon increasing the reversal rate, D sim C. Non-monotonic transport behavior: effective diffusion Most prominently, the MSDs indicate that the longtime diffusivities D sim eff of run-reverse polymers display a non-monotonic behavior with respect to the reversal rate λ. We further introduce the scaled path length Λ = L c,max / run = L c,max λ/v c , which characterizes the tradeoff between the maximal pore length of the medium and the run length of the polymer. We find that the nonmonotonic behavior persists for a broad range of Péclet Hence, we propose that equation (2) serves as geometric criterion for optimal transport of active agents in porous media, which occurs when the run length run is comparable to the maximal pore length L c,max characteristic of the porous the environment. Further, we found that the non-monotonic transport behavior persists for porous media with fewer obstacles, N = 500, 750, but becomes rather weak for dilute environments, N 500, where it approaches the monotonic behavior of active polymers in an unconfined environment, see Fig. 3b. The chord-length distributions, ϕ Lc ( ) [ Fig. 1e], show a strong decay for densely packed environments. Most importantly, our results demonstrate a data collapse for N 500 upon rescaling with the longest pore length available in the environment, L c,max , see Fig. 3c. This result emphasizes that indeed the longest pore length of the environment is the characteristic length scale dictating the transport behavior of active polymers and furthermore confirms our geometric criterion.
In contrast to previous work [40], we have also found non-monotonic spreading of active agents in porous environments with concave pore shapes [ Fig. 5 in Methods]. This highlights that solely the interplay of length scales dictates this intricate behavior and not the pore shape.
Our findings offer an interpretation for the observations of the seminal experimental work by Wolfe and Berg [1], who used engineered strains of non-chemotactic  E. coli. The tumbling rate and thus the run length of these mutants can be tuned by varying the concentration of an external inducer (isopropyl β-D-thiogalactoside (IPTG)) in the medium. Monitoring the flagella bundles of cells tethered to a surface suggested an increase of the tumbling rate upon increasing the concentration of IPTG. Qualitatively, the experiments of Wolfe and Berg revealed that bacterial swarms that tumble at an intermediate rate spread more than either incessantly tumbling cells or smooth swimming, i.e., non-tumbling, strains. However, the dynamical behavior of individual cells has not been quantified. Therefore, we anticipate that our simulations, which qualitatively confirm these experimental findings, enable us to elucidate the physics of the non-monotonic behavior and provide predictions for the optimal spreading.

D. Individual trajectories
To elucidate this non-monotonic behavior of the diffusivities, we investigate individual trajectories for Pe= 50 [ Fig. 4a]. For rarely reversing agents, λτ 0 = 5 × 10 −4 , the trajectories in the xy plane indicate highly confined motion, whereas the trajectories of a polymer reversing at higher rates, λτ 0 = 0.5 and 5, cover significantly larger areas including many pores [ Fig. 4a(inset)]. By inspecting the associated time evolution of the displacements, we observe that the motion of rarely reversing polymers λτ 0 = 5 · 10 −4 , i.e., the purple curve in Fig. 4, is characterized by fast hopping events, at which the polymer moves through the pore space, and long, extended phases of trapping, where the polymer is trapped inside a pore. This motility pattern agrees with experimental findings of wild-type E. coli moving in a porous structure [4,5].
We further observe that upon increasing the reversal rate the trapping events become shorter, see orange curve for λτ 0 = 0.5 in Fig. 4a. In fact, they disappear for large reversal rates, λτ 0 = 5, where the trajectory is dominated by hopping events, i.e., the green curve in Fig. 4a. Now, we explain the mechanism for the trapping and hopping, which identifies the main features necessary for an active agent to explore a porous medium. Suppose the active agent is exploring the porous environment with an intrinsic run time τ run ∼ 1/λ and after a while it enters a dead-end pore at a time t. Consequently, the agent will be trapped within the dead-end pore for the remaining time interval [t, τ run ]. This implies that the trapping time of a self-propelled agent can potentially be reduced by increasing the reversal rate λ. The exemplary trajectory for λτ 0 = 0.5, which corresponds to the optimal reversal rate, exhibits the largest displacements with few trapping events [ Fig. 4a]. Moreover, we observe that despite the sparse and short trapping events of frequently reversing polymers, λτ 0 = 5, their hops become significantly shorter.
We quantify this behavior by extracting the distributions (1) ϕ T (t) for the trapping time, i.e., the time the polymer spends trapped inside a pore, and (2) ϕ H ( ) for the hop length, i.e., the distance the polymer moves from one trapping event to another or from one trapping event to the next reversal, from the individual trajectories (see Methods).

E. Trapping time distribution
We find that the trapping time distributions for active polymers with Pe= 50 exhibit a power-law scaling at long times with an exponent that increases from 2 to 7/2 for increasing reversal rates λ [ Fig. 4b]. To rationalize this behavior, we extend the concept of an 'entropic trap' model introduced in the study of diffusing long polymers escaping small outlets [41] and recently for E. coli bacteria moving in a porous medium [4,5]. In this model, the escape of an active polymer from an obstruction in the porous medium is determined by the number of orientations (Ω t ) keeping it trapped inside the pore and the number of orientations (Ω e ) allowing the polymer to leave it. In short (see Methods), our model predicts a power-law trapping time distribution ϕ T (t) ∼ t −(1+β) at long times. Here, the exponent β = X/C 0 is related to the active energy X, which has dimensions of F p L, and the average depth of the entropic trap C 0 , which can be quantified by the average free energy difference between the two states: C 0 = k B T ln(Ω t /Ω e ) , where the brackets correspond to the ensemble average over all pores.
In accord with this phenomenological model, we find that the exponent β increases for increasing reversal rate λ [ Fig. 4b]. Specifically, at a fixed Pe (corresponding to X = const.), a larger reversal rate increases the probability Ω e of leaving the (e.g., dead-end) pores and leads to a relatively lower trap depth C 0 , and thus a larger exponent β = X/C 0 . Further, we note that for β ≤ 1, corresponding to the case where the active energy becomes smaller than the trap depth X ≤ C 0 , the entropic trap model predicts a divergence of the mean trapping duration τ T ∼ (β − 1) −1 . This becomes evident in the monotonic increase of the mean trapping duration by orders of magnitude with decreasing reversal rates [ Fig. 4b(inset) open circles].

F. Hopping length distribution
We further extract the distributions of the hopping length of the agents ϕ H ( ) [ Fig. 4c for Pe= 50] and rationalize that the hop lengths are determined by the interplay of the intrinsic run length of the active polymers run ≡ v c /λ and the pore geometry. In particular, we find that the chord length distribution [ Fig. 4c black solid line] agrees with the hopping length distribution of rarely reversing polymers, λτ 0 = 5 · 10 −4 (and 5 · 10 −3 ). These agents have an intrinsic run length of run 10 3 σ (and 10 2 σ), which is much longer than the longest available straight pathways, L c,max 20σ, and therefore their hopping motion is fully determined by the chord lengths of the porous medium [ Fig. 4c]. Further, the average hop length is comparable to the average chord length for small λ, H ∼ L c [ Fig. 4c(inset)]. In contrast, for large reversal rates λτ 0 5 the distribution decays even faster and the average hop length approaches the run length of the polymer, H → run .
Our findings further show that polymers with an intermediate reversal rate, e.g., λτ 0 = 0.5, can follow the straight path available in the porous structure and continue to explore another successive pore without getting trapped, which leads to hopping lengths longer than the longest chord length. Furthermore, we find that the probability for longer hops becomes larger at an intermediate reversal rate than the chord-length distribution, which indicates that the polymer explores these more often than the shorter pores of the medium.
Most significantly, our results demonstrate that the probability for long hops and also the average hop length H [ Fig. 4c(inset)] vary non-monotonically with λτ 0 : they are lowest for λτ 0 = 5 · 10 −4 and λτ 0 = 15 and highest for λτ 0 = 0.5. This strong non-monotonic variation of the average hop length could explain the optimal spreading and thus the maximal long-time effective diffusivity of run-reverse polymers in porous media [ Fig. 3a]. In particular, agents with Péclet number Pe= 50 and reversal rate λτ 0 = 0.5 have an intrinsic run length of run = 20σ, comparable to the longest pore length L c,max . Our findings reveal that the active polymer continuously explores the pore space without getting trapped too often [ Fig. 4a], which suggests the qualitative picture: the agent moves through the pores until it reaches an obstruction, where subsequent reversals allow it to continue to explore the pore space.

G. Coarse-grained dynamics
To quantitatively characterize the long-time effective diffusivities, we develop a coarse-grained model for the 3D dynamics of the alternating hopping and trapping phases. During the hopping phase the agent, modeled as a point particle, moves straight at effective velocity v and reverses its swimming direction at rate λ with propagator P H (∆r, t), which denotes the probability that the particle displaces ∆r during time t while hopping. During the trapping phase the particle cannot move and after the trapping event it hops along a new, random direction. The time the particle spends in a hopping or trapping phase is determined by the hopping and trapping time distributions extracted from the simulations [ Fig. 4b-c]. As the swim speed of the polymer during the hopping events remains roughly constant, the hopping time, i.e., the duration of a hopping phase of length , t ∼ /v, also follows an exponential distribution ϕ H (t) = exp(−t/τ H )/τ H with mean duration τ H . The trapping time distribution obeys a power-law behavior ϕ T (t) = β(1 + t/τ ) −1−β /τ with mean duration τ T = τ (β − 1) −1 , where τ is a characteristic time scale for trapping. This power-law behavior suggests that the dynamics are non-Markovian and therefore we use a renewal theory [42,43].
The probability density P (∆r, t) for the particle to have displaced ∆r during lag time t is the sum of the probability densities for the particle to be in a hopping or a trapping phase: P (∆r, t) = P H (∆r, t) + P T (∆r, t). We further introduce the probability densities (per time) for the particle to start a hopping phase and a trapping phase at displacement ∆r and at lag time t by H(∆r, t) and T (∆r, t), respectively. Then the probabil-ity P T (∆r, t) that the particle is trapped at ∆r for lag time t is obtained as the sum of the probability of never having hopped before P (0) T (∆r, t) and the probability of having hopped at least once before getting trapped at ∆r and at an arbitrary earlier time t − t : Here, ϕ where T (1) (∆r, t) denotes the probability for the first trapping event. The second term corresponds to the sum over the probabilities that the agent has escaped a trap and started to hop at ∆r − at an earlier time, t − t , until it gets trapped again at ∆r and t. Similar equations hold for P H (∆r, t) and H(∆r, t) (see equations (A6)-(A7) in Methods). We can derive an analytic solution for the probability density in Fourier-Laplace [42]. By isotropy, it can be expanded for small wave numbers k up to O(k 4 ) via P (k, s) ∞ 0 dt e −st 1−k 2 |∆r(t)| 2 /3! = s −1 − k 2 |∆r(s)| 2 /3!, which allows for an analytical derivation of the Laplace transform of the mean-square displacement |∆r(s)| 2 .
Finally, we obtain the long-time transport behavior by taking the limit of |∆r(s)| 2 for s → 0. We find lim s→0 |∆r(s)| 2 6D theo eff s −2 , which corresponds to a long-time diffusive behavior |∆r(t)| 2 6D theo eff t for t → ∞ with effective diffusivity: It is determined by the fraction of time spent hopping between traps, τ H /(τ H + τ T ), the hopping time, τ H , the effective velocity, v, and the tumbling rate λ. The tumbling angle for run-reverse motion, as employed here, is ϑ 0 = −π [9]. The expression [equation (5)] depends on the exponent β of the trapping time distribution via the average trapping time τ T . To compare the simulation data with our coarsegrained theory, we extract the average hopping and trapping times, τ H and τ T , from the individual trajectories and obtain the effective diffusivity from the hop-andtrap model [equation (5)]. We observe that our model captures the non-monotonic trend of the simulation results, D sim eff ∼ D theo eff [ Fig. 3a dashed lines]. It is interesting that the hopping and trapping times, used as inputs for the model, can semi-quantitatively predict the nonmonotonic trend of the effective diffusivities of a complex system. The deviations at small reversal rates, corresponding to the scaled path length Λ 10 −1 , are expected because at such small reversal rates the polymers remain trapped most of the time. Therefore, the extracted mean trapping times are underestimated and even longer simulations would be required.
We rationalize the non-monotonic behavior by inspecting the fraction of time spent hopping, p = τ H /(τ H + τ T ) [ Fig. 4b inset]. As expected, for rare reversal events the fraction of time spent hopping is small p 1/2. This corresponds to τ H τ T and therefore equation (5) reduces to D eff = v 2 τ 2 H /(3τ T ) → 0, as large trapping times suppress transport. In contrast, the trapping times for frequently reversing polymers are negligible compared to the hopping times, τ T τ H as p 1/2. Therefore, the behavior of the effective diffusivity is dominated by the reversal rate and equation (5) simplifies to D eff ∼ v 2 /λ → 0, which vanishes as the reversal rate increases. Hence, maximal transport occurs at an intermediate reversal rate λ, in agreement with our simulation results for the long-time effective diffusivities and the average hopping lengths.
Our results indicate that the non-monotonic behavior of the effective diffusivities is a clear signature of hopand-trap dynamics. This behavior vanishes for dilute environments, where the dynamics are solely determined by the run-reverse motion of the polymers. The transition between both motility modes is shown in Fig. 3b.

II. SUMMARY AND CONCLUSION
To rationalize the seminal experiments of Wolfe and Berg [1], which characterized spreading in a porous environment of bacteria with different tumbling rates, we have performed Brownian dynamics simulations of runreverse stiff polymers in a porous environment and find a non-monotonic transport behavior as a function of the reversal rate. We introduce, for the first time, a geometric criterion for the optimal spreading of active agents in a porous environment, which occurs when the intrinsic run length agrees with the longest straight path available in the environment. We further show that individual polymer trajectories exhibit a hop-and-trap mechanism, in accord with recent experiments of E. coli [4,5]. Our results suggest that optimal transport at an intermediate reversal rate is characterized by a maximal average hop length. In contrast, the motion of rarely and frequently reversing polymers is set by the pore length or the runlength, respectively. We corroborate these findings using a renewal theory for the coarse-grained hop-and-trap dynamics.
Our results demonstrate that this non-monotonic transport behavior of active agents in a porous environment persists irrespective of the shape of the pores, in contrast to earlier predictions [40]. These findings indi-cate that the 'size of the pores, not their shape, matters' for this large-scale spreading and we therefore anticipate that this behavior is universal for densely packed environments.
In the future, it will be interesting to study the effect of swimmer shape anisotropy on active transport in a porous environment. In particular, recent work has shown that the non-monotonic behavior persists even for point particles on a 2D lattice with obstacles [32]. However, the interplay between pore geometry and run length remains to be explored. Most importantly, our study and recent experiments [4,5] predict power-law distributions for the trapping times, which could be amplified by particle shape. Taking into account this memory dependence in a coarse-grained description goes beyond earlier theoretical predictions [32,35].
We anticipate that our theoretical findings can be tested in various biological systems, such as bacteria [3,[8][9][10], algae [11], or archaea [12]. This could shed light on fundamental microbiological processes, which include the adaption of the tumbling behavior under varying environmental conditions. Experimental observations of Bacillus subtilis [44] have shown that these cells can reverse their swimming direction once they encounter an obstacle and, similarly, Spirochetes [45] display an increase of reversal rate while invading a heterogeneous fibrous medium. It would be interesting to elucidate if these spatially varying reversal mechanisms allow them to optimize their motion in a 3D porous medium. Similarly, our results could guide the design of future synthetic swimmers with, e.g., specific magnetic properties [46]. Their spreading could be enhanced by reorienting them via externally applied magnetic fields. On the macroscale, our findings could find application in instructing robots during search and rescue operations in disaster zones [40].
Beyond porous media, our results lay the foundation for studying transport of active semiflexible polymers in dynamically re-arranging environments composed of other polymers, such as the interior of cells [39,47]. The strongly interacting, crowded environment could lead to a relaxation of the extended trapping events of the active agents and entail complex entanglement effects of the individual constituents.

DATA AVAILABILITY
The data used for this paper are available from the corresponding authors upon request.
(see later section on details about the porous environment).
The elasticity of the chain is characterized by the wellestablished wormlike chain (WLC) model pioneered by Kratky and Porod [48]. The discretized interaction energy is where we have introduced the persistence length p of a semiflexible polymer. It is a measure of the decay length of its tangent-tangent correlations and allows distinguishing between flexible, p /L 1, semiflexible, p /L 1, and stiff polymers, p /L 1.
Interactions between neighboring beads are modeled by the finitely extensible non-linear elastic (FENE) potential, which ensures a finite, maximal distance between two monomers, L/N p + δ. For monomer-monomer distances larger than L/N p + δ the interaction potential diverges, U FENE = ∞. Consequently, the entire polymer chain has a maximal length of L + N p δ. The interactions between non-neighboring polymer beads and the interactions between the polymer and the surrounding, porous environment are modeled using the Weeks-Chandler-Anderson (WCA) potential [49], Here, d = σ for the interaction of two monomers and d = (σ + σ s )/2 for the interaction between a monomer and an obstacle of the porous environment with diameter σ s (see later section on the Porous environment). The full contribution is the sum over all monomers and beads of the porous environment U WCA = i,j;i =j U WCA ij . The total interaction energy is then U = U WLC + U FENE + U WCA , which includes (amongst others) three important dimensionless parameters: the persistence length of the polymer chain relative to the contour length p /L, and FENE and WCA , which measure the relative importance of the interaction energies with respect to thermal energy. In the present study, we simulate stiff polymers with p /L = 50. To assure minimal polymer extension and prevent overlap of monomers and the environment, we set FENE = 10 3 and WCA = 5, respectively. Furthermore, we use a Brownian time step of τ B = 10 −6 τ 0 with τ 0 the diffusive time scale of a monomer and wait for 10 9 τ B before taking measurements.
In addition to self-propulsion, the polymer performs a run-reverse motion, as sketched in Fig. 1b. The reversal events occur at randomly drawn times, t, which follow an exponential distribution λ exp(−λt) with reversal rate λ. At the run-reverse event, the active forces, acting along the polymer chain, change sign, F p . Thus, the polymer instantaneously moves along the new, opposite swimming direction. In our framework the reversal events occur instantaneously, so that the polymer does not stop before moving into the opposite direction.

Porous environment
The porous structure of the environment is generated by randomly distributed obstacles of diameter σ s , which are allowed to overlap. By varying the diameter of the obstacles, we tune the average pore diameter, σ p , of the medium, inspired by the experimental set-up of Refs. [4,5]. The porous environment is characterized by the average pore diameter and the chord length. For ensemble averaging, we use 20 statistically independent structures.
We extract the average pore diameter of the medium by monitoring the transport behavior of Brownian particles with radius σ B in the porous environment. Therefore, we measure the mean-square displacement, (∆r B (t)) 2 with ∆r B (t) = r B (t) − r B (0), and calculate the local exponent, α(t) = d ln |∆r B (t)| 2 /d ln t, as a function of time. The local exponent displays a minimum at an intermediate time τ min , which allows introducing the pore diameter by σ p = σ B + |∆r B (τ min )| 2 1/2 . For simplicity, we choose σ B = σ.
To measure the chord length distribution ϕ Lc ( ), we image several two dimensional planes randomly passing through the simulated porous medium [Fig. 1d]. These images provide a map of the pore space. We then binarize these images where two phases represent solid obstacles and open pores, respectively. We calculate the distribution of chords of length , which fit within each binarized pore space image. This protocol yields a direct measurement of straight pathways available in the pore space.

Effect of pore shape
We have addressed the effect of pore shape on the large-scale spreading of active agents. In particular, we have replaced the WCA potential [equation (A4)], corresponding to pores with convex boundaries, and modeled the interaction between the polymer and concave pores by the interaction potential: where R i is the distance between monomer i of the polymer and the nearest obstacle. We choose C = 100 and a = 1.75, corresponding to concave voids of diameter ∼ 3.5σ, as illustrated in Fig. 5c. We find that also for pores with concave walls the long-time behavior of the mean-square displacement is diffusive, see Fig. 5a. Extracting the long-time effective diffusivities shows that the non-monotonic behavior persists as a function of the reversal rate [ Fig. 5b]. The maximal diffusivity occurs at a reversal rate of λτ 0 = 0.5 and a run length of run = 20σ. This optimal behavior agrees with the findings for a porous environments with convex pore shapes, as discussed in the main text. Overall, the effective diffusivities are smaller than for a convex environment, which may indicate that the polymer explores the individual pores for a longer time and should depend on the energy depth C .

Data analysis of the individual trajectories
We extract the distributions for the trapping time, hopping time, and the hopping length from the individual trajectories. We note that the hopping time is defined as the time between hopping from one trap to the next or from one trapping event to the next reversal and the hopping length is the length the polymer moved during this time. To extract these quantities, we follow the approach from Refs. [4,5] and first measure the average velocity of a non-tumbling polymer in a free environment, v . Then we calculate the instantaneous velocities of the center monomer, v i = |r c (t i+1 )−r c (t i )|/(t i+1 −t i ), where t i corresponds to the i-th time step, of individual particle trajectories. Thus, we can classify a hopping phase by v i ≥ v /2 and a trapping phase by v i < v /2, which allows extracting the trapping and hopping time distributions, ϕ T (t) and ϕ H (t) with mean durations, τ T and τ H . In addition, we keep track of the reversal times, which are input to our simulations, and compare it with the hopping duration. This provides the hopping length distribution, ϕ H (t) with mean hopping length H . In particular, if the reversal time since the last trapping event t λ is shorter than the hopping phase t hop , the hopping length corresponds to the length displaced until the reversal event. For t λ > t hop the hopping length is the displacement from one trapping event to the next. For the comparison of the effective diffusivities extracted from simulations, D sim eff , to the predictions of the hop-and-trap model [equation (5)], we further use as effective velocity the cut-off velocity, v = v /2. We can fully recover the non-monotonic behavior and explain the simulation data up to a constant pre-factor.
To test our approach, we have varied the cut-off velocity between v /3 to 2 v /3 and found that it does not change our conclusions: the power-law behavior of the trapping time distributions, the behavior of the hopping length distributions, and the semi-quantitative agreement of the effective diffusivities of the coarse-grained model and the simulations remain preserved.

Entropic trap model for the trapping time distribution
Motivated by other disordered media [4,5,41], the probability for the entropic trap C is assumed to follow P (C) = C −1 0 exp(−C/C 0 ) with average trap depth C 0 . By analogy to equilibrium physics, the probability of an active polymer to escape a trap of depth C is assumed to obey an Arrhenius-like relation. Then the trapping duration is given by t = τ (exp(C/X) − 1), where X characterizes the active energy due to its swimming motion and τ corresponds to a characteristic time scale for trapping. In particular, the trapping duration vanishes, t → 0, for small entropic traps, i.e., C/X → 0. For a passive polymer the active energy X is replaced by the thermal energy k B T . The probability distribution of the trapping times can be obtained as ϕ T (t) = P (C)(∂t/∂C) −1 = β(1 + t/τ ) −1−β /τ with β = X/C 0 .

Renewal theory for the hop-and-trap dynamics
The probability density for a particle to be in a hopping phase follows: Here, P H (∆r, t) denotes the probability that the particle has never been trapped before and the second term corresponds to the sum over all hopping phases, which started after at least one trapping event. Further, ϕ (0) is the probability that the hopping time exceeds t. The probability density (per time) that a new hopping phase starts obeys the equation of motion H(∆r, t) =H (1) (∆r, t) where H (1) (∆r, t) is the probability that the particle starts the first hop. After a Fourier transform of the probability densities, ∆r → k, and by the convolution theorem, the renewal equations [equations (3),(4),(A6),(A7)] simplify to P T (k, t) = P with p = τ H /(τ H + τ T ), which account for the fact that the system starts in a stationary state [42]. In particular, the probability to have never hopped before, P T (k, t), depends on the probability that the particle is in a trapped state, 1 − p, and on the time integral, which represents the probability that the trapping phase exceeds time t. It can be rationalized as follows: The probability density that the trapping phase is of length t is given by t ϕ T (t )/τ T . The probability that after lag time t the particle is still trapped is (t − t)Θ(t − t)/t , where Θ(·) denotes the Heaviside function. Then the probability that the particle has not yet started a hopping phase at time t is obtained by integrating over all durations t : ∞ t dt ϕ T (t )(t − t)/τ T . Similarly, we can derive P Moreover, the probability densities (per time) for the first trapping and hopping event are Here, the probability for the first trapping event T (1) (k, t) depends on the probability that the particle is in a hopping state p and has hopped for a time t with propagator P H (k, t). Further, the probability density for a hopping phase to be of length t is given by t ϕ H (t )/τ H . The probability for the lag time t = 0 to be in the same interval is uniformly distributed, 1/t . Then the probability density that the trapping phase starts at t given a hopping phase of length t is obtained by the integral over all possible t , ∞ t dt ϕ H (t )/τ H . Similar considerations hold for H (1) (k, t).
Finally, we specify the propagator in the hopping phase P H (k, t). Since we are interested in terms up to O(k 2 ), we expand it in k: P H (k, t) 1 − k 2 |∆r(t)| 2 RR /3!, where the mean-square displacement of a run-reverse particle is |∆r(t)| 2 RR = 2v 2 /λ 2 exp(−λt) +λt − 1 . The effective rate depends on the turning angle ϑ 0 viã λ = λ(1−cos ϑ 0 ) [9]. Subsequently, we perform a Laplace transform, t → s, of equations (A8a)-(A8d) and use the convolution theorem to derive an analytical solution for the renewal equations in Fourier-Laplace space [42]. The formal solution has been presented elsewhere [42]. We insert the expansion of P H and the hop and trapping time distributions, ϕ H and ϕ T , into the theoretical predictions in Fourier-Laplace space and keep only terms up to O(k 4 ). Using the expansion from the main text, we derive an analytical solution for the mean-square displacement in Laplace space, |∆r(s)| 2 . An analytical backtransform to time space is not possible, however, we can extract the long-time (corresponding to s → 0) behavior, see main text.
We further note that the long-time effective diffusivity can also be derived analytically for a truncated powerlaw distribution, ϕ T (t) = exp(−γt)(1 + t/τ ) −1−β f (γτ, β) with f (γτ, β) = [e γτ β ExpIntE(1 + β, γτ )] −1 , where ExpIntE(·, ·) denotes the generalized exponential integral function [50]. It assumes the same form as equation (5) with average trapping time τ T = We note that the average trapping time of the truncated powerlaw distribution reduces for γ = 0 to that of the powerlaw distribution used in our manuscript. Details of the distribution become apparent in the short-time behavior of the mean-square displacement, but do not affect our data analysis for the long-time effective diffusivities.