Analytical Solution of Oscillatory Stokes Flow in a Porous Pipe with Spatiotemporally Periodic Suction/Injection

Abstract

The study presents a generalized analytical solution of the laminar, oscillatory, creeping flow of an incompressible Newtonian fluid in a porous circular pipe with spatiotemporally periodic suction/injection at the wall. The analytical solution is examined for a variety of values of the dimensionless parameters, namely the Womersley number and the dimensionless suction/injection number.

1. Introduction

Flow in pipes with porous walls is of considerable importance for biological and industrial filtration processes. The mathematical modelling of steady flow in circular porous pipes was initially addressed by Yuan and Finkelstein assuming constant seepage velocity at the wall [1]. The Yuan and Finkelstein exact solution of the incompressible Navier–Stokes equations was obtained via a perturbation method for small and large suction/injection flow. Later, Terrill [2] obtained an exact solution for steady axisymmetric flow in a porous pipe with spatially variable suction/injection. Tsangaris et al. [3] extended the work of Terrill introducing unsteady wall injection and/or suction. The pulsatile flow in a porous tube with increasing seepage rate along the walls was studied analytically by Chang [4], using a linear approximation of the Navier–Stokes equations. Recently, Sindawi et al. [5] applied an approximate analytical method to solve the pulsatile flow in a porous tube with prescribed leakage flow rate which was considered to be smaller than axial flow rate.

Blood flow in arterioles and tubular fluid flow in the renal tubules and during hemodialysis are prominent examples of biological flows that are usually modelled as Stokes flows in porous pipes. Macey [6] modelled the tubular fluid flow devising a porous pipe with a reabsorption rate that varied linearly with distance, which is a better assumption about the renal filtration process than constant reabsorption rate, especially at the distal part of the tubules. Two years later, Macey [7] improved the renal flow modelling considering the exponential decay of the reabsorption rate in the proximal tubule, following the study of Kelman [8]. Kozinski et al. [9] extended the work of Macey for arbitrary variation in the seepage velocity using Fourier analysis. Bhatti et al. [10] derived an analytical solution of the flow in an infinitely long circular porous pipe with temporally periodic suction/injection. Siddiqui et al. [11] developed an exact analysis of steady Stokes flow through tubes with periodic permeable wall.

The Stokes flow in permeable annulus seems to fit in the theory of the recently discovered cerebrospinal fluid flow in the perivascular spaces of the brain, which is part of the brain’s glymphatic system [12]. The glymphatic system serves an important function in the clearance of waste, accomplished through the exchange of fluid between the perivascular spaces and the interstitial spaces [13]. Recent experimental results based on particle tracking velocimetry unveiled the basic fluid mechanics aspects of the glymphatic system [14,15]. Briefly, the cerebrospinal fluid flow in the perivascular spaces was found to be pulsatile, driven by the heartbeat and the arterial wall pulsations, termed also as perivascular pumping. Moreover, the cerebrospinal fluid flow is best described by Navier–Stokes equations as the perivascular system is an unobstructed fluid-filled space. The experimental results showed also that the arterial wall wave is not sinusoidal but rather follows the cardiac cycle dynamics with a fast expansion followed by a small contraction in accordance with the cardiac systole and diastole. A last notice is that the perivascular space is not best described as a refined circular annulus but shall be more eccentric and flattened. The fluid mechanics problem considered in the current study may be regarded as a model of the cerebrospinal fluid flow in perivascular spaces, as they share many aspects. Worth mentioning is that the fluid mechanics modelling of peristaltic flow in the glymphatic system is a mostly undiscovered area of research.

Our group has recently presented an analytic solution of oscillatory Stokes flow between porous plates, accounting for magnetohydrodynamic effects [16]. The purpose of the current study is to expand the analytical framework of Stokes flow problems in porous pipes, by considering unsteady oscillatory state and periodic suction/injection at the permeable wall. Previous solutions in the literature on the steady state case prove to be a limit of the derived, generalized solution that is introduced here.

2. Methods

2.1. Mathematical Formulation

A pipe with permeable wall is considered, as depicted in Figure 1. The suction/injection is considered to be periodic, with a period of $L=2\pi /\alpha$. The flow is assumed to be slow, laminar, and axisymmetric governed by the continuity and the Stokes equations, written in the cylindrical coordinate reference system. Furthermore, the fluid is assumed to be incompressible, with a constant density ρ, and Newtonian, with a constant dynamic viscosity μ. Introducing the stream function $\mathsf{\Psi }\left(r,z, t\right)$, the axial velocity ${\upsilon }_z$ and the radial velocity ${\upsilon }_r$ are defined as follows,

$${\upsilon }_z=\frac{1}{r}\frac{\partial \mathsf{\Psi }}{\partial r},{\upsilon }_r=-\frac{1}{r}\frac{\partial \mathsf{\Psi }}{\partial z}$$

(1)

Subsequently, the continuity equation is satisfied identically and the Stokes equations are written as follows,

$$\frac{\partial p}{\partial z}=-\rho \frac{1}{r}\frac{{\partial }^2\mathsf{\Psi }}{\partial r\partial t}+\mu \frac{1}{r}\frac{\partial }{\partial r}\left(E^2\mathsf{\Psi }\right)$$

(2)

$$\frac{\partial p}{\partial r}=\rho \frac{1}{r}\frac{{\partial }^2\mathsf{\Psi }}{\partial z\partial t}-\mu \frac{1}{r}\frac{\partial }{\partial z}(E^2\mathsf{\Psi })$$

(3)

is a second order differential operator, defined by [17],

$$E^2\equiv \frac{{\partial }^2}{\partial r^2}-\frac{1}{r}\frac{\partial }{\partial r}+\frac{{\partial }^2}{\partial z^2}$$

(4)

Eliminating pressure $p$ from Equations (2) and (3), we obtain a partial differential equation with respect to $\mathsf{\Psi }$ [18],

$$\rho \frac{\partial }{\partial t}\left(E^2\mathsf{\Psi }\right)=\mu E^2{E}^2\mathsf{\Psi }$$

(5)

The corresponding boundary conditions for the velocity components are

$${\upsilon }_z\left(R,z,t\right)=0, {\upsilon }_r\left(R,z,t\right)={\upsilon }_{0}sin\left(az\right)sin\left(\omega t\right) or {\upsilon }_{0}sin\left(az\right)con\left(\omega t\right)$$

(6a,b,c)

A constant pressure condition $p_o$ is assumed along the $z$ axis of the pipe,

$$\begin{array}{c} p\left(0,z,t\right)=p_0\end{array}$$

(7)

2.2. Problem Solution

We seek the analytical solution of the above problem, using the following Ansatz for the stream function Ψ,

$$\mathsf{\Psi }\left(r,z,t\right)=f\left(r\right) cos\left(az\right)e^{i\omega t}$$

(8)

where $i=\sqrt{-1}$.

According to the above expression, the velocity components can be written as

$${\upsilon }_z=\frac{1}{r}\frac{df}{dr}con\left(az\right)e^{i\omega t},{\upsilon }_r=a\frac{f}{r}sin\left(az\right)e^{i\omega t}$$

(9)

The $E^2$ and $E^2{E}^2$ operators in Equations (2), (3) and (5) reduce to the following expressions

$$E^2\mathsf{\Psi }\equiv \frac{{\partial }^2\mathsf{\Psi }}{\partial r^2}-\frac{1}{r}\frac{\partial \mathsf{\Psi }}{\partial r}+\frac{{\partial }^2\mathsf{\Psi }}{\partial z^2}=\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)con\left(az\right)e^{i\omega t}$$

(10)

$$E^2{E}^2\mathsf{\Psi }\equiv \frac{d^2}{d{r}^2}\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)-\frac{1}{r}\frac{d}{dr}\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)-a^2\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)con\left(az\right)e^{i\omega t}$$

(11)

So, the partial differential Equation (5) of the problem reduces to the following fourth order linear ordinary differential equation

$$\rho i\omega \left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)=\mu \left[\frac{d^2}{d{r}^2}\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)-\frac{1}{r}\frac{d}{dr}\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)-a^2\left(\frac{d^{2}f}{d{r}^2}-\frac{1}{r}\frac{df}{dr}-a^{2}f\right)\right]$$

(12)

By the introduction of the dimensionless variables,

$${\beta }^2=i\frac{\rho \omega }{\mu a^2},\sigma =ra$$

(13)

the differential Equation (12) takes the following form

$${\beta }^2\left(\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f\right)=\frac{d^2}{d{\sigma }^2}\left(\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f\right)-\frac{1}{\sigma }\frac{d}{d\sigma }\left(\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f\right)-\left(\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f\right)$$

(14)

and by substituting

$$F\equiv \frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f$$

(15)

it reduces to a linear Bessel equation of the form,

$$\frac{d^{2}F}{d{\sigma }^2}-\frac{1}{\sigma }\frac{dF}{d\sigma }-F\left({\beta }^2+1\right)=0$$

(16)

Using the variable transformation,

$$s=\sigma \sqrt{{\beta }^2+1}$$

(17)

the differential Equation (16) takes the form of a Bessel equation,

$$\frac{d^{2}F}{d{s}^2}-\frac{1}{s}\frac{dF}{ds}-F=0$$

(18)

The general solution of the above Bessel equation is [19],

$$F\left(s\right)=c_{3}s{I}_1\left(s\right)+c_{4}s{K}_1\left(s\right)=c_3\sigma \sqrt{{\beta }^2+1} I_1\left(\sigma \sqrt{{\beta }^2+1}\right)+c_4\sigma \sqrt{{\beta }^2+1} K_1\left(\sigma \sqrt{{\beta }^2+1}\right)$$

(19)

$I_1\left(s\right)$, $K_1\left(s\right)$ are the first-order-modified Bessel functions of the first and second kind, respectively. Substituting the expression (19) in (15), an inhomogeneous Bessel equation is derived,

$$\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f=c_3\sigma \sqrt{{\beta }^2+1} I_1\left(\sigma \sqrt{{\beta }^2+1}\right)+c_4\sigma \sqrt{{\beta }^2+1} K_1\left(\sigma \sqrt{{\beta }^2+1}\right)$$

(20)

The Bessel function of the second kind is rejected as it obtains infinite values at the pipe centerline, and Equation (20) reduces to,

$$\frac{d^{2}f}{d{\sigma }^2}-\frac{1}{\sigma }\frac{df}{d\sigma }-f=c_3\sigma \sqrt{{\beta }^2+1} I_1\left(\sigma \sqrt{{\beta }^2+1}\right)$$

(21)

The inhomogeneous Bessel Equation (21) has the following solution [19],

$$f\left(\sigma \right)=c_1\sigma I_1\left(\sigma \right)+c_2\sigma I_1\left(\sigma \sqrt{{\beta }^2+1}\right)$$

(22)

The above expression is the addition of the solution of the corresponding homogeneous equation and a partial solution of the inhomogeneous equation.

Rewriting the boundary conditions (6) for the velocity at the wall of the porous pipe in imaginary formulation,

$${\upsilon }_z\left(R,z,t\right)=0,{\upsilon }_r\left(R,z,t\right)={\upsilon }_{0}sin\left(az\right)e^{i\omega t}$$

(23a,b)

and accordingly with respect to f,

$$\frac{df\left(R\right)}{dr}=0,f\left(R\right)=\frac{{\upsilon }_0 R}{a}$$

(24a,b)

The constants $c_1$ and $c_2$ are determined by introducing the boundary conditions (24a,b) in the solution (22),

$$c_1=\frac{{\upsilon }_0}{a^2}\frac{k{I}_0\left(kRa\right)}{k{I}_0\left(kRa\right)I_1\left(Ra\right)-I_0\left(Ra\right)I_1\left(kRa\right)}$$

(25)

$$c_2=-\frac{{\upsilon }_0}{a^2}\frac{I_0\left(Ra\right)}{k{I}_0\left(kRa\right)I_1\left(Ra\right)-I_0\left(Ra\right)I_1\left(kRa\right)}$$

(26)

Finally, the analytical expressions for the velocity components are,

$$k=\sqrt{{\beta }^2+1}$$

(27)

$${\upsilon }_z={\upsilon }_{0}k\frac{I_0\left(kRa\right)I_0\left(\sigma \right)-I_0\left(Ra\right)I_0\left(k\sigma \right)}{k{I}_0\left(kRa\right)I_1\left(Ra\right)-I_0\left(Ra\right)I_1\left(kRa\right)}con\left(az\right)e^{i\omega t}$$

(28)

$${\upsilon }_r={\upsilon }_0\frac{k{I}_0\left(kRa\right)I_1\left(\sigma \right)-I_0\left(Ra\right)I_1\left(k\sigma \right)}{k{I}_0\left(kRa\right)I_1\left(Ra\right)-I_0\left(Ra\right)I_1\left(kRa\right)}sin\left(az\right)e^{i\omega t}$$

(29)

and in dimensionless form,

$${\overline{\upsilon }}_{\overline{z}}=\sqrt{{\epsilon }^2+i{\gamma }^2}\frac{I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)I_0\left(\overline{r}\epsilon \right)-I_0\left(\epsilon \right)I_0\left(\overline{r}\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}{\sqrt{{\epsilon }^2+i{\gamma }^2} I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)I_1\left(\epsilon \right)-\epsilon I_0\left(\epsilon \right)I_1\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}con\left(\epsilon \overline{z}\right)e^{i\tau }.$$

(30)

$${\overline{\upsilon }}_{\overline{r}}=\frac{\sqrt{{\epsilon }^2+i{\gamma }^2} I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)I_1\left(\overline{r}\epsilon \right)-\epsilon I_0\left(\epsilon \right)I_1\left(\overline{r}\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}{\sqrt{{\epsilon }^2+i{\gamma }^2} I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)I_1\left(\epsilon \right)-\epsilon I_0\left(\epsilon \right)I_1\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}sin\left(\epsilon \overline{z}\right)e^{i\tau }$$

(31)

where the dimensionless coordinates and parameters are defined as,

$$\overline{r}=\frac{r}{R},\overline{z}=\frac{z}{R},\tau =t\omega ,{\overline{\upsilon }}_{\overline{z}}=\frac{{\upsilon }_z}{{\upsilon }_0},{\overline{\upsilon }}_{\overline{r}}=\frac{{\upsilon }_r}{{\upsilon }_0},\epsilon =aR \left(suction/injection parameter\right),\gamma =R\sqrt{\frac{\rho \omega }{\mu }} \left(Womersley number\right)$$

(32)

The mean dimensionless velocity over the cross-section is calculated by the integration of the axial velocity component, resulting in the following relation,

$${\overline{\upsilon }}_{\overline{z}m}=\frac{{\upsilon }_{zm}}{{\upsilon }_o}=\frac{2}{\epsilon }con\left(\epsilon \overline{z}\right)e^{i\tau }$$

(33)

The pressure distribution is calculated by integrating Equations (2) and (3). The boundary condition and the calculated pressure distribution in dimensionless form are given by the following expressions assuming $p_0=0$,

$$\overline{p}\left(0,\overline{z},\tau \right)=0$$

(34)

$$\overline{p}\left(\overline{r},\overline{z},\tau \right)=-i\frac{{\gamma }^2}{\epsilon }\frac{\sqrt{{\epsilon }^2+i{\gamma }^2} I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}{\sqrt{{\epsilon }^2+i{\gamma }^2} I_0\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)I_1\left(\epsilon \right)-\epsilon I_0\left(\epsilon \right)I_1\left(\sqrt{{\epsilon }^2+i{\gamma }^2}\right)}\left[I_0\left(\epsilon \overline{r}\right)-1\right]sin\left(\epsilon \overline{z}\right)e^{i\tau }$$

(35)

where $\overline{p}$ is the dimensionless pressure defined as,

$$\overline{p}=\frac{p}{P},P=\frac{\mu {\upsilon }_0}{R}$$

(36)

Rewriting Equations (30), (31) and (35) as,

$${\overline{\upsilon }}_{\overline{z}}=\left({\overline{V}}_{\overline{z}R}+i{\overline{V}}_{\overline{z}I}\right)con\left(\epsilon \overline{z}\right)e^{i\tau }$$

(37)

$${\overline{\upsilon }}_{\overline{r}}=\left({\overline{V}}_{\overline{r}R}+i{\overline{V}}_{\overline{r}I}\right)sin\left(\epsilon \overline{z}\right)e^{i\tau }$$

(38)

$$\overline{p}=\left({\overline{P}}_{aR}+i{\overline{P}}_{aI}\right)sin\left(\epsilon \overline{z}\right)e^{i\tau }$$

(39)

and taking into consideration the imaginary part of the relations, which correspond to the boundary conditions (6a,b), we arrive at the following relations for the velocity components and the pressure,

$${\overline{\upsilon }}_{\overline{z}}={\overline{V}}_{\overline{z}amp}sin\left(\tau +{\phi }_{\overline{z}}\right)con\left(\epsilon \overline{z}\right), {\overline{V}}_{\overline{z}amp}=\sqrt{{\overline{V}}_{\overline{z}R}{}^2+{\overline{V}}_{\overline{z}I}{}^2}, tan\left({\phi }_{\overline{z}}\right)=\frac{{\overline{V}}_{\overline{z}I}}{{\overline{V}}_{\overline{z}R}}$$

(40)

$${\overline{\upsilon }}_{\overline{r}}={\overline{V}}_{\overline{r}amp}sin\left(\tau +{\phi }_{\overline{r}}\right)sin\left(\epsilon \overline{z}\right), {\overline{V}}_{\overline{r}amp}=\sqrt{{\overline{V}}_{\overline{r}R}{}^2+{\overline{V}}_{\overline{r}I}{}^2}, tan\left({\phi }_{\overline{r}}\right)=\frac{{\overline{V}}_{\overline{r}I}}{{\overline{V}}_{\overline{r}R}}$$

(41)

$$\overline{p}={\overline{P}}_{amp}sin\left(\tau +{\phi }_{\overline{p}}\right)sin\left(\epsilon \overline{z}\right), {\overline{P}}_{amp}=\sqrt{{\overline{P}}_{aR}{}^2+{\overline{P}}_{aI}{}^2}, tan\left({\phi }_{\overline{p}}\right)=\frac{{\overline{P}}_{aI}}{{\overline{P}}_{aR}}$$

(42)

$${\overline{\upsilon }}_{\overline{z}m}=\frac{2}{\epsilon }sin\left(\tau \right)con\left(\epsilon \overline{z}\right), {\overline{\upsilon }}_{\overline{z}m, amp}=\frac{2}{\epsilon }, {\phi }_{\overline{z}m}=0$$

(43)

The obtained solution, due to the linear form of the creeping flow equations, can be straightforwardly added to the Poiseuille solution of steady flow and/or to the oscillatory flow in a circular pipe. As a limit of the present solution for very small Womersley number ($\gamma \to 0$), the steady state solution of Siddiqui et al. [11] can be obtained.

3. Results and Discussion

Calculations were carried out for different values of the parameters $\gamma$ and $\epsilon$ entering the problem and the results are depicted graphically. Figure 2 presents the velocity field at various time instances for three values of the Womersley number. The periodicity of the velocity is evident both spatially and temporally. Peak values migrate towards the pipe wall with increasing oscillation frequency, while the velocity acquires a flat profile around the pipe axis. A dominant characteristic of the flow, clear in the stream density plots of Figure 2, is the recurrent recirculation zones with alternating direction every half a time period.

Similar observations are made in Figure 3, which depicts the amplitude profile of axial and radial velocity components, and pressure, for three different values of the Womersley parameter (γ = 1, 5, 10) and the suction/injection parameter (ε = 1, 2, 3). For small values of Womersley parameter (γ = 1) and suction injection parameter (ε = 1), the axial velocity acquires an oscillating, parabolic-like profile (quasi steady Poiseuille flow), with mean value ${\overline{\upsilon }}_{\overline{z}m}=2sin\left(\tau \right)con\left(\overline{z}\right)$ and maximum value approaching ${\overline{\upsilon }}_{\overline{z}m}=4sin\left(\tau \right)con\left(\overline{z}\right)$. It can be noticed that the amplitude of the axial velocity reduces with a factor of 2/ε. For larger values of the Womersley parameter, the axial velocity profile flattens around the pipe axis, while higher values appear near the pipe wall. Such behavior is typical for oscillating pipe flows at high frequencies and is due to the boundary layer character of the flow near the pipe wall [20,21]. Regarding the radial velocity, its peak value shifts closer to the pipe wall with increasing values of the Womersley parameter.

The pressure amplitude is an increasing function of the Womersley parameter (Figure 3). A physical explanation of the pressure amplitude increase under increasing Womersley number lies on the balance of inertial, pressure and viscous forces dictated by the law of momentum conservation. An oscillatory flow with high frequency leads to the increase in unsteady inertial forces that is compensated mainly by the increase in pressure. It is noticeable that, for γ = 1, the pressure amplitude increases for increasing ε values, while the reverse behavior is observed for higher γ values. The reversal of the pressure amplitude profiles under increasing suction/injection parameter takes place gradually within the interval 3 < γ < 4. After that, any increase in the frequency leads to bigger differences between the pressure amplitudes for the various ε values. Physically, the latter effect is a result of high suction/injection frequency both in space and time. Τhe reversed pattern of pressure (but not of velocity), predicted by the solution, is an interesting finding that is worth further investigation.

The phase angle profiles of the axial and radial velocity components, and pressure are shown in Figure 4. For small values of γ, the phase angle is close to zero, for all ε values. The phase angles of ${\phi }_{\overline{z}}$ and ${\phi }_{\overline{r}}$ have the same value at the axis of the pipe and different values at the pipe wall, where ${\phi }_{\overline{r}}=0$, according to the applied boundary condition. For large values of Womersley number γ, the phase angle ${\phi }_{\overline{z}}$ near the pipe wall approaches to π/4, for all values of ε. The phase angle of pressure is constant along the $r$ axis, because the expression (42) is independent of $\overline{r}$. When γ tends to very small values, ${\phi }_{\overline{p}}\to 0$, and when γ tends to very high values, the pressure phase angle approaches a value ${\phi }_{\overline{p}}\to \pi$. Overall, the phase angle is also an increasing function of the Womersley number like the pressure amplitude.

4. Conclusions

The current study presents an analytical solution that generalizes the literature on flows driven by the suction/injection occurring at the wall of permeable circular pipes. Introducing the stream function, the momentum equation reduces to a fourth order ordinary differential equation and after several transformations, to an inhomogeneous Bessel equation that is solved analytically. Here, the suction/injection is considered to be periodic both in space and time generating an oscillatory flow with recurrent recirculation zones and points of flow separation. The increase in suction/injection temporal frequency, controlled by the Womersley number, results in the shift in peak velocity towards the permeable wall and in the increase in pressure. Concerning the spatial frequency, the increase in the suction/injection parameter leads to a decrease in velocity for all Womersley numbers. For small Womersley numbers, an increase in pressure is noticed. The combined increase in both parameters was found to be causing a reduction in both velocity and pressure. The problem description mostly fits in biofluidic systems, such as the renal tubules and the glymphatic system in the brain, where fluid exchange seems to be taking place between permeable interfaces. The mathematical model along with its analytic solution could serve as a basis for the study of the aforementioned biofluidic systems.