On the Asymptotic Boundary Condition at the Free-Fluid/ Porous-Medium Interface in Periciliary Layer due to the Ciliary Movement

/e purpose of this study is to study the movement of the periciliary (PCL) fluid due to the ciliary locomotion. In this research, because the bundle of cilia is considered instead of a single cilium, Stokes–Brinkman equations in a macroscopic scale are employed to find the velocity of the PCL fluid. When the cilia are perpendicular to the horizontal plane, the PCL consists of only the cilia./e inclination of cilia (cilia make an angle θ(θ< 90∘) to the horizontal plane) results in two different domains in the PCL, the regions comprising and not comprising cilia. /e main objective of this study is to determine the appropriate boundary conditions of the velocity between these two regions where the PCL fluid is moved by self-propelled cilia rather than the pressure gradient. A matched asymptotic expansion method is applied to the Stokes–Brinkman equations to determine the constraints. Two boundary conditions at the interface are obtained and the velocity of the PCL fluid at the top of PCL can be used to be the boundary conditions at the bottom of the mucus layer to determine the velocity of mucus. /is model can help engineers to build devices to treat patients who have problem with the respiratory system. Applications include modeling fluid flow through filters such as engine filter and rice fields.


Introduction
Breathing is one of the most important bodily functions that keep us alive. When we breathe air in through the nose it is not just oxygen but also small substances containing particles and pathogens that pass through lungs. However, the respiratory system has filters to trap and remove the strange pathogens, thus allowing us to breathe without irritation. e filters lining the bronchus within the respiratory system are called cilia, hair-like structure, scattered throughout the ciliated cells, which are vital components of the immune system. Beside cilia, there are goblet cells scattering among the ciliated cells in the respiratory system. One of the primary functions of the goblet cells is to secrete mucus to trap the particles or microorganisms that pass through our respiratory system. en mucus forms a layer at the tips of cilia. e cilia aid in transporting these foreign particles out of the body by sweeping the mucus upwards towards the throat. e tips of the cilia make contact with the mucus layer on the forward stroke and bend sideways and backwards on the reverse stroke. Essentially, this causes the mucus to be propelled only in the forward direction on the forward stroke [1]. e excreted mucus is ejected through the vocal cords and into the pharynx. e primary innate defense mechanism is called mucociliary clearance [2]. Under optimal lung conditions, cilia beat at 12 to 15 Hz in coordinated waves propelling mucus at 4 to 20 mm/min [1]. Figure 1 shows a small portion of the trachea cross section. Inside the segment there exist three main layers, air, mucus, and the periciliary layers (PCL). e PCL is a moist layer composed of cilia covering the ciliated cells. e fluid in this layer is called PCL fluid. Below the PCL is epithelial cells. Goblet cells secret mucus to trap the strange particles entering the body to protect the lungs from contaminants found in the epithelium. Among the goblet cells, there exist ciliated cells, where cilia are found on the top and wave in a rhythmic or pulsating motion and use that motion to keep passage ways free of mucus or foreign particles.
In this work, we determine the velocity of the PCL fluid due to the movement of cilia. We consider that the movement of cilia results in two different patterns in the PCL as shown in Figure 2. e left of Figure 2 shows the cilia perpendicular to the horizontal plane. e PCL has only one domain, which contains only cilia. When the cilia move forward and form an angle θ, θ < 90 ∘ with the horizontal plane, the PCL can be divided into two domains, the layer containing cilia and the layer above the cilia as shown on the right of Figure 2. Since the PCL layer consists of the PCL fluid, assisting to treat the cilia to operate normally, and the solid phases, cilia, this part in PCL is considered as a porous medium. e layer above the cilia in the PCL is called freefluid domain, as shown on the right of Figure 2. e PCL domain is denoted by ω � ω 1 ∪ ω 2 where ω 1 is the free-fluid region above the porous layer ω 2 when the cilia form an angle θ, θ < 90 ∘ with the horizontal plane. With the perpendicular case, we have only ω 2 . e variable y Stoke is the height of the domain ω 2 . Notice that y Stoke changes according to the angle θ. In this study, the first aim is to figure out the velocity of the PCL fluid in both domains to determine the boundary conditions between the domains ω 1 and ω 2 .
e equation applied in the free-fluid region was described by the Stokes equation or Navier-Stokes equations, while flow in the porous-medium region was usually modeled by using Darcy's Law or Brinkman equations [6,[8][9][10][11][12][13]. In general, the Brinkman model is used for the porous media with high porosity [14]. Examples of flow in free-fluid region are diverse. For instance, Abbas et al. [15] studied the effect of the flow for an incompressible liquid with nonlinear differential equations. For flow in a porous medium, Ado-Elkhair et al. [16] presented the influence of the magnetic field and Darcy medium on oscillatory wavy walled (Peristalsis) flow due to the cilia motion by using a Adomian decomposition method (ADM) with a system of nonlinear partial differential equations. Zaher et al. [17] studied the effect of the ciliated walls in the uterine tube to strengthen the sperm to reach the egg. ey used Darcy model with an oscillating wall and non-Newtonian magnetized fluid. An Adomian analysis method was employed to solve the equation. e illustrations of simultaneous flow in both free-fluid and porous regions are various; for example, Valdes-Parada and Lasseux [10] provided a new approach for modeling flow in porous and free-fluid domains based on Beaver and Joseph (BJ) configuration by using Darcy's Law and Navier-Stokes equation. Carraro et al. [12] figured out a high-precision simulation of the slow viscous flow over a porous bed by applying a finite element method to Stokes equation in free-fluid domain and Darcy's Law in porous bed. Naqvi and Bottaro [13] studied boundary conditions between an isotropic porous-medium and a free-fluid region by using the Stokes-like equation and homogenization theory. Du and Zu [18] proposed a local and parallel finite element method to the mixed Navier-Stokes/Darcy model with the Beaver-Joseph interface condition for an incompressible fluid flow. Wuttanachamsri [6] used one-dimensional Stokes-Brinkman equations to find the shape of free interface between a porous-medium and a free-fluid region. e two-layer configuration that uses different equations in different domains requires proper boundary conditions at the interface between the free-fluid and porous-medium regions. e difficulty in determining boundary conditions results from the fact that the orders of the corresponding differential operators are often different [14]. e boundary conditions at the interface between the porous-medium and the free-fluid regions had been studied by several authors [10,[19][20][21][22][23]. For example, Beavers and Joseph [22] proposed a slip boundary condition based on repeated experiments with a Poiseuille flow over a permeable porous block. Sahraoui and Kaviany [23] determined an appropriate hydrodynamic boundary condition at the interface between a porous-medium and a plain domain using the pressure correction method. Valdes-Parada and Lasseux [10] studied three different boundary conditions: creeping flow under no-slip conditions, inertial flow no-slip conditions, and creeping flow with slip conditions with Navier-Stokes-Darcy equations. e solutions were calculated by using Direct Numerical Simulation (DNS).
Although the boundary conditions had been investigated extensively by many researchers, the studies of boundary conditions at the interface between the free-fluid region and the porous medium in the human lungs are limited. is study aims to determine the boundary conditions at the interface between the free-fluid region and the porous medium in PCL by using the method of matched asymptotic expansion, which has been used widely by considerable researchers. For example, Chandesris and Jamet [9] studied the Poiseuille flow over a permeable block by using the method of matched asymptotic expansion and determined the boundary conditions at the top of the permeable block. Poopra and Wuttanachamsri [11] calculated the velocity of the PCL fluid with the BJ boundary condition when cilia formed angles with the horizontal plane by using the asymptotic expansion method. e results in the reviewed literature, however, are not appropriated for our research problem. In most of the studies, solid phases in the porous domain are stationary. In our study, on the contrary, cilia in the porous domain are mobile. Although Poopra and Wuttanachamsri [11] studied the PCL with self-propelled solid phases, they did not seek for matching conditions to match the solutions at the transition zone at the free-fluid/porous-medium interface and they did not provide a proper boundary condition at the interface, which was practical to the problem. us, this research has been conducted to figure out the boundary conditions suiting for the problem.
We first provide our mathematical model and boundary conditions at the top and bottom of our domain in dimensionless forms in Section 2. Second, the velocity of the PCL fluid is calculated by using the asymptotic expansion method in Sections 3.1 and 3.2. ird, because unknown constants occurred during finding the asymptotic solutions, matching conditions presented in Section 3.3 are needed to match the outer and inner results presented in Sections 3.1 and 3.2. Fourth, the zeroth-order solution after applying the matching conditions to the inner and outer results is shown in Section 3.4. Next, to verify our solutions, comparisons of our results with experimental data are described in Section 4. en, the boundary conditions at the free-fluid/porous-medium interface for this problem are illustrated in Section 5. Finally, the conclusion is drawn in the last section.

Mathematical Model and Boundary Conditions
We provide the governing equations, the Stokes-Brinkman equations, and their normalized models, in a one-dimensional domain and boundary conditions in this section. Because we consider the movement of the carpet of cilia (not just the movement of a single cilium), the equations used are in a macroscopic scale obtained from the Hybrid Mixture eory (HMT), which is an upscaling method of changing a microscale equation to be a macroscale equation [24].

Stokes-Brinkman Equations.
In this section, we provide the Stokes-Brinkman equations derived in [7] in a macroscopic scale in the one-dimensional domain ω. at is, where μ is the dynamic viscosity, k −1 11 is the inverse of the permeability tensor, ε l is the porosity, v l and v S are the velocities of the liquid and solid phases, respectively, p is pressure, g is gravity, and ρ is the fluid density. We assume that the fluid flows along the x-axis. With the macroscopic scale, we project the average variables on the y-axis. at is, other variables depend on only the y-direction and in the y-direction the pressure gradient is zero and dP/dx is a constant. e Stokes equation is equation (1) without the terms consisting of the permeability k 11 and ε l � 1. at is, and the Brinkman equation is equation (1). erefore (1) is called Stokes-Brinkman equations. e Brinkman equation used in ω 2 has been normalized in [11], which is and the normalized Stokes equation is where h is the characteristic length, U 0 is a volumetric average velocity in the porous medium, K 0 is characteristic permeability, g 0 is characteristic gravity, and B(y + ) � h 2 ε l v + /K 0 K + where the normalized variables For the Stokes equation, the variable U + � v l /U 0 because ε l � 1. In this study, we assume that the velocity of cilia is linear, that is, where a 0 and a 1 are constants. Next, we provide the boundary conditions at the top and the bottom of the domains used in this research.

Boundary Conditions.
Because cilia are immobilized at the bottom of the domain, we assume that the PCL fluid has a no-slip condition at y � 0. at is, For the boundary condition at the top of the domain, we assume that the velocity is bounded by a constant. at is, where D 0 is a constant as illustrated in Figure 3.
For the special cases of the boundary condition, Vanaki et al. [25] used and U + is periodic on the left and right sides of the boundaries for the study of the movement of the PCL fluid in case of diseased and healthy conditions. Liron and Meyer [26] studied the PCL flow at y + � 1 by using the boundary condition zU + /zy � 0. In this study, we use bounded boundary at y + � 1 so that one may apply this study to other related problems. We now have the normalized Stokes-Brinkman equations and boundary conditions. Next, we find the solution U + by using the method of matched asymptotic expansion in Section 3.

Asymptotic Solutions
In this section we employ the method of asymptotic expansion to find the solutions in both domains ω 1 and ω 2 . Because there exists a transition layer at the free-fluid/porous-medium interface, the approximation of this part of the domain is called inner solution. e approximation of other parts of the domain ω is named outer solution. en the solutions are combined through matching conditions provided in Section 3.3. e inner and outer solutions, derived in [11], are briefly introduced in Sections 3.1 and 3.2. e zeroth-order results after applying the matching conditions are provided in Section 3.4.

Outer Solutions.
In this section, we briefly provide the outer solutions found in [11] in two different domains, ω 1 and ω 2 , by using the method of asymptotic expansion. e solutions in domain ω 2 and ω 1 are obtained by substituting into (3) and (4), respectively. us, we have the zeroth-order and first-order solutions for the outer region ω 2 with the boundary condition (7) are respectively, where w 1 and w 5 are constants, and the zerothorder and first-order solutions for the outer region ω 1 with the boundary condition (8) are respectively, where w 3 and w 7 are constants. Next, we provide the solution in the inner region.

Inner Solutions.
In this section, we provide the inner solutions in the transition zone between the porous medium and the free-fluid domain. In the inner region, for any variable f, we use f instead of a physical variable f. For example, y + in the inner region is y + and U + in the inner region is U + . Define Applying (15) to (3) with the properties, obtained from (15), and using the same process of finding the outer solution, we have the zeroth-, first-, and second-order solutions as follows: where m 1 , m 2 , m 3 , and m 4 are constants. e outer and inner solutions will be matched by using the matching conditions provided in the next section.

Matching Conditions.
In this section, we derive matching conditions to match the variables defined in the outer and inner regions. e matching condition used in this work is developed from Van Dyke's rule [9,27]. With the outer and inner variables, and applying Van Dyke's rule matching principle, we have We find the matching condition of f (0) by substituting (20) and (21) into (22). erefore, We suppose that lim y + ⟶ y ± Stoke f (i) (y + ) and lim y + ⟶ ± ∞ f (i) (y + ) exist, i � 0, 1, 2, ... en, considering the zeroth-order term of ε, we have Next, we find another matching condition by applying Taylor's series expansion to (20) around the point y + Stoke or y − Stoke , where the superscript plus and minus signs mean the right and left sides of y Stoke and we assume that f (i) (y + ), i � 0, 1, 2, 3, . . . are continuous at y � y Stoke so that en, we have Mathematical Problems in Engineering 5 Substituting (15) into (26), we obtain Differentiating (27) with respect to y + , taking y + ⟶ y ± Stoke on both sides, and omitting the big O, then, we obtain where in this work we assume that the function f is continuous enough so that the differentiation can be switched with the limit. To find the matching condition for the first derivative of the functions f and f, we differentiate (21) with respect to y + and taking limit y + ⟶ ± ∞ on both sides, we obtain lim y + ⟶ ± ∞ df y + , ε dy + � lim Since we assume that the differentiation can be switched with the limit, from (22), we have is fact demonstrates that (28) and (29) are equal. en, we get ε lim Collecting the same order of ε in (31), we obtain Next, we determine the zeroth-order matching solutions for the couple Stokes-Brinkman equations.
To illustrate our results, the solution (45) is plotted in Figure 4, where the angles between the cilia and the horizontal plane considered in this work are 40 ∘ , 50 ∘ , 60 ∘ , 70 ∘ , 80 ∘ , and 90 ∘ , respectively. In this research, we assume that the cilia have the highest velocity at θ � 90 ∘ and stop beating at θ � 40 ∘ . e values of the variables used in (45) are presented in Tables 1 and 2. e values in Table 1 are used for all angles θ and the values of the variables in Table 2 are different for each angle θ, which are obtained from [8,28], where they have obtained the data from a biological experiment [29]. We concentrate on only one viscosity because we intend to find the solutions that are close to the real problem and the solutions do not change with the different viscosities. e results for the different angles are interesting because we determine the PCL velocities due to the movement of cilia, which means the angle changes.
e Earth gravity, solid velocity, and permeability variables in (5) are fixed according to our designed problem. Pressure can be modified but for this problem the most effect on the velocity of the PCL fluid is the solid phases, rather than pressure gradient. e variables x and y in (5) are independent variables, which give the position of the fluid velocity on the domain. Figure 4 shows the velocity of the PCL fluid with different angles θ in domain ω 2 . e highest velocity occurred at θ � 90 ∘ as expected.
is is because, in this study, we assume that the maximum velocity of cilia occurs at θ � 90 ∘ . Notice that the velocities decrease when the angles between the cilia and the horizontal plane decrease.
In the next section, we compare our solution with an experimental data where the velocity of the solid phases in our solution is set to be zero.

Comparison with an Experimental Study
In this section, we compare the mass flow rate Φ of our result with the study of Beavers and Joseph [22], where where M f is the mass flow rate in the free-fluid region and M p is the mass flow rate in the porous medium. In their study, BJ carried out the experiment to measure the flow rate of a long porous block and a small uniform gap immediately Mathematical Problems in Engineering above the block, while in this study, we analytically find the mass flow rate of the fluid flow in the PCL in human lung by using the asymptotic expansion method. e porous medium in the experimental study is stationary porous blocks, while the porous medium in this study is the ciliary layer, which is a motion medium. erefore, we compare our solution with the experimental data at only one porosity and the solid velocity is set to be zero. Since Beavers and Joseph provide the data, the mass flow rate, only on the free-fluid domain, as the result, only the solution on the free-fluid domain is compared.
By applying the zeroth-order solutions (45), (46) provided in Section 3.4 to (47), the mass flow rate obtained from both equations is    [27]. e gravity g + g 0 is 9.8[m/s 2 ]. e length of the free-fluid domain ranges from 3.37046 × 10 − 4 to 1, where y Stoke is assumed to be 3.37046 × 10 − 4 . e porosity is 0.78 and the characteristic permeability K 0 is 7.1 × 10 −9 . e values of other variables used in (48) are K + � k 11 /K 0 , k 11 � 0.4, w 1 � −2.5, U 0 � 1, D 0 � 1 and h � 7.5, and Φ is normalized by 10 4 . Figure 5 shows the comparison of the mass flow rate in this study (solid line) through the free-fluid domains as a function of σ, with the mass flow rate in the experimental study of Beavers and Joseph [22]. As can be seen from Figure 5, the solution of our study (solid line) is in a good agreement with the experimental data of Beavers and Joseph [22].

Boundary Conditions at the Free-Fluid/ Porous-Medium Interface
In this section we find the boundary condition at the freefluid/porous-medium interface by using the two solutions in the two domains ω 1 and ω 2 . Based on solutions (45) and (46), we can now obtain the boundary condition at the freefluid/porous-medium interface, which is Next, we provide another boundary condition at the free-fluid/porous-medium interface.
Differentiating the zeroth-order solutions (45) and (46) with respect to y + both sides, we have respectively. erefore we obtain another boundary condition at the free-fluid/porous-medium interface, which is us, we obtain two boundary conditions at the freefluid/porous-medium interface, which are proper for the fluid flow in the PCL due to the movement of cilia.

Summary and Conclusions
e movement of the cilia in the PCL layer plays an important role in the clearance of mucus and in preventing foreign materials from entering the lower respiratory tract. Determining the velocity of the PCL fluid flow in motile cilia is necessary for a better understanding of the clearance mechanism. e velocity of PCL fluid can be determined using an ordinary differential equation (ODE), which requires boundary conditions. is research quantifies the boundary conditions at the interface between the porousmedium and free-fluid regions using Stokes-Brinkman equations. To achieve our aim: (i) We apply the Stokes equation to the free-fluid domain and the Brinkman equation to the porousmedium region. (ii) We also assumed two boundary conditions: e first boundary condition is no-slip boundary condition, applied at the bottom of the porous-medium domain; the second boundary condition is that the velocity of the PCL fluid is bounded when y tends to a constant, applied at the top of the free-fluid region. (iii) e velocity of the PCL fluid due to the movement of cilia is obtained by using the matched asymptotic expansion. (iv) e results are illustrated in Figure 4 with different angles, θ � 40 ∘ , 50 ∘ , 60 ∘ , 70 ∘ , 80 ∘ , 90 ∘ , that the cilia make with the horizontal plane. e variables employed to plot the graph of the solutions are shown in Tables 1 and 2. (v) We find that the velocity of the PCL fluid reaches its maximum value when θ � 90 ∘ and decreases when θ decreases, which satisfies our assumption that the highest velocity of cilia occurs when the cilia are perpendicular to the horizontal plane and the cilia velocities decrease with decreasing angle θ. (vi) We also found that the solution in this study is consistent with the experimental data studied by Beavers and Joseph, which is shown in Figure 5.
As mentioned above, not only are there two domains, porous domain and free-fluid domain, but also there is the interface between the free-fluid domain and porous medium.
e boundary conditions at the interface between the porous-medium and the free-fluid regions are determined from the zeroth-order solutions provided in Section 3.4.
(i) e first boundary condition is the relationship between the two zeroth-ordered solutions in the two domains while the second boundary condition is the derivative of those two solutions. (ii) It should be noted that, in this study, the solid phases in the porous domain are self-propelled motion and we provide proper boundary conditions for the problem. On the contrary, in previous studies, they are stationary. Although they are moved, available literature has not provided appropriate boundary conditions for the fluid flow in PCL. (iii) e limitation of this study is that the asymptotic expansion method applied in this study is good for a linear equation, but it is not easy to find the solutions for nonlinear equations. (iv) e disadvantage of this study is that we calculate the zeroth-order approximation solution, where the other orders are cut off. (v) In future work, the more realistic solution such as every angle is considered in one calculation rather than one angle. Applications of the boundary conditions to some other problems are, for example, flows in rice field [30][31][32][33][34].

Conflicts of Interest
All authors declare that they have no conflicts of interest.