Pulsatile pressure enhanced rapid water transport through flexible graphene nano/Angstrom-size channels: a continuum modeling approach using the micro-structure of nanoconfined water

Several researchers observed a significant increase in water flow through graphene-based nanocapillaries. As graphene sheets are flexible (Wang and Shi 2015 Energy Environ. Sci. 8 790–823), we represent nanocapillaries with a deformable channel-wall model by using the small displacement structural-mechanics and perturbation theory presented by Gervais et al (2006 Lab Chip 6 500–7), and Christov et al (2018 J. Fluid Mech. 841 267–86), respectively. We assume the lubrication assumption in the shallow nanochannels, and using the microstructure of confined water along with slip at the capillary boundaries and disjoining pressure (Neek-Amal et al 2018 Appl. Phys. Lett. 113 083101), we derive the model for deformable nanochannels. Our derived model also facilitates the flow dynamics of Newtonian fluids under different conditions as its limiting cases, which have been previously reported in literature (Neek-Amal et al 2018 Appl. Phys. Lett. 113 083101; Gervais et al 2006 Lab Chip 6 500–7; Christov et al 2018 J. Fluid Mech. 841 267–86 ; White 1990 Fluid Mechanics; Keith Batchelor 1967 An Introduction to Fluid Dynamics ; Kirby 2010 Micro-and Nanoscale Fluid Mechanics: Transport in Microfluidic Devices). We compare the experimental observations by Radha et al (2016 Nature 538 222–5) and MD simulation results by Neek-Amal et al (2018 Appl. Phys. Lett. 113 083101) with our deformable-wall model. We find that for channel-height Ho<4 Å, the flow-rate prediction by the deformable-wall model is 5%–7% more compared to Neek-Amal et al (2018 Appl. Phys. Lett. 113 083101) well-fitted rigid-wall model. These predictions are within the errorbar of the experimental data as shown by Radha et al (2016 Nature 538 222–5), which indicates that the derived deformable-wall model could be more accurate to model Radha et al (2016 Nature 538 222–5) experiments as compared to the rigid-wall model. Using the model, we study the effect of the flexibility of graphene sheets on the flow rate. As the flexibility α increases (or corresponding thickness T and elastic modulus E of the wall decreases), the flow rate also increases. We find that the flow rate scales as m˙flexible∼α0 for (αΔpW/EHo)≪1 ; m˙flexible∼α for (αΔpW/EHo)∼O(10−1) ; and m˙flexible∼α3 for (αΔpW/EHo)∼O(1) , respectively. We also find that, for a given thickness T , the percentage change in flow rate in the smaller height of the channel is more than the larger height of the channels. As the channel height decreases for the given reservoir pressure and thickness, the Δm˙/m˙ increases with Ho−1 followed by Ho−3 after a height-threshold. Further, we investigate how the applied pulsating pressure influences the flow rate. We find that due to the oscillatory pressure field, there is no change in the averaged mass flow rate in the rigid-wall channel, whereas the flow rate increases in the flexible channels with the increasing magnitude of the oscillatory pressure field. Also, in flexible channels, depending on the magnitude of the pressure field, either of the steady or oscillatory or both kinds of pressure field, the averaged mass flow rate dependence varies from Δp to Δp4 as the pressure field increases. The flow rate in the rigid-wall channel scales as m˙rigid∼Δp , whereas for the deformable-wall channel it scale as m˙flexible∼Δp for (αΔpW/EHo)≊0 , m˙flexible∼Δp2 for (αΔpW/EHo)∼O(10−1) , and m˙flexible∼Δp4 for (αΔpW/EHo)∼O(1) . We find that both the flexibility of the graphene sheet and the pulsating pressure fields to these flexible channels intensify the rapid flow rate through nano/Angstrom-size graphene capillaries.


Introduction
Water transport through nanocapillaries, such as in biological systems and nanoscale materials, is ubiquitous in the natural world and has many applications in technology [1][2][3][4][5][6][7][8]. In biological systems, capillary action in nano capillaries enables efficient water transport in plants, allowing them to absorb water from the soil and distribute it throughout their tissues [9,10]. It also facilitates water movement in the vascular systems of humans and animals, such as in blood vessels and the microcapillaries of the human body [11]. From a technological standpoint, flexible nanocapillary flow has been harnessed for various applications. Some notable examples include flexible nanochannels, which are utilized in lab-on-a-chip devices for precise manipulation and analysis of small volumes of fluids [12,13]. These flexible nanochannels enable miniaturized systems for tasks such as chemical and biological sensing, DNA sequencing, and drug delivery [14][15][16]. Flexible graphene nanochannels are employed in water purification, desalination, filtration, and separation of biomolecules processes due to their ability to selectively allow specific molecules or particles to pass through while blocking others [17,18]. Due to the high surface area and confined flow of water, the flexible nanochannels are incorporated into heat exchange systems for efficient cooling [19][20][21][22]. Inkjet printers use oscillating pressure pulses to eject the jet through nano-micro channels to control ink flow onto paper or other surfaces. The precise manipulation of tiny droplets enables high-resolution printing [23][24][25]. Even these nano-micro channels are also used to generate shock-wave by using the piezoelectric effect and controlling the pressure field in various ways [26][27][28]. Radha et al [29] has observed fast water flow within graphene-based nanocapillaries.
Neek-Amal et al [30] used a combination of the Hagen-Poiseuille theory, with confined water properties under continuum modeling, and the molecular dynamic simulations to study the capillaries. They represented the capillaries as rigid wall channels and presented a model for the flow rate. They attributed the fast water flow to the high density and viscosity of water molecules confined within these tiny nanocapillaries. Although the Graphene sheets are flexible [31]. The flexibility of a shallow channel plays a crucial role in influencing both the effective pressure drop across the channel and the resulting flow profile [32,33]. This is primarily because the flow rate is susceptible to the cross-section length scale, exhibiting a fourth power dependence [33]. Consequently, even minor deformations in the channel's geometry can significantly alter the pressure drop and flow behavior within the system.
Gervais et al [32] provided a satisfactory model of the deformation-induced change in the flow rate by empirically connecting a Hookean elastic response with the lubrication approximation for Stokes flow. Nonetheless, their model involves a fitting parameter that needs to be determined through experimentation for every channel shape. Christov et al [33] relate the fitting parameter in Gervais et al [32] using a perturbation technique for the flow. Therefore, using the small displacement structural mechanics and perturbation theory presented by Gervais et al [32] and Christov et al [33], respectively, in this paper, we represent these graphene-based nanochannels with a deformable channel-wall model. Further, Akram et al [34] investigate the Stokes second problem for a Burgers' fluid flowing over a plain wall, oscillating between two perpendicular side walls. They derived analytical solutions for the velocity distribution and the corresponding shear stress. However, their formulation did not use the microstructure properties or deformability. Unlike in our case, where pulsatile pressure is considered, they studied the flow induced by the vibrating wall. Akram et al [35][36][37][38][39] studied the convection phenomenon involving double diffusion in the context of peristaltic transport using Sisko, Jeffreys, Prandtl, fourth-grade and Johnson-Segalman nanofluids, respectively, within a presumed sine/cosine-shaped channel exposed to a magnetic field. However, their formulation did without incorporating the pulsatile pressure field and the deformability.
Under the lubrication assumption in the shallow nanochannels (specifically, the ratio of the channel's height to its width and the ratio of the channel's height to its length are both assumed small), and using the microstructure of confined water along with slip at the capillary boundaries and disjoining pressure [30], we study the effect of flexibility of graphene sheet to the flow rate. Additionally, we examine the influence of applied pulsating pressure on the water flow rate within these flexible nano capillaries. The application of pulsating pressure fields or vibrations is a well-known technique for enhancing the flow rate of complex fluids through porous mediums, including channels and tube capillaries [23][24][25][40][41][42][43]. However, for a Newtonian fluid (such as water), the pulsating pressure has no effect on the flow rate in a rigid-wall channel/tube [44]. We conduct a comparison between the predicted flow rates considering both the flexibility of the channel wall and the application of pulsating pressure with the experimental findings reported by Radha et al [29], as well as from the molecular dynamic simulations predictions by Neek-Amal et al [30], where their results were well-fitted using the rigid-wall model.
We consider a shallow rectangular channel with a length L, width W, and height H, satisfying the condition H ≪ W ≪ L as shown in figure 1. The channel's upper wall/lower wall is a flexible graphene sheet sealed at the edges of the vertical channel wall, which can undergo deformation. The time-dependent flow In (a), we show the complete view of the channel without deformability. In (b), we show the front view of the channel with deformability (solid-line) and without (dashed-line) deformability. In (c), we show the side-view of the channel, where inside pressure is p(x) at any longitudinal axis x and outside pressure is assumed to be zero. Q(t) because of the applied pulsating pressure field occurs along the x-direction. As a result of the normal stresses exerted by the flow on the walls, the soft top wall of the channel deforms in the positive z direction, away from the x − y plane. This deformation defines the steady shape of the channel's top wall as is the deformation along the vertical direction as depicted in figure 1. A pulsatile pressure field p(t) is applied at the reservoir at x = 0, and the exit pressure is assumed to be zero for reference purposes. Currently, we do not assume any particular magnitude for the displacement. Still, as the pressure is pulsatile, δ can be greater and less than zero, and the magnitude |δ| ≪ W is in our problem.
The structure of the paper is outlined as follows. Section nomenclature describes the nomenclature with symbols, units, and description. Section 2 describes the governing equations, while section 3 presents the model's derivation. In the section 4, we present and analyze the results. Finally we present conclusions in section 5.

Cauchy equations
The Cauchy's equation [45] and the continuity equation for an incompressible fluid are given by where v = [u v w] is the fluid velocity, p is the fluid pressure, ρ is the fluid density, g is the gravitational body force, and τ is the total deviatoric stress tensor.

Boundary conditions
Boundary conditions play an essential role in determining the solution. We assume that the fluid cannot penetrate the channel wall; therefore on the boundary Γ, it gives v · n wall = 0, where n wall is the unit outward normal vector on the wall. Generally, the no-slip boundary condition at the fluid-solid interface is a fundamental notion in fluid mechanics. However, fluid flow at the nano and Angstrom scale require a certain degree of tangential velocity (Navier slip) in order to match experimental observations [46][47][48][49]. This leads to where m wall is the tangential unit vector along the channel wall. Also, the arbitrary parameter θ meets 0 ⩽ θ ⩽ 1. Here, θ = 0 and 1 corresponds to pure-slip and no-slip boundary conditions. The symmetry boundary condition at the centreline of the channel z = 0 demands the velocity normal to the centreline, and the Cauchy traction vector t tangential to the centreline is zero. These two conditions can be expressed as v · n centreline = 0, and t · m centreline = 0, respectively, where n centreline and m centreline are the unit normal and unit tangent vector to the symmetry boundary, respectively. The traction on the boundary, which is equivalent to a Neumann boundary condition, is expressed as 3.
where p(x) is the pressure at any longitudinal direction x and 0 < α < 2/3. Young's modulus of a graphene membrane is E = 1 TPa [31,50]. H o is the initial height of the channel when δ = 0. However, the α is a fitting parameter that varies with the different geometry of the channel and needs to be calculated explicitly from the experiments. To overcome this issue, Christov et al [33] performed the perturbation analysis using the isotropic quasi-static plate bending and the Stokes equations and figured out by comparing with Gervais et al [32] model that for rectangular cross-section α = (1/60)(W/T ) 3 (1 − ν 2 ), where T is the thickness of the upper horizontal wall and µ is the Poisson's ratio of the material (for incompressible material ν = 0.5 [51]).

2D planar model
We consider fully developed, time-periodic, incompressible laminar flow in a rectangular channel of height H and width W as shown in the schematic diagram (figure 1). The channel is assumed to be sufficiently long and wide in comparison to the height (i.e. H/W ≪ 1, and H/L ≪ 1) to use a two-dimensional planar model [33,52,53]. The external forcing vibrations are applied using a pulsatory pressure gradient in the water reservoirs similar to the pistons placed at infinity from the channel wall [44]. We also exclude any hydrodynamic instability caused by pulsatory pressure in the transience flow field. We further assume a minimal expansion or contraction due to deformability in comparison to the height of the channel, δ/H ≪ 1, which is caused by the pressure difference between the fluid and the atmospheric conditions in the deformable channel. We assume the Cartesian velocity components u and w along longitudinal and vertical directions x and z, respectively. The z coordinate is measured from the channel's mid-plane. Therefore, using the lubrication assumptions in the shallow cross-section of the channel, we retain the leading order terms as described in [33]. Using the impermeable solid-wall boundary condition, we get w(z = −H(x)/2) = w(z = H(x)/2) = 0. In the leading order terms, using the impermeable solid-wall boundary condition, the normal velocity vanishes everywhere, i.e.
Also, we neglect gravitational forces over capillary forces. Under these assumptions for H/W ≪ 1, and H/L ≪ 1, the Cauchy's equation (1) can be written as where η and ρ are the viscosity and density of confined water, respectively. Using the Navier slip, the tangential solid-wall boundary condition gives where

Fourier representation
We use a similar mathematical solution methodology as described by O'Brien [55]. Due to the linear nature of equation (9), it is possible to generate any time-varying flow with a period T within a rectangular duct using the Fourier representation. First, the pulsating pressure wave is decomposed into its fundamental frequency and harmonic components. Similarly, the velocity is decomposed into its steady state and time-periodic harmonic components. The overall solution is obtained by combining all the steady state and harmonic solutions, each scaled by its respective Fourier coefficient.
In the standard lubrication assumption, we assume the pressure depends only on the longitudinal direction. We consider the width-averaged pressure gradient to be where p o ⩾ 0 is the steady pressure gradient. Also, p s,n , and p c,n represent the sine and cosine amplitudes of the pressure gradient of nth harmonics of the oscillatory pressure field. In deformable channels under small displacement approximations, all coefficients p o , p s,n ,and p c,n might weakly depend on 'x' . In compact form, the above function can be written as where p n = p c,n − i p s,n , only the real part of the pressure gradient represents the solution. Similarly, the velocity can be described as (u s,n sin (ωnt) + u c,n cos (ωnt)) , where u o , and u n are the width-averaged velocities. Also, the u n = u c,n − i u s,n and only the real part of the flow profile represent the solution.

Analytical solution
Substituting equations (12) and (13) in the model equation (9) along x direction, we get To satisfy equation (14), the steady and the time-dependent terms for each harmonic requires and respectively. After solving equations (15) and (16), we get respectively. From equation (17), the velocity gradients are given by respectively. Using the symmetry boundary condition at z = 0 for all time t, i.e. ( ∂u ∂z ) z=0,t = 0 in equation (18), we get Similarly, using the tangential solid-wall boundary condition with Navier slip for all time t from equation (10), we obtain which yields and respectively. Here sgn(c 4,n ) is the 'sign function of argument (c 4,n )' . As we know, the pressure is decreasing along the channel length for a steady part. Therefore, we assumed p o ⩾ 0. For a meaningful solution we require sgn(c 4,n ) = 1. Hence substituting equations (17), (19), (21) and (22) in equation (13), the generalized velocity profile under pulsating field can be written as We define k as a non-dimensional frequency or kinetic Reynolds number as Equation (23) can be further written as, In equation (25), the real part of the pressure gradient p n e (iωnt) is the valid physical condition. Therefore, the 'Term A' mentioned using underbrace must also be real for a valid physical flow profile. We are using properties of confined water in nano and angstrom-scale Pa s (see section 3.5 on the microstructure of confined water), the kinetic Reynolds number k is approximately O(10 −12 ω) ≈ 0. Assuming β = √ ikn, therefore, even for very large values of applied ω, and nth harmonic of applied pulsating pressure, the β tends to zero. Hence, in equation (25), the 'Term A' can be written to where β = √ ikn. We find that irrespective of the sign of β sinh(β), the limit of Term B will not be affected and is real. Also, Term C has to be real for a physically valid flow profile. The real part of the limit of Term C depends on the sign of the real part of β sinh(β). Also with in the experimental range, the 0.5H ρωn/η ≈ O(10 −6 √ ωn) ≪ 10, which upon substitution we get the real part of Hence for the real solution, equation (26) can be simplified to Substituting equation (27) back in equation (25), we get The volume flow rate in a deformable nanochannel at any time t is given by As we know from equation (7), that Using equation (12), and substituting equation (30) in equation (29), we get We assume p(x = 0) to be the relative oscillatory inlet pressure with respect to the pressure at the outlet of the channel, i.e. p(x = L) = 0. In nano and Angstrom scale channels, the viscosity η(H) and density ρ(H) is a function of the height of the channel [30]. Even in many of the literature we found that at the nano-Angstrom scale, the slip-length λ(H) is also a function of the height of the channel [56,57]. Under small displacements of the wall, we assume the viscosity, density and slip-length to be η( , respectively. We integrate equation (31) along the channel length L for a given applied oscillatory pressure field at the inlet of the channel, which yieldŝ The flow rate Q(t) is not a function of longitudinal direction. Therefore substituting x = 0 in equation (33) gives the volume flow rate in the channel as where p(x = 0) = ∆p contains both the steady and the oscillatory pressure field at x = 0. Therefore the mass flow rate can be written asṁ In order to calculate the mean volumetric flow in a time-periodic quasi-steady state over a complete oscillatory pressure field cycle in the channel, we can integrate it as Similarly the mean mass flow rate across one channel of dimensions L × H × W is given bȳ

Microstructure of confined water
Radha et al [29], and Neek-Amal et al [30], observed that the density ρ and viscosity η of the confined water inside a nanochannel is significantly large (for H o ≲ 12Å) and varies with the size of the nanocapillary. In hydrophobic nanocapillaries, the presence of excluded volume near the confining walls reduces the available volume for a given mass of water, leading to increased density. Additionally, the confining walls induce a structuring effect that enhances the viscosity of water [29,30,58,59]. Therefore, for a given channel height, Neek-Amal et al [30] assumed a uniform density in the channel under a continuum approach and proposed a model for density and viscosity of the confined water as a function of the height of the channel, which are and   [29] and Neek-Amal et al [30], respectively, we choose ζ(H o ) = 1 for all channel height.
In nano and angstrom-scale capillaries, the hydrostatic pressure is likely to have little effect on the dynamics of flow, which is many orders of magnitude smaller than the capillary pressure ∆p c , and disjoining  [30]. The solid black line is a fit using f(Ho) = 1 + ae −Ho/δ1 . In the inset, we show the nondimensionalized viscosity of the confined water. The green symbols and solid black line represent the results using the molecular dynamic simulation, and a fit using g(Ho) = 1 + be −Ho/δ2 [30], respectively. Reprinted from [30], with the permission of AIP Publishing. Figure 3. Variation of the entropic pressure ∆pe (kbar) and capillary pressure ∆pc (bar) (in the inset) with width-averaged channel height. The data is taken from [29,30]. Reprinted from [30], with the permission of AIP Publishing.
pressure ∆p d . The capillary pressure is due to interfacial tension, which is given by ∆p c = 2σ cos(ϕ)/H, where σ ≊ 70 mN m −1 is the surface tension and ϕ is the contact angle between water and surface (graphene/graphite). In different studies, the contact angle varies between 55 • and127 • [30,60,61]. Further, the disjoining pressure is due to the enhanced ordering of water structure in nanocapillaries. It is given by  [29]). Using these experimental parameters, we show the variation of the entropic pressure ∆p e (kbar) and capillary pressure ∆p c (bar) (in the inset) with width-averaged channel height in figure 3 [29,30].

Results and discussion
The mass flow rateṁ without oscillatory pressure field (i.e. p n = 0) in the rigid-wall N nanochannels of dimensions L × H × W using equations (7), (29) and (35) can be calculated aṡ where Q is the volume flow rate in a single nanochannel. Also, ∆p = p(x = 0) is the inlet pressure. Therefore, using the microstructure properties of confined water, we can writė where ∆p = ∆p vdW + ∆p e + ∆p c = as presented by [30].
Using equations (34) and (35), we can write the mass flow rate for the flexible wall N graphene channels using the microstructure properties of confined water and under an oscillatory pressure fielḋ where ∆p = p(x = 0) is the applied pressure at the inlet of the channel, which contains both the steady and oscillatory pressure field, i.e. p o ̸ = 0, and p n ̸ = 0. Using equation (43), we can easily identify the following limits: I. The mass flow rate in the rigid wall N nanochannels under pulsating pressure field, i.e. α = 0: using equation (43), the mass flow rate for the 'N' rigid channel under an oscillatory pressure field can be written asṁ where ∆p = p(x = 0) is the applied pressure at the inlet of the channel, which contain both the steady and oscillatory pressure field, i.e. p o ̸ = 0, and p n ̸ = 0. To the best of our knowledge, we have not seen the above derived equation (44) in the literature so far. II. The mass flow rate in the flexible wall N nanochannels without pulsating pressure field, without slip, and without confined water properties, i.e. λ = 0, p n = 0, f(H o ) = g(H o ) = 1: using equation (43), we can write the mass flow rate for the 'N' flexible channel without pulsating pressure field, without slip, and without confined water propertieṡ Equation (45) is the same as the analytical model expression under lubrication assumptions derived by Gervais et al [32]. III. The mass flow rate in the rigid wall N nanochannels without pulsating pressure field, i.e. α = 0, p n = 0: under these limits, we obtain equation (42) as presented by Neek-Amal et al [30].
IV. The mass flow rate in the rigid wall N nanochannels without pulsating pressure field, without slip and without confined water properties, i. e. α = 0, λ = 0, p n = 0, f(H o ) = g(H o ) = 1: under these limits, we obtainṁ which is a classical result of Hagen-Poiseuille flow in channels [45,62,63]. Further, using equations (43) and (44), we see that the difference between the mass flow rate from the flexible and rigid nanochannel can be written as  [29] and from the molecular dynamic simulations by Neek-Amal et al [30] using their well-fitted rigid-wall model (shown with the black line) to the flow rates when flexibility is considered.
In equation (47), we observe that the additional non-linear terms α∆pW/EH o because of flexibility increases the mass flow rate. Also, in the following section, we will find that the mean mass flow rate in a complete pure oscillatory cycle of pressure field increases the mass flow rate from the channel, whereas in the case of a rigid channel, the increment mass flow rate due to half oscillatory cycle of pressure gets nullified from the other half cycle of pressure-field. In the following section, we will discuss the effect of flexibility and pulsating pressure field on the flow rate of these nanochannels and compare the flow rate observed in the experiments by Radha et al [29] and the predictions from molecular simulations by Neek-Amal et al [30] using the well-fitted rigid-wall model described by Neek-Amal et al [30].

Effect of thickness of the top capping graphene sheet wall on mass flow rate
To solve the above analytical equations (41)-(47), we wrote MATLAB code. Using the experimental parameters, pressure [29], and the fitted confined water properties [30], we show the mass flow rate as a function of the height of the nano/Angstrom-sized channels in figure 4(a). In figure 4(b), we show the percentage change in flow rate when the flexibility of the graphene upper capping wall is considered. In both figures 4(a) and (b), the black line shows the mass flow rate in the experimental graphene channels (where we choose α = 0, i.e. T = ∞ because α ∝ (W/T ) 3 ) to include the experimental and molecular dynamic simulation's predictions which were described and shown well-fitted using the rigid channel model [29,30]. The black arrow indicates the increasing values of upper capping wall thickness T , which varies from 35 nm (green color) to 100 nm (magenta color) in intervals of 5 nm. We find that the flow-rate prediction for the magenta line (when T = 100 nm) matches well with the experimental predictions for the rigid channel as shown by a black line (T = ∞) for all channel heights H o > 4 Å. This indicates that the pressure fields due to disjoining pressure and capillary interfacial tension were not enough to deform the upper capping of the graphene channel in experiments. Although for H o < 4 Å, we see the deviation between the magenta and black line, where upto 5-7% increment of flow rate is predicted as shown in figure 4(b). This deviation is within the error bar of the experimental values [29]. This shows that the derived flexible channel wall model could be more accurate in modeling flexible graphene nanochannel than the rigid wall model by Neek-Amal et al [30]. Assuming graphene as the incompressible material (ν = 0.5), we calculated the corresponding increasing values of α with decreasing T from the perturbation theory, where α = (1/60)(W/T ) 3 (1 − ν 2 ) as shown by the red arrow, which varies between (0.027 ⩽ α ⩽ 0.64).  We also observe as the thickness of the upper capping wall decreases, the mass flow rate inside the channel increases. For T = 35 nm, the percentage change in flow rate due to the flexibility of the wall increased by more than 80%, as shown by the green line. We noticed, for channel heights H o > 10 Å, that the change in flow rate due to the flexibility of the wall is negligible. It is because the pressure field due to disjoining pressure and capillary interfacial tension decreases with channel height (as shown in figure 3) and becomes relatively low enough to deform the upper wall of the channel.
By applying an additional constant pressure at the reservoir (i.e. at x = 0), which is ten times the disjoining pressure and capillary pressure at channel height H o = 1 Å, we show the mass flow rate prediction in figures 5(a) and (b). Here, the black and red arrows indicate the varying values of T , and α where the varying interval remains the same as in 4(a) and (b).
In 5(a), we find that the flow rate increases with channel height, unlike in figure 4(a). This is because we applied the additional pressure supply to all channel heights, whereas in Radha et al [29], the disjoining and interfacial pressure inside the channel decreased drastically with the channel height (shown in figure 3).
We find that the flow rate in the rigid channel scales asṁ rigid ∼ ∆p, whereas for the flexible channels, for negligible wall-displacement perturbation, the flow rate scales asṁ flexible ∼ ∆p for (α∆pW/EH o ) ≪ 1, and for large perturbationṁ flexible ∼ ∆p 4 . We also find that, for a given thickness T , the percentage change in flow rate in the smaller height of the channel is much larger than the larger height of the channels. For the larger height of the channel, where (α∆pW/EH o ) ≪ 1, the change in flow rate scales as (ṁ flexible −ṁ rigid )/ṁ rigid = ∆ṁ/ṁ ∼ (∆p/H o ), whereas for the smaller height channel, the wall perturbation due to deformity is significant and the change in flow rate scales as ∼(∆p/H o ) 3 . Hence, as the channel height decreases for the given reservoir pressure, ∆ṁ/ṁ increases H −1 o followed by H −3 o after a threshold with the channel height as shown in 5(b).

Effect of varying pressure on the mass flow rate
In this section, we take a channel height H o = 10 Å and the confined water properties at the H o = 10 Å. We keep the other experimental parameters and show the mass flow rate as a function of varying reservoir pressure in figure 6(a). In figure 6(b), we show the percentage change in flow rate when the flexibility of the graphene upper capping wall is considered a function of varying reservoir pressure.
In both figures 6(a) and (b), the data shown from the green line to the purple line indicate α is increasing (shown with red arrow) from 0 to 0.1, respectively. Corresponding wall thickness T is decreasing (shown with black arrow). For all data values, the nonlinear term due to flexibility is between 0 ⩽ (α∆pW/EH o ) ⩽ 0.9, where the maximum value 0.9 occurs for the purple data for α = 0.1, and ∆p/∆p Ho=1Å = 175. We find that for the rigid channel, theṁ rigid ∼ ∆p as shown with green data in figure 6(a), whereas for the flexible channels (shown between red to purple colors), the flow rate scales aṡ (1). We find that the nonlinear increment of mass flow rate in figure 6(a) is due to (α∆pW/EH o ) ∼ O(10 −1 ).  On the other hand, in 6(b), we find that the percentage change in mass flow rate for the rigid channel, and the flexible channels are ∆ṁ rigid /ṁ rigid ∼ ∆p 0 , ∆ṁ flexible /ṁ rigid ∼ ∆p for (α∆pW/EH o ) ≪ 1 and ∆ṁ flexible /ṁ rigid ∼ ∆p 3 for (α∆pW/EH o ) ∼ O(1), respectively. We also find that as the flexibility parameter α increases for the given reservoir pressure, the mass flow increases in the channel. We notice for α = 0.02 and α = 0.1, at ∆p/∆p Ho=1Å = 175, the percentage change in mass flow rate (∆ṁ flexible −ṁ rigid )/ṁ rigid is 20% (shown in red circle), and 120% (shown in purple left-side triangle), respectively. This states that, due to the flexibility of graphene sheets, the channel deformations effect becomes significantly important and substantially increases the mass flow rate.

Effect of pulsating/oscillatory pressure in the reservoir on the mass flow rate
We take a channel height H = 10 Å and the confined water properties at the H = 10 Å in this section. We use α = 0 for rigid channels and α = 0.1 for flexible channels. We give a pulsating reservoir pressure ∆p = p 1 L sin(ωt), where p 1 L (amplitude of the applied oscillatory pressure field at the reservoir, i.e. at x = 0) is varying between 0 to p 1 L = 325∆p Ho=1Å indicated by the data lines from black to pastel green. We keep the other experimental parameters as it is and show the mass flow rate in one oscillatory cycle for rigid (α = 0) and flexible channel (α = 0.1) as a function of non dimensionalize time in figures 7(a) and (b), respectively. Also, please note that we choose the ∆p has only one sinusoidal pressure field here to show the pulsation effect simply, but our derivation and model (i.e. equation (43)) for ∆p is not restricted to only one oscillatory field, infact, it can have any number of oscillatory pressure field with different amplitude and frequencies.
We find that the flow rate increases with pressure in both the rigid and flexible channels for the first half cycle (when the sin(ωt) ⩾ 0). We find that the ratio of maximum flow rate (when p 1 L = 325∆p Ho=1Å and ωt = π/2 ) isṁ flexible /ṁ rigid = 3.42. This means that the pulsating/oscillatory pressure in the reservoir Figure 8. We show the averaged mass flow ratem in one oscillatory cycle with and without oscillatory pressure field for the rigid channel (α = 0) in black and for the flexible channel (α = 0.1) with and without oscillatory pressure field in red and green, respectively as a function of applied pressure.
intensifies the mass flow rate in a flexible channel. We also noticed that, in the next half cycle (when the sin(ωt) ⩽ 0), the flow rate is in the opposite direction for both channels. In the case of a rigid wall channel, the amount of flow reversal is the same as the amount of flow in the first half cycle; this gives the net averaged flow over a complete cycle zero for all pure oscillatory pressure fields. On the other hand, in the case of a flexible channel, the amount of flow reversal is negligibly small in the next half cycle (when the sin(ωt) ⩽ 0) as compared to the amount of flow in the first half cycle. This states that the pulsating flow field shows a net positive averaged mass flow rate in a flexible channel.

Effect of pulsating pressure and flexibility of graphene sheet on the mass flow rate
In this section, we take a channel height H o = 10 Å and the confined water properties at the H o = 10 Å. We use α = 0 for rigid channels and α = 0.1 for flexible channels. We give a pulsating reservoir pressure ∆p = p o L + p 1 L sin(ωt), where (p o L and p 1 L are the amplitudes of the applied steady and oscillatory pressure field at the reservoir, i.e. at x = 0). The p o L is varying between 0 to p o L = 175∆p Ho=1Å , and p 1 = 2p o , respectively. We keep the other experimental parameters as it is and show the averaged mass flow rate in one oscillatory cycle for the flexible channel (α = 0.1) without oscillatory pressure field (i.e. p 1 = 0) in green color and with oscillatory pressure field (i.e. p 1 = 2p o ) in red color in figure 8 as a function of varying pressure, respectively. We also show the averaged mass flow rate for the rigid channel (α = 0) with the black line in the same figure 8.
We find that the flexibility of the wall increases the averaged mass flow rate for all varying pressure fields in the case of without an oscillatory pressure field compared to the rigid wall channel, as shown with green and black data lines. We find that for p o L = 162.5∆p Ho=1Å , the averaged mass flow rate in rigid channel wall is 1.7 × 10 −3 g s −1 , whereas for flexible channel (α = 0.1) in the absence of the oscillatory pressure field, it is 3.6 × 10 −3 g s −1 . Therefore, slight deformability in the channel increases the averaged mass flow rate by 2.12 times, which is significant. We further find that, as we introduced the oscillation field in the pressure in addition to the constant pressure field, the averaged mass flow rate remains unchanged in the rigid wall channel, as shown with a black line. In contrast, the oscillatory field intensifies the flow rate further in the flexible channel (shown in red). We find that for p o L = 162.5∆p Ho=1Å , the averaged mass flow rate in rigid channel wall with oscillations is 1.7 × 10 −3 g s −1 , whereas, for flexible channel (α = 0.1) with oscillatory pressure field, it is 8.8 × 10 −3 g s −1 , which is 5.2 times higher in flow rate concerning rigid channel and 2.44 times higher when there is no oscillatory pressure field in flexible channel.
We also find that as the magnitude of the pulsating pressure increases in a flexible channel, the time-averaged flow rate also increases with it. The reason is that in the case of a rigid channel, the mass flow rate is linearly proportional to ∆p; therefore, when we time-averaged the flow rate, the integration of the oscillatory part of the pressure field becomes zero in the complete cycle. On the other hand, in the case of the flexible channel, the flow rate (from equation (43)) consists of non-linear terms of ∆p; due to the nonlinearity, the time-averaging of the oscillatory part of the pressure field in a complete cycle is non-zero and increases with the amplitude of the oscillatory pressure field (shown in red).

Conclusion
In this paper, we derived a model for the mass flow rate in the rigid and deformable nanochannels for the nanoconfined water transport by using the small displacement structural mechanics and perturbation theory presented by Gervais et al [32] and Christov et al [33], respectively under the lubrication approximation. For the validation purpose, we show that the newly derived model also facilitates the flow dynamics of Newtonian fluids under different conditions as its limiting cases, which have been previously reported in the literature [30,32,33,45,62,63]. Thorough validation tests have revealed that the newly derived model produce mathematically and physically sensible flow rate in diverse situations of applied pressure, deformability, shallow channel geometry, and boundary conditions in the nanochannels. In our study, we also compare the predictions by our deformable-wall model with the experimental results by Radha et al [29] and the MD simulation predictions by Neek-Amal et al [30], which were well-predicted by the rigid-wall model. We find that for H o < 4 Å, the flow rate prediction by the deformable-wall model is 5%-7% more than the rigid-wall model. These predictions are within the errorbar of the experimental values as shown by Radha et al [29], which indicate that the derived flexible channel wall model could be more accurate to model flexible graphene nanochannel in comparison to the rigid wall model by Neek-Amal et al [30].
In our study, we focus on investigating the impact of two key factors: the flexibility of the graphene sheet and the application of a pulsating pressure field on the flow rate. We find that as the flexibility α increases (or corresponding thickness T and elastic modulus E of the wall decreases), the flow rate inside the channel increases. We find that the flow rate in flexible channels scales asṁ flexible ∼ α 0 for (α∆pW/EH o ) ≪ 1, m flexible ∼ α for (α∆pW/EH o ) ∼ O(10 −1 ) andṁ flexible ∼ α 3 for (α∆pW/EH o ) ∼ O (1). We also find that, for a given thickness T , the percentage change in flow rate in the smaller height of the channel is much larger than the larger height of the channels. As the channel height decreases for the given reservoir pressure and thickness, the ∆ṁ/ṁ increases with H −1 o followed by H −3 o after a threshold with the channel height as shown in 5(b).
Additionally, to gain insights into the mechanisms and dynamics of fluid flow in nanocapillaries and to understand how they can be manipulated or controlled for various applications, we analyze the effects of applying a pulsating pressure field, which involves periodically varying the pressure exerted on the reservoir. Due to the oscillatory pressure field, there is no change in the averaged mass flow rate in the rigid channel. On the otherhand, the flow rate increases in the flexible channels with the increasing magnitude of the oscillatory pressure field. Also, in flexible channels, depending on the magnitude of the pressure field, either of the steady or oscillatory or both kinds of pressure field, the averaged mass flow rate dependence varies from ∆p to ∆p 4 as the pressure field increases. The flow rate in the rigid-wall channel scales asṁ rigid ∼ ∆p, whereas for the deformable-wall channel it scale asṁ flexible ∼ ∆p for (α∆pW/EH o ) ≊ 0,ṁ flexible ∼ ∆p 2 for (α∆pW/EH o ) ∼ O(10 −1 ), andṁ flexible ∼ ∆p 4 for (α∆pW/EH o ) ∼ O (1).
We find that although the confined water properties increase the density upto 1.6 times in the channel and the slip also plays a significant role in increasing the mass flow rate in the nanocapillaries, the flexibility of these channels and oscillatory pressure field at the reservoir intensify the flow rate through these channels several times.

Data availability statement
All data that support the findings of this study are included within the article (and any supplementary files).