Flow Rate Prediction for a Semi-permeable Membrane at Low Reynolds Number in a Circular Pipe

A concise and accurate prediction method is required for membrane permeability in chemical engineering and biological fields. As a preliminary study on this topic, we propose the concentration polarization model (CPM) of the permeate flux and flow rate under dominant effects of viscosity and solute diffusion. In this model, concentration polarization is incorporated for the solution flow through a semi-permeable membrane (i.e., permeable for solvent but not for solute) in a circular pipe. The effect of the concentration polarization on the flow field in a circular pipe under a viscous-dominant condition (i.e., at a low Reynolds number) is discussed by comparing the CPM with the numerical simulation results and infinitesimal Péclet number model (IPM) for the membrane permeability, strength of the osmotic pressure, and Péclet number. The CPM and IPM are confirmed to be a reasonable extension of the model for a pure fluid, which was proposed previously. The application range of the IPM is narrow because the advection of the solute concentration is not considered, whereas the CPM demonstrates superior applicability in a wide range of parameters, including the permeability coefficient, strength of the osmotic pressure, and Péclet number. This suggests the necessity for considering concentration polarization. Although the mathematical expression of the CPM is more complex than that of the IPM, the CPM exhibits a potential to accurately predict the permeability parameters for a condition in which a large permeate flux and osmotic pressure occur.


Introduction
Biological membranes consist of phospholipid bilayers with integral and peripheral membrane proteins, and they have a variety of functions based on the type of a cell. For example, a glycocalyx in endothelial cell membranes functions as the mechanotransductor and maintains homeostasis (Weinbaum et al. 2011;Jiang et al. 2020). The functions and structures of the biological membranes are extremely complex, and the transport mechanism of the biological membranes has not been understood completely. Conversely, artificial membranes are relatively simple when compared to biological membranes, and they are mainly used for chemical engineering applications such as filtration. The filtration technology is categorized into microfiltration, ultrafiltration, nanofiltration, and reverse osmosis based on the size of the substance to be separated. In microfiltration, the pores of the membrane can be confirmed by electron microscopy (Holdich et al. 2006), while for reverse osmosis membranes, which filter sea salts of small molecular weight, the pores are too small for observation. To understand the transport mechanism of substances through membranes, a variety of mathematical models have been proposed, to date.
The mathematical models of membrane transport are classified into two models (Mazid 1984;Wang et al. 2014): the mechanistic model, which is based on the viewpoint of membrane structures (Wijmans and Baker 1995), and the phenomenological model, which is based on non-equilibrium dynamics (Katchalsky and Curran 1965;Schultz 1980;De Groot and Mazur 1984;Demirel 2007;Friedman 2010). The solution-diffusion and pore-flow models are well-known in the mechanistic model. In the solution-diffusion model, the pressure inside the membrane is assumed as constant, and the substance moves from the surface on the feed side to that on the permeate side by diffusion inside the membrane. The solution-diffusion model can describe transport phenomena for nanofiltration and reverse osmosis. On the other hand, in the pore-flow model, the substance moves through the pores of a straight circular tube by pressure-driven advection, and the concentration does not affect the permeate flux. The pore-flow model is applied to membranes with relatively large pores that are observed in microfiltration and ultrafiltration. In the phenomenological model, the Kedem-Katchalsky model (Katchalsky and Curran 1965) and Speigler-Kedem model (Spiegler and Kedem 1966) are widely used. The phenomenological models do not specify membrane structure, and the microscopic membrane structure is not considered. In these models, the flux through the membrane is generally driven by the gradient in the chemical potential, and the flux of the solution exhibits a linear relationship with hydrostatic and osmotic pressures. A major difference in the mechanistic and phenomenological models is that the phenomenological model can deal with various types of permeabilities, such as selective permeability and impermeability, by adjusting the value of a parameter termed as the reflection coefficient. However, mechanistic models cannot deal with various types of permeabilities because of the absence of the reflection coefficient. The physical interpretation of the permeability parameters in the Kedem-Katchalsky model is detailed in an extant study (Kedem and Katchalsky 1961).
To calculate the filtration efficiency using the aforementioned models, the permeability parameters in the models should be determined. To measure the permeability parameters, the concentration on the feed side of the membrane is required. However, it is difficult to directly measure the accurate concentration value on the membrane due to the presence of an unstirred layer (i.e., the heterogeneous concentration field) in the vicinity of the membrane. The elevation in the concentration on the feed side of the membrane surface due to the flow toward the membrane is termed as concentration polarization, which significantly decreases the filtration efficiency. Instead of the direct measurement of the concentration on the membrane, the boundary film theory (Sablani et al. 2001;Porter 1972), used for concentration polarization, with non-dimensional correlation formulae (Berg et al. 1989;Lévêque 1928) or velocity variation method (Nakao and Kimura 1981) is used to predict the concentration value on the membrane surface. In the theoretical studies on the membrane transport, the Kedem-Katchalsky model has been modified to include the effect of the concentration boundary layer (Kargol 1996(Kargol , 2000Bryll and Ślȩzak 2017). However, the prediction method has not been established because it is difficult to obtain an accurate value of the mass transfer coefficient for the boundary layer thickness in the boundary film theory.
We previously developed a relationship between the flow rate and permeability coefficient of the membrane, which is the proportional coefficient between the permeate volume flux and pressure difference on the membrane, for a pure fluid through a circular pipe under low Reynolds number and proposed a formula for predicting the permeability coefficient from the flow rate (Takeuchi et al. 2019). This formula is based on the phenomenological model and does not assume the types of materials. Therefore, this formula is universal as long as the permeate flux is proportional to the pressure difference on the membrane. In this study, we propose a formula that describes the effects of the concentration polarization on the solution permeation through a semi-permeable membrane, as a preliminary study, to establish a concise and accurate prediction method to determine membrane permeability parameters including the reflection coefficient. With respect to concentration polarization, we develop a model for the permeate flux and flow rate from the cross-sectional average steady advection-diffusion equation with an assumption of a unidirectional flow along the boundary pressure gradient. The effect of the concentration polarization on the flow field in the circular pipe is discussed by comparing the numerical simulation results and the model, in which the infinitesimal concentration boundary layer (i.e., the infinitesimal Péclet number condition) is assumed.

Flow Rate Through a Membrane in a Circular Pipe
The flow rate of a pressure-driven flow through a circular pipe impeded by a membrane is derived. The fluid domain is setup, as shown in Fig. 1, and the flow is driven by a constant pressure difference, p L − p R (≡ Δp) , across distance in the longitudinal direction x. The membrane is placed at x = � in parallel with the r axis (i.e., perpendicular to the centerline). The deflection of the membrane is not considered.
In the case of no membrane, the fully developed velocity profile u P corresponds to the axisymmetric Hagen-Poiseuille flow: where r is the radial distance from the center axis, R the pipe radius, the viscous coefficient, and p the prescribed pressure gradient −Δp∕ . Hereafter, the flow rate of the Hagen-Poiseuille flow is denoted as Q P (= − p R 4 ∕8 ).
By placing a permeable membrane, the above parabolic velocity profile is no longer valid. However, with a membrane of large permeability L p , we can reasonably expect that the parabolic profile is asymptotically recovered (Takeuchi et al. 2019). Therefore, in the following model, we assume a large permeability case, and an approximate velocity profile can be obtained.
To reflect the effect of the permeable membrane on the local velocity profile, p in Eq. (1) is temporarily replaced with dp∕dx , and the equation is integrated with respect to x in the ranges 0 ≤ x ≤ ′ and ′ ≤ x ≤ separately: where u L and u R denote the fluid velocities in the left-and right-hand sides of the membrane, respectively, and p ∓ denotes the limiting pressure values on the interface as shown in Fig. 1. These velocities are equated with the on-membrane velocity ũ(r) , and ′ is eliminated from the above equations by adding Eqs. (2a), (2b): where p denotes the pressure difference at the membrane defined as p ≡ p − − p + . Given that the permeate flux through the membrane corresponds to the local flow rate per unit area (Katchalsky and Curran 1965), ũ is also termed as the permeate flux in this paper.
Based on the Kedem-Katchalsky model (Katchalsky and Curran 1965), the volumetric flux is proportional to the difference between the hydrostatic and osmotic pressure differences across the semi-permeable membrane with a permeability coefficient L p . Furthermore, by equating the volumetric flux through the membrane with ũ and using the van't Hoff law, the following equation is derived: where [[c]] denotes the concentration difference on the membrane defined as [[c]] ≡ c − − c + , R g denotes the gas constant, and T the temperature. Substituting the above p into Eq. (3) and solving the equation with respect to ũ , the following equation is obtained: In the following Sects. 2.1 and 2.2, we predict the concentration difference with and without consideration of the concentration boundary layer and propose two novel models that relate the permeate flux and flow rate.

Infinitesimal Péclet Number Model (IPM)
In this subsection, we propose the model for the permeate flux and flow rate without considering the concentration boundary layer. The Péclet number Pe is defined as V max R∕D , where D denotes the diffusion coefficient, V max the maximum fluid velocity in the no-membrane case. In the condition with infinitesimal Péclet number ( Pe → 0 ), convective transport is negligible, and the following approximations can be made under steady state: where c L and c R denote the boundary concentrations (as shown in Fig. 1). Substituting [[c]] into Eq. (5), the permeate flux ũ and flow rate Q are obtained as follows: where M denotes the non-dimensional parameter defined as follows: where L(= L p ∕R) denotes the non-dimensional permeability coefficient. Then, the normalized form is as follows: where denotes the ratio of the osmotic pressure drop to the prescribed hydrostatic pressure difference between the inlet and outlet of the circular pipe as follows: Under → 0 (i.e., no concentration), Eq. (7) reduces to the formula proposed by Takeuchi et al. (2019). This indicates that the above flow rate is a reasonable extension of the pure Δp .

Concentration Polarization Model (CPM)
In this subsection, we propose the model for the permeate flux and flow rate by considering the concentration boundary layer. For the analysis of the concentration polarization, it is important to model the increase and decrease in the solute concentration due to advection, diffusion, and permeation. In this model, the concentration difference on the membrane is derived from the steady advection-diffusion equation in the cylindrical coordinate system as follows: with the following boundary conditions: where u denotes the fluid velocity. We assume that the variation in the concentration in the circumferential direction and the advection of the concentration in the radius direction can be ignored. These assumptions along with the boundary condition (12d) and replacement of the convective velocity u with the cross-sectional average velocity ū in the longitudinal direction, yield the following one-dimensional differential equation: where ū and c are defined as ū = 2 R 2 ∫ R 0 rudr and c = 2 R 2 ∫ R 0 rcdr , respectively. Thus, the concentration difference on the membrane is derived from Eqs. (12a), (12b), (12c) and (13) as follows: in Eq. (5) as an approximate value, permeate flux ũ and flow rate Q and its normalized form Q∕Q P are as follows: Δp , Using the relation of Q = R 2ū and Maclaurin series for the exponential functions in Eq. (14), ū is expressed as follows: The permeate flux and flow rate in Eqs. (15) and (17) are calculated using Eq. (14) with Eq. (18). As the mathematical expression Q∕Q P = F(M) is complex, the inverse formula M = F −1 (Q∕Q p ) cannot be obtained in a straightforward manner.

Computational Conditions
The applicability of the models explained above is assessed by comparing the fully validated numerical simulation, which directly solves the basic equations of fluid motion and mass transport. The numerical method is based on the discrete-forcing immersed boundary method (Sato et al. 2016;Takeuchi et al. 2018), which undertakes (i) the consistency between the incompressible velocity and pressure fields and (ii) direct discretization even at the grid points adjacent to the interface, and thereby enables the sharp representation of the pressure discontinuity across the membrane (Takeuchi et al. 2019). The details of the numerical method and its validity for mass transport problems, including the semipermeable membrane, are described in the Supplementary material. A msore generalized numerical method considering the coupling with solute permeation is provided in Yamada et al. (2021). All the governing equations in the simulation are non-dimensionalized by the maximum velocity V max in the no membrane case, pipe radius R, reference pressure V 2 max (with fluid density ), reference time R∕V max , and boundary concentration difference Δc . The simulation is conducted in an axisymmetric two-dimensional domain. The domain size is and R in the longitudinal x and radial r directions, respectively, and = 5 throughout this study. The semi-permeable membrane is placed at x = 2.5 . For the boundary condition, p = 2, c = 2 at the inlet ( x = 0 ), p = 1, c = 1 at the outlet ( x = 5 ), u = 0, c∕ r = 0 at the impermeable wall ( r = 1 ), and a gradient-free boundary condition at the centerline ( r = 0 ) are imposed. The Reynolds number Re is defined as V max R∕ , and it is set to 1. To investigate the applicability of the models in the wide ranges of L, , and Pe , L = 10 −5 , 10 −4 , 10 −3 , 10 −2 , 10 −1 , 10 0 and = 0, 0.3, 0.5, 0.9 are employed. The Péclet number is set as Pe = 1 unless specified otherwise. The width of the grid points and time increment for the numerical simulation are 0.025(= R∕40) and 10 −4 , respectively.

Results and Discussion
The radial profiles of the permeate flux for the combination of L = 10 −2 , 10 0 , and = 0.3, 0.9 are shown in Fig. 2. In the cases with the small permeability coefficient L, the magnitude of the permeate flux is small with an almost flat profile near the centerline ( r = 0 ) when compared to those with large L in the numerical simulation results and analytical models. A similar trend is observed in the previous study without the osmotic pressure (Takeuchi et al. 2019). The magnitude of the permeate flux decreases as increases due to the elevation of the osmotic pressure. In all the results, the CPM is closer to the numerical simulation results than the IPM near the center, particularly for cases with large values of L. For cases with small values of L, the CPM and IPM overestimate and underestimate the permeate flux near the center and in the vicinity of the wall ( r = 1 ), respectively.
The pressure and concentration distributions along the centerline obtained from the numerical simulation are shown in Fig. 3. When comparing Fig. 3a, b, the magnitude of the pressure gradient in the pipe for the case with large L (Fig. 3b) is higher than that in the case with small L (Fig. 3a) because of the higher flow velocity. A large magnitude of the concentration gradient is also observed in Fig. 3b because of the advection effect. In the case of the lowest flow velocity, as shown in Fig. 3c, the pressure and concentration profiles are flat. In Fig. 3d, a slight pressure gradient can be observed, whereas the concentration distribution appears to be flat. This nearly flat profile (i.e., [[c]] ≃ Δc ) can appear as inconsistent with the significant deviation of the flux in the CPM when compared to that in the IPM in Fig. 2d. This is explained in the following. From Eqs. (6) and (15) Fig. 2d. The normalized flow rates Q∕Q P of the IPM and CPM are compared with that from the numerical simulation in Fig. 4. The difference between the IPM and numerical simulation increases as the non-dimensional parameters M and increase. Conversely, the curves of the CPM are in good agreement with the numerical simulation results in all the cases of over the entire M range, which demonstrates that the CPM can deal with the solution flow through the semi-permeable membrane. To investigate more details, the errors in the flow rates of the IPM and CPM, based on the numerical simulation results, are plotted in Fig. 5. Given that the osmotic pressure on the membrane does not affect the velocity field when = 0 , the symbols denoting IPM and CPM overlap in Fig. 5. The error of the CPM is lower than that of the IPM at large permeability coefficient values of L = 10 −1 , 10 0 . However, at small permeability coefficient values of L = 10 −2 , 10 −3 , the error of the CPM is not necessarily lower than that of the IPM. In cases with small values of permeability coefficient, both models overestimate and underestimate the permeate flux near the center and wall, respectively, as shown in Fig. 2, and these errors cancel each other, resulting in the lower errors in the IPM than those in the CPM, coincidentally. In the IPM, the errors tend to increase as the permeability coefficient increases, whereas the CPM succeeds in suppressing the increases in error even in the conditions with large L and by considering concentration polarization. The errors tend to increase with for L = 10 −1 , 10 0 .
The change in the flow field and the applicable range of the models for the Péclet number are investigated at L = 10 0 and = 0.9 , in which the highest osmotic pressure occurs. The radial profiles of the permeate flux for Pe = 1, 10, 100 are shown in Fig. 6. The result of Pe = 1 is identical to the result "Numerical" (in red) in Fig. 2d. The magnitude of the flux of Pe = 1 is more than twice as large as those of Pe = 10, 100, whereas the magnitudes of the flux with Pe = 10, 100 are comparable. The profile becomes flat near the center as the Péclet number increases and the position at which the permeate flux is maximum changes from r = 0 to r ≈ 0.5 for Pe = 100. Figure 7 displays the distributions of the ratio of the velocity component in the radial direction to that in the longitudinal direction. In all the cases, the flows toward the wall and center are observed in the vicinity on the feed side and permeate side of the membrane, respectively. The magnitude of the ratio increases as the Péclet number increases. The advection in the radial direction near the membrane affects the osmotic pressure and has a significant effect on the profile of the permeate flux as shown in Fig. 6. For low values of Pe, the assumption of negligible advection in the radial direction in the CPM can be justified. The difference between the analytical flow rate model [(Eq. (7) or Eq. (16)] and the numerical simulation result is plotted in Fig. 8a for a wide range of Péclet numbers, and Fig. 8b shows the error defined as the absolute value of the difference divided by the numerical value. In Fig. 8a, the differences in the IPM increase sharply in the range 0 < Pe ≤ 20 , Fig. 6 Radial profiles of the permeate flux for Pe = 1, 10, 100 at L = 10 0 and = 0.9 and then decrease as Pe increases. On the other hand, the differences in the CPM in the range 0 < Pe ≤ 10 are almost the same, and then decrease. In Fig. 8b, the tendency of the error in the IPM is similar to that of the difference. In the CPM, a non-smooth variation of the error is observed at Pe ≈ 10 because the sign of the difference changes from positive to negative. In the range higher than Pe = 10, the errors of the CPM are in the comparable order of magnitude to that in the IPM.
Based on the results explained above, the applicable range of the IPM is confined to the low values of and Pe . However, the CPM is effective for the solution flow through the semi-permeable membrane in a wide range of L, , and Pe when compared to IPM. This indicates that it is necessary to consider concentration polarization even for almost homogeneous concentration fields.

Conclusion
In this study, we proposed a concentration polarization model (CPM) for two-component solutions across a semi-permeable membrane. The model describes the effects of concentration polarization on solvent permeation through a semi-permeable membrane by using the cross-sectional average steady advection-diffusion equation. In the CPM, an unidirectional flow is assumed along the boundary pressure gradient. Although the model developed by Takeuchi et al. (2019) for a pure fluid is extended to the infinitesimal Péclet number model (IPM) in a straightforward manner, the application range is narrow because the concentration boundary layer is not considered. Conversely, the CPM demonstrates superior applicability in a wide range of permeability coefficients, strength of the osmotic pressure, and Péclet numbers. Hence, the results indicate that it is necessary to consider concentration polarization even for almost homogeneous concentration fields.
In our previous study (Takeuchi et al. 2019), we formulated the relationship between a permeability coefficient and flow rate for pure fluids, which allows one to accurately determine the membrane permeability coefficient from the measured flow rate of the simple experiment. Recalling that the aim of the present study is to establish a concise and accurate method for predicting membrane permeability parameters, the CPM can potentially predict the permeability parameters accurately for cases involving large permeation flux and osmotic pressure, although the mathematical expression of the CPM is more complex  2,4,6,8,10,20,40,60,80, 100 at L = 10 0 and = 0.9