Abstract

For over a century, scientists have sought to understand how fish orient against an incoming flow, even without visual and flow cues. Here, we elucidate a potential hydrodynamic mechanism of rheotaxis through the study of the bidirectional coupling between fish and the surrounding fluid. By modeling a fish as a vortex dipole in an infinite channel with an imposed background flow, we establish a planar dynamical system for the cross-stream coordinate and orientation. The system dynamics captures the existence of a critical flow speed for fish to successfully orient while performing cross-stream, periodic sweeping movements. Model predictions are examined in the context of experimental observations in the literature on the rheotactic behavior of fish deprived of visual and lateral line cues. The crucial role of bidirectional hydrodynamic interactions unveiled by this model points at an overlooked limitation of existing experimental paradigms to study rheotaxis in the laboratory.

Editor's evaluation

The authors present a simple model of fish swimming in a channel and reacting to the surrounding flow with their lateral line and no other sensory system. They demonstrate that the fish stably orients upstream in certain conditions. Particularly, rheotaxis can emerge even in the absence of sensory feedback, purely as a consequence of passive hydrodynamic interactions in the presence of the walls.

https://doi.org/10.7554/eLife.75225.sa0

eLife digest

One fascinating and perplexing fact about fish is that they tend to orient themselves and swim against the flow, rather than with it. This phenomenon is called rheotaxis, and it has countless examples, from salmon migrating upstream to lay their eggs to trout drift-foraging in a current. Yet, despite over a century of experimental studies, the mechanisms underlying rheotaxis remain poorly understood. There is general consensus that fish rely on water- and body-motion cues to vision, vestibular, tactile, and other senses. However, several questions remain unanswered, including how blind fish can perform rheotaxis or whether a passive hydrodynamic mechanism can support the phenomenon. One aspect that has been overlooked in studies of rheotaxis is the bidirectional hydrodynamic interaction between the fish and the surrounding flow, that is, how the presence of the fish alters the flow, which, in turn, affects the fish.

To address these open questions about rheotaxis, Porfiri, Zhang and Peterson wanted to develop a mathematical model of fish swimming, one that could help understand the passive hydrodynamic pathway that leads to swimming against a flow. Unlike experiments on live animals, a mathematical model offers the ability to remove cues to certain senses without interfering with animal behavior.

Porfiri, Zhang and Peterson modeled a fish as a pair of vortices located infinitely close to each other, rotating in opposite directions with the same strength. The vortex pair could freely move through an infinitely long channel with an imposed background flow, devoid of all sensory information expect of that accessed through the lateral line. Analyzing the resulting system revealed that there is a critical speed for the background flow above which the fish successfully orients itself against the flow, resulting in rheotaxis. This critical speed depends on the width of the channel the fish is swimming in. Depriving the fish of sensory information received through the lateral line does not preclude rheotaxis, indicating that rheotaxis could emerge in a completely passive manner.

The finding that the critical speed for rheotaxis depends on channel width could improve the design of experiments studying the phenomenon, since this effect could confound experiments where fish are confined in narrow channels. In this vein, Porfiri, Zhang and Peterson’s model could assist biologists in designing experiments detailing the multisensory nature of rheotaxis. Evidence of the importance of bidirectional hydrodynamic interactions on fish orientation may also inform modeling research on fish behavior.

Introduction

Swimming animals display a complex behavioral repertoire in response to flows (Chapman et al., 2011). Particularly fascinating is the ability of several fish species to orient and swim against an incoming flow, a behavior known as rheotaxis. While intuition may suggest that vision is necessary for fish to determine the direction of the flow, several experimental studies of midwater species swimming in a channel have documented rheotaxis in the dark above a critical flow speed (Coombs et al., 2020). When deprived of vision, fish lose the ability to hold station and they may perform sweeping, cross-stream movements from one side of the channel to other (Bak-Coleman et al., 2013; Bak-Coleman and Coombs, 2014; Elder and Coombs, 2015; Figure 1).

Figure 1

Download asset Open asset

Fish rheotaxis.

(a) Illustration of the problem with notation, showing a fish swimming in a background flow described by Equation 4. (b) Schematic of the cross-stream sweeping movement of some fish species swimming without visual cues; snapshots of fish at earlier time instants are illustrated by lighter shading.

In addition to vision, fish may rely on an array of compensatory sensory modalities to navigate the flow, which utilizes tactile, proprioceptive, olfactory, electric, kinematic, and hydrodynamic signals (Montgomery et al., 2000; von der Emde, 1999). For example, fish could sense and actively respond to linear accelerations caused by the surrounding flow using their vestibular system (Pavlov and Tjurjukov, 1995). Similarly, with the help of tactile sensors on their body surface, fish could maintain their orientation against a current through momentary contacts with their surroundings (Arnold, 1969; Lyon, 1904). Several modern studies have unveiled the critical role of the lateral line system, an array of mechanosensory receptors located on the surface of fish body (Montgomery and Baker, 2020), in their ability to orient against a current (Baker and Montgomery, 1999; Montgomery et al., 1997), hinting at a hydrodynamics-based rheotactic mechanism that has not been fully elucidated. When deprived of vision, can fish rely only on lateral line feedback to perform rheotaxis? Is there a possibility for rheotaxis to be achieved through a purely passive hydrodynamic mechanism that does not need any sensing?

Through experiments on zebrafish larvae swimming in a laminar flow in a straight tube, Oteiza et al., 2017 have recently unveiled an elegant hydrodynamic mechanism for fish to actively perform rheotaxis. Utilizing their mechanosensory lateral line, fish can sense the flow along different parts of their body, which is sufficient for them to deduce local velocity gradients in the flow and adjust their movements accordingly. As further elaborated upon by Dabiri, 2017, the insight offered by Oteiza et al., 2017 is grounded in the fundamental relationship between vorticity and circulation given by the Kelvin-Stokes’ theorem, so that fish movements will be informed by local sampling of the vorticity field. While offering an elegant pathway to explain rheotaxis, the framework of Oteiza et al., 2017 does not include a way for rheotaxis to be performed in the absence of information about the local vorticity field. Several experimental studies have shown that fish can perform rheotaxis even when their lateral line is partially or completed ablated, provided that the flow speed is sufficiently large (Bak-Coleman et al., 2013; Bak-Coleman and Coombs, 2014; Baker and Montgomery, 1999; Elder and Coombs, 2015; Montgomery et al., 1997; Oteiza et al., 2017; Van Trump and McHenry, 2013).

Mathematical modeling efforts seeking to clarify the mechanisms underlying rheotaxis are scant (Burbano-L and Porfiri, 2021; Chicoli et al., 2015; Colvert and Kanso, 2016; Oteiza et al., 2017), despite experiments on rheotaxis dating back more than a century (Lyon, 1904). A common hypothesis of existing mathematical models is that the presence of the fish does not alter the flow physics with respect to the background flow, thereby neglecting interactions between the fish and the walls of the channel. For example, the model by Oteiza et al., 2017 implements a random walk in a virtual flow, matching experimental measurements of the background flow in the absence of the animal through particle image velocimetry. A similar line of approach was pursued by Burbano-L and Porfiri, 2021 for the study of multisensory feedback control of adult zebrafish.

Thus, according to these models, the fish acts as a perfectly non-invasive sensor that probes and reacts to the local flow environment without perturbing it. There are countless examples in fluid mechanics that could question the validity of such an approximation, from coupled interactions between a fluid and a solid in vortex-induced vibrations (Williamson and Govardhan, 2004) to laminar boundary layer response to environmental disturbances that range from simple decay of the perturbation to bypass transition (Saric et al., 2002). We expect that accounting for bidirectional coupling between the fluid flow and the fish will help clarify many of the puzzling aspects of rheotaxis.

To shed light on the physics of rheotaxis, we formulate a mathematical model based on the paradigm of the finite-dipole, originally proposed by Tchieu et al., 2012. Within this paradigm, a fish is viewed as a pair of point vortices of equal and opposite strength separated by a constant distance in a two-dimensional plane. The application of the finite-dipole has bestowed important theoretical advancements in the study of hydrodynamic interactions between swimming animals (Gazzola et al., 2016; Filella et al., 2018; Kanso and Cheng Hou Tsang, 2014; Kanso and Michelin, 2019; Porfiri et al., 2021), although numerical validation of the framework against full solution of Navier-Stokes equations is lacking – conducting such a validation is also part of this study. Upon validating the dipole model, we investigate the bidirectional coupling between a fish and the surrounding fluid flow in a channel. Our work contributes to the recent literature on minimal models of fish swimming (Gazzola et al., 2014; Gazzola et al., 2015; Sánchez-Rodríguez et al., 2020) that builds on seminal work by Lighthill, 1975, Taylor, 1997 and Wu, 2006 to elucidate the fundamental physical underpinnings of locomotion and inform the design of engineering systems.

We focus on an ideal condition, where fish are deprived of all sensing systems, other than the lateral line that gives them access to information about the flow. Such flow information is coupled, however, to the motion of the fish itself, which acts as an invasive sensor and perturbs the background flow. Just as fish motion influences the local flow field, so too does the local flow field alter fish motion through advection. Predictions from the proposed model are compared against existing empirical observations on fish rheotaxis, compiled through a comprehensive literature review of published work since 1900. Data presented in the literature are used to offer context to the predicted dependence of rheotaxis performance on local flow characteristics, individual fish traits, and lateral line feedback.

Results

Model of the fluid flow

Consider a single fish swimming in an infinitely long two-dimensional channel of width $h$ (Figure 1(a)). Let one wall of the channel be at $y = 0$ and the other at $y = h$ , with $x$ pointing along the channel. The fish position at time $t$ is given by ${\vec{r}}_{f} (t) = x_{f} (t) \hat{i} + y_{f} (t) \hat{j}$ , where $\hat{i}$ and $\hat{j}$ are the unit vectors in the $x$ and $y$ directions, respectively. The orientation of the fish with respect to the $x$ axis is given by $θ_{f} (t)$ (positive counter-clockwise) and its self-propulsion velocity is ${\vec{v}}_{f} = v_{0} (\cos θ_{f} \hat{i} + \sin θ_{f} \hat{j}) = v_{0} {\hat{v}}_{f}$ , where v₀ is the constant speed of the fish and ${\hat{v}}_{f}$ is a unit vector in the swimming direction.

The flow is modeled as a potential flow, which is a close approximation of the realistic flow field around a fish. This simple linear fluid model is intended to capture the mean flow physics, thereby averaging any turbulence contribution. The fish is modeled as a dipole, the potential field of which at some location $\vec{r} = x \hat{i} + y \hat{j}$ is given by

ϕ_{f} (\vec{r}, {\vec{r}}_{f}, θ_{f}) = - r_{0}^{2} (\frac{(\vec{r} - {\vec{r}}_{f}) \cdot {\vec{v}}_{f}}{{|| \vec{r} - {\vec{r}}_{f} ||}^{2}}),

where r₀ is the characteristic dipole length-scale (on the order of the amplitude of the fish tail beating), so that the circulation of each vortex is $2 π r_{0} v_{0}$ . This potential field is constructed assuming a far-field view of the dipole (Filella et al., 2018), wherein r₀ is small in comparison with the characteristic flow length scale, which is satisfied for $ρ = r_{0} / h ≪ 1$ . The velocity field at $\vec{r}$ due to the dipole (fish) is ${\vec{u}}_{f} = \nabla ϕ_{f}$ .

A major contribution of the proposed model is the treatment of the fish as an invasive sensor that both reacts to and influences the background flow, thereby establishing a coupled interaction between the fish and the surrounding environment. A fish swimming in the vicinity of a wall will induce rotational flow near the boundary. In the inviscid limit, this boundary layer is infinitesimally thin and can be considered as wall-bounded vorticity (Batchelor, 2000). Employing the classical method of images (Newton, 2011), the influence of the wall-bounded vorticity on the flow field is equivalent to that of a fictitious fish (dipole) mirrored about the wall plane. For the case of a fish in a channel, this results in an infinite number of image fish (dipoles) (Figure 2), the position vectors for which are

{\vec{r}}_{<, n}^{+} = x_{f} \hat{i} + (y_{f} - 2 (n + 1) h) \hat{j},

{\vec{r}}_{<, n}^{-} = x_{f} \hat{i} + (- y_{f} - 2 n h) \hat{j},

{\vec{r}}_{>, n}^{+} = x_{f} \hat{i} + (y_{f} + 2 (n + 1) h) \hat{j},

{\vec{r}}_{>, n}^{-} = x_{f} \hat{i} + (- y_{f} + 2 (n + 1) h) \hat{j},

Figure 2

Download asset Open asset

Method of images.

Schematic of the fish (black) in the channel (thick lines) and the set of images (gray) needed to generate the channel. The streamlines generated by the fish in an otherwise quiescent fluid are shown in the channel colored by local velocity magnitude (red: high; blue: low). Dashed and solid lines are mirroring planes for the method of images, the pattern for which continues *ad infinitum*.

where $n$ is a non-negative integer representing the $n$ -th set of images. Subscripts “<” and “>” correspond to position vectors of the images at $y < 0$ and $y > h$ , respectively. Likewise, superscript “±” denotes the orientation of the image dipole as $\pm θ_{f}$ ; that is, a position vector with superscript “+” indicates that the associated image has the same orientation as the fish.

The potential function for a given image is found by replacing ${\vec{r}}_{f}$ in Equation 1 with its position vector from Equation 2d and adjusting the sign of $θ_{f}$ in Equation 1 to match the superscript of its vector. The potential field at $\vec{r}$ due to the image dipoles is

ϕ_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) = \sum_{n = 0}^{\infty} (ϕ_{f} (\vec{r}, {\vec{r}}_{<, n}^{+}, θ_{f}) + ϕ_{f} (\vec{r}, {\vec{r}}_{<, n}^{-}, - θ_{f}) + ϕ_{f} (\vec{r}, {\vec{r}}_{>, n}^{+}, θ_{f}) + ϕ_{f} (\vec{r}, {\vec{r}}_{>, n}^{-}, - θ_{f})) .

Thus, the velocity field due to the wall is computed as ${\vec{u}}_{w} = \nabla ϕ_{w}$ , and the overall velocity field induced by the fish is ${\vec{u}}_{f} + {\vec{u}}_{w}$ . (A closed-form expression for the series in terms of trigonometric and hyperbolic functions is presented in Appendix 1) Overall, the presence of the walls distorts the flow generated by the dipole, both compressing the streamlines between the fish and the walls in its proximity and creating long-range swirling patterns in the channel (Figure 2).

The presence of a background flow in the channel is modelled by superimposing a weakly rotational flow,

{\vec{u}}_{b} (\vec{r}) = U_{0} (1 - 4 ϵ {(\frac{y}{h} - \frac{1}{2})}^{2}) \hat{i},

which has speed U₀ at the channel centerline and $U_{0} (1 - ϵ)$ at the walls, $ϵ$ being a small positive parameter. As $ϵ \to 0$ , a uniform (irrotational) background flow is recovered: such a flow is indistinguishable from the one in Figure 2, provided that the observer is moving with the background flow.

For $ϵ ≪ 1$ , the imposed velocity profile approximates that of a turbulent channel flow, wherein a modest degree of velocity profile curvature is present near the channel centerline. We note that this velocity profile does not satisfy the no-slip boundary condition (zero velocity on the walls), and the flow is entirely described by only two parameters (U₀ and $ϵ$ ). For $ϵ ≃ 1$ , the profile approaches that of a laminar flow with parabolic dependence on the cross-stream coordinate. The overall fluid flow in the channel is ultimately computed as $\vec{u} = {\vec{u}}_{f} + {\vec{u}}_{w} + {\vec{u}}_{b}$ .

The circulation in a region $ℛ$ in the flow field centered at some location $y$ is $Γ = \int_{ℛ} ω d A$ , where $ω = (\nabla \times \vec{u}) \cdot \hat{k}$ is the local fluid vorticity ( $\hat{k} = \hat{i} \times \hat{j}$ ). For the considered flow field, we determine

ω (\vec{r}) = \frac{8 U_{0} ϵ}{h} (\frac{y}{h} - \frac{1}{2}),

whereby the irrotational component of the flow field does not contribute to the circulation, and the circulation at a point (per unit area) is equivalent to the local vorticity.

Numerical validation of the dipole model

Despite the success of dipole-based models in the study of fish swimming (Filella et al., 2018; Gazzola et al., 2016; Porfiri et al., 2021; Tchieu et al., 2012), their accuracy against complete Navier-Stokes simulations remains elusive. The potential flow framework in which these models are grounded neglects boundary layers and the resulting wakes that emerge from viscous effects. Quantifying the extent to which these effects influence the flow field generated by the fish is part of this study.

Specifically, computational fluid dynamics (CFD) simulations were conducted to detail the flow field around a fish during steady swimming. The simulation setup was based upon a giant danio of body length $l = 7.3 cm$ in a channel of width $h = 15.0 cm$ . The length of the simulation domain was $L = 50.0 cm$ ( $\sim 6.8 l$ ) with the fish model placed at the channel centerline $20 cm$ ( $\sim 2.7 l$ ) downstream of the inlet. The body undulation of the giant danio was imposed a priori in the simulation, based on data from Najafi and Abtahi, 2022. The time-resolved flow field around the fish was quantified by solving the incompressible Navier-Stokes equations. Details on the setup of the numerical framework and convergence analysis supporting its accuracy are included in Appendix 2.

The mean velocity field averaged over a tail beating cycle is displayed in Figure 3(a). The predominant flow feature observed in the mean field is a flow circulation from the head to the tail of the fish with left-right symmetry and compression of the streamlines near the channel walls. The highest velocity is found at the head of the fish with a thrust wake and recirculation region downsteam of the animal – both at a lower velocity than the anterior flow. The largest velocity in the wake is less than 20% of the peak values recorded ahead of the fish. Additional simulation results are included in Appendix 2. The flow field predicted by the dipole model is in good agreement with numerical simulations, as shown in Figure 3. The dipole model is successful in capturing the circulation from the head to the tail and the compression of the streamlines near the walls. These features are expected to be the main drivers of the interaction between the fish and the channel walls, thereby supporting the value of a dipole model for a first-order analysis of the hydrodynamics of rheotaxis.

Figure 3

Download asset Open asset

Velocity field around a swimming fish from computational fluid dynamics.

(a) Mean velocity field around the steady swimming giant danio relative to the background flow. (b) Velocity field predicted by a dipole with $θ = π$ located at $0.315 l$ from the fish head along its centerline relative to the background flow. The selection of the dipole location and strength is detailed in Appendix 2.

Model of fish dynamics

From knowledge of the fluid flow in the channel, we compute the advective velocity $\vec{U} ({\vec{r}}_{f}, θ_{f})$ and hydrodynamic turn rate $Ω ({\vec{r}}_{f}, θ_{f})$ at the fish location, which encode the influence of the confining walls and background flow on the translational and rotational motion of the fish, respectively. Neglecting the inertia of the fish so that it instantaneously responds to changes in the fluid flow, we determine (Filella et al., 2018)

{\dot{\vec{r}}}_{f} (t) = \vec{U} ({\vec{r}}_{f} (t), θ_{f} (t)) + {\vec{v}}_{f} (θ_{f} (t)),

{\dot{θ}}_{f} (t) = Ω ({\vec{r}}_{f} (t), θ_{f} (t)) + λ ({\vec{r}}_{f} (t), θ_{f} (t)),

where λ is the feedback mechanism based on the circulation measurement through the lateral line.

The advective velocity is found by de-singularizing the total velocity field $\vec{u}$ at $\vec{r} = {\vec{r}}_{f}$ , which is equivalent to calculating the sum of the velocity due to the walls and the background flow in correspondence of the fish (Milne-Thomson, 1996)

\begin{matrix} \vec{U} ({\vec{r}}_{f}, θ_{f}) = {({\vec{u}}_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) + {\vec{u}}_{b} (\vec{r})) |}_{\vec{r} = {\vec{r}}_{f}} \\ = - \frac{π^{2} v_{0} ρ^{2}}{12} [(1 + 3 \csc^{2} (\frac{π y_{f}}{h})) \cos θ_{f} \hat{i} - (1 - 3 \csc^{2} (\frac{π y_{f}}{h})) \sin θ_{f} \hat{j}] \\ + U_{0} (1 - 4 ϵ {(\frac{y_{f}}{h} - \frac{1}{2})}^{2}) \hat{i} . \end{matrix}

Equation (7) indicates that the walls have a retarding effect on the swimming speed of the fish that increases in magnitude the closer the fish gets to either wall of the channel. A fish swimming with orientation $θ_{f} = 0$ at the center of the channel, for example, will swim with velocity ${\dot{\vec{r}}}_{f} (t) = v_{0} (1 - (π^{2} / 3) ρ^{2}) \hat{i} + U_{0} \hat{i}$ . This effect should not be mistaken as traditional viscous drag, which is not included in potential flow theory; rather, it should be intended as the impact of nearby solid boundaries.

Hydrodynamic turn rate is incorporated by considering the difference in velocity experienced by the two constituent vortices comprising the dipole, namely,

\begin{array}{ll} Ω ({\vec{r}}_{f}, θ_{f}) & = - {\hat{v}}_{f} \cdot [\nabla ({\vec{u}}_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) + {\vec{u}}_{b} (\vec{r})) |_{\vec{r} = {\vec{r}}_{f}}] {\hat{v}}_{f}^{⊥} \\ = - \frac{π^{3} ρ^{2} v_{0}}{4 h} \cot (\frac{π y_{f}}{h}) \csc^{2} (\frac{π y_{f}}{h}) \cos θ_{f} + \frac{8 U_{0} ϵ}{h} (\frac{y_{f}}{h} - \frac{1}{2}) \cos^{2} θ_{f}, \end{array}

where ${\hat{v}}_{f}^{⟂} = \hat{k} \times {\hat{v}}_{f}$ ; see Materials and methods section for the mathematical derivation. Equation (8) indicates that interaction with the walls causes the fish to turn towards the nearest wall; for example, a fish at $y_{f} = 3 / 4 h$ , will experience a turn rate due to the wall of $(π^{3} ρ^{2} v_{0}) / (2 h) \cos θ_{f}$ , such that it will be rotated counter-clockwise if swimming downstream and clockwise if swimming upstream. On the other hand, the turning direction imposed by the background flow is always positive (counter-clockwise) in the right half of the channel and negative (clockwise) in the left half, irrespective of fish orientation, so that a fish at $y_{f} = 3 / 4 h$ will always be rotated counter-clockwise. As a result, the fish may turn towards or away from a wall, depending on model parameters and orientation.

Based on experimental observations and theoretical insight (Burbano-L and Porfiri, 2021; Oteiza et al., 2017), we hypothesize that hydrodynamic feedback, that is, lateral line measurements of the surrounding fluid that fish can employ to navigate the flow, is related to the measurement of the circulation in a region surrounding the fish. This hypothesis is supported by experimental evidence presented in Oteiza et al., 2017, which indicated the ability of fish to sense variations in the local velocity gradients to perform rheotaxis. Therein, the authors also found that partial ablation of the lateral line leads to the loss of the rheotactic behavior, hinting at the importance of the flow information on both sides of the fish body for the estimation of the flow circulation.¹¹1 Temporal fluctuations, as in turbulent flows, are neglected in our model. As such, we assume time-averaged circulation sensing from the lateral line. We consider a rectangular region $ℛ$ of width r₀ along the fish body length $l$ . For simplicity, we assume a linear feedback mechanism, $λ ({\vec{r}}_{f}, θ_{f}) = K Γ ({\vec{r}}_{f}, θ_{f}),$ where we made evident that circulation is computed about the fish location and $K$ is a non-negative feedback gain. Assuming that the fish size is smaller than the characteristic length scale of the flow, we linearize the vorticity along the fish in Equation 5 as $ω (\vec{r}) \approx ω ({\vec{r}}_{f}) + \nabla ω ({\vec{r}}_{f}) \cdot {\hat{v}}_{f} Δ l$ . By computing the integral from $Δ l = - l / 2$ to $l / 2$ , we obtain

λ ({\vec{r}}_{f}, θ_{f}) = K r_{0} l \frac{8 U_{0} ϵ}{h} (\frac{y_{f}}{h} - \frac{1}{2}) .

Compared to established practice for modeling fish behavior in response to visual stimuli (Calovi et al., 2014; Couzin et al., 2005; Gautrais et al., 2009; Zienkiewicz et al., 2015), the proposed model introduces nonlinear dynamics arising from the bidirectional coupling between the motion of the fish and the flow physics in its surroundings. We note that the employed feedback in Equation 9 neglects additional potential sensing mechanisms, including vision (Lyon, 1904), acceleration sensing through the vestibular system (Pavlov and Tjurjukov, 1995), and pressure sensing through sensory afferents in the fins (Hardy et al., 2016), which might enhance the ability of fish to navigate the flow.

Analysis of the planar dynamical system

Given that the right hand side of equation set Equation 6a and Equation 6b is independent of the streamwise position of the fish, the equations for the cross-streamwise motion and the swimming direction can be separately studied, leading to an elegant nonlinear planar dynamical system. We center the cross-stream coordinate about the center of the channel and non-dimensionalize it with respect to $h$ , introducing $ξ = y_{f} / h - 1 / 2$ . The governing equations become

\dot{ξ} = [1 - \frac{π^{2} ρ^{2}}{12} (3 \csc^{2} (π (ξ + \frac{1}{2})) - 1)] \sin θ_{f},

{\dot{θ}}_{f} = - \frac{π^{3} ρ^{2}}{4} \cot (π (ξ + \frac{1}{2})) \csc^{2} (π (ξ + \frac{1}{2})) \cos θ_{f} + 8 α ξ (\cos^{2} θ_{f} + κ),

where we non-dimensionalized by the time needed for the fish to traverse the channel in the absence of a background flow, that is, $h / v_{0}$ , and introduced $α = U_{0} ϵ / v_{0}$ and $κ = K r_{0} l$ (see Materials and methods section for estimation of these parameters from experimental observations).

In search of the equilibria of the dynamical system, we note that swimming downstream or upstream ( $θ_{f} = 0$ and π, respectively) solves Equation 10a for any choice of the cross-stream coordinate, the value of which is determined from the solution of Equation 10b for the corresponding orientation $θ_{f}$ . In the case of downstream swimming, the only solution of the resulting transcendental equation is $ξ = 0$ . For upstream swimming, depending on the value of the parameter $β = (α (1 + κ)) / ρ^{2}$ , we have one or three solutions: if $β < β^{*} = π^{4} / 32$ , the only solution is $ξ = 0$ , otherwise, in addition to $ξ = 0$ , there are two solutions symmetrically located with respect to the centerline that approach the walls as $β \to \infty$ (Figure 4(a), see Materials and methoods section for mathematical derivations).

Figure 4

Download asset Open asset

Qualitative dynamics of equation set Equation 10a.

(a) Cross-stream equilibria for upstream swimming as a function of β. (**b,c**) Phase plot for downstream and upstream swimming in the case $α = 0.1$ , $ρ = 0.1$ , and $κ = 1$ , so that $β = 20$ . In all panels, red refers to unstable equilibria and green to stable equilibria.

Local stability of these equilibria is determined by studying the eigenvalues of the state matrix of the corresponding linearized dynamics. For all the considered dynamics, the trace of the state matrix is zero, so that the equilibria can be saddle points (unstable) or neutral centers (stable), if the determinant is negative or positive, respectively (Bakker, 1991; see Materials and methods section for mathematical derivations). In the case of downstream swimming, the determinant is always negative, such that the equilibrium $(θ_{f} = 0, ξ = 0)$ is a saddle point (Figure 4(b)). For upstream swimming, the equilibrium $(θ_{f} = π, ξ = 0)$ is stable if $β > β^{*}$ , leading to periodic oscillations similar to experimental observations (Bak-Coleman et al., 2013; Bak-Coleman and Coombs, 2014; Elder and Coombs, 2015; Figure 1(b)); the other two equilibria located away from the centerline are always unstable (Figure 4(b and c)). Oscillations about the centerline during rheotaxis have a radian frequency $ω_{0} ≃ (π^{2} / 2) ρ \sqrt{β / β^{*} - 1}$ , such that the frequency increases with the square root of β and is zero at $β^{*}$ (see Materials and methods section for the mathematical derivation).

Discussion

There is overwhelming evidence that fish can negotiate complex flow environments by responding to even small flow perturbations (Liao, 2007). However, seldom are these perturbations included in mathematical models of fish behavior, which largely rely on vision cues (Calovi et al., 2014; Couzin et al., 2005; Gautrais et al., 2009; Zienkiewicz et al., 2015). In this paper, we proposed a hydrodynamic model for the bidirectional coupling between fish swimming and fluid flow in the absence of any sensory input but lateral line feedback – encapsulated by a simple linear feedback mechanism. The model reduces to a nonlinear planar dynamical system for the cross-stream coordinate and orientation, of the kind that are featured in nonlinear dynamics textbooks for their elegance, analytical tractability, and broad physical interest (Sastry, 2013).

The planar system anticipates several of the surprising features of rheotaxis. In particular, this study provides some potential answers to the question raised by Coombs et al., 2020: “…what role, if any, do passive (e.g. wind vane) mechanisms play in rheotaxis and how are these influenced by fish factors (e.g. body shape) and flow dynamics?” Through the mathematical analysis of the model, we uncovered an equilibrium at the channel centerline for upstream swimming whose stability is controlled by a single non-dimensional parameter that summarizes flow speed, lateral line feedback, flow gradient, channel width, and fish size. Above a critical value of this parameter, the model predicts that rheotaxis is stable and fish will begin periodic cross-stream sweeping movements whose amplitude can be as large as the channel width. Interestingly, the model anticipates rheotaxis even without sensory feedback, through only passive hydrodynamic mechanisms.

Our mathematical proof of the existence of a nontrivial threshold for β above which upstream swimming becomes stable finds partial support in experimental observations on a number of species in the absence of visual cues (see Appendix 3, where we have performed a bibliographical survey on experimental studies about rheotaxis). Several of these experiments have indicated the existence of a threshold in the flow speed or flow gradient above which fish successfully perform rheotaxis. Importantly, we predict that the presence of channel walls is necessary for the emergence of such a threshold, since for $ρ \to 0$ , $β \to \infty$ , thereby automatically guaranteeing the stability of upstream swimming. Based on our estimation of α and ρ from available data, β can be as small as 10^-1 and exceed 10², thereby encompassing the critical value $β^{*} ≃ 3$ (see Materials and methods section for estimation of model parameters). We should exercise care in drawing comparisons with experiments, which only control for visual feedback, in contrast with the model where we block all sensory modalities except of the lateral line. As reviewed by Coombs et al., 2020, water-motion cues can also be accessible to tactile or other cutaneous senses, beyond the lateral line that is included in our model. In addition, body-motion cues are not limited to visual senses, whereby they can be accessed by tactile and vestibular senses. Hence, a one-to-one comparison between experiments and theory is presently not possible.

The model predicts the emergence of rheotaxis in the absence of any sensory information. Setting $κ = 0$ in our model eliminates hydrodynamic feedback, yet, the fish is able to perform rheotaxis at sufficiently large flow speeds and steep flow gradients. This finding would support the possibility of a completely passive mechanism for rheotaxis. To date, there is no experiment on live animals that can be used to support this claim, owing to the necessity to eliminate all sensory modalities without compromising fish ability to swim. In practice, this may be unfeasible to do. As discussed in Coombs et al., 2020, existing approaches for disabling senses suffer from potential pitfalls, including: (i) unintended effects on the overall fish behavior, which are likely to occur in an effort to block at once vestibular, tactile, and lateral line senses, and (ii) difficulty in guaranteeing complete blockage of a sensory modality, which, like the lateral line, can be distributed throughout the whole body.

A potential line of approach to explore the possibility of a complete passive form of rheotaxis is through experiments with robotic fish ( Duraisamy et al., 2019; Wang et al., 2020; Zhang et al., 2016), mimicking locomotory patterns of live animals and allowing to precisely control sensory input. In this vein, we foresee experiments with robotic fish in a complete open-loop operation that does not utilize any sensory input. The robotic fish developed in our previous study (Kopman and Porfiri, 2013; Kopman et al., 2015) could offer a versatile platform to conduct such an experiment. Such a robot is actuated by a built-in step motor to undertake a periodic tail beating with a predetermined frequency. All its electronics is encased in the frontal section of the robot, so that its size and shape can be readily adjusted though rapid prototyping.

Although free swimming experiments on the robotic fish would be ideal, practicality may suggest to constraint the streamwise and vertical location of the robot while allowing cross-stream motions and heading changes. Such a setup shall also include a load cell to measure the drag on the robot, providing an independent measure to set the tail beat frequency – similar to CFD simulations (see Appendix 2). Also, we would recommend measuring energy expenditure by the robotic fish to gain insight into the hydrodynamic costs and benefits of rheotaxis. Experimental parameters encapsulated in β, including inlet flow speed, channel width, and robot length can be all controlled, and the flow curvature at the centerline can be measured through velocimetry techniques. In the experiments, one shall track the motion of the robotic fish to score conditions in which stable rheotaxis is observed and extract other salient information, such as the frequency of cross-stream sweeping, if present.

The model predicts that increasing κ broadens the stable region, leading to more robust rheotaxis, which is in qualitative agreement with experimental observations – including experiments on animals with intact versus compromised lateral lines (see Appendix 3). The model prediction on the influence of the environment on rheotaxis, including the flow gradient and flow channel size, also parallels the literature on rheotaxis (see Appendix 3). For example, observations by Oteiza et al., 2017 suggest that increasing the flow gradient $ϵ$ enhances hydrodynamic feedback in zebrafish, resulting in improved rheotaxis.

Our analysis also indicates that wider channels should promote rheotaxis by lowering the critical speed above which swimming against the flow becomes a stable equilibrium. This mathematical finding is indirectly supported by experimental observations (see Appendix 3), and bears relevance in the design of experimental protocols for the study of rheotaxis. Confining the subject in a narrow channel will promote bidirectional hydrodynamic interactions with the walls, so that small movements of the animal will reverberate into sizeable changes in the flow physics that will mask the gradient of the background flow. Similarly, in partial alignment with experimental observations, the model predicts a lower threshold for longer fish, owing to a magnification of the hydrodynamic feedback received by a longer body (see Appendix 3). Again, we warn care in drawing comparisons due to the presence of other senses in real experiments, which are not modelled in our work.

Finally, the model anticipates the onset of periodic cross-stream sweeping, which has been studied in some experiments on fish swimming in channels without vision (Coombs et al., 2020). While there is not conclusive experimental evidence regarding the dependence of the frequency of oscillations on flow conditions, the model is in qualitative agreement with experiments by Elder and Coombs, 2015, showing a sublinear dependence on the flow speed. Therein, it is shown that the radian frequency has a weak positive tendency with respect to the flow speed for Mexican tetra swimming with or without cues from the lateral line. Above $2 cm s^{- 1}$ , the animals can successfully perform rheotaxis and display sweeping oscillations at about three cycles per minute and increase to about four cycles per minute at $12 cm s^{- 1}$ . These correspond to a radian frequency on the order of $0.1 rad s^{- 1}$ , which is similar to what we would predict for β ranging from 10⁰ to 10¹ and ρ of the order of $10^{- 1}$ (recall that the time is scaled with respect to time required by the animal to traverse the channel from wall to wall in the absence of a background flow). We acknowledge that the current model does not describe contact and impact with the walls of the channel, which could be important in further detailing the onset of cross-sweeping motions that could involve stick-and-slip at the bottom of the channel (Van Trump and McHenry, 2013).

Just as other minimal models of fish swimming have helped resolve open questions on scaling laws (Gazzola et al., 2014), gait (Gazzola et al., 2015), and drag (Sánchez-Rodríguez et al., 2020), the proposed effort addresses some of the baffling aspects of rheotaxis through a transparent and intuitive treatment of bidirectional hydrodynamic interactions between fish and their surroundings. The crucial role of these bidirectional interactions hints that active manipulation of their surroundings by fish offers them a pathway to overcome sensory deprivation and sustain stable rheotaxis.

The proposed model is not free of limitations, which should be addressed in future research. The current model neglects the elasticity and inertia of the fish, which might reduce the accuracy in the prediction of rheotaxis, especially transient phenomena. Future research should refine the dipole paradigm toward a dynamic, unsteady model that accounts for added mass effects and distributed elasticity, similar to those used in the study of swimming robots (Colgate and Lynch, 2004; Sfakiotakis et al., 1999). The model could also be expanded to account for additional sensory modalities, such as vision, vestibular system, and tactile sensors on the fish body surface. We argue that pursuing any of these extensions shall require detailed experimental data, beyond what the literature can currently offer. Experiments are also needed to refine the linear hydrodynamic feedback mechanism that we hypothesized for the lateral line; in this vein, future experiments could be designed to parametrically vary the flow speed and quantify the activity level of lateral line nerve fibers through neurophysiological recordings (Mogdans, 2019). Beyond the inclusion and refinement of individual sensory modalities, we envision research toward the incorporation of a multisensory framework, as the one introduced by Coombs et al., 2020. Such a framework identifies several motorsensory integration sites in the central nervous system that could contribute to rheotaxis, thereby calling for modeling efforts at the interface of neuroscience and fluid mechanics.

Despite its limitations, the proposed minimalistic model is successful in anticipating some of the puzzling aspects of rheotaxis and points at the possibility of attaining rheotaxis in a purely passive manner, without any sensory input. Most importantly, the model brings forward a potential methodological oversight of laboratory practice in the study of rheotaxis, caused by bidirectional hydrodynamic interactions between the swimming fish and the fluid flow. To date, there is no gold standard for the selection of the size of the swimming domain, which is ultimately chosen on the basis of practical considerations, such as facilitating behavioral scoring and creating a laminar background flow. The model demonstrates that the width of the channel has a modulatory effect on the threshold speed for rheotaxis and the cross-stream swimming frequency, which challenges the comparison of different experimental studies and confounds the precise quantification of the role of individual sensory modalities on rheotaxis. Overall, our effort warrants reconsidering the behavioral phenotype of rheotaxis, by viewing fish as an invasive sensor that modifies the encompassing flow and hydrodynamically responds to it.

Materials and methods

Derivation of the turn rate equation for the fish dynamics

Request a detailed protocol

The expression for the turn rate in Equation 8 is obtained from the original finite-dipole model by Tchieu et al., 2012, in the limit of small distances between the vortices in the pair ( $r_{0} \to 0$ ).

Specifically, eqution (2.11) from Tchieu et al., 2012, adapted to the case of a single dipole reads

{\dot{θ}}_{f} = R e [\frac{(U ({\vec{r}}_{f, r}) - i V ({\vec{r}}_{f, r})) - (U ({\vec{r}}_{f, l}) - i V ({\vec{r}}_{f, l}))}{r_{0}} e^{i θ_{f}}],

where subscript $l$ and $r$ refer to the left and right vortices forming the pair and $\vec{U} = U \hat{i} + V \hat{j}$ is the advective velocity field acting on the dipole. The advective field consists of the interactions with the walls and the background flow, so that $\vec{U} (\vec{r}) = {\vec{u}}_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) + {\vec{u}}_{b} (\vec{r})$ ; in the case of Tchieu et al., 2012, such a field encompasses the velocity field induced by any other dipole in the plane. Left and right vortices are defined so that ${\vec{r}}_{f, l} = {\vec{r}}_{f} + r_{0} {\hat{v}}_{f}^{⟂} / 2$ and ${\vec{r}}_{f, r} = {\vec{r}}_{f} - r_{0} {\hat{v}}_{f}^{⟂} / 2$ , which yields ${\vec{r}}_{f, l} - {\vec{r}}_{f, r} = {\hat{v}}_{f}^{⟂} r_{0}$ .

By carrying out the complex algebra in Equation 11, we determine

{\dot{θ}}_{f} = (\frac{\vec{- U} ({\vec{r}}_{f, l}) + \vec{U} ({\vec{r}}_{f, r})}{r_{0}}) \cdot {\hat{v}}_{f},

which supports the intuition that the dipole will turn counter-clockwise if the right vortex would experience a stronger velocity along the swimming direction. Upon linearizing the term in parenthesis in the neighborhood of ${\vec{r}}_{f}$ , this expression becomes

{\dot{θ}}_{f} = - \nabla \vec{U} (\vec{r_{f}}) {\hat{v}}_{f}^{⊥} \cdot {\hat{v}}_{f} .

The chosen approach is consistent from the standpoint of vortex dynamics, by which each vortex in the pair advects in response to local fluid velocity. In this vein, the fish is interpreted as a bluff body, which rotates according to a difference in the drag experienced by its left and right sides. Such a difference is amplified by the pectoral fins, which enhance the effect of any left-to-right asymmetry in the surrounding fluid flow. In the literature, this description is termed T-dipole, in opposition to the so-called A-dipole that introduces two fiducial points along the direction of motion of the dipole that govern its turning (Kanso and Cheng Hou Tsang, 2014). Whether one representation is superior to the other in terms of accuracy is yet to be clarified; our choice of using a T-dipole is based on its theoretical consistency and intuition on the underlying flow physics. Potential avenues for resolution include detailed CFD simulations of free swimming fish or experiments with robotic fish. For completeness, in Appendix 4, we report model predictions based on the A-dipole.

Determination of the equilibria of the planar dynamical system

Request a detailed protocol

By setting $θ_{f} = 0$ or $θ_{f} = π$ in equation set (Equation 10a, Equation 10b), we determine that ξ should be equal to some constant, which is a root of the following transcendental equation:

\frac{π^{3}}{32} \cot (π (ξ + \frac{1}{2})) \csc^{2} (π (ξ + \frac{1}{2})) = \pm β ξ,

where the positive sign corresponds to $θ_{f} = 0$ and the negative sign to $θ_{f} = π$ . Here, $β = α (1 + κ) / ρ^{2}$ as introduced from the main text.

As shown in Figure 5, for $θ_{f} = 0$ , there is only one root of the equation ( $ξ = 0$ ; see the intersection between the solid red line and the solid black curve), while up to three roots can rise for $θ_{f} = π$ depending on the value of β. For β smaller than a critical value $β^{*}$ , only $ξ = 0$ is a solution (see the intersection between the solid blue line and the solid black curve), while for $β > β^{*}$ two additional solutions, symmetrically located with respect to the origin emerge (see the intersections between the dashed blue line and the solid black curve). The critical value $β^{*}$ is identified by matching the slope of the black curve at $ξ = 0$ , so that $β^{*} = π^{4} / 32$ . Notably, the two solutions symmetrically located with respect to the centerline approach the walls as $β \to \infty$ .

Figure 5

Download asset Open asset

Visual illustration of the process of determining the roots of Equation 14.

(a) Plot of the function $\frac{π^{3}}{32} \cot (π (ξ + \frac{1}{2})) \csc^{2} (π (ξ + \frac{1}{2}))$ (black), superimposed with three lines of different slope: 200 (red), -200 (dashed blue), and -2 (solid blue). (b) Zoomed-in view of the curves in (a) showing that the blue line can only intersect the black curve at the origin.

Local stability analysis of the planar dynamical system

Request a detailed protocol

To examine the local stability of the equilibria of the planar dynamical system, we linearize equation set (Equation 10a, Equation 10b). The state matrix of the linearized dynamics, $A$ , describes the local behavior of the nonlinear system when perturbed in the vicinity of the equilibrium, that is,

\dot{δ q} (t) = A δ q (t),

where $δ q = [δ ξ, δ θ_{f}]^{T}$ is the variation about the equilibrium. The eigenvalues of the $A$ are indicative of local stability about each equilibrium.

For $θ_{f} = 0$ and $ξ = 0$ , the state matrix is given by

A = [\begin{matrix} 0 & 1 - \frac{π^{2} ρ^{2}}{6} \\ 8 (1 + κ) α + \frac{π^{4} ρ^{2}}{4} & 0 \end{matrix}] .

Given that the trace of the matrix is zero ( $tr A = 0$ ), the analysis of the stability of the equilibrium resorts to ensuring the sign of the determinant to be positive ( $det A > 0$ ). Specifically, if the determinant is positive, the eigenvalues are imaginary and the equilibrium is a neutral center (stable, although not asymptotically stable), otherwise one of the eigenvalues is positive and the equilibrium is a saddle point (unstable) (Bakker, 1991). Hence, stability requires that

\frac{1}{24} (- 6 + π^{2} ρ^{2}) (32 α (1 + κ) + π^{4} ρ^{2}) > 0.

Since the first factor is always negative ( $ρ ≪ 1$ ) and the second is positive, the inequality is never fulfilled and the equilibrium is a saddle point (unstable) (Figure 4(a and b)).

For $θ_{f} = π$ and $ξ = 0$ , the state matrix is given by

A = [\begin{matrix} 0 & - 1 + \frac{π^{2} ρ^{2}}{6} \\ 8 (1 + κ) α - \frac{π^{4} ρ^{2}}{4} & 0 \end{matrix}] .

Similar to the previous case, stability requires that $det A > 0$ , that is,

\frac{1}{24} (- 6 + π^{2} ρ^{2}) (- 32 α (1 + κ) + π^{4} ρ^{2}) > 0.

Due to the sign change in the first summand appearing in the second factor with respect to the previous case, stability becomes possible. Specifically, the equilibrium is a neutral center (stable) for $β > β^{*} = π^{4} / 32$ , which is also the necessary condition for the existence of the two equilibria symmetrically located with respect to the channel centerline (Figure 4(a and c)).

When $β > β^{*}$ , we register the presence of two more equilibria at $\pm ξ \neq 0$ . The state matrix takes the form

A = [\begin{matrix} 0 & - 1 - \frac{π^{2} ρ^{2}}{12} + \frac{1}{4} π^{2} ρ^{2} \sec^{2} (π ξ) \\ 8 (1 + κ) α - \frac{1}{4} π^{4} ρ^{2} (2 - \cos (2 π ξ)) \sec^{4} (π ξ) & 0 \end{matrix}],

Also in this case, stability requires that $det A > 0$ , that is,

\frac{1}{48} (- 12 + 3 π^{2} ρ^{2} \sec^{2} (π ξ) - π^{2} ρ^{2}) (- 32 α (1 + κ) + π^{4} ρ^{2} (2 - \cos (2 π ξ)) \sec^{4} (π ξ)) > 0

Once again, for $ρ ≪ 1$ , we can assume that the first factor in parenthesis is negative. (This assumption is grounded upon Equation 14, which yields that $(ξ \pm 1 / 2) = O (ρ^{2 / 3})$ ; since $\cos (π ξ)^{2} = O ((ξ \pm 1 / 2)^{2})$ , we have that $ρ^{2} \sec^{2} (π ξ) \to 0$ as $ρ \to 0$ .) Hence, we obtain

β > \frac{π^{4}}{32} (2 - \cos (2 π ξ)) \sec^{4} (π ξ),

which is not satisfied for any choice of $β > β^{*}$ . Thus, the two equilibria away from the channel centerline, close to the walls are always saddle points (unstable) (Figure 4(a and c)).

We comment that the local stability analysis requires only knowledge of the curvature of the flow field at the centerline of the channel. Hence, should one contemplate alternative profiles for the background flow, linear stability results shall not change. Higher-order parameterizations for the flow profile will result into nonlinear dependencies on ξ that do not affect the linear analysis. Likewise, while we considered a linear feedback mechanism to integrate lateral information via a simple gain, one may explore nonlinear relationships between λ and $Γ$ . The linear stability analysis shall not change, whereby these nonlinear forms will result into dependencies on higher powers of ξ.

Frequency of cross-stream sweeping

Request a detailed protocol

The linearized planar system about the stable focus in Equation 18 is equivalent to a classical second-order system in terms of the cross-stream coordinate, similar to a mass-spring model. Hence, the radian resonance frequency of the system is

ω_{0} = \sqrt{det A} ≃ \frac{π^{2}}{2} ρ \sqrt{\frac{β}{β^{*}} - 1} .

where the last approximation holds for $ρ ≪ 1$ . Equation (23) shows that, close to the threshold, the frequency of oscillations is small and it increases with β and ρ.

Estimation of model parameters

Request a detailed protocol

In a typical experimental setup on rheotaxis, the width of the channel, $h$ , is on the order of three to ten times the body length of the animal, $l$ . For example, experiments from Elder and Coombs, 2015 on Mexican tetras of $l = 8.3 cm$ were conducted in a channel with $h = 25 cm$ . Similarly, in the experiments on adult zebrafish from Burbano-L and Porfiri, 2021, $l = 3.6 cm$ and $h = 13.8 cm$ , and in the experiments on zebrafish larvae from Oteiza et al., 2017, $l = 4.2 mm$ (inferred from the animals’ age) and $h = 1.27$ – $4.76 cm$ . The distance between the vortices simulating a fish, r₀, should be on the order of a tail beat, which has a typical value of $0.2 l$ (Gazzola et al., 2014). As a result, it is tenable to assume that $ρ^{2}$ is between $10^{- 4}$ and $10^{- 2}$ .

A safe estimation of the velocity of the animal in the absence of the background flow, v₀, would be on the order of few body lengths per second (Gazzola et al., 2014). The speed used for the background flow across experiments, U₀, tend to be of the same order as the magnitude of v₀, leaning toward values close to one body length per second (Coombs et al., 2020). For instance, data on zebrafish from Burbano-L and Porfiri, 2021 suggest $v_{0} = 5.7 cm s^{- 1}$ and $U_{0} = 3.2 cm s^{- 1}$ . The estimation of the non-dimensional parameter $ϵ$ associated with the shear in the flow is more difficult, since data on the velocity profiles are seldom reported. That being said, for channel flow of sufficiently high Reynolds number, the velocity profile in the channel is expected to be blunt, approximating a uniform flow profile near the channel center (White, 1974). Thus, it is tenable to treat $ϵ$ as a small parameter, between 10^-2 and 10^-1. For flow of low Reynolds number (Oteiza et al., 2017) ( $R e < 100$ ), the velocity gradient in the channel has been observed to be large, corresponding to $ϵ$ values in the range of 10^-1 and 1. By combining these estimations, we propose that α ranges between 0 and 1.

An estimation of κ is difficult to offer, whereby feedback from the lateral line has only been included in few studies (Burbano-L and Porfiri, 2021; Chicoli et al., 2015; Colvert and Kanso, 2016; Oteiza et al., 2017). Using the data-driven model from Burbano-L and Porfiri, 2021, it is tenable to assume values on the order of 10¹ for individuals showing high rheotactic performance. This gain can also be estimated by comparing the threshold speeds of fish, $U_{c}$ , with and without the lateral line, through $\frac{U_{c} (LL -)}{U_{c} (LL +)} = 1 + κ$ , according to Equation 29 in Appendix 3. The significant increase in the threshold speed following lateral line ablation in Baker and Montgomery, 1999 indicates that $κ \in [2, 7]$ , while the indistinguishable threshold speed between LL+ and LL- fish in a few other studies (Bak-Coleman and Coombs, 2014; Elder and Coombs, 2015; Van Trump and McHenry, 2013) may suggest that $κ \sim 0$ . In Table 1, we summarize the model parameters identified from data in the experimental studies detailed in Appendix 3.

Table 1

Estimation of model parameters from data in the literature.

Reference	$ρ$	$ϵ$	$α$	$κ$	$β$
Bak-Coleman et al., 2013	$\sim 0.05$	$[10^{- 2}, 10^{- 1}]$	$[0, 0.17]$	—	—
Bak-Coleman and Coombs, 2014	$\sim 0.04$	$[10^{- 2}, 10^{- 1}]$	$[0, 0.16]$	$\sim 0$	$[0, 100]$
Baker and Montgomery, 1999 and Montgomery et al., 1997	$\sim 0.1$	$[10^{- 2}, 10^{- 1}]$	$^{*} [0, 0.32]$	$[2, 7]$	$[0, 256]$
Elder and Coombs, 2015	$\sim 0.066$	$[10^{- 2}, 10^{- 1}]$	$^{*} [0, 0.24]$	$\sim 0$	$[0, 55]$
Kulpa et al., 2015	$\sim 0.04$	$\sim 1$ near center of jet	$\sim 1.3$	—	—
Oteiza et al., 2017	$[0.018, 0.066]$	$[0.20, 0.82]$	—	—	—
Peimani et al., 2017	$\sim 0.044$	$\sim 1$	—	—	—
Suli et al., 2012	$\sim 0.018$	$[0.1, 1]$	—	—	—
Van Trump and McHenry, 2013	$[0.055, 0.127]$	$[10^{- 2}, 10^{- 1}]$	$^{*} [0, 0.32]$	$\sim 0$	$[0, 106]$

LL+ cavefish swimming speed $v_{0} \sim 5 c m / s$ in zero background flow in Bak-Coleman and Coombs, 2014 is used to estimate α.

Appendix 1

Complete expression for the velocity field caused by image dipoles

The velocity field ${\vec{u}}_{f} = u_{f} \hat{i} + v_{f} \hat{j}$ at $\vec{r}$ induced by the single dipole at ${\vec{r}}_{f}$ , given by the potential function in Equation 1, is

u_{f} (\vec{r}, {\vec{r}}_{f}, θ_{f}) = r_{0}^{2} v_{0} (\frac{({(x - x_{f})}^{2} - {(y - y_{f})}^{2}) \cos θ_{f} + 2 (x - x_{f}) (y - y_{f}) \sin θ_{f}}{{({(x - x_{f})}^{2} + {(y - y_{f})}^{2})}^{2}}),

v_{f} (\vec{r}, {\vec{r}}_{f}, θ_{f}) = r_{0}^{2} v_{0} (\frac{- ({(x - x_{f})}^{2} - {(y - y_{f})}^{2}) \sin θ_{f} + 2 (x - x_{f}) (y - y_{f}) \cos θ_{f}}{{({(x - x_{f})}^{2} + {(y - y_{f})}^{2})}^{2}}) .

The potential function describing the image vortex system for a dipole in a channel presented in Equation 3 can be simplified using Mathematica, yielding

\begin{array}{ll} ϕ_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) & = \frac{r_{0}^{2} v_{0}}{4} [4 \frac{(x - x_{f}) \cos θ_{f} + (y - y_{f}) \sin θ_{f}}{(x - x_{f})^{2} + (y - y_{f})^{2}} \\ - \frac{π e^{- i θ_{f}}}{h} (e^{2 i θ_{f}} (\coth (π A) + \coth (π B^{*})) + \coth (π A^{*}) + \coth (π B))], \end{array}

where $A = ((x - x_{f}) + i (y - y_{f})) / (2 h)$ , $B = ((x - x_{f}) + i (y + y_{f})) / (2 h)$ , $i = \sqrt{- 1}$ , and a superscript * indicates complex conjugate. The velocity field at due to the walls, is

\begin{array}{ll} v_{w} = \frac{r_{0}^{2} v_{0}}{4} [\frac{π^{2} e^{- i θ_{f}}}{2 h^{2}} (e^{2 i θ_{f}} (c s c h^{2} π A + c s c h^{2} π B^{*}) + (c s c h^{2} π A^{*} + c s c h^{2} π B)) \\ + \frac{4 \cos θ_{f}}{(x - x_{f})^{2} + (y - y_{f})^{2}} - \frac{8 (x - x_{f}) ((x - x_{f}) \cos θ_{f} + (y - y_{f}) \sin θ_{f})}{((x - x_{f})^{2} + (y - y_{f})^{2})^{2}}] . \end{array}

\begin{array}{ll} v_{w} = \frac{r_{0}^{2} v_{0}}{4} [\frac{i π^{2} e^{- i θ_{f}}}{2 h^{2}} (e^{2 i θ_{f}} (c s c h^{2} π A - c s c h^{2} π B^{*}) - (c s c h^{2} π A^{*} + c s c h^{2} π B)) \\ + \frac{4 \sin θ_{f}}{(x - x_{f})^{2} + (y - y_{f})^{2}} - \frac{8 (y - y_{f}) ((x - x_{f}) \cos θ_{f} + (y - y_{f}) \sin θ_{f})}{((x - x_{f})^{2} + (y - y_{f})^{2})^{2}}] . \end{array}

Superimposing the velocity fields from the dipole and its images and setting $y = 0$ (or $y = h$ ) yields $v_{f} + v_{w} = 0$ , thereby confirming that the walls of the channel are streamlines.

Appendix 2

Computational fluid dynamics

Framework

To quantify the flow around the swimming fish, we solved the incompressible Navier-Stokes equations numerically in the commercial software COMSOL Multiphysics (version 5.6). We focused on the steady swimming of a giant danio exhibiting a carangiform swimming pattern, consisting of large body undulations in the posterior of the animal and minimal lateral movement in the anterior portion. The lateral movement of a giant danio was mathematically described through a local coordinate system, $x^{'}$ - $y^{'}$ , such that the undeformed fish aligns its centerline along the $x^{'}$ -axis with head at $x^{'} = 0$ and tail at $x^{'} = l$ (see Appendix 2—figure 1a). The lateral movement were described as a traveling wave from the head to the tail as (Najafi and Abtahi, 2022).

y^{'} = [c_{1} (x^{'} - 0.3 l) + c_{2} {(x^{'} - 0.3 l)}^{2}] \sin [k_{L} (x^{'} - 0.3 l) - 2 π f t]

where c₁ and c₂ are two constants that describe the shape of the undulation envelope, $k_{L}$ is the wavenumber, and $f$ is the tail beating frequency. The values of $l$ , c₁, c₂, and $k_{L}$ are taken from Najafi and Abtahi, 2022, and the value of $f$ is taken from Zhang et al., 2019; these parameters are summarized in Appendix 2—table 1. For the chosen tail beat frequency, the oscillatory Reynolds number (Gazzola et al., 2014) is $S w = \frac{2 π C f l}{ν} = 8, 100$ , where $C$ is the tail beat amplitude, given as $C = 0.7 \times c_{1} l + 0.49 \times c_{2} l^{2}$ , and ν is the kinematic viscosity of water at room temperature ( $ν = 1.0 \times 10^{- 6} m^{2} s^{- 1}$ ).

The undeformed giant danio body shape was approximated as a NACA0013 airfoil, which has a maximum thickness-to-chord length matching that of the animal in Zhang et al., 2019. A NACA airfoil offers a reasonable approximation of the fish body in terms of its streamlined shape and its ability to emulate thrust production of a swimming fish (Lucas et al., 2020). The movement of the fish body was imposed in the numerical simulations using a moving boundary. No-slip boundary conditions was set on the fish body, whereas the channel walls satisfied only the no-penetration boundary condition, that is, they were treated as slip walls. In the simulations, we set a uniform flow with speed U₀ at the inlet, imposed a zero pressure boundary condition at the outlet, and fixed the axial location of the fish at the channel centerline.

Appendix 2—figure 1

Download asset Open asset

Details about the implementation of the computational fluid dynamics simulations.

(a) Mesh implemented in the simulations, with definitions of coordinate systems and a zoomed-in view of the refined mesh around the fish. (b) Mesh convergence analysis, showing the mean drag as a function of number of elements in the simulation.

The fluid domain was discretized using triangular elements, as illustrated in Appendix 2—figure 1a, and was allowed to deform in time to accommodate the body undulations. A grid convergence study was performed to ensure the solution was independent of the mesh. As shown in Appendix 2—figure 2b, a total number of $76 k$ elements were sufficient to guarantee the mesh independence of the simulations, such that implementing a higher number of $122 k$ elements introduced negligible variation in the mean drag prediction on the fish. As a result, all simulations were conducted using $76 k$ elements. A total of 10 tail beating cycles were simulated in each simulation, which was sufficiently long for the flow in the channel to fully develop. All simulations were conducted on 12 Intel Xeon Platinum 8268 CPUs with a base frequency of $2.90 GHz$ and a total of $48 GB$ memory. A typical simulation required approximately $6 h$ computational time to complete.

To represent the realistic swimming condition, the imposed inlet flow speed, U₀, should match the fish swimming speed, v₀, which depends on the body undulations described by (27). To identify the value of v₀, we conducted a series of preliminary simulations by varying U₀ until the time-averaged total drag on the fish, ${\bar{F}}_{D}$ was zero. That is, the drag experienced by the fish exactly balanced the thrust generated by the body undulations. The resulting swimming speed was determined to be $v_{0} = 19.74 cm / s$ . The resulting Reynolds number based upon channel width is $Re = \frac{h U_{0}}{ν} = 29, 600$ .

Results

The instantaneous velocity fields in the vicinity of the giant danio relative to the background flow are presented at five time instants during half of a tail beating cycle in Appendix 2—figure 2. In comparison with the mean flow shown in Figure 3, we observe distortions of the streamlines, with a left-right asymmetry caused by the body undulations. We also identify a series of alternating vortices generated through tail beating that form the wake. The magnitude of the wake flow decays as it is advected downstream.

CFD results are also useful to quantify the thickness of the boundary layer along the animal and offer support in favor of the proposed inviscid model for the study of the interactions between a fish and the channel walls. The thickness of the boundary layer can be quantified from the flow velocity component along the $x$ -direction, $u_{t}$ , in the presence of the background flow. As shown in Appendix 2—figure 3, the value of $u_{t}$ is zero at the body surface due to the no-slip boundary condition. A large gradient is observed within a thin layer at the boundary, in which $u_{t}$ rapidly increases to v₀. The boundary layer thickness can be estimated by identifying the location at which $u_{t}$ reaches 99% of its asymptotic value away from the fish body. During a tail beating cycle, the boundary layer thickness ranged between 2.8% and 15.7% of the fish body length, thereby supporting the feasibility of neglecting the influence of the viscous boundary layers as a first-order approximation.

Appendix 2—figure 2

Download asset Open asset

Instantaneous velocity fields around a swimming fish relative to the background flow from computational fluid dynamics.

(a) – (e) correspond to $t = 0$ , $T / 8$ , $T / 4$ , $3 T / 8$ , and $T / 2$ , respectively, where $T$ is the period of a tail beat. White curves are streamlines with arrows indicating flow directions.

Appendix 2—figure 3

Download asset Open asset

Analysis of the boundary layer thickness along a swimming fish.

x-component of the flow velocity, $u_{t}$ , extracted across the half-length of the fish body. The values of $u_{t}$ are measured in a coordinate system moving at the speed of v₀ along the fish swimming direction.

Comparison with the dipole model

The velocity field predicted by the dipole model is validated against the mean velocity quantified through CFD (see Figure 3(a)). Consistent with the fish location and heading direction, we set $y_{f} = h / 2$ and $θ = π$ for the dipole. The axial location, x_f, and the characteristic length scale, r₀, of the dipole are treated as fitting parameters. The velocity field associated with the dipole model in Figure 3(b) corresponds to a set of fitting parameters that minimizes the discrepancy between the model prediction and the simulated simulation. Denoting the model-predicted velocity field as ${\vec{u}}^{dipole}$ and the simulated one as ${\vec{u}}^{CFD}$ , the discrepancy between them is quantified through

E = \sqrt{\frac{\int_{D} {| {\vec{u}}^{d i p o l e} - {\vec{u}}^{C F D} |}^{2} d A}{\int_{D} d A}},

where $D$ is the computational domain. The optimal values of the fitting parameters that minimize $E$ are determined through an exhaustive search of the parameter space, leading to $(x_{f} - x_{s}) = 0.315 l$ and $r_{0} = 0.062 l$ , where x_s is the axial location of the fish head (the origin of the $x^{'}$ - $y^{'}$ coordinate system).

Appendix 2—table 1

Parameters employed in (27) to describe the locomotory pattern of a giant danio.

Parameters	$l$	$c_{1}$	$c_{2}$	$k_{L}$	$f$
Values	$7.3 cm$	0.004	$- 2.33 m^{- 1}$	$78.5 m^{- 1}$	$3 Hz$

Appendix 3

Bibliographical survey

We surveyed over three hundred publications cited by Arnold, 1974 and Coombs et al., 2020 – two review papers on rheotaxis, with the former focusing on early investigations from 1900 s to 1970 s, and the latter highlighting more recent works conducted between 1970 s and 2020. Publications were selected through the following inclusion and exclusion criteria.

Criteria

Inclusion criteria

We selected studies where: (i) the subject animals were fish; (ii) fish demonstrated rheotactic behavior; (iii) no unsteady flow events were present in the swimming domain, such as the wake structure of obstacles; (iv) the sensory cues available to fish could be identified with some confidence; (v) fish behavior was not influenced by social interactions; (vi) fish swam without visual cues; and (vii) the publication was written in English. Within criterion (iii), we focused on experiments with steady flows where the flow gradient is consistent over time, thus excluding swimming in random flow events. Criterion (vi) was introduced to direct our search toward the effects of hydrodynamic cues and lateral line sensing, which limited our search to experiments using blind fish or experiments in the dark. We acknowledge that this condition is not reflective of the model hypotheses, which, in fact, block any sense except for the lateral lines (such as vestibular and tactile senses).

Exclusion criteria

Among studies identified through the selection criteria, we excluded experiments on pleuronectiform flatfishes, which swim on their side and generate propulsive undulations in a vertical plane (Webb, 2002). This locomotory pattern differs fundamentally from the current model, derived on the assumption the fish align their bodies vertically and undulate on a horizontal plane, which is the swimming strategy of the majority of fishes.

Dataset

Appendix 3—table 1 presents data extracted from the selected studies, including the fish species, size of the swimming domain, flow conditions, sensory cues available to the fish, and the measured rheotaxis threshold speed. Swimming domains with rectangular cross-sections are defined by their length ( $L$ ), width ( $h$ ), and depth ( $W$ ), while swimming domains with circular cross-sections by their length and diameter ( $D$ ). Flow conditions are quantified through the flow speed and flow gradient. If information about the flow gradient was not available in a study, we qualitatively estimated its value through the Reynolds number of the flow, defined based on the width (diameter) of the channel and the background flow speed as $Re = \frac{h U_{0}}{ν}$ ( $Re = \frac{D U_{0}}{ν}$ ). For a sufficiently high $Re$ , the flow gradient near the center of the channel is expected to be low.

Comparison against model predictions

The studies identified in Appendix 3—table 1 are utilized to offer some context to the proposed theoretical framework with respect to the rheotaxis stability threshold. We first express the stability threshold $β = β^{*}$ in dimensional form in terms of the rheotaxis threshold speed

U_{c} = \frac{π^{4} r_{0}^{2} v_{0}}{32 h^{2} (1 + K r_{0} l) ϵ},

such that $U_{0} > U_{c}$ corresponds to the stable condition $β > β^{*}$ , and vice versa. From most studies, the values of $U_{c}$ can be identified and its confidence level can be inferred (see Appendix 3—table 1).

As evidenced in Equation 29, a series of parameters could influence the rheotaxis threshold speed, including the lateral line feedback, flow gradient, swimming domain size, and fish body length. Specifically, Equation 29 predicts that increasing the lateral line feedback, flow gradient, and/or width of swimming channel promotes rheotaxis at lower flow speeds, whereas increasing fish size will require higher flow speeds to elicit rheotactic behavior. The effects of these parameters are validated independently in Appendix 3—table 2, where we include experimental evidence garnered within each study and, when possible, carry out a comparison, across them. We compare each of these empirical observations to model predictions, and assess if they support the model, contradict the model, or are inconclusive. An observation is considered supportive of (contradictive to) our model if the measured $U_{c}$ exhibits with statistical significance the same (opposite) dependence on a certain parameter. Data that lack statistical significance are considered inconclusive.

The confidence intervals of the measured $U_{c}$ values were estimated to determine if $U_{c}$ were significantly different across studies. When the mean and standard error of the mean (s.e.m.) of $U_{c}$ were provided, we estimated its 95% confidence interval as $(mean - 1.96 s . e . m ., mean + 1.96 s . e . m .)$ . If the confidence intervals of two $U_{c}$ values did not overlap, we considered them significantly different. For instance, in Bak-Coleman and Coombs, 2014 and Elder and Coombs, 2015, the confidence intervals of $U_{c}$ were determined to be $(0.63, 1.17) cm / s$ and $(1.27, 2.64) cm / s$ , respectively, and thus the $U_{c}$ values were considered significantly different. In several studies, such as Baker and Montgomery, 1999 and Van Trump and McHenry, 2013, the threshold speeds were only estimated as intervals, where fish swimming below a lower bound did not perform rheotaxis, while they exhibited rheotaxis above an upper bound. We treated this speed interval as the confidence interval for $U_{c}$ in our statistical analysis.

Effect of lateral line feedback

Several studies provide some support in favor of the prediction of our model of the beneficial role of lateral line feedback, showing a significant reduction in rheotactic performance when the lateral line is compromised (Kulpa et al., 2015; Oteiza et al., 2017; Suli et al., 2012), see Appendix 3—table 2. In these studies, fish locomotion was measured in steady background flows, so that a fish holding station would experience minimal linear acceleration and marginally engage the vestibular system. Throughout these studies, fish were not observed to make contact with the swim channel, indicating that tactile senses played a negligible role in rheotaxis.

Effect of flow gradient

We identified two studies (Lyon, 1904; Oteiza et al., 2017) that could back the predicted effect of the flow velocity gradient on rheotaxis, as summarized in Appendix 3—table 2. In both studies, fish locomotion was recorded in flows with varying velocity gradients. In qualitative agreement with the proposed model, the rheotaxis performance of zebrafish larvae significantly improved with increasing gradient magnitudes (Oteiza et al., 2017). Similar observations were obtained by Lyon on blind Fundulus (Lyon, 1904), where in a flow with a small gradient, fish performed rheotaxis only when tactile cues were available, while in a jet flow with a large flow gradient, rheotaxis could be elicited solely by the flow. Although qualitatively in line with our predictions, we conservatively considered this study as inconclusive due to a lack of quantitative data for statistical tests.

Effect of size of swimming domain

To elucidate the role of the swimming domain size on rheotaxis threshold speed, we conducted cross-study comparisons as shown in Appendix 3—table 2. As evidenced through comparisons between two experiments on zebrafish larvae (Oteiza et al., 2017; Peimani et al., 2017) in swim channels of drastically different sizes, rheotaxis was elicited at a higher threshold speed in a smaller flow channel, which supports our model prediction. Our model is also qualitatively supported by comparisons between Bak-Coleman and Coombs, 2014 and Baker and Montgomery, 1999, or Bak-Coleman and Coombs, 2014 and Van Trump and McHenry, 2013, where experiments on blind cavefish of comparable body sizes uncovered higher threshold speeds in smaller flow channels. In the experiments of Bak-Coleman and Coombs, 2014, blind cavefish were observed to receive transient tactile senses while swimming, which could have contributed to its lower rheotaxis threshold. As a result, experimental data on blind cavefish were conservatively deemed to be inconclusive.

Effect of fish size

The relationship between the threshold speed and fish body size is less straightforward, as the body size not only determines the value of $l$ , but also influences r₀, which is on the order of the fish tail beat amplitude. We assumed $r_{0} = 0.2 l$ , which is a typical tail-beat-amplitude-to-body-length ratio (Gazzola et al., 2014). For fish with functional lateral lines that produce positive feedback, $K > 0$ , we obtain $U_{c} \sim 1 - \frac{1}{1 + 0.2 K l^{2}}$ . For fish with disabled lateral line, $K = 0$ , we find $U_{c} \sim l^{2}$ . In both cases, the model predicts that $U_{c}$ is larger for fish with larger body length. Some evidence can be garnered by contrasting a pair of studies by Bak-Coleman and Coombs, 2014 and Elder and Coombs, 2015, where experiments were conducted on fish of the same species (Astyanax mexicanus) in swim tunnels of the same size, and tested in flows at a similar range of speeds. The high flow speeds in both studies suggest that the flow gradients in these experiments were small. We assume that the lateral line feedback were equivalent in both studies, as the subjects were conspecific. Although the tactile cues present in the experiments by Bak-Coleman and Coombs, 2014 hinder our ability to reach a definitive conclusion on the effect of fish body size, the higher threshold speed observed in larger fish qualitatively supports our model prediction. The paucity of data for validation of the effect of fish size is a result of a lack of studies with matching experimental conditions, including dimensions of the flow facilities, flow conditions, and functionality of the lateral line.

In summary, we identified a total of five sets of experiments in support of our model, and nine sets of studies that offer inconclusive evidence. None of the data contradicted predictions from the proposed model.

Most experiments used in our comparison listed in Appendix 3—table 1 were conducted in the past 25 years, and only two studies date back to before 1970. This disparity is attributed to an evolution of the methodologies for the study of rheotaxis over time. Among the earlier efforts, a large portion relied on observations of fish behavior in the field (Arnold, 1974). Although these studies minimized the introduction of external stimuli stemming from human presence and unfamiliar environments that could alter the behavior of fish in the wild, a lack of flexibility in the design of controlled experiments in the field, together with an insufficient measurement resolution, has led to only a limited number of works that could distinguish the impact of one sensory cue from another. As a result, a large number of earlier efforts do not meet our inclusion criteria.

Likely, the interest in fish rheotaxis was recently reignited owning to the advancements in technologies that could selectively deactivate specific fish sensory organs, thereby allowing for the targeted investigation of the role of each sensory cue in rheotaxis. For instance, pharmacological methods that could disable the lateral line led to studies (Baker and Montgomery, 1999; Montgomery et al., 1997) challenging the long-standing perception that the lateral line could not mediate rheotaxis. In addition, the development of high speed cameras with infra-red sensing capabilities enabled precise measurements of fish behavior in the dark, allowing for the elimination of visual cues from the study of lateral line functionality in rheotaxis.

Some early experiments that have been considered in the past as evidence against the role of the lateral line are not listed in Appendix 3—table 1 due to a lack of a controlled experimental design in the field setting. For example, some species of salmonids, including salmon and trout, were observed to swim against the current in the day and rest on the bottom of a stream at night (Davidson, 1949; Edmundson et al., 1968; Gibson, 1966), leading to a conclusion that the lateral line played a minimal role in rheotaxis (Arnold, 1974). However, we did not include these experiments for our model validation due to confounding factors posed by the field settings, such as variations in water temperature (Edmundson et al., 1968; Fraser et al., 1993; Needham and Jones, 1959) and current speeds at different hours of the day. Daily fluctuations in the availability of food (Waters, 1962; Elliott, 1965) is another factor that could influence the activity levels of fish at night, as observed in white bass (McNaught and Hasler, 1961) and trout (Elliott, 1965). Another class of experiments that led to the previous rejection of lateral line was the demonstration of the imperative role of vision in rheotaxis. Experiments on salmon (Hoar, 1954) and herring (Brawn, 1960) showed a reduction in rheotaxis when vision was obscured in the dark or in muddy water, conflating the role of visual cues in rheotaxis. Again, these observations do not directly contradict the proposed model, which suggests that in the absence of visual cues, rheotaxis could still manifest provided that the flow speed is sufficiently high.

Appendix 3—table 1

Relevant publications on fish rheotaxis in the absence of visual cues, identified through literature review.

Reference	Fish		Swimming domain	Flow properties		*Sensory cues	Rheotaxis threshold speed
Reference	Species	Length	Swimming domain	Flow speed	Flow gradient	*Sensory cues	Rheotaxis threshold speed
Bak-Coleman et al., 2013	Giant danio (Devario aequipinnatus)	6.0 –7.3 cm	Flow tank of $25 \times 25 \times 25 c m$ ( $L \times h \times W$ )	0, 3, and $7 cm / s$	$Re \sim 7500$ at LL+ threshold speed; flow gradient expected to be small near center of tank	LL+/LL-	$\leq 3 cm / s$
Bak-Coleman and Coombs, 2014	blind cavefish (Astyanax mexicanus)	4.2 –5.0 cm	Flow tank of $25 \times 25 \times 25 c m$ ( $L \times h \times W$ )	0, 1, 2, 3, 4, 7 and $8 cm / s$	$Re \sim 2000$ at LL+ threshold speed; flow gradient expected to be small near center of tank	LL+/LL-; fish made transient contacts with substrate	LL+: $0.90 c m / s$ ; LL-: $0.54 c m / s$
^†Baker and Montgomery, 1999	blind cavefish (Astyanax fasciatus)	4 –7 cm	Flow tank of $51 \times 9 \times 20 c m$ ( $L \times h \times W$ )	0, 2, 3, 5, 9 and $16 cm / s$	$Re \sim 2000$ at LL+ threshold speed; flow gradient expected to be small near center of tank	LL+/LL-; tactile senses	LL+: 2– $3 cm / s$ ; LL-: 9– $16 cm / s$
Elder and Coombs, 2015	Mexican tetras (Astyanax mexicanus)	$8.3 cm$	Flow tank of $25 \times 25 \times 25 c m$ ( $L \times h \times W$ )	0, 1, 2, 4, 7, and $12 cm / s$	$Re \sim 5000$ at threshold speed; flow gradient expected to be small near center of tank	LL+/LL-	$\sim 2 cm / s$ for LL+ and LL-
Kulpa et al., 2015	blind cavefish (Astyanax mexicanus)	4.4 –5.3 cm	Flow tank of $25 \times 25 \times 10 c m$ ( $L \times h \times W$ )	Maximum speed of $8 cm / s$	Jet flow across center of tank; flow gradient expected to be large	LL+/LL-	$\leq 8 cm / s$
^‡Lyon, 1904	blind Fundulus	unspecified	Trough with unspecified dimensions; tideway leading to pond	“not too strong current” in trough and current with “more or less eddy and irregularity” in tideway	Flow gradient expected to be small	LL+; some fish gained tactile senses	Not measured; rheotaxis elicited only by tactile cues
^‡Lyon, 1904	blind Fundulus	unspecified	Trough with unspecified dimensions	flow “gushing rather violently”	Jet flow; flow gradient expected to be large	LL+	Not measured; rheotaxis elicited by flow
^†Montgomery et al., 1997	blind cavefish (Astyanax fasciatus)	4 –7 cm	^§Flow tank of $51 \times 9 \times 20 c m$ ( $L \times h \times W$ )	0, 2, 3, 5, 9, and $16 cm / s$	$Re \sim 2000$ at LL+ threshold speed; flow gradient expected to be small near center of tank	LL+/LL-; tactile senses	LL+: 2–3 cm/s; LL-: 9–16 cm/s
Oteiza et al., 2017	zebrafish (Danio rerio) larva 5–7 days post fertilization (dpf)	unspecified	13 cm-long circular tube with diameter 1.27– $4.76 c m$	0.2–0.8 cm/s	Low to high flow gradients identified through particle image velocimetry	LL+/LL-	LL+: rheotaxis observed as low as $0.2 cm / s$
Peimani et al., 2017	zebrafish (Danio rerio) larva 5–7 dpf	estimated $\sim 0.35 c m$	Flow channel of $63.3 \times 1.6 \times 0.55 m m$ ( $L \times h \times W$ )	0.95–3.8 cm/s	$Re \sim 10$ at threshold speed; flow gradient expected to be large	LL+	$0.95 cm / s$
Suli et al., 2012	zebrafish (Danio rerio) larva 5 dpf	$\sim 0.33 cm$	Flume of $110 \times 3.7 \times 2.8 c m$ ( $L \times h \times W$ )	0.075, 0.15, $0.2 cm / s$	$R e < 75$ ; flow gradient expected to be large	LL+/LL-	Not quantified
Van Trump and McHenry, 2013	blind Mexican cavefish (Astyanax fasciatus)	3 –7 cm	Cylindrical channel of $150 \times 11 c m$ ( $L \times D$ )	0, 1, 2, 4, 6, 8, 10, 13, $16 cm / s$	$R e > 2000$ at threshold speed; flow gradient expected to be small near center of tank	LL+/LL-	2–4 cm/s

*

LL+: lateral line enabled; LL−: lateral line disabled
†

Data are extracted from the same set of experiments
‡

Two experiments are considered from the same paper
§

Data are from Baker and Montgomery, 1999.

Appendix 3—table 2

Results of the bibliographical research on fish rheotaxis in the absence of visual cues, used to validate the proposed model.

Reference	Fish species	*Evidence	Comparison with model
Reference	Fish species	*Evidence	Supportive	Inconclusive
Within studies
Effect of lateral line
Bak-Coleman et al., 2013	Giant danio (Devario aequipinnatus)	No significant difference in fish heading angle against current was detected between LL+ and LL-		×
Bak-Coleman and Coombs, 2014	blind cavefish (Astyanax mexicanus)	Rheotaxis threshold speed was slightly (but not significantly) lower in LL- condition		×
Baker and Montgomery, 1999 and Montgomery et al., 1997	blind cavefish (Astyanax fasciatus)	Rheotaxis threshold speed was significantly higher in LL- condition; fish received intermittent tactile senses		×
Elder and Coombs, 2015	Mexican tetras (Astyanax mexicanus)	No significant influence of LL condition was detected on rheotactic performance		×
Kulpa et al., 2015	blind cavefish (Astyanax mexicanus)	Significantly higher rheotaxis index in LL+ fish than LL- fish in jet stream	×
Oteiza et al., 2017	zebrafish (Danio rerio) larva 5–7 dpf	Posterior lateral line ablation or chemical neuromast ablation severely reduced rheotaxis	×
Suli et al., 2012	zebrafish (Danio rerio) larva 5 dpf	LL hair cell damage led to a significant decrease in rheotaxis; regeneration of LL hair cells restored rheotaxis	×
Van Trump and McHenry, 2013	blind Mexican cavefish (Astyanax fasciatus)	In LL+ and LL-, fish exhibited statistically indistinguishable rheotaxis behavior		×
Effect of flow gradient
Lyon, 1904	blind Fundulus	In a flow with small gradient, rheotaxis was elicited only when fish received tactile cues; in jet flow with large gradient, rheotaxis was elicited by flow without tactile cues. Lack of data on statistical significance		×
Oteiza et al., 2017	zebrafish (Danio rerio) larva 5–7 dpf	Rheotaxis of fish improved with increasing gradient magnitudes	×
Across studies
Effect of channel width
Bak-Coleman and Coombs, 2014; Baker and Montgomery, 1999	blind cavefish (Astyanax mexicanus); blind cavefish (Astyanax fasciatus)	Significantly different threshold speed for LL+ fish: $0.90 \pm 0.137 cm / s$ (mean ±s.e.m.) in $25 cm$ wide tunnel; between $2 cm / s$ and $3 cm / s$ in $9 cm$ wide tunnel. Tactile cues available to fish in Bak-Coleman and Coombs, 2014		×
Bak-Coleman and Coombs, 2014; Van Trump and McHenry, 2013	blind cavefish (Astyanax mexicanus); blind cavefish (Astyanax fasciatus)	Significantly different threshold speed for LL+ fish: $0.90 \pm 0.137 cm / s$ (mean ±s.e.m.) in $25 cm$ wide tunnel; between $2 cm / s$ and $4 cm / s$ in $\sim 11 cm$ diameter tunnel. Tactile cues available to fish in Bak-Coleman and Coombs, 2014		×
Oteiza et al., 2017; Peimani et al., 2017	zebrafish (Danio rerio) larva 5–7 dpf	Onset of rheotaxis in LL+ fish observed at flow speed $0.95 cm / s$ in $1.6 mm$ wide tunnel; rheotaxis observed in LL+ fish at flow speed $0.2 cm / s$ in $2.22 cm$ diameter tunnel	×
Effect of body length
Bak-Coleman and Coombs, 2014; Elder and Coombs, 2015	blind cavefish (Astyanax mexicanus); Mexican tetras (Astyanax mexicanus)	Significantly different threshold speed for LL+ fish: $0.90 \pm 0.137 cm / s$ (mean ±s.e.m.) for 4.2– $5.0 cm$ long fish; $1.96 \pm 0.350 cm / s$ (mean ±s.e.m.) for $8.3 cm$ long fish. Tactile cues available to fish in Bak-Coleman and Coombs, 2014		×
Total			5	9

*

LL+: lateral line enabled; LL−: lateral line disabled

Appendix 4

Stability analysis of the dynamical system based on the A-dipole

Derivation of the turn rate equation

Different from Equation 8, the A-dipole model treats fish as a slender body, which rotates in response to a gradient in the cross flow along its body. As such, under the A-dipole paradigm, the hydrodynamic turn rate becomes

\begin{matrix} Ω^{(A)} ({\vec{r}}_{f}, θ_{f}) = {\hat{v}}_{f}^{⟂} \cdot [{\nabla ({\vec{u}}_{w} (\vec{r}, {\vec{r}}_{f}, θ_{f}) + {\vec{u}}_{b} (\vec{r})) |}_{\vec{r} = {\vec{r}}_{f}}] {\hat{v}}_{f} \\ = \frac{π^{3} ρ^{2} v_{0}}{4 h} \cot (\frac{π y_{f}}{h}) \csc^{2} (\frac{π y_{f}}{h}) \cos θ_{f} + \frac{8 U_{0} ϵ}{h} (\frac{y_{f}}{h} - \frac{1}{2}) \sin^{2} θ_{f} . \end{matrix}

While the equation for the cross-stream motion of the A-dipole remains the same as Equation 10a, the governing equation for its orientation requires modification. Substituting Equation 30 into Equation 6b yields

{\dot{θ}}_{f} = \frac{π^{3} ρ^{2}}{4} \cot (π (ξ + \frac{1}{2})) \csc^{2} (π (ξ + \frac{1}{2})) \cos θ_{f} + 8 α ξ (\sin^{2} θ_{f} + κ),

which is analogous to Equation 10b for a T-dipole. Equation 10a and Equation 31 constitute the governing equations of the nonlinear planar dynamical system for an A-dipole fish model.

Determination of the equilibria of the dynamical system

Similar to our analysis of the T-dipole, the equilibria of the dynamical system for the A-dipole can be determined by setting $θ_{f} = 0$ and $θ_{f} = π$ in Equation 10a and Equation 31. Analogous to Equation 14, the cross-stream location ξ of the equilibria can be determined through a transcendental equation,

\frac{π^{3}}{32} \cot (π (ξ + \frac{1}{2})) \csc^{2} (π (ξ + \frac{1}{2})) = \mp β^{(A)} ξ,

where the negative sign corresponds to $θ_{f} = 0$ and the positive sign to $θ_{f} = π$ . Similar to the analysis of the T-dipole, here, we introduce a nondimensional parameter $β^{(A)} = α κ / ρ^{2}$ .

For downstream swimming ( $θ_{f} = 0$ ), there exists one or three roots to Equation 32 depending on the value of $β^{(A)}$ . For $β^{(A)}$ smaller than the critical value of $β^{*} = π^{4} / 32$ , $ξ = 0$ is the only solution, while for $β^{(A)} > β^{*}$ , two additional roots symmetrically located with respect to the centerline emerge. On the other hand, for upstream swimming ( $θ_{f} = π$ ), there is only one root at $ξ = 0$ .

Local stability analysis of the dynamical system

Local stability of the equilibria is studied by linearizing Equation 10a and Equation 31 about the equilibria. Specifically, for $θ_{f} = 0$ and $ξ = 0$ , the state matrix reads

A = [\begin{matrix} 0 & 1 - \frac{π^{2} ρ^{2}}{6} \\ 8 κ α - \frac{π^{4} ρ^{2}}{4} & 0 \end{matrix}] .

Stability requires that the determinant of $A$ is positive ( $det A > 0$ ), that is,

\frac{1}{24} (6 - π^{2} ρ^{2}) (π^{4} ρ^{2} - 32 α κ) > 0.

Since the first factor is always positive ( $ρ ≪ 1$ ), the inequality is fulfilled when $β^{(A)} < β^{*}$ . Hence, the system is stable if $β^{(A)} < β^{*}$ and unstable if $β^{(A)} > β^{*}$ .

For $θ_{f} = 0$ and $β^{(A)} > β^{*}$ , we identify two additional equilibria symmetrically located with respect to the centerline. The state matrix for these equilibria reads

A = [\begin{matrix} 0 & 1 + \frac{π^{2} ρ^{2}}{12} - \frac{1}{4} π^{2} ρ^{2} \sec^{2} (π ξ) \\ 8 κ α - \frac{1}{4} π^{4} ρ^{2} (2 - \cos (2 π ξ)) \sec^{4} (π ξ) & 0 \end{matrix}] .

Stability requires that $det A > 0$ , that is,

\frac{1}{48} (12 - 3 π^{2} ρ^{2} \sec^{2} (π ξ) + π^{2} ρ^{2}) (- 32 α κ + π^{4} ρ^{2} (2 - \cos (2 π ξ)) \sec^{4} (π ξ)) > 0.

As demonstrated in Methods and Materials Section, the first factor in Equation 36 can be assumed to be positive. Thus, the inequality in Equation 36 is equivalent to

β^{(A)} < \frac{π^{4}}{32} (2 - \cos (2 π ξ)) \sec^{4} (π ξ),

which is always satisfied for any choice of $β^{(A)}$ . Therefore, the additional equilibria emerging for $β^{(A)} > β^{*}$ at $ξ \neq 0$ are always stable.

For upstream swimming along the centerline ( $θ_{f} = π$ and $ξ = 0$ ), the state matrix is

A = [\begin{matrix} 0 & \frac{π^{2} ρ^{2}}{6} - 1 \\ 8 κ α + \frac{π^{4} ρ^{2}}{4} & 0 \end{matrix}] .

Similarly, stability requires that $det A > 0$ , that is,

\frac{1}{24} (6 - π^{2} ρ^{2}) (32 α κ + π^{4} ρ^{2}) > 0.

Since both the first and the second factors are positive, the condition in Equation 39 is always satisfied and the equilibrium is stable.

Comparing predictions with the T-dipole, results are similar in that there is a critical threshold that modulates the stability of upstream/downstream swimming, with a preference for rheotactic behavior above the threshold. For the T-dipole, below the critical threshold, neither swimming downstream nor upstream along the centerline of the channel are stable, and above the threshold, the latter becomes stable. Put simply, at low flow speed it is difficult for fish to orient in the flow, but they can successfully orient against the flow at high speed. For the A-dipole, below the critical threshold, both swimming downstream and upstream along the centerline of the channel are stable, and above the threshold, only the latter remains stable. In other words, fish can orient in the flow at low speed, without a preference for downstream or upstream swimming, but at high speed only swimming upstream becomes feasible.

Data availability

The authors declare that the data supporting the findings of this study are available within the paper. The Mathematica notebook used to derive the governing equations, study the planar dynamical, and generate associated figures, together with the CFD data discussed in the text, are also available at https://github.com/dynamicalsystemslaboratory/Rheotaxis; (copy archived at swh:1:rev:0c2828a439cb0df9971cdf12dc7dd93d16ce3450).

References

1. Arnold GP
(1969) The reactions of the plaice (Pleuronectes Platessa L.) to water currents
Journal of Experimental Biology 51:681–697.

https://doi.org/10.1242/jeb.51.3.681
- Google Scholar
1. Arnold GP
(1974) Rheotropism in fishes
Biological Reviews of the Cambridge Philosophical Society 49:515–576.

https://doi.org/10.1111/j.1469-185x.1974.tb01173.x
- PubMed
- Google Scholar
(2013) The spatiotemporal dynamics of rheotactic behavior depends on flow speed and available sensory information
The Journal of Experimental Biology 216:4011–4024.

https://doi.org/10.1242/jeb.090480
- PubMed
- Google Scholar
1. Bak-Coleman J.
2. Coombs S
(2014) Sedentary behavior as a factor in determining lateral line contributions to rheotaxis
The Journal of Experimental Biology 217:2338–2347.

https://doi.org/10.1242/jeb.102574
- PubMed
- Google Scholar
1. Baker CF
2. Montgomery JC
(1999) The sensory basis of rheotaxis in the blind Mexican cave fish, Astyanax fasciatus
Journal of Comparative Physiology A 184:519–527.

https://doi.org/10.1007/s003590050351
- Google Scholar
Book
1. Bakker P
(1991)
Some Elements of the Qualitative Theory of Differential Equations

Heidelberg, Germany: Springer.
- Google Scholar
Book
1. Batchelor GK
(2000) An Introduction to Fluid Dynamics
Cambridge Mathematical Library. Cambridge University Press.

https://doi.org/10.1017/CBO9780511800955
- Google Scholar
1. Brawn VM
(1960) Underwater television observations of the swimming speed and behaviour of captive herring
Journal of the Fisheries Research Board of Canada 17:689–698.

https://doi.org/10.1139/f60-054
- Google Scholar
1. Burbano-L DA
2. Porfiri M
(2021) Modeling multi-sensory feedback control of zebrafish in a flow
PLOS Computational Biology 17:e1008644.

https://doi.org/10.1371/journal.pcbi.1008644
- PubMed
- Google Scholar
1. Calovi DS
2. Lopez U
3. Ngo S
4. Sire C
5. Chaté H
6. Theraulaz G
(2014) Swarming, schooling, milling: phase diagram of a data-driven fish school model
New Journal of Physics 16:015026.

https://doi.org/10.1088/1367-2630/16/1/015026
- Google Scholar
1. Chapman JW
2. Klaassen RHG
3. Drake VA
4. Fossette S
5. Hays GC
6. Metcalfe JD
7. Reynolds AM
8. Reynolds DR
9. Alerstam T
(2011) Animal orientation strategies for movement in flows
Current Biology 21:R861–R870.

https://doi.org/10.1016/j.cub.2011.08.014
- PubMed
- Google Scholar
(2015) Rheotaxis performance increases with group size in a coupled phase model with sensory noise
The European Physical Journal Special Topics 224:3233–3244.

https://doi.org/10.1140/epjst/e2015-50080-x
- Google Scholar
1. Colgate JE
2. Lynch KM
(2004) Mechanics and control of swimming: A review
IEEE Journal of Oceanic Engineering 29:660–673.

https://doi.org/10.1109/JOE.2004.833208
- Google Scholar
1. Colvert B
2. Kanso E
(2016) Fishlike rheotaxis
Journal of Fluid Mechanics 793:656–666.

https://doi.org/10.1017/jfm.2016.141
- Google Scholar
(2020) Rheotaxis revisited: a multi-behavioral and multisensory perspective on how fish orient to flow
The Journal of Experimental Biology 223:23.

https://doi.org/10.1242/jeb.223008
- PubMed
- Google Scholar
1. Couzin ID
2. Krause J
3. Franks NR
4. Levin SA
(2005) Effective leadership and decision-making in animal groups on the move
Nature 433:513–516.

https://doi.org/10.1038/nature03236
- PubMed
- Google Scholar
1. Dabiri JO
(2017) How fish feel the flow
Nature 547:406–407.

https://doi.org/10.1038/nature23096
- PubMed
- Google Scholar
1. Davidson VM
(1949) Salmon and eel movement in constant circular current
Journal of the Fisheries Research Board of Canada 7c:432–448.

https://doi.org/10.1139/f47-035
- Google Scholar
(2019) Design, modeling, and control of biomimetic fish robot: A review
Journal of Bionic Engineering 16:967–993.

https://doi.org/10.1007/s42235-019-0111-7
- Google Scholar
(1968) Permanence of station in juvenile chinook salmon and steelhead trout
Journal of the Fisheries Research Board of Canada 25:1453–1464.

https://doi.org/10.1139/f68-126
- Google Scholar
1. Elder J
2. Coombs S
(2015) The influence of turbulence on the sensory basis of rheotaxis
Journal of Comparative Physiology. A, Neuroethology, Sensory, Neural, and Behavioral Physiology 201:667–680.

https://doi.org/10.1007/s00359-015-1014-7
- PubMed
- Google Scholar
1. Elliott JM
(1965) Daily fluctuations of drift invertebrates in a dartmoor stream
Nature 205:1127–1129.

https://doi.org/10.1038/2051127c0
- Google Scholar
1. Filella A
2. Nadal F
3. Sire C
4. Kanso E
5. Eloy C
(2018) Model of collective fish behavior with hydrodynamic interactions
Physical Review Letters 120:198101.

https://doi.org/10.1103/PhysRevLett.120.198101
- PubMed
- Google Scholar
(1993) Temperature-dependent switch between diurnal and nocturnal foraging in salmon
Proceedings of the Royal Society of London. Series B 252:135–139.

https://doi.org/10.1098/rspb.1993.0057
- Google Scholar
1. Gautrais J
2. Jost C
3. Soria M
4. Campo A
5. Motsch S
6. Fournier R
7. Blanco S
8. Theraulaz G
(2009) Analyzing fish movement as a persistent turning walker
Journal of Mathematical Biology 58:429–445.

https://doi.org/10.1007/s00285-008-0198-7
- PubMed
- Google Scholar
(2014) Scaling macroscopic aquatic locomotion
Nature Physics 10:758–761.

https://doi.org/10.1038/nphys3078
- Google Scholar
(2015) Gait and speed selection in slender inertial swimmers
PNAS 112:3874–3879.

https://doi.org/10.1073/pnas.1419335112
- PubMed
- Google Scholar
(2016) Learning to school in the presence of hydrodynamic interactions
Journal of Fluid Mechanics 789:726–749.

https://doi.org/10.1017/jfm.2015.686
- Google Scholar
1. Gibson RJ
(1966) Some factors influencing the distributions of brook trout and young atlantic salmon
Journal of the Fisheries Research Board of Canada 23:1977–1980.

https://doi.org/10.1139/f66-185
- Google Scholar
(2016) Touch sensation by pectoral fins of the catfish Pimelodus pictus
Proceedings. Biological Sciences 283:1824.

https://doi.org/10.1098/rspb.2015.2652
- PubMed
- Google Scholar
1. Hoar WS
(1954) The behaviour of juvenile pacific salmon, with particular reference to the sockeye (Oncorhynchus nerka
Journal of the Fisheries Research Board of Canada 11:69–97.

https://doi.org/10.1139/f54-009
- Google Scholar
1. Kanso E
2. Cheng Hou Tsang A
(2014) Dipole models of self-propelled bodies
Fluid Dynamics Research 46:061407.

https://doi.org/10.1088/0169-5983/46/6/061407
- Google Scholar
1. Kanso E.
2. Michelin S
(2019) Phoretic and hydrodynamic interactions of weakly confined autophoretic particles
The Journal of Chemical Physics 150:044902.

https://doi.org/10.1063/1.5065656
- PubMed
- Google Scholar
1. Kopman V
2. Porfiri M
(2013) Design, modeling, and characterization of a miniature robotic fish for research and education in biomimetics and bioinspiration
IEEE/ASME Transactions on Mechatronics 18:471–483.

https://doi.org/10.1109/TMECH.2012.2222431
- Google Scholar
1. Kopman V
2. Laut J
3. Acquaviva F
4. Rizzo A
5. Porfiri M
(2015) Dynamic modeling of a robotic fish propelled by a compliant tail
IEEE Journal of Oceanic Engineering 40:209–221.

https://doi.org/10.1109/JOE.2013.2294891
- Google Scholar
(2015) The lateral line is necessary for blind cavefish rheotaxis in non-uniform flow
The Journal of Experimental Biology 218:1603–1612.

https://doi.org/10.1242/jeb.119537
- PubMed
- Google Scholar
1. Liao JC
(2007) A review of fish swimming mechanics and behaviour in altered flows
Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences 362:1973–1993.

https://doi.org/10.1098/rstb.2007.2082
- PubMed
- Google Scholar
Book
1. Lighthill SJ
(1975) Mathematical Biofluiddynamics
Philadelphia, United states: SIAM.

https://doi.org/10.1137/1.9781611970517
- Google Scholar
(2020) Airfoil-like mechanics generate thrust on the anterior body of swimming fishes
PNAS 117:10585–10592.

https://doi.org/10.1073/pnas.1919055117
- PubMed
- Google Scholar
1. Lyon EP
(1904) On rheotropism. i. — rheotropism in fishes
American Journal of Physiology-Legacy Content 12:149–161.

https://doi.org/10.1152/ajplegacy.1904.12.2.149
- Google Scholar
1. McNaught DC
2. Hasler AD
(1961) Surface schooling and feeding behavior in the white bass, ROCCUS CHRYSOPS (RAFINESQUE), in lake mendota1
Limnology and Oceanography 6:53–60.

https://doi.org/10.4319/lo.1961.6.1.0053
- Google Scholar
Book
1. Milne-Thomson LM
(1996)
Theoretical Hydrodynamics

Courier Corporation.
- Google Scholar
1. Mogdans J
(2019) Sensory ecology of the fish lateral-line system: Morphological and physiological adaptations for the perception of hydrodynamic stimuli
Journal of Fish Biology 95:53–72.

https://doi.org/10.1111/jfb.13966
- PubMed
- Google Scholar
(1997) The lateral line can mediate rheotaxis in fish
Nature 389:960–963.

https://doi.org/10.1038/40135
- Google Scholar
1. Montgomery J
2. Carton G
3. Voigt R
4. Baker C
5. Diebel C
(2000) Sensory processing of water currents by fishes
Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences 355:1325–1327.

https://doi.org/10.1098/rstb.2000.0693
- PubMed
- Google Scholar
Book
1. Montgomery JC
2. Baker CF
(2020)
Lateral line and fish behavior

In: Fritzsch B, editors. The Senses: A Comprehensive Reference (second edition). Oxford: Elsevier. pp. 133–142.
- Google Scholar
1. Najafi S
2. Abtahi B
(2022) Numerical simulation of aquatic animals locomotion using unsteady panel method
Ocean Engineering 244:110380.

https://doi.org/10.1016/j.oceaneng.2021.110380
- Google Scholar
1. Needham PR
2. Jones AC
(1959) Flow, temperature, solar radiation, and ice in relation to activities of fishes in sagehen creek, california
Ecology 40:465–474.

https://doi.org/10.2307/1929764
- Google Scholar
Book
1. Newton PK
(2011)
The N-Vortex Problem: Analytical Techniques

Heidelberg, Germany: Springer Science & Business Media.
- Google Scholar
1. Oteiza P
2. Odstrcil I
3. Lauder G
4. Portugues R
5. Engert F
(2017) A novel mechanism for mechanosensory-based rheotaxis in larval zebrafish
Nature 547:445–448.

https://doi.org/10.1038/nature23014
- PubMed
- Google Scholar
1. Pavlov DS
2. Tjurjukov SN
(1995) Reactions of dace to linear accelerations
Journal of Fish Biology 46:768–774.

https://doi.org/10.1111/j.1095-8649.1995.tb01600.x
- Google Scholar
(2017) A microfluidic device for quantitative investigation of zebrafish larvae’s rheotaxis
Biomedical Microdevices 19:1–6.

https://doi.org/10.1007/s10544-017-0240-x
- PubMed
- Google Scholar
Book
(2021)
Emergence of In-Line Swimming Patterns in Zebrafish Pairs

Cambridge, United Kingdom: cambridge university.
- Google Scholar
(2020) A minimal model of self propelled locomotion
Journal of Fluids and Structures 97:103071.

https://doi.org/10.1016/j.jfluidstructs.2020.103071
- Google Scholar
(2002) Boundary-layer receptivity to freestream disturbances
Annual Review of Fluid Mechanics 34:291–319.

https://doi.org/10.1146/annurev.fluid.34.082701.161921
- Google Scholar
Book
1. Sastry S
(2013)
Nonlinear Systems: Analysis, Stability, and Control

Springer Science & Business Media.
- Google Scholar
(1999) Review of fish swimming modes for aquatic locomotion
IEEE Journal of Oceanic Engineering 24:237–252.

https://doi.org/10.1109/48.757275
- Google Scholar
1. Suli A
2. Watson GM
3. Rubel EW
4. Raible DW
5. Riley B
(2012) Rheotaxis in larval zebrafish is mediated by lateral line mechanosensory hair cells
PLOS ONE 7:e29727.

https://doi.org/10.1371/journal.pone.0029727
- Google Scholar
1. Taylor GI
(1997) Analysis of the swimming of long and narrow animals
Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences 214:158–183.

https://doi.org/10.1098/rspa.1952.0159
- Google Scholar
(2012) The finite-dipole dynamical system
Proceedings of the Royal Society A 468:3006–3026.

https://doi.org/10.1098/rspa.2012.0119
- Google Scholar
1. Van Trump WJ
2. McHenry MJ
(2013) The lateral line system is not necessary for rheotaxis in the Mexican blind cavefish (Astyanax fasciatus)
Integrative and Comparative Biology 53:799–809.

https://doi.org/10.1093/icb/ict064
- PubMed
- Google Scholar
1. von der Emde G
(1999) Active electrolocation of objects in weakly electric fish
The Journal of Experimental Biology 202:1205–1215.

https://doi.org/10.1242/jeb.202.10.1205
- PubMed
- Google Scholar
Book
1. Wang R
2. Wang S
3. Wang Y
4. Cheng L
5. Tan M
(2020)
Development and Motion Control of Biomimetic Underwater Robots: A Survey

IEEE Transactions on Systems.
- Google Scholar
1. Waters TF
(1962) Diurnal periodicity in the drift of stream invertebrates
Ecology 43:316.

https://doi.org/10.2307/1931988
- Google Scholar
1. Webb PW
(2002) Kinematics of plaice, Pleuronectes platessa, and cod, Gadus morhua, swimming near the bottom
The Journal of Experimental Biology 205:2125–2134.

https://doi.org/10.1242/jeb.205.14.2125
- PubMed
- Google Scholar
Book
1. White FM
(1974)
Viscous Fluid Flow

New York, United States: McGraw-Hill.
- Google Scholar
1. Williamson C
2. Govardhan R
(2004) Vortex-induced vibrations
Annual Review of Fluid Mechanics 36:413–455.

https://doi.org/10.1146/annurev.fluid.36.050802.122128
- Google Scholar
1. Wu TYT
(2006) Hydromechanics of swimming propulsion. Part 1. Swimming of a two-dimensional flexible plate at variable forward speeds in an inviscid fluid
Journal of Fluid Mechanics 46:337–355.

https://doi.org/10.1017/S0022112071000570
- Google Scholar
1. Zhang F
2. Lagor FD
3. Lei H
4. Tan X
5. Paley DA
(2016) Robotic Fish
Mechanical Engineering 138:S2–S5.

https://doi.org/10.1115/1.2016-Mar-6
- Google Scholar
(2019) An information-theoretic study of fish swimming in the wake of a pitching airfoil
Physica D 396:35–46.

https://doi.org/10.1016/j.physd.2019.02.014
- Google Scholar
(2015) Data-driven stochastic modelling of zebrafish locomotion
Journal of Mathematical Biology 71:1081–1105.

https://doi.org/10.1007/s00285-014-0843-2
- PubMed
- Google Scholar

Article and author information

Author details

Maurizio Porfiri
1. Department of Biomedical Engineering, New York University, New York, United States
2. Department of Mechanical and Aerospace Engineering, New York University, New York, United States
3. Center for Urban Science and Progress, New York University Tandon School of Engineering, New York, United States
Contribution
Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, M.P. conceived the study, developed the theoretical model and performed data analysis, wrote a first draft of the manuscript, which was consolidated in its present form by all the authors., Methodology, Project administration, Resources, Software, Supervision, Validation, Visualization, Writing – original draft, Writing – review and editing

For correspondence
mporfiri@nyu.edu

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0002-1480-3539
Peng Zhang
1. Department of Mechanical and Aerospace Engineering, New York University, New York, United States
2. Center for Urban Science and Progress, New York University Tandon School of Engineering, New York, United States
Contribution
Data curation, Investigation, Methodology, P.Z. conducted the literature review and performed the model validation. The manuscript was consolidated in its present form by all the authors., Validation, Writing – original draft

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0001-8237-1259
Sean D Peterson

Mechanical and Mechatronics Engineering Department, University of Waterloo, Waterloo, Canada

Contribution
Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, S.D.P. conceived the study, developed the theoretical model and performed data analysis, wrote a first draft of the manuscript, which was consolidated in its present form by all the authors., Software, Supervision, Validation, Visualization, Writing – original draft, Writing – review and editing

For correspondence
peterson@uwaterloo.ca

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0001-8746-2491

Funding

National Science Foundation (CMMI-1901697)

Maurizio Porfiri

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

The authors acknowledge financial support from the National Science Foundation under Grant No. CMMI-1901697. The authors acknowledge Alain Boldini and Simone Macrì for useful discussions.

Copyright

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.