2.1. Problem formulation

The model developed in Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) is based on lubrication theory. More specifically, the following assumptions have been made.

1. (i) Planar geometry: it is assumed that ureteral peristalsis is planar and can be qualitatively modelled by two-dimensional geometry.

2. (ii) Long wavelength: it is assumed that the peristaltic wavelength is much larger than the lumen height, a necessary condition for applicability of lubrication theory. This is justified by the anatomy of the ureter in a typical adult person as noted in table 1.

3. (iii) Inertia-free flow: the Reynolds number is assumed to be sufficiently small.

4. (iv) Fluid properties: the fluid is assumed to be Newtonian with constant density and viscosity.

5. (v) Infinite sinusoidal wavetrain: the flow is periodic and the peristaltic wave is assumed to be an infinite progressive train of sinusoidal waves, moving longitudinally at a constant velocity. The more realistic peristaltic waveform described by Lykoudis & Roos (Reference Lykoudis and Roos1970) is eschewed in favour of simplicity and comparison to the majority of previous studies. Also, the periodic boundary condition assumed here is challenged for a small number of waves as pointed out by Li & Brasseur (Reference Li and Brasseur1993).

Table 1. Normal anatomical information of a typical ureter and peristalsis in an adult person (Zheng et al. Reference Zheng, Carugo, Mosayyebi, Turney, Burkhard, Lange, Obrist, Waters and Clavica2021).

Using standard non-dimensionalization, we set the wave velocity $c$, fluid viscosity $\mu$ and lumen half-height $a$ to unity (see Appendix C for more detail). In the laboratory frame the problem is unsteady and the wall transverse coordinate $H$ is given by

(2.1)\begin{equation} H = 1 - b \cos(k X - \omega t), \end{equation}

where $b$ is the peristaltic amplitude, $X$ is the longitudinal coordinate in the laboratory frame, $t$ is time, $k = 2 {\rm \pi}/ \lambda$ is the wavenumber, $\omega = 2 {\rm \pi}/ T$ is the angular frequency and $T$ is the waveperiod. The lubrication approximation is valid for small normalized wavenumbers, $k \ll 1$. Attaching the frame of reference to the travelling wave and assuming the pressure gradient between the ‘ends’ (i.e. the bladder and kidney) remains constant (i.e. assuming periodic boundary conditions), the flow becomes steady, and the wall transverse coordinate becomes

(2.2)\begin{equation} h = 1 - b \cos kx. \end{equation}

The transformation between the two frames is given by

(2.3)\begin{equation} \left.\begin{array}{c@{}} x = X - t, \quad y = Y, \\ u(x, y) = U(X - t, Y) - 1, \quad v(x, y) = V(X - t, Y), \end{array}\right\} \end{equation}

where $(x, y)$ and $(X, Y)$ are parameterizing longitudinal and transverse coordinates in the wave and laboratory frames of reference, respectively. Similarly, $(u, v)$ and $(U, V)$ are longitudinal and transverse velocity components in the wave and laboratory frames of reference. A schematic of a peristalsis wave, as described above, in the steady frame of reference with non-dimensional variables is shown in figure 1.

Figure 1. A schematic of a peristalsis wave in the wave frame (moving left to right with the peristalsis wave) with non-dimensional variables. The dashed lines represent peristalsis of zero amplitude ($b = 0$). The dashed-dotted line is the symmetry plane.

With the lubrication assumptions, the Navier–Stokes equation in the longitudinal and transverse directions becomes

(2.4a,b)\begin{equation} \frac{\partial p}{\partial x} = \frac{\partial^2 u}{\partial{y}^2}, \quad \frac{\partial p}{\partial y} = 0. \end{equation}

To account for wall longitudinal motion, we set $u(y = h) = -1 + \beta (x)$, where $\beta (x)$ is a $2 {\rm \pi}$-periodic function in space and time (in the Lagrangian sense). Note that $\beta = 0$ recovers the Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) solution. Therefore, the boundary conditions become

(2.5)\begin{equation} \left.\begin{array}{c@{}} \displaystyle \dfrac{\partial u}{\partial y} = 0 \quad\text{ at } y = 0\ (\text{symmetry condition}), \\ u ={-}1 + \beta(x) \quad\text{ at } y= h \ (\text{imposed wall velocity}). \end{array}\right\} \end{equation}

While the correct form of $\beta (x)$ should be determined based on clinical experiments, these quantitative data are not yet available, and we defer investigating a more realistic and clinical form of $\beta (x)$ for future studies. Here, we examine general effects of wall longitudinal motion on peristalsis, to do so we assume $\beta$ is in the form of periodic Gaussian as follows:

(2.6)\begin{equation} \beta = \beta_0 + \sum_{m ={-}\infty}^{\infty} \beta_1 \exp{\frac{-(kx - m k\lambda - \phi)^2}{2 \beta_2^2}}. \end{equation}

Here $\phi, \beta _1$ and $\beta _2$ are parameters controlling phase, amplitude and width of $\beta$, respectively. Parameter $\beta _0$ is a constraint ensuring that the motion of the wall in the wave frame remains $2{\rm \pi}$-periodic in time (in the Lagrangian sense). Parameter $\beta _0 (b, \beta _1, \beta _2$) must be obtained such that the wall in the wave frame travels one wavelength in one waveperiod, i.e. to satisfy

(2.7)\begin{equation} \int_0^{\omega T = 2{\rm \pi}} u(x, y = h)\, {\rm d}t ={-}k \lambda ={-}2{\rm \pi}. \end{equation}

The summation in (2.6) is required to make $\beta$ periodic. We find that for the range of $\beta _2$ considered in this study, it is sufficient to truncate the summation from $m = -5$ to $5$. The form of $\beta$ given in (2.6) is far from general but allows consideration of localized wall motion near the lumen extremas with phase lags.

2.2. Flow kinematics

The problem given by (2.4a,b) can be rewritten in terms of a stream function $\psi$ (satisfying $u = \partial \psi / \partial y, v = -\partial \psi / \partial x$) as

(2.8)\begin{equation} \frac{\partial^4 \psi}{\partial y^4} = 0. \end{equation}

Boundary conditions given by (2.5) can also be expressed in terms of the stream function as

(2.9)\begin{equation} \left.\begin{array}{c@{}} \displaystyle \dfrac{\partial^2 \psi}{\partial{y}^2} = 0\quad \text{ at } y = 0,\\ \displaystyle \dfrac{\partial \psi}{\partial y} ={-}1 + \beta(x)\quad \text{ at } y = h. \end{array}\right\} \end{equation}

The two remaining conditions needed to solve (2.8) can be obtained by noting that in the wave frame both planes $y = 0$ and $y = h$ are streamlines and the volume flow rate $q$ is constant at all cross-sections:

(2.10)\begin{equation} \left.\begin{array}{c@{}} \psi = 0 \quad\text{ at } y = 0, \\ \psi = q \quad\text{ at } y = h. \end{array}\right\} \end{equation}

Solving (2.8)–(2.10) gives $u$ and $v$ as

(2.11)\begin{gather} u = \frac{3}{2} \frac{h (1 - \beta) + q}{h} \left[1 - \left(\frac{y}{h}\right)^2 \right] - (1 - \beta), \end{gather}

(2.12)\begin{gather}v =\frac{h\beta' - h'(1 - \beta)}{2} \frac{y}{h} \left[ 3 - \left(\frac{y}{h}\right)^2 \right] - \frac{3}{2} \frac{h(1 - \beta) + q}{h} \frac{y}{h} h' \left(\frac{y}{h} - 1\right)^2 - y \beta ', \end{gather}

where the primes denote derivative with respect to $x$. An example of longitudinal velocity $u$ for $\beta = 0$ – Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) case – is shown in Appendix D, figure 13. The parabolic, Poiseuille type flow can be seen in the figure.

2.3. Pressure-flow characteristics

Considering that the peristaltic wave acts as a pump (i.e. working on the fluid to increase its pressure), two quantities of practical interest are the time-mean volume flow rate at each cross-section $\bar {Q}$, which measures the net volume flow or mean discharge rate (per lumen width), see Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) for derivation, and the pressure rise per wavelength, ${\rm \Delta} p_\lambda$:

(2.13)\begin{gather} \bar{Q} = q + 1, \end{gather}

(2.14)\begin{gather}{\rm \Delta} p_{\lambda} = \int_{0}^{k\lambda = 2 {\rm \pi}} {\frac{{\rm d} p}{{\rm d}\kern 0.06em x}}\, {{\rm d}\kern 0.06em x}, \end{gather}

where $\textrm {d}p/{\textrm {d}\kern0.06em x} = \partial p /\partial x$ can be obtained from (2.4a,b) and (2.11) as

(2.15)\begin{equation} \frac{{\rm d} p}{{\rm d}\kern 0.06em x} = \frac{-3}{h^3} [q + (1 - \beta ) h ]. \end{equation}

For complete lumen closure (i.e. $b = 1$), $\bar {Q} = 1$ and the peristaltic wave acts as a positive-displacement pump while, for no net flow, $\bar {Q} = 0$.

As is evident from (2.15) and consistent with Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969), due to the low-Reynolds-number assumption employed in this study, ${\rm \Delta} p_\lambda$ has a linear relationship with $\bar {Q}$; see figure 14(a) in Appendix D. Also, as noted by Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969), this linear relationship can be characterized by horizontal and vertical intercepts of these curves with the following physical significance.

1. (i) Vertical intercepts with $\bar {Q} = 0$ are at ${\rm \Delta} p^0_\lambda$ and indicate a pressure rise per wavelength for zero time-mean flow that is similar to when there is a blockage in the ureter or if the pressure difference between the bladder and the kidney are such that the net flow is zero.

2. (ii) Horizontal intercepts with ${\rm \Delta} p_\lambda = 0$ are at $\bar {Q}_0$ and indicate a time-mean flow for a zero pressure rise per wavelength. At this flow rate the peristaltic wave is not doing useful work and flow is not controlled by the peristalsis. Under such conditions, depending on differential pressure between the kidneys and bladder, we might have no flow, forward flow or backward flow. In these situations, all the work done by the peristaltic wave will be converted to heat by viscous dissipation.

The above intercepts mark the range where the peristalsis acts as a pump, defined as the range of volume flow rates where mean flow is in the direction of the pressure rise per wavelength – non-negative $\bar {Q}$ and non-negative ${\rm \Delta} p_\lambda$ in figure 14(a) in Appendix D – i.e. $0 \le \bar {Q} \le \bar {Q}_0$. Consistent with similar works in the literature (Shapiro et al. Reference Shapiro, Jaffrin and Weinberg1969; Liron Reference Liron1976), the peristalsis activity in this range is the focus here. For convenience, $\bar {Q}$ is normalized by its pumping range limit $\bar {Q}_0$ (see figure 14b in Appendix D):

(2.16)\begin{equation} \varGamma = \frac{\bar{Q}}{\bar{Q}_0}. \end{equation}

2.4. Peristaltic reflux

Peristaltic reflux is a topic of controversy in the literature. The debate is whether to use a Lagrangian description, advocated by Shapiro (Reference Shapiro1967), or an Eulerian description, promoted by Fung & Yih (Reference Fung and Yih1968), who argue that the Eulerian description gives the same result as the Lagrangian. Despite Shapiro & Jaffrin (Reference Shapiro and Jaffrin1971) definitively showing the difference between the two descriptions and concluding that the Lagrangian description is the only correct one, the Eulerian description is still being commonly misused (Abd Elnaby & Haroun Reference Abd Elnaby and Haroun2008; Jiménez-Lozano, Sen & Dunn Reference Jiménez-Lozano, Sen and Dunn2009; Vahidi et al. Reference Vahidi, Fatouraee, Imanparast and Moghadam2011; Gad Reference Gad2014; Hosseini et al. Reference Hosseini, Ji, Xu, Rezaienia, Avital, Munjiza, Williams and Green2018). In such analysis, a time-mean velocity of less than zero is incorrectly identified as the criteria for peristaltic reflux. Furthermore, they predict reflux to happen in the centre of the lumen, which is not consistent with experimental results of Weinberg et al. (Reference Weinberg, Eckstein and Shapiro1971), showing reflux occurs near the wall. As shown previously (Shapiro et al. Reference Shapiro, Jaffrin and Weinberg1969; Takabatake & Ayukawa Reference Takabatake and Ayukawa1982) and in this study, peristaltic reflux can happen even when $\bar {Q} > 0$. Therefore, consistent with Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969), Shapiro & Jaffrin (Reference Shapiro and Jaffrin1971), Takabatake & Ayukawa (Reference Takabatake and Ayukawa1982), we define peristaltic reflux as retrograde motion of fluid particles, resulting in a net backward displacement of more than one wavelength over one waveperiod. The quantity of interest becomes the trajectories of fluid particles over one waveperiod, which can be obtained by solving an initial value problem given by

(2.17a,b)\begin{equation} \frac{{\rm D}x}{{\rm D}t} = u, \quad \frac{{\rm D}y}{{\rm D}t} = v, \end{equation}

where $\textrm {D}/\textrm {D}t$ is the total derivative and $u$ and $v$ are given by (2.11) and (2.12), respectively. The particle trajectories in the laboratory frame can be simply obtained by the coordinate transformation given in (2.3), i.e. $X = x + t$, $Y = y$. Examples of these trajectories showing refluxing and non-refluxing fluid particles, based on the above definition, are shown in figures 2(a) and 2(b) in the wave and laboratory frames, respectively. Note that the fluid particle starting closer to the wall refluxes while the fluid particle closer to the lumen centre does not. The corresponding animation of particle trajectories of figure 2(b) is shown in supplementary movie 1 available at https://doi.org/10.1017/jfm.2023.363. As evident in figure 2, the wave and particle periods are not the same. In the above figures, parameters $b = 0.8$ and $q = -0.7$ are chosen as an example case and are based on clinical expertise and anatomical measurements (see § 4 for more detail) and, unless stated otherwise, the same parameters are used for the remainder of the paper.

Figure 2. Example peristaltic reflux illustration. (a) Fluid pathlines are shown starting at $kx/2 {\rm \pi}= 0$ over one waveperiod ($0 < \omega t /2 {\rm \pi}< 1$) in the wave frame. The refluxing pathlines are shown as red dashed lines. The timeline originally vertically spanning the small gap at $kx/2 {\rm \pi}= 0$ is shown with a green dotted line at $\omega t / 2 {\rm \pi}= 1$. (b) Two examples of pathlines in the laboratory frame starting at $kX/2{\rm \pi} = 0, Y = 0.15$ (with net negative displacement, i.e. refluxing) and $kX/2{\rm \pi} = 0, Y = 0.05$ (with net positive displacement, i.e. not refluxing). The associated particles in (b) are marked in (a) as well. The empty circles show initial locations while the filled circles are showing locations at the end of one waveperiod. The results are shown for $\beta = 0$, i.e. no wall longitudinal motion. The wall in (b), i.e. in the lab frame, translates and is shown at $\omega t / 2 {\rm \pi}= 1$. The corresponding animation of (b) is shown in supplementary movie 1.

In light of the above discussion, a more rigorous criteria for peristaltic reflux can be defined based on the mean velocity of the fluid particles over one waveperiod ($\omega T = 2 {\rm \pi}$). The mean velocity of fluid particles $\bar {u}$ (or $\bar {U}$ in the laboratory frame) can be obtained from the horizontal motion of particles over one waveperiod ${\rm \Delta} x_\lambda$ as

(2.18a,b)\begin{equation} \bar{u} = \frac{{\rm \Delta} x_\lambda}{2 {\rm \pi}}, \quad \bar{U} = \bar{u} + 1. \end{equation}

Mean velocities of $\bar {u} < -1$ in the wave frame (or $\bar {U} < 0$ in the laboratory frame) indicate fluid particles with a net backward motion of more than one wavelength over one waveperiod, i.e. refluxing condition, this is shown in figure 3. Therefore, the range of refluxing flow rates $\varGamma _r$ can be defined as

(2.19)\begin{equation} \varGamma_{r} = \{\varGamma \in (0, 1) \mid \bar{U} (\varGamma) < 0 \}. \end{equation}

Figure 3. Reflux criteria. Average longitudinal velocity of fluid particles at $kX / 2{\rm \pi} = 0$ over one waveperiod ($0 < \omega t / 2 {\rm \pi}< 1$) in the laboratory frame. The dashed red line is showing refluxing fluid particles. The forward and backward (reflux) time-mean volume flow rate can be calculated by integrating $\bar {U}$ with respect to $Y$ for $\bar {U} > 0$ and $\bar {U} < 0$, respectively. This figure is identical to the timeline of figure 2(a) at $\omega t / 2 {\rm \pi}= 1$ (green dotted line), albeit at a different scale.

2.5. Reflux quantification

Consistent with Takabatake & Ayukawa (Reference Takabatake and Ayukawa1982), since we use the waveperiod not the particle period for calculating the mean velocity of figure 3, we can integrate $\bar {U}$ with respect to $Y$ in the refluxing range to calculate the time-mean backward (i.e. refluxing) volume flow rate per unit width in the out-of-plane direction:

(2.20)\begin{equation} \bar{Q}_{r} = \int_{\bar{U} < 0}{\bar{U}\, {\rm d}Y}. \end{equation}

Similarly, the time-mean forward (non-refluxing) volume flow rate $\bar {Q}_{nr}$ can be calculated as

(2.21)\begin{equation} \bar{Q}_{nr} = \int_{\bar{U} > 0}{\bar{U} \,{\rm d}Y}. \end{equation}

The flow rates $Q_r$ and $Q_{nr}$ are represented by the areas enclosed by the red dashed line bounded by $\bar {U} = 0$, and by the solid curve bounded by $Y = 0$ and $\bar {U} = 0$, respectively. The algebraic sum of the refluxing and non-refluxing flow rates gives the net time-mean flow rate $\bar {Q} = \bar {Q}_r + \bar {Q}_{nr}$. Lastly, the reflux fraction ${\mathcal{X}}$ can be defined as the ratio of magnitudes of the refluxing to net volume flow rates, i.e.

(2.22)\begin{equation} {\mathcal{X}} = \frac{| \bar{Q}_r |}{\bar{Q}}. \end{equation}

2.6. Efficiency

In the peristaltic wave pumping range ($0 \le \varGamma \le 1$), the mechanical efficiency $E$ can be obtained based on its common definition in pumping machinery as

(2.23)\begin{equation} E = \frac{\bar{Q} {\rm \Delta} p_\lambda}{\bar{W}}, \end{equation}

where $\bar {Q} {\rm \Delta} p_\lambda$ is the useful pumping power and $\bar {W}$ is the mean rate of wall mechanical work, both per wavelength and per unit width (Shapiro et al. Reference Shapiro, Jaffrin and Weinberg1969; Liron Reference Liron1976). In the above equation $\bar {W}$ can be calculated by noting that it is the sum of the rate of useful pumping power ($\bar {Q} {\rm \Delta} p_\lambda$) and viscous dissipation $D$:

(2.24)\begin{equation} \bar{W} = \bar{Q} {\rm \Delta} p_\lambda + D, \end{equation}

where $D$, based on the small slope approximation, is

(2.25)\begin{equation} D = \int_0^{2 {\rm \pi}} \int_0^h \left(\frac{{\rm d} p}{{\rm d}\kern 0.06em x} y\right)^2 \, {\rm d}y \,{\rm d}\kern 0.06em x. \end{equation}

We note that efficiencies outside this pumping range (where $\bar {W}$ changes sign) have some anomalies that we do not yet understand. Peristaltic efficiency and reflux fraction in the pumping range, calculated numerically with (2.23)–(2.25), and (2.20)–(2.22), are plotted in figures 4(a) and 4(b), respectively for $\beta = 0$, i.e. the Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) case. For comparison purposes, the Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) results – their (17) and (51) for efficiency and reflux fraction, respectively – are overlaid on top of the results obtained here.

Figure 4. Peristaltic efficiency and reflux. (a) Pump efficiency and (b) reflux fraction plotted against normalized volume flow rate for two different peristaltic amplitudes $b = 0.5$ and $0.8$. The results are obtained with $\beta = 0$ (without wall longitudinal motion). The lines show the results from the analytical model (Shapiro et al. Reference Shapiro, Jaffrin and Weinberg1969) – their (17) and (51) – while the discrete points are obtained with our numerical model. The refluxing conditions in (a) are shown with empty circles and dotted lines. The short vertical bar in (a) is the limit of reflux as predicted by the Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) model (the minimum needed normalized volume flow rate to prevent reflux). The numerical and analytical solutions match for the case of $\beta = 0$ – the discrepancy in (b) is due to the expansion solution approximation employed by Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969). For clarity, the discrete points with zero reflux fractions are omitted in (b). In (b) the example critical reflux fraction ${\mathcal{X}}_c = 0.02$, described in the text, is shown with a dashed horizontal line as an example.

The spatial resolution for numerical integrations involved in (2.23) is ${\rm \Delta} x = {\rm \Delta} y/h = 0.01$. The ordinary differential equation (2.17a,b) for reflux quantification (2.20), on the other hand, is solved with a relative tolerance of $10^{-5}$ and spatial resolution of ${\rm \Delta} y/h = 0.01$. The average error for efficiency compared with the closed-form solution of Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969), independent of the peristaltic amplitude $b$, is of the order of ${\sim }10^{-6}$. The closed-form solution of Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) for the reflux fraction is an expansion for small $b$ (see their figure 13, where they compare the closed-form solution to the more accurate numerical integration results). Therefore, for the reflux fraction, the numerical error is estimated with the solution for finer temporal and spatial resolutions (relative tolerance of $10^{-6}$ and ${\rm \Delta} y/h = 0.005$). This numerical error is of the order of ${\sim }10^{-5}$ for the largest $b$ considered ($b = 0.8$) and it decreases for smaller peristaltic amplitudes (${\sim } 10^{-8}$ for $b = 0.3, 0.5$).

In both figures the lines show Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) results while the discrete points indicate the results obtained here. In figure 4(a) (i.e. efficiency curves) the range of refluxing $\varGamma$, obtained from (2.19), are marked with empty circles. Similarly, the range of refluxing $\varGamma$ predicted by Shapiro et al. (Reference Shapiro, Jaffrin and Weinberg1969) – their (47) – are indicated with dotted lines ending with a vertical bar. In figure 4(b) the additional horizontal dashed line, indicating an example for a critical reflux fraction ${\mathcal{X}}_c$, is also plotted. This is because, although not well understood, small amounts of peristaltic reflux can be harmless, especially if they are not accompanied by other abnormal conditions such as high bladder pressure, poor contractility or a UTI (Boyarsky & Labay Reference Boyarsky and Labay1981). Therefore, the more important quantity for clinical applications would be the reflux volume fraction. The critical reflux fraction, beyond which the peristaltic reflux becomes harmful, depends on clinical scenarios and varies on a case-by-case basis. In this study and as shown in figure 4(b), we choose the critical reflux fraction of ${\mathcal{X}}_c = 0.02$ as an example. The following important observations can be made from figure 4.

1. (i) In the absence of ureteral wall longitudinal motion, peristaltic reflux is common for flow rates within the pumping range.

2. (ii) Depending on the peristaltic amplitude, there is a minimum flow rate required for preventing peristaltic reflux. This minimum required flow rate increases with increasing peristaltic amplitude.

3. (iii) The reflux fraction increases with peristaltic amplitude.

4. (iv) As the net volume flow rate increases, the reflux fraction decreases and, depending on the peristaltic amplitude, can go to zero.

5. (v) As the net volume flow rate decreases, the reflux fraction increases and at the limit $\varGamma \rightarrow 0$ it goes to infinity. This is because while the net flow goes to zero, the backward flow remains finite.

6. (vi) As the peristaltic amplitude increases ($b \rightarrow 1$), the maximum peristalsis efficiency increases but at the cost of an increase in the refluxing range and refluxing volume fractions.