Large-amplitude, short-wave peristalsis and its implications for transport

Valveless, tubular pumps are widespread in the animal kingdom, but the mechanism by which these pumps generate fluid flow is often in dispute. Where the pumping mechanism of many organs was once described as peristalsis, other mechanisms, such as dynamic suction pumping, have been suggested as possible alternative mechanisms. Peristalsis is often evaluated using criteria established in a technical definition for mechanical pumps, but this definition is based on a small-amplitude, long-wave approximation which biological pumps often violate. In this study, we use a direct numerical simulation of large-amplitude, short-wave peristalsis to investigate the relationships between fluid flow, compression frequency, compression wave speed, and tube occlusion. We also explore how the flows produced differ from the criteria outlined in the technical definition of peristalsis. We find that many of the technical criteria are violated by our model: Fluid flow speeds produced by peristalsis are greater than the speeds of the compression wave; fluid flow is pulsatile; and flow speed have a nonlinear relationship with compression frequency when compression wave speed is held constant. We suggest that the technical definition is inappropriate for evaluating peristalsis as a pumping mechanism for biological pumps because they too frequently violate the assumptions inherent in these criteria. Instead, we recommend that a simpler, more inclusive definition be used for assessing peristalsis as a pumping mechanism based on the presence of non-stationary compression sites that propagate unidirectionally along a tube without the need for a structurally fixed flow direction.


Introduction
Tubular pumps are found across the animal kingdom in all stages of development (e.g., Griffiths et al. 1987;Gashev 2002;Xavier-Neto et al. 2007; Lee and Socha 2009;Glenn et al. 2010;Xavier-Neto et al. 2010;Krenn 2010;Greenlee et al. 2013;Harrison et al. 2013a). In all vertebrate embryos, the heart first forms as a valveless tubular pump (e.g., Taber 2001). Similarly, the hearts of many non-vertebrate chordates such as tunicates are also valveless, tubular pumps (Anderson 1968;Kalk 1970;Glenn et al. 2010;Waldrop and Miller 2015). In the last 10 years, the mechanism of pumping in tubular hearts and other tubular pumps has been described as either peristalsis or dynamic suction pumping (Forouhar et al. 2006;Taber et al. 2007;Männer et al. 2010;Maes et al. 2011).
For the case of peristalsis, a wave of active contraction propagates down the length of the tube. Since the advent of mechanical peristaltic pumps (i.e., roller pumps), the technical definition of a peristaltic pump has been refined to include these characteristics (summarized in Männer et al. 2010): 1. Peristaltic pumps are positive-displacement pumps, which displace a fixed volume to create fluid flow. 2. Peristaltic pumps have non-stationary compression sites, i.e., waves of compression that unidirectionally propagate down a flexible tube. 3. Peristaltic pumps produce continuous flow. 4. There are no structurally fixed directions of flow (e.g., there are no one-way valves), and the direction of flow is determined by the direction of the compression wave.
5. Flow velocity is equivalent to the speed of the compression wave. Peak flow velocity does not exceed the speed of the compression wave. 6. There is a linear relationship between the frequency of the compression wave and flow rate it produces.
The principles of the technical definition are based upon a body of analytical work using small-amplitude and/or longwave approximations of peristaltic pumping where nonlinear effects, such as inertia and large flows in the radial direction, are small (Jaffrin and Shapiro 1971;Hanin 1968;Shapiro et al. 1969;Fung and Yih 1968). These results may not apply to cases when inertia is non-negligible or where the compression amplitude is not small relative to the diameter or the wavelength. These nonlinear cases do, however, characterize many biological pumps (Santhanakrishnan and Miller 2011).
There are, however, some numerical and analytic results available for large-amplitude and/or short-wavelength peristalsis. Childress (2009) showed analytically for large occlusion ratios and long wavelengths that the peak flow velocity can exceed the peristaltic wave speed. Pozrikidis (1987) modeled peristalsis for relatively small wavelengths in a Stokes fluid. He found that peak flow velocities are nearly twice the wave speed for an 80 % occlusion of the tube. Large-amplitude, short-wave peristalsis has also been studied numerically for Stokesian viscous and viscoelastic fluids. Teran et al. (2008) simulated up to 50 % occlusion and reported only mean flow rates for viscous and viscoelastic fluids. Aranda et al. (2011) simulated high-amplitude peristaltic pumping in a three-dimensional tube with a phase-shifted asymmetry. The mean flow rate approached the wave speed for full occlusion, implying that the peak speeds are higher than the wave speed given the large spatial and temporal variations in the flow. Ceniceros and Fisher (2012) also simulated viscoelastic fluids at high occlusion ratios but also only reported mean flow rates. For Newtonian fluids in nearly occluded tubes, mean flow rates approached the wave speed. Given the large spatial variations in flow as evidenced by the given vorticity plots, it is likely that peak flow speeds were much higher than the wave speed.
In contrast, dynamic suction pumping is defined by an isolated region of active contraction that is asymmetrically located in a section of flexible tube connected to relatively stiffer inflow and outflow tracts (Liebau 1954(Liebau , 1955(Liebau , 1957. Passive elastic traveling waves emanate from the active contraction site, and these waves drive flow. Analytical models (Auerbach et al. 2004;Bringley et al. 2008), physical experiments (Hickerson et al. 2005a;Bringley et al. 2008), and numerical simulations (Jung and Peskin 2000;Jung et al. 2008;Baird et al. 2014) support that this pumping mechanism can effectively transport fluid under certain conditions. Furthermore, these pumps are characterized by a nonlinear frequency-flow relationship, reversals in flow directions, and flow speeds higher than the wave speed.
Based on the technical definition of peristalsis, Forouhar et al. (2006) make the case that the zebrafish embryonic heart is not a peristaltic pump by observing the kinematics of heart compression and blood flow in vivo. They found these observations violated some of the principles that underpin peristalsis, specifically: They observed bidirectional compression and reflection waves of the heart muscle and imply that activation of the muscle was limited to one site (#2 and #3), the flow velocity inside the heart exceeded the compression wave speed (#5), and there was a nonlinear relationship between blood flow speed and the frequency of heart compressions (#6). As a result of these observations, Forouhar et al. (2006) reject peristalsis and suggest that the fluid dynamics of the embryonic vertebrate heart share more of the properties of dynamic suction pumping than peristalsis (a localized site of active compression is placed off-center, and the resulting traveling waves are passive elastic).
Since the publication of this paper nearly 8 years ago, some researchers have speculated that other biological pumps drive blood and other fluids using dynamic suction pumping (Davidson 2007;Vogel 2007;Harrison et al. 2013b). Other researchers continue to describe tubular pumping in the embryonic heart as peristalsis (Christoffels and Moorman 2009;Postma et al. 2008;Taber et al. 2007). Some of these tubular pumps, like tunicate hearts and many embryonic hearts, violate the low-amplitude, long-wave assumptions of the technical definition of peristalsis. In many of these cases, the compression wave almost completely occludes the tube, action potentials are known to propagate the length of the entire tube by activating helically wound muscle fibers (Anderson 1968;Kalk 1970), and the width of the tube is not significantly longer than its length. Furthermore, these pumps are often at a scale where neither inertial nor viscous effects can be neglected.
In this paper, we use direct numerical simulation of the fully coupled fluid-structure interaction problem to quantify the flows generated by large-amplitude, short-wave peristalsis. We then evaluate which aspects of the technical definition of peristalsis a large-amplitude, short-wave peristaltic pump can fulfill. The simulations were done using the immersed boundary method (Peskin 2002). The peristaltic wave was prescribed as a traveling sine wave. The length of the tube was four times the diameter, and the region of prescribed motion was centered in the middle of the tube. The ends of the elastic section were allowed to bend as determined by the coupling between the fluid and the elastic tube. We quantify peak and average fluid speeds downstream of the compression region, the pressure at the inflow and outflow points of the compression tube, and flow speeds at a point within the compression region.

The immersed boundary method
The modeling and numerical simulations of peristalsis were conducted with the immersed boundary method (Peskin 2002) using the library IBAMR (Griffith 2014). The IB method can be used either to simulate a boundary moving with a prescribed motion or for solving for the motion based on the interaction of the fluid with an elastic boundary. The following outline describes the two-dimensional formulation of the immersed boundary method, but the three-dimensional extension is straightforward. For a full review of the method, see Peskin (2002). The equations of fluid motion are given by the Navier-Stokes equations: where u(x, t) is the fluid velocity, p(x, t) is the pressure, F(x, t) is the force per unit area applied to the fluid by the immersed boundary, ρ is the density of the fluid, and μ is the dynamic viscosity of the fluid. The independent variables are the time t and the position x.
The interaction equations between the fluid and the boundary are given by: where f(r, t) is the force per unit length applied by the boundary to the fluid as a function of Lagrangian position and time, δ(x) is a two-dimensional delta function, X(r, t) gives the Cartesian coordinates at time t of the material point labeled by the Lagrangian parameter r . Equation 3 applies force from the boundary to the fluid grid, and Eq. 4 evaluates the local fluid velocity at the boundary. The boundary is then moved at the local fluid velocity, and this enforces the noslip boundary condition. Each of these equations involves a two-dimensional Dirac delta function, δ, which acts as the kernel of an integral transformation. These equations convert Lagrangian variables to Eulerian variables and vice versa. The force equations are specific to the application. In a simple case where a preferred motion is enforced, boundary points are tethered to target points via springs. The equation describing the force applied to the fluid by the boundary in Lagrangian coordinates is given by: where f(r, t) is the force per unit length, k targ is a stiffness coefficient of the tethering springs, and Y(r, t) is the prescribed position of the target boundary.

Dimensionless numbers
To compare the flows of different pumps over a range of length scales and velocities, it is a useful exercise to nondimensionalize the terms in the Navier-Stokes equations as follows, where L, U , and 1/ f * are characteristic flow length, velocity, and time scales, respectively. In this paper, we choose L to be the diameter of the elastic section of the tube, f * to be the beat frequency of our base case, and the characteristic velocity as U = L f * . Note the U describes the velocity in terms of heart tube lengths per beat.
The dimensionless bending stiffness of the boundary may then be calculated as where k bend is the flexural stiffness of the boundary. The dimensionless spring stiffness may be written as where k is the spring stiffness which describes the resistance to stretching. The dimensionless pressure can then be defined as where p is the dimensional pressure and U is the characteristic velocity.
Of particular relevance to internal biological flows is the Womersley number (Wo). The Wo describes to what extent unsteady effects matter in pulsatile flows. It is given by the following equation, where ω is the frequency of the pulse, d is the diameter of the tube, and ν is the kinematic viscosity of the fluid. Within the context of a blood vessel, when the value of Wo is high, the velocity profile is nearly flat over most of the cross section and there is a small region near the vessel wall known as the boundary layer where viscous effects are important. Also, the flow at the center of the tube is inertial and pulsatile. At the other low extreme end of Wo, the velocity profile over the vessel cross section is parabolic, and the flow is quasi-steady and viscous dominated. The transient effects can be ignored when Wo is sufficiently small, generally when Wo < 1, and this is common in the case of microcirculation. Unless otherwise noted, we vary the Wo by changing the viscosity of the fluid only.

Numerical method
We used an adaptive and parallelized version of the immersed boundary method, IBAMR (Griffith 2014). IBAMR is a C++ framework that provides discretization and solver infrastructure for partial differential equations on block-structured locally refined Eulerian grids (Berger and Oliger 1984;Berger and Colella 1989) and on Lagrangian (structural) meshes. IBAMR also includes infrastructure for coupling Eulerian and Lagrangian representations. The Eulerian grid on which the Navier-Stokes equations were solved was locally refined near the immersed boundaries and regions of vorticity with a threshold of |ω| > 0.1. This Cartesian grid was organized as a hierarchy of four nested grid levels, and the finest grid was assigned a spatial step size of dx = D/512, where D is the length of the domain. The ratio of the spatial step size on each grid relative to the next coarsest grid was 1:4. The numerical parameters used for the simulations are given in Table 1.

Model of peristalsis
A numerical model of an elastic heart tube connected to a rigid racetrack was constructed to study high-amplitude peristaltic flows. The racetrack design was used for easy comparison with previous models of tubular heart pumping (e.g., Jung and Peskin 2000;Hickerson et al. 2005b;Baird et al. 2014;Avrahami and Gharib 2008;Lee et al. 2012). The racetrack section was constructed by connecting two sections of Compression ratio 0.40-0.95 straight tube (one of which represents the elastic heart) to curved sections. The resting diameter of the racetrack was constant throughout its length. Dimensions and elastic properties of the racetrack are given in Table 2, and the model set up is shown in Fig. 1a.
The bottom portion of the racetrack consisted of an elastic section of dimensionless length L = 4.0 and dimensionless diameter d = 1.0. The middle 3/4 of the tube was tethered to target points that drove the peristaltic motion. Note that the tube was initialized at rest in a straight configuration. The amplitude of the peristaltic wave increased from zero to its maximum amplitude after one pulse. To ensure conservation of volume inside the racetrack during pumping, the ends of the elastic section were allowed to deform as a result of the coupling between the elastic boundary and the fluid. Note that the sides and top of the racetrack were also tethered to target points that did not move. The target point stiffness, k targ was set equal to the spring stiffness k. This parameter can be adjusted to minimize the distance between the actual boundary point and its target position. The position of the target points was determined by the following equation, where y top,bot is the y-coordinate of the top or bottom of the tube, R top,bot is the distance of the top or bottom of the tube from the horizontal centerline of the racetrack, A is the amplitude of the contraction, f is the frequency of contraction, c is the wave speed, and x t is the horizontal distance from the beginning of the prescribed motion. The compression ratio gives the percent occlusion, and is equal to 2 A.
In addition to the racetrack design, we constructed an open model of peristalsis which consisted of only the elastic heart tube on the same domain (Fig. 2). The heart tube was driven in an identical way to the racetrack design. This design Fig. 1 Panels showing racetrack circulatory system and compression region with speed across the simulation area (background color) and light blue flow marker particles. a Initial simulation conditions (t = 0.0) showing racetrack circulatory system (black lines), tube diameter d (white vertical line), compression region, L tube (white horizontal line between white boxes 3 and 4), tethered region (yellow bars) and marker particles (light blue dots). Areas over which data were collected are labeled as white boxes and a white/black point (see text for details on calculations). b Simulation at t = 0.4. c Simulation at t = 0.8. d Simulation at t = 1.2. e Simulation at t = 1.6. f Simulation at t = 2.0. Color scale units non-dimensional speed allowed us to eliminate the possibility of Liebau pumping in our closed, racetrack model by removing the flexible regions of heart tube that connected the contracting region to the rigid racetrack.

Parameter sweeps
We varied the values of four parameters in five separate sets of simulations. (1) Wo was changed by altering the dynamic viscosity μ of the fluid within the tube, which alters ν in Eq. 9 since ν = μ/ρ, where ρ is fluid density. Wo ranged from 0.1 to 10. (2) Tube occlusion was changed by altering the amplitude of contraction A in Eq. 10 by some factor (compression ratio), which ranged from 40 % tube occlusion (compressionratio = 0.4) to nearly complete occlusion of the tube (compressionratio = 0.95). (3) The speed of compression wave. (4) The frequency of compressions ( f in Eq. 10) was changed, ranging from 0.5 to 2.0, while allowing the speed of compression wave, c, to increase linearly. (5) The frequency of compressions ( f in Eq. 10) was changed while holding the speed of the compression wave, c, constant. Note that the no-racetrack design was used only for sweeps of compression wave frequency.
The default parameter values used for simulation sweeps other than the parameter of interest are: W o = 1, f * = 1, c = 3.0, and compression ratio = 0.8.    Fig. 1a. These mean speeds were then temporally averaged over the entire simulation to find the average speed per simulation U avg , presented in Fig. 5a-c (black circles) and Fig. 8a (black and gray circles).
U max,t , U max : The non-dimensional maximum magnitude of velocity across the same section of tube (white box labeled 1 in Fig. 1a) was also calculated for each time point of the simulation, U max,t . Fig. 3 shows U max,t versus dimensionless time. These numbers were temporally averaged over the entire simulation to find the average maximum speed per simulation, U max , presented in Fig. 5a-c (white diamonds) and Fig. 8a (white and gray diamonds).
U peak : Using the non-dimensional, maximum speeds across the upper section of tube (white box labeled 1 in Fig. 1a) for each simulation time step, peak speeds (max-imum of maximum speed) were found for each pulse. For each simulation, the peak speeds over each pulse were averaged to find the average peak speed for each simulation, U peak . u m , u x,t , u x : To characterize some of the dynamics in the compression region, the instantaneous, non-dimensional magnitude of velocity (|u |, speed) and x-component of fluid flow velocity, u x,t , were sampled at each time at one point in the center of the compression tube, indicated by the white/black dot labeled 2 in Fig. 1a. For the x-component of velocity, the instantaneous component of velocity for each time step u x,t versus t are presented in Fig. 4. These speeds were then temporally averaged across the entire simulation to find the average speed or magnitude of velocity (u m ) and average x-component of velocity (u x ) for each simulation. These speeds are presented in Fig. 5d-f as black squares (u m ) and triangles (u x ) and in Fig. 8b as black and gray squares (u m ) and black and gray triangles (u x ).
p in , p out , Δp: Non-dimensional pressure, p , at each time step in the simulation was spatially averaged across two cross sections of the racetrack near the compression region: the inflow region (Fig. 1a, white box 3) and outflow region (Fig. 1a, white box 4). The inflow pressures were subtracted from the outflow pressures to find the pressure difference at each time step and averaged across the entire simulation to find the average inflow pressure p in , the average outflow pressure p out , and average pressure difference (Δp) for each simulation. These pressures are presented in Fig. 7a-d as connected circles ( p in ), connected diamonds ( p out ), and dotted triangles (Δp).
For the no-racetrack design, the average volume flow rate v avg was calculated using MATLAB. The average xcomponent of velocity across the opening of the heart tube was calculated and then divided by the diameter of the opening at that time point. These values were then temporally averaged across simulation time to find v avg .

Results
Results discussed below are for the racetrack design unless otherwise noted.

Womersley number
Flow speeds in the racetrack U max,t are graphed against simulation time for several Wo in Fig. 3a. Note that the Wo was varied by changing the viscosity of the fluid rather than the pulse frequency and/or the wave speed. The dashed line shows the constant wave speed for these simulations.  likely due to inertial effects and the formation of strong vortices in the tube. In all cases, the maximum dimensionless flow speed is greater than the wave speed for part of the cycle.
U avg , U max , and U peak are graphed against Womersley number (Wo) in Fig. 5a. Nonlinear relationships exist between each of the dimensionless flow speed metrics and Wo, where the Wo was varied by changing the viscosity. U avg tend to steadily increase with increasing Wo, while average peak flow speeds are greatest around Wo = 1 and 50. U peak are greater than the speed of the compression wave, c, for all Wo and show variability with simulation time for values of Wo > 10 (Fig. 3a). Flow speeds within the compression tube u m and u x are plotted against Wo in Fig. 5d. u x,t is plotted against time for several simulations in of u x decreases and becomes negative for Wo > 1. This indicates that flow is moving on average in clockwise direction, or opposite the direction of the propagating compression wave. This phenomenon is due to significant backflow that occurs when the tube rapidly expands behind the compression wave.
Pressure For values Wo < 2, there are large differences in pressure (Δp). Since Wo is varied by changing viscosity, the Wo < 1 cases are viscous dominated, and this results in a high resistance of the fluid to move through the tube. Δp are greatest at Wo = 0.3, with p > 4000 (Fig. 7a). Of note is the case where Wo = 0.1 and the pressure difference decreases relative to Wo = 0.3. Note in Fig. 5a that the flow speed at this Wo drops to almost zero. For a compression ratio of 0.8, the peristaltic pump is not able to drive highly viscous flow around the racetrack, and the pressure difference between the inflow and outflow tracts drops. For values W ≥ 2, pressure differences between the two regions approach zero as the resistance to flow decreases with decreasing viscosity. Nondimensional pressure p is reported versus time in Fig. 6a for varying values of Wo.

Tube occlusion
Flow speeds in the racetrack U max,t is graphed against dimensionless time for several compression ratios in Fig. 3b. Note that U max,t does not exceed c for compression ratios set to 0.6 and less. For large compression ratios of 0.8-0.9, U max,t reaches speeds of nearly double the value of c. In all cases considered, U max,t is positive, and the flow moves in the counterclockwise direction at the top of the racetrack. U avg , U max , and U peak is graphed against compression ratio in Fig. 5b. Each metric of the flow speed has a nonlinear, increasing relationship with increasing compression ratios. For compression ratios under 0.7, corresponding to 70 % tube occlusion, U peak does not exceed c which is consistent with a low-amplitude approximation of peristalsis (Fig. 3b). However, for values 0.7 and above, the average peak flow speeds exceed c. Flow speeds within the compression tube Figure 4b reports u x,t measured at the center of the compression region as a function of dimensionless time for different compression ratios. Note that in all cases, significant backflow occurs as the tube rapidly expands behind the compression wave. Backflow is minimized at the highest compression ration of 0.9. Figure 5e reports the temporally averaged speeds, u m and u x , versus compression ratio. Higher compression ratios lead to higher forward flow speeds, as evidenced by the increase in u x relative to overall magnitude in Fig. 5b. u x exceeds the speed of the compression wave for the highest values (≥ 0.9). u m starts to exceed c as low as compression ratio = 0.6 (Fig. 4b).
Pressure Pressure differences between the inflow region ( Fig. 1a at 3) and outflow region (Fig. 1a at 4) are small for lower compression ratios and grow nonlinearly with increasing tube occlusion (Fig. 7b). The higher pressure differences correspond to the larger flow speeds generated by the higher compression ratios. Non-dimensional pressure p is reported versus time in Fig. 6b for varying values of compression ratio.

Speed of compression wave, c
Flow speeds in the racetrack In this set of simulations, the pulsing frequency is held constant and c is varied. U max,t is graphed against dimensionless time for several wave speeds in Fig. 3c. Note that changing c while holding the frequency constant has the effect of changing the wave form (see Fig. 2). The base case wave speed c = 3.0 was halved and doubled. As expected, lower values of c resulted in slow flows and higher values of c resulted in faster flows. Figure 5c presents U avg , U max , and U peak versus c for the upper section of the racetrack. All three speeds have a nonlinear and increasing relationship with increasing wave speed, with average peak speeds being consistently higher than c. There is a small dip in fluid speeds around value c = 3.33 which is consistent across all three measures of fluid speed. Flow speeds within the compression tube Figure 4c reports u x,t as a function of dimensionless time for different wave speeds. Note that there is significant backflow for the higher wave speed cases. In Fig. 5f, the temporally averaged speeds u m and u x measured at a point within the compression tube are plotted against the wave speed. All three speed metrics have a similar nonlinear relationship with c. Speeds increase overall with increasing c. Pressure Pressures across both the inflow region ( Fig. 1a at 3) and outflow region (Fig. 1a at 4) shares a nonlinear relationship with c, similar to the relationship seen between fluid flow speed and wave speed (Fig. 7c). Generally, increasing wave speed also increases the pressure difference between the two points. Note that the pressure difference dips around c = 3.33, just as U avg , U max , and U peak also drop at this value of c. Non-dimensional pressure p is reported versus time in Fig. 6c for varying values of c.

Compression wave frequency
Here, we consider the following two cases: (1) Frequency and wave speed are varied proportionally such that the waveform along the compression tube is unchanged, and (2) frequency is varied while the wave speed is held constant. This allows us to consider the cases when wave speed is coupled to and decoupled from changes in the pulse frequency.

Variable wave speed
Flow speed Figure 8a shows U avg , U max , and U peak as a function of the pulse frequency, f * , for the racetrack design. The gray data connected by lines show the case when the compression wave speed is allowed to change proportionally with f * . In this case, all measures of speed show a linear relationships with f * . The dashed line shows the wave speed, and U max is close to the wave speed across all frequencies considered. U peak also varies linearly with f * and is higher than the wave speed. Figure.  and u x (triangles) measured at a point within the compression region of the tube tionship is observed. The linear relationships between flow speeds (U avg and U peak ) and f for flows generated in the noracetrack model as well (Fig 9a). Additionally, U peak exceed c for the no-racetrack design.
Volume flow rate For the no-racetrack design, the volume flow rate, v avg , has a linear relationship with compression frequency, f , when c increases with increasing f (Fig. 9c).  Figure 7d shows the relationship between pressure, p , and compression frequency, f , for cases when the wave speed, c, changes proportionally with frequency (black items). p in shows a nonlinear but regular increase with increasing f , and p out shows a nonlinear decrease with increasing frequency. The pressure difference Δp between the inflow and outflow tract increases in a regular and nonlinear way with f .

Constant wave speed
Flow speeds in the racetrack When compression wave speed is held constant while the pumping frequency increases, all measures of fluid speed have markedly nonlinear relationships with compression frequency. Note that this has the effect of changing the wave form along the length of the compression region. Within the upper section of the racetrack, U avg , U max , and U peak show an increasing, oscillatory pattern with distinct peaks near f = 0.8, 1.6, and 2.5 (Fig. 8a). This pattern is more pronounced for the U peak . U peak greatly exceed the compression wave speed, c. U peak is triple the wave speed near f = 2.5. Flow speeds within the compression tube for the racetrack design u m and u x also share a nonlinear relationship with the contraction frequency, f . This graph is characterized by oscillations, but peak flow speeds within the compression tube do not correspond to the pattern seen within the racetrack (Fig. 8b).
Flow speeds for the no-racetrack design U avg and U peak show a similar nonlinear relationship with f when c is held constant (the racetrack and no-racetrack designs are compared in Fig. 9b). Peaks in both speeds are seen around f = 0.8, 1.6, and 2.5. Values of U peak are consistently above c. Volume flow rate For the no-racetrack design, the volume flow rate, v avg , has a nonlinear relationship with compression frequency, f , when c is held at 3.0 (Fig. 9d).
Pressure Figure 7d shows the relationship between pressures, p in , p out and Δp, and frequency, f , for constant wave speed (gray items). Pressures show a nonlinear relationship with frequency, as in the case of fluid flow speed. Männer et al. (2010) summarizes six qualities of technical peristaltic pumps often used to evaluate biological pumps based on small-amplitude and/or long-wave approximation of peristalsis. Our direct, numerical simulation of peristalsis adheres to a simpler, more inclusive definition of peristalsis: having a non-stationary compression wave that travels unidirectionally down a tube with no structurally fixed direction of flow. The geometry considered in this paper violates the long-wave approximation (the full wavelength is four times the diameter for the base case) and the small-amplitude approximation (compressions up to 95 % are considered). Furthermore, the Wo was varied from 0.1 to 50, spanning cases where both inertia and viscosity are significant. Note that this model was not developed to capture the specific dynamics of any particular tubular heart but rather to consider highly nonlinear dynamics.

Evaluating large-amplitude, short-wave peristaltic pumps against the technical definition of peristalsis
Our models incorporates the most basic features of peristalsis and produces flow characteristics that run counter to those listed in the technical definition. Our model generates pulsatile fluid flow for Wo of 5 and less at high compression ratios (Fig. 3a, b); flow speeds that show nonlinear relationships when either the compression wave speed or the pulsing frequency is fixed (Fig. 5c); and peak flow speeds that exceed compression wave speeds in most cases . These results are similar to the results of other valveless, peristaltic models with cardiac cushions (Taber et al. 2007). When flow speeds are measured in the compression section, significant backflow can be observed. Additionally, the model shows that flow speeds can have either a linear or nonlinear relationship with compression frequency, depending on how the speed of the compression wave is handled. A linear relationship between flow speed and compression frequency results when the speed of the compression wave increases linearly with compression frequency (Figs. 8a, 9a, c), demonstrating that flow is being driven by peristalsis and not another pumping mechanism. In this case, the wave form along the tube is the same, and the pump just operates faster. This is a feature of most mechanical pumps that does not necessarily translate to biological pumps. In the case of a constant wave speed, the relationship is strikingly nonlinear between fluid speed and compression frequency (Figs. 8a,9b,d). For these cases, the wave form along the compression region changes (demonstrated in Fig. 2). This scenario corresponds to situations when the pacemaker activity is altered while the propagation speed of an action potential is constant.

Implications of large-amplitude, short-wave peristaltic pumps in biological fluid transport systems
Since many of the principles of the technical definition derived from mechanical pumps are violated by our simple model, we suggest that the technical definition of peristalsis is inappropriate for evaluating the pumping mechanism of biological pumps. While the definition appropriately describes some peristaltic pumps, many biological pumps may fail simply because they violate the low-amplitude, long-wave approximation used to establish the criteria in the technical definition. Simply decoupling compression wave speed and frequency can also cause a pump to violate the technical definition even though it retains the essential features of peristalsis (traveling wave of compression and no structurally fixed flow direction) that mathematicians and biologists use. There are many biological pumps that fail criteria in the technical definition of peristalsis but fit the more inclusive definition. In tunicate hearts, heart beat frequency increases with increasing ambient temperature but the conduction velocity that passes down the heart tube remains constant (Kriebal 1967). In mosquito hearts, hemolymph flow is not continuous and flow speed is greater than compression wave speed (Glenn et al. 2010). For many tubular pumps, the contraction wave may nearly or even completely occlude the tube. Other tubular pumps, such as the embryonic hearts of tunicates, are not much longer than they are wide. Tubes with diameters on the order of their length violate the long-wave assumption.

Significance of the pumping mechanism to the embryonic heart
In the past 8 years, both peristalsis and dynamic suction pumping have been proposed as the mechanism by which the early embryonic heart tube drives the flow of blood. The important point of this difference for the case of understanding the evolution and development of the heart is the distinction between which regions of the heart actively contract. If tubular hearts pump using dynamic suction pumping, then only one region of the heart actively contracts, and the waves observed are passive and elastic. On the other hand, if the pumping mechanism is peristalsis, then the heart actively contracts along its entire length. This distinction has consequences for how the cardiac conduction system and the musculature develops. In Forouhar et al. (2006) work on embryonic zebrafish hearts, they present observations that demonstrate this valveless, tubular heart violates the technical definition of peristalsis. As we have shown with this simple model of flow being driven by peristalsis, fluid flow can mimic each observation made on the embryonic heart, including pulsatile flow both within the compressing tube (Fig. 4a-c) and far away from the compression tube (Fig. 3a-c); peak flow speeds that exceed compression wave speeds ( Fig. 5a-f, 9a, b); and a nonlinear relationship between fluid flow and compression frequency that is very similar to the relationship observed in their study (Fig. 8a, 9b). The results of our study suggest that peristalsis cannot be ruled out as the pumping mechanism of the vertebrate embryonic heart and that the exact pumping mechanism requires further study since both mechanisms can result in similar features for cases with nearly complete occlusion and diameter to length ratios on the order of 1:4.
The only observation that supports rejection of peristalsis as the pumping mechanism for the zebrafish heart is that no active contraction down the length of tube was observed (Forouhar et al. 2006). However, this observation has been disputed by others (Männer et al. 2010;Maes et al. 2011;Goenezen et al. 2012) and was not supported by direct measurement of muscle activation or morphological features that would suggest electrical activation of cardiac muscle was limited to one area of the heart. Many studies show that the architecture required to propagate such a signal are present early during heart development and that conduction throughout the myocardium is itself a critical component of normal cardiogenesis (Paff 1938;Paff et al. 1964;Rottbauer et al. 2001;Chi et al. 2010).