Transport across thin membranes: Effective solute flux jump

A model to describe the transport across membranes of chemical species dissolved in an incompressible ﬂow is developed via homogenization. The asymptotic matching between the microscopic and macroscopic solute concentration ﬁelds leads to a solute ﬂux jump across the membrane, quantiﬁed through the solution of diffusion problems at the microscale. The predictive model, written in a closed form, covers a wide range of membrane behaviors, in the limit of negligible Reynolds and P (cid:2) eclet numbers inside the membrane. The closure problem at the microscale, found via homogenization, allows one to link the membrane microstructure to its effective macroscopic properties, such as solvent permeability and solute diffusivity. After a validation of the model through comparison with the corresponding full-scale solution, an immediate application is provided, where the membrane behavior is a priori predicted through an analysis of its microscopic properties. The introduced tools and considerations may ﬁnd applications in the design of thin microstructured membranes


I. INTRODUCTION
Transport phenomena across porous membranes, defined as thin microstructured permeable surfaces, are massively exploited in industry.These flows play a crucial role in a wide range of industrial procedures, such as desalination, sterile filtration, food processing, petroleum refining, and other medical and environmental applications (cf.Mohanty and Purkait 1 for a review).To give an example, fog water harvesting systems are largely employed in arid climates and are usually composed of nets 2 or harps. 3The interplay between the aerodynamic flow of such permeable structures 4 and the capture of water drops is a crucial issue in these systems. 5Membranes are ubiquitous in the biological world.Some unicellular organisms use thin permeable structures in their displacement and feeding strategies, 6 or plants use them to spread their seeds. 7,8Solvent (water) and solute (sugar) translocation across aquaporin porous channels constellating the cellular membranes is a primary process for the good performance of organisms, spanning from plants 9,10 to animals. 11Serious genetic diseases are directly linked to mutations of aquaporin channels. 12Each abovementioned filtration process and the related transport phenomenon are characterized by at least two length scales, the macroscopic size of the membrane itself and the characteristic length of the pores of the membrane, which spans from the nanometer 13 to the millimeter [14][15][16] scale.The interplay between these two intrinsically different length scales renders the description of this phenomenon extremely complex, while their deep understanding and predictive modeling are of great interest for engineering applications and medical progress.In particular, from the industrial viewpoint, separation processes constitute 10%-15% of the world energy consumption, 17,18 thus highlighting the necessity of a better understanding of these mechanisms.Several attempts to address this issue are proposed in the literature in the case of classical porous media whose dimension along the filtration direction is comparable to the macroscopic length at play (called here thick membranes).For thick membranes, two principal approaches used to analyze the solute-solvent flow across membranes are identified.The first approach relies on the employment of macroscopic models to mimic the presence of the solid skeleton of the membrane through a set of equivalent equations modeling the behavior of the coupled solute-solvent flow through the membrane.These formulations are based on the so-called Darcy law 19 or its Brinkmann extension, 20 coupled with a solute transport equation.The advective and diffusion terms are properly modified to account for the presence of the microstructure of the pores. 21A critical issue of these models stems from their weak predictive power since, in the macroscopic equations, some local properties of the porous medium are treated as free parameters and, thus, ad hoc calibrated for the physical configuration studied or supported by experimental measurements. 22n alternative to the Darcy or Brinkman approaches consists of the massive and consuming task of resolving the full-scale physics through continuous first principles, e.g., the Navier-Stokes equations for the fluid solvent and the advection-diffusion equation for the transported solute. 23In this case, despite the high level of detail in the solution at the pore-scale, the global physics underlying the phenomenon remains hidden by the particular geometry considered.
Recent patents on ultra-thin carbon-based membranes have increased the interest also in the quasi-two-dimensional, so-called thin, membranes. 44While these classical approaches used to analyze porous media are suitable to model fluid flows across thick membranes, this is not the case for thin membranes, whose literature is more limited [45][46][47] and suffers of little generality in terms of flow configurations and membrane geometry.Macroscopic models for thin membranes semipermeable to the solute are commonly based on the so-called Kedem and Kaltchasky empirical law, 48 which qualitatively describes the average solute and solvent flux due to osmosis.On the other extreme of the spectrum, molecular dynamics simulations across Angstrom-sized pores are widely employed, 49 but still incur computational and conceptual limitations, not allowing to link them to a continuum-mechanics approach.Links between these two approaches are still in a developmental phase, e.g., in Bacchin 50 or Cardoso and Cartwright, 51 where the microscopic characterization of the solvent-solute behavior was included in a continuum macroscopic model via the introduction of a parameter found via kinetic theory.
Homogenization showed great potential in describing the flow at the interface between a free-fluid region and a porous medium 35 or the flow in the vicinity of a rough surface, 52 through a Navier slip condition which recovers the stress-jump at the interface, similar to the pioneering studies of Beavers and Joseph. 53Zampogna and Gallaire 54 developed a homogenization-based model to describe solvent flows across microstructured thin membranes, which provides a predictive and explicit link between the microscopic geometry of the membrane and their permeability and slip properties.The so-called effective stress jump model expresses the proportionality between the velocity and the stress-jump across the membrane and already showed potential in the optimization-based design of membranes to build full-scale geometries satisfying a given macroscopic hydrodynamic objective. 55n the present paper, we study the effects of the microstructured surface on the concentration field for thin membranes fully permeable to both solute and solvent.While, in Ref. 54, only the hydrodynamics of the fluid across the membrane was considered, here, we extend this solvent model to the passive transport of a dilute solute through a microstructured permeable surface characterized by different physical-chemical behaviors.The paper is structured as follows.In Sec.II, we introduce the first principles governing the solvent and solute flow.Sections III and IV are devoted to the derivation of the equivalent model and the effective equations describing the solute concentration flow, for different behaviors of the membrane.The microscopic problems associated with the macroscopic properties of the membrane are solved in Sec.V for a particular geometry.Once the microscopic solution is known, the macroscopic model is validated in Sec.VI through comparisons with the full-scale solution in some representative configurations.In Sec.VII, we conclude by proposing an analysis of the microscopic results for a given fractal-like micro-structure, with the aim of understanding the macroscopic physical-chemical properties of the membrane and in the perspective of limiting the trial-and-error procedure which often rules these tasks.

II. EFFECTIVE CONCENTRATION JUMP MODEL
We consider a solute of molecular diffusivity D transported by a flow of an incompressible Newtonian fluid, denoted as the solvent, of constant density q and viscosity l.We neglect variations of the solvent properties due to the concentration of the solute.The solute-solvent flows across a microstructured permeable surface (from now on a membrane, cf.Fig. 1) formed by a periodic repetition of a given solid inclusion so that a microscopic elementary cell [highlighted by the dashed rectangles in Figs.1(b) and 1(c)] can be identified.A separation of scales related to this transport phenomenon subsists, i.e., l ( L; where l denotes the characteristic thickness of the membrane and L the characteristic length of the flow, thus leading to the definition of an infinitesimally small separation of scales parameter The membrane is fully permeable to the solute.In the case of negligible inertia, the solvent velocity and pressure ðû i ; pÞ and solute concentration ĉ are governed by the Stokes and advection-diffusion equations, 56 valid in the whole fluid domain (cf.Fig. 1).Exploiting the Einstein summation convention, these equations can be written as and where R ij and F i denote the stress tensor and concentration flux, i.e., The differential equations describing the solute-solvent flow past the membrane M (cf.Fig. 1) are closed by specifying the boundary conditions on @M, i.e., the walls of the microscopic structure forming the membrane.The solvent velocity satisfies the no-slip condition ûi ¼ 0 on @M; (5) while the solute concentration satisfies a generic boundary condition of the form In Eq. ( 6), ni denotes the outward normal to the membrane walls @M [cf.Fig. 1(c)], ĵ the adsorption/desorption coefficient, and Ĉw an eventual constant value of concentration imposed at the membrane.
The coefficients a, b, and c allow one to change the boundary condition imposed on @M, i.e., homogeneous and non-homogeneous Neumann, Dirichlet, and Robin conditions.Some combination of a, b, and c corresponds to precise chemical-physical interactions between the solute molecules and the membrane walls, thus describing a variety of operative conditions of the membrane walls and their effect on the solute concentration, spanning from non-reactive walls (homogeneous Neumann, a ¼ 1 and b ¼ c ¼ 0), 40 to chemostats (Dirichlet, a ¼ 0 and b ¼ c ¼ 1), 57 through partially absorbing/desorbing materials (Robin, a ¼ b ¼ 1 and c ¼ 0). 58he interplay between the characteristic length of the pores l and the macroscopic length scale L defined by Eq. (1) indicates the multiscale nature of the phenomenon.Here, we propose a simple macroscopic model which faithfully reproduces the flow behavior without the requirement of representing the full-scale geometry.A homogenization-based procedure has been developed in Zampogna and Gallaire 54 and applied to Eq. (2), leading to the following solvent velocity conditions: where the bar over the variables represents the spatial average which, for a generic function f, is defined as In the above equation, jC F j is the surface area of the fluid part of C and ð t ;ŝÞ are the tangential-to-the-surface coordinates, as indicated in Fig. 1.All variables in Eqs. ( 7) and ( 8) are purely macroscopic, i.e., averaged exploiting Eq. ( 9).While Eq. ( 7) imposes the continuity of the solvent velocity across the membrane, Eq. ( 8) describes a jump in the fluid stress due to the presence of the membrane, which depends on the macroscopic solvent velocity and on the tensors M ijk and N ijk .The stress-jump, thus, states that the average velocity embeds the macroscopic effect of the microscopic viscous stresses at the solid walls, similar to the flow through a porous medium, 32 at the interface of a porous medium 35 or on a rough wall. 52These proportionality tensors are calculated via Stokes problems within the microscopic domain F deduced via homogenization.They are, respectively, called upward and downward Navier tensors.Their non-zero components represent the ability of the flow to permeate across or slip over the membrane.We refer to Zampogna and Gallaire and Ledda et al. 54,55 for further insights on the solvent model.In the following, we develop a macroscopic boundary condition for the solute flow, in analogy with Eqs. ( 7) and (8).

III. PREAMBLE
Before applying the homogenization technique to the advection-diffusion problem, we perform a preliminary step.We focus on the region of space in the vicinity of the membrane, and we consider a control volume whose size corresponds to that of the microscopic elementary cell defined in Fig. 1.In the absence of unsteady nearmembrane phenomena, upon employment of Gauss' theorem and enforcing periodicity in the membrane directions, the integration of Eq. (3) over the control volume gives where dA represents the difference of fluxes between the upward and downward side of the control volume, U and D. Substituting Eq. ( 6) in the RHS, we obtain where hÁi denotes the surface average over @M, i.e., If the diffusive regime is dominant in the vicinity of the membrane, i.e., when variations of ĉ in F are negligible, the following relation holds hĉi % j@Mjĉ : In this case, Eq. ( 11) is rewritten as with ĵeff ¼ ĵj@Mj.Equation ( 14) can be used together with model ( 7) and ( 8) to calculate a reference value of concentration Ĉ0 at the macroscopic membrane, which depends on the values of a, b, and c.While in the case of Neumann and Robin condition the value of Ĉ0 is not known a priori, for the case of Dirichlet Ĉ0 assumes the constant value of Ĉw on the membrane.
Figure 2 shows a comparison between the solution of Eq. ( 11) and full-scale simulations for a test flow configuration 73 where a linear membrane formed by ten solid square inclusions splits a twodimensional channel of height L ¼ 10 À3 m.The channel ideally represents a portion of a microfluidic device, freely inspired by highperformance computer refrigeration systems, and is filled by water at 20 C (q ¼ 10 3 kg=m 3 and l ¼ 10 À6 m 2 =s) and ethylene glycol (D ¼ 10 À9 m 2 =s), resulting in one of the most common coolant solutions. 59The membrane therein represents a micro-structured baffle, i.e., a solid structure which favors heat exchange. 60The velocity profile and solute concentration are constant at the inlet L (u 1 ¼ 10 À6 m=s and c ¼ 1 mol=m 3 ).With these values, the P eclet number associated with the channel, defined as u 1 L=D, is unitary, while viscous diffusion is faster than mass diffusivity since the Schmidt number Sc ¼ =D is equal to 10 3 .The top and bottom walls of the channel, T and B, are made of aluminum alloy 3003 (Ref.61) and react with ethylene glycol.The macroscopic observable of this reaction is an adsorption condition @ n ĉ ¼ Àĵĉ with ĵ ¼ 0:01 mol=m 2 .At the outlet R, a zero-normal solvent stress and solute flux have been imposed.Panels (b) and (c) show the values of solute concentration over the centerline of the membrane, C, for the chemostat-like (a In panel (b), the solid inclusions of the membrane have side equal to 10 À6 m, while in panel (c), their side is 9:9 Â 10 À5 .In both cases, the membrane thickness corresponds to the side of the solid inclusions.The discrepancies between the macroscopic model and the averaged full-scale results can be relatively important, depending on the case.Interestingly, the case of constant concentration on the solid inclusions does not give a constant concentration profile, but such value of concentration varies along the membrane.
To improve this first estimate, a formal interface condition for the passive transport of species across micro-structured membranes is developed via a homogenization-based procedure.This procedure is applied to c Ã ¼ ĉ À Ĉ0 , which still satisfies the advection-diffusion problem (3), while the boundary condition ( 6) on @M is rewritten as

IV. HOMOGENIZATION PROCEDURE
In the vicinity of the membrane, variations of the quantities occur at the scale l ( L. We refer to this problem as the microscopic or inner problem, defined within the microscopic elementary cell, in opposition to the macroscopic or outer problem, far from the membrane, where the effects of the microscale are not present.The flow equations are normalized by employing the typical inner and outer scales.A multiple scale decomposition of the normalized inner problem is then performed taking into account the far-field conditions on U and D (see Fig. 1 for their definition), which ensure the asymptotic matching between the inner and outer problems.The successive step consists of the spatial average of the microscopic quantities to obtain a homogenized macroscopic model that describes the microscopic behavior of the membrane.The result of this procedure is a set of purely macroscopic interface conditions for the outer problem defined over the macroscopic homogeneous membrane C, where the fluid region F and the solid region M are indistinct.

A. Outer problem
We focus on the problem far from the membrane, in the outer upward and downward regions, whose variables are denoted with the apex Á O .The reference scales for the outer quantities are where DC O ¼ ĈO À Ĉ0 and ĈO is a reference value for the concentration, defined by the macroscopic, far from the membrane, problem.
The outer Peclet number is defined as The non-dimensional governing equation in the outer upward and downward region is

B. Inner problem
The inner problem is valid in the microscopic domain F identified in Figs.1(b) and 1(c) by the dashed rectangles.Variables within the inner domain are denoted with the superscript Á I .The boundaries of the microscopic domain S 1 ; S 2 and T 1 ; T 2 are characterized by periodicity owing to the repeated microstructure of the membrane.We non-dimensionalize Eq. (3) using the following relations: where DC I is set by the balance of diffusive fluxes between the inner and outer domain, i.e., while U I is the inner solvent velocity and l is the microscopic length scale defined in Fig. 1(b).With these normalizations, Eq. ( 3) is dimensionalized as follows: with the non-dimensional solute flux defined as The non-dimensional version of the boundary condition on @M introduced in Eq. ( 6) is To formally develop the macroscopic model, we assume that j ¼ OðÞ, i.e., the membrane can absorb/desorb a Oð1Þ quantity of solute over a length equal to L [cf.Dalwadi et al. 37 and Chernyavsky et al. 62 where it is shown that, for bulk porous media, a value of j larger than OðÞ implies an order 1= effect on the solute flux with a subsequent breakdown of the asymptotic expansion adopted to develop the model].

C. Matching conditions between the inner and outer problem
The problem introduced in Sec.IV B is well defined only upon the specification of the boundary conditions on U and D. These boundary conditions are deduced by the matching between the inner and outer problems, given by the following dimensional equations: and where the superscripts Á þ and Á À denote the inner and outer side of U and D, according to Fig. 1(c).Equations ( 24)-( 26) represent the continuity of solvent velocity, solute concentration, and normal-tothe-membrane solute flux, respectively.Using the reference scales introduced in Secs.IV B and IV A, the non-dimensional version of interface conditions (25) and ( 26) reads and In Eq. ( 28), the continuity of the advective part of the flux u O i c O n i is ensured by the continuity of the solute concentration and solvent velocity fields across the membrane.Hence, Eq. (28) reduces actually to While always considering that the actual contribution to the flux jump is given by the diffusive contribution, we pursue our analysis with the complete flux interface condition (28) to preserve the conservative structure of the governing equations.

D. Multiple scale decomposition
Since the structure of the membrane allows one to uniquely identify the separation of scales parameter ( 1 introduced in Eq. (1), we decompose the inner unknown fields with a multiple scale expansion.We introduce the fast (microscopic) and slow (macroscopic) variables together with the expansion n¼0 n c I;ðnÞ ðx; X; tÞ: The spatial derivatives are transformed following the rule where capital indices denote derivation with respect to X i .Assuming the inner Peclet number , we obtain the following leading order set of equations: with The flux interface condition at the outer boundaries of the microscopic domain U and D assumes the form The boundary condition on @M can be finally written as where the terms multiplied by j disappear at the leading order since j $ , and the RHS has been accordingly modified to treat the three types of boundary conditions in a compact form.

E. Solution of the leading-order problem
We now proceed to the solution of the leading-order problem.Since F O i j U and F O i j D are source terms for the set of linear equations ( 33)- (36) and do not depend on the integration variable x i , they form a basis for the space of the solutions.As a consequence, the solution is formally written as Substituting Eq. ( 37) in the set of Eqs. ( 33)- (36), the quantities T i and Y i satisfy the microscopic problems and Each problem above is a diffusion problem enforced by nonhomogeneous Neumann boundary conditions on U; D and with periodic boundary conditions over S 1 ; S 2 and T 1 ; T 2 .Note that the two problems are formally analogous in the bulk domain since they differ only for the inversion of the boundary conditions between U and D.

F. Retrieval of the final macroscopic solution
Applying the spatial average defined in Eq. ( 9) to the solution (37) of the leading-order problem, we write the macroscopic concentration field at the membrane as where, in the outer fluxes F O i j U;D , only the diffusive part counts because of continuity of velocity and concentration fields.Rearranging Eq. ( 40) for the original, dimensional concentration field ĉ, we obtain Equation ( 41) states the existence of a jump in the diffusive flux across the membrane, and thus, the presence of the membrane changes the slope of the solute concentration profile.Condition ( 41) is also reminiscent of the stress-jump condition for the solvent flow derived in Zampogna and Gallaire, 54 in which the velocity was proportional to the stress-jump (which represents the flux for the Navier-Stokes equations) weighted by so-called upward and downward Navier slip tensors.The vectors T i and Y i can be interpreted as effective diffusion vectors of the equivalent membrane, which measure the modifications of the membrane-normal diffusion due to the presence of the solid inclusions.
Both solvent and solute interface conditions are hence perfectly in line with previous homogenization-based results, such as the case of the Navier slip condition on a rough wall, 52 where the velocity at the wall was proportional to the stress via a slip tensor.In the particular case of an impermeable rough surface, condition (41) reads ĉj C À Ĉ0 ¼ T n @ N c O j C À .The vector component T n is, therefore, interpreted as a chemical slip, which quantifies the deviations of the solute concentration from the case in which the membrane corresponds to a smooth, single-scale impermeable surface.
The macroscopic interface condition is based on the average of the solution of the microscopic problems (38) and (39).After this step, their averages are used in Eq. ( 41), together with the continuity of the concentration field over the interface C across the membrane, i.e.,

V. EFFECTIVE DIFFUSION VECTORS
In this section, problems (38) and (39) are numerically solved to deduce the components of the effective diffusion vectors for the different values of a and c.The values assumed by these parameters have been discussed after their introduction in Eq. (6).Since during the homogenization procedure the parameter b is dropped [cf.Eq. ( 36)], only two different sets of values for a and c are meaningful; the couple ða ¼ 1; c ¼ 0Þ represents insulating and adsorbing membranes, while ða ¼ 0; c ¼ 1Þ corresponds to chemostat-like membranes.The numerical solution relies on a weak form implementation in the finite-element solver COMSOL Multiphysics.The spatial discretization is based on cubic triangular elements. 74s a test case, two-dimensional square parametric inclusions are considered.The microscopic problems (38) and (39) are solved in the frame of reference of the membrane denoted with (s, n), where spatial variables are normalized with respect to the microscopic length l.Figures 3(a For this reason, we restrict our analysis to T n , while similar considerations apply for Y n using the geometrical argument represented by Eq. (43).To apply the macroscopic model described by Eq. ( 40), the spatial average defined in ( 9) is computed.
The size of the microscopic domain along the normal-to-themembrane direction, jU À Dj, is not specified in the inner problem since, in principle, these boundaries are located at an infinite distance from the membrane to correctly match the inner and outer solutions.Table I shows that the values of T n and Y n reach an asymptotic value for jU À Dj tending to infinity, assessing the well-posedness of the problem.
In Fig. 3(c), the values of T n are shown in the case of square inclusions as a function of the fluid-to-total ratio along C, i.e., # ¼ jC F j=jC F [ C M j, in the range of 0:01 < # < 0:99.The blue curve corresponds to the insulating/adsorbing case (a ¼ 1; c ¼ 0), while the orange one to the chemostat case (a ¼ 0; c ¼ 1).In the case of a ¼ 1; c ¼ 0; T n has a diverging behavior approaching # ¼ 0, and hence, it can be interpreted as a measure of the deviation of the concentration field with respect to the case when the inclusions are not present.In the absence of inclusions, one can write lim meaning that the concentration field and its normal-to-the-membrane flux are continuous or, in other words, the membrane is not present anymore.On the contrary, the normal component of the vector T for a ¼ 0; c ¼ 1, i.e., when the membrane acts as a chemostat, tends to 0 for # tending to zero.In this case, T n measures the deviation of the concentration from the values imposed on the solid inclusions, which increases with the fluid-to-total ratio, confirming the analogy with the Navier slip condition 52 outlined in Sec.IV.The coefficient j eff =j ¼ j@Mj in Fig. 3(d) is equal to the perimeter of the inclusion M and represents the main contribution of the effective adsorption rate of the membrane.As shown in Sec.VII, j eff is of fundamental importance in the design of efficient adsorbing/ desorbing membranes.
In this section, we showed the results of the microscopic calculations for a particular geometry of the inclusions, i.e., squares.The validation of the homogenized model requires its comparison with full-scale simulations.In the following, we, thus, investigate different macroscopic configurations to validate the microscopic model proposed for the membrane.

VI. MACROSCOPIC VALIDATION OF THE EQUIVALENT MODEL
In Sec.VI, we proceed with the validation of the equivalent macroscopic model by means of comparisons with full-scale simulations of the coupled solvent-solute flow.For the sake of conciseness, we will focus on the solute flux since the solvent flow has been already extensively validated in Zampogna and Gallaire and Ledda et al. 54,55 The solvent-solute flow configurations chosen as testing benches are sketched in Fig. 4, while in Table II, conditions on top (T), bottom (B), left (L), and right (R) boundaries of the macroscopic computational domain are listed.The results presented in this section have been rendered non-dimensional using the macroscopic outer reference scales for the solvent velocity and pressure, the solute concentration, and the lengths.With such normalization, the knowledge of the inner representative quantities introduced in Sec.IV is not needed to apply the macroscopic model.In particular, the outer macroscopic reference concentration has been defined, for both cases considered, as where ĉj L is the value of concentration imposed at the inlet L of the macroscopic domain, while ĉmin is the minimum value of concentration which tends to zero as x 1 !þ1 since the boundaries T and B adsorb the concentration at a given constant rate j.The approximation of the solute transport up to OðÞ requires the knowledge of both C 0 and c I;ð0Þ .The mathematical solution of the macroscopic problem, thus, relies on two steps.The first step consists of the solution of the non-dimensional version of Eq. ( 3) together with interface conditions ( 14) and ( 42).In the second step, Eq. ( 46) is solved together with the interface conditions ( 40) and (42)  where the solution of the first step is employed.The Stokes equations satisfied by the solvent flow, once non-dimensionalized, assume the form together with the continuity of velocities across C and the normalized stress jump condition All equations mentioned above are numerically implemented via their weak formulation in the finite element solver COMSOL Multiphysics, using a domain decomposition method 63 to couple the upward and downward solvent flow and solute fluxes.In this framework, macroscopic models ( 42) and ( 40) and ( 7) and ( 8) are interface conditions between two domains, respectively.To exchange information from the upward to the downward domain, the stress jump and concentration flux conditions are implemented by exploiting the interface integral emerging from the weak formulation of the corresponding governing equations, while, to exchange information from the downward to the upward domain, the continuity of velocity and concentration is imposed via a Dirichlet boundary condition.The spatial discretization is based on the Taylor-Hood (P2-P1) triangular elements for the solvent flow solution and cubic triangular elements for the solute concentration solution. 75 Configuration C1 In Sec.VI A, we verify the accuracy of the model (40) in the case of configuration C1.A channel flow is considered.The channel is split by a membrane placed on the line x 1 ¼ 0, whose microscopic structure (represented in the red circular inset of Fig. 4) is made of square inclusions such that the porosity is # ¼ 0: and adsorbing behaviors (a ¼ b ¼ 1; c ¼ 0), respectively.In the present case, for the chemostat-like membrane, the value of the concentration field C 0 ¼ C w on @M has been set equal to 0.2 while, for the absorbing case, the value of the absorption rate j on @M is .In each frame, the colored and grayscale isocontours are referred to as the fullscale solutions, while the red corresponding isolevels refer to the solution of macroscopic model (40).We notice that the grayscale in frames (a) and (e) has been adopted to appreciate the small variations of concentration in the downward region, which is rendered almost constant by the presence of the membrane.In the lower half-part of each frame, the solvent flow streamlines have been sketched to have a global idea of the coupled solvent-solute flow behavior.The comparison between isocontours and isolevels is satisfactory in all cases shown.A quantitative comparative analysis between the full-scale solution and the macroscopic model depicted in Figs.5(b), 5(d), and 5(f) confirms the good agreement.In Fig. 5(b), the concentration field is sampled over C, corresponding to the vertical line x 1 ¼ 0 in the fluid domain.For each membrane behavior indicated in the graph, black curves represent the full-scale concentration profile, the colored curves correspond to the macroscopic solution, and the blue stars the full-scale profile averaged using definition (9).In the remaining two frames, Figs.5(d) and 5(f), the concentration c and its gradient @ 1 c along the normal-to-the membrane direction have been sampled on the horizontal line x 3 ¼ 0 depicted by a dashed black line in the left frames of Fig. 5.The macroscopic solution in color well reproduces the full-scale profiles in black.Besides, the simulations confirm that the membrane produces a jump in the normal-to-the-membrane direction that can be estimated via the effective flux jump condition.We finally focus on the correction FIG. 5. Overview of the solute concentration behavior in configuration C1 and comparison between the full-scale and the macroscopic solution [Eq. ( 40)].Panels ðaÞ; ðcÞ, and (e): isocontours of the full-scale concentration field (in colors and gray scale) and isolevels of the equivalent, macroscopic fields in red for the membrane behaviors indicated in the title.The grayscale in frame (a) and (e) has been adopted to appreciate the small variations of concentration in the downward region, which is rendered almost constant by the presence of the membrane.Solvent flow streamlines are represented in the lower half-part of each frame.Panels ðbÞ; ðdÞ, and (f): quantitative comparison between macroscopic (colored continuous lines) and full-scale quantities (black lines).In panels (b) and (d), the concentration is sampled along C (line x 1 ¼ 0 and x 3 ¼ 0, respectively), while in panel (f), the concentration horizontal flux is sampled over the line x 3 ¼ 0. Blue stars represent the average of the full-scale solution using definition (9).
given by the homogenization-based condition (40) to the value of concentration C 0 estimated via Eq.( 14).For the insulating case, the solution C 0 has a precision similar to C 0 þ c I;ð0Þ , with a very small difference with the full-scale solution.In the other cases, we observe a larger relative difference between C 0 and C 0 þ c I;ð0Þ , which appears to be of the same order of .A remarkable conceptual improvement given by Eq. ( 40) can be noticed for chemostat-like membranes (a ¼ 0; b ¼ c ¼ 1) since it is able to reproduce the macroscopic variability of the concentration profile along C, poorly predicted to a constant value by C 0 (Fig. 6).

B. Configuration C2
A free uniform flow of a solvent past a vertical membrane positioned at x 1 ¼ 0 and for À0:5 < x 3 < 0:5 is considered, carrying a solute whose concentration is maintained constantly equal to 1 on L for À0:5 < x 3 < 0:5, and it is adsorbed at a rate j ¼ on T and B (cf. Table II).The membrane, whose structure is sketched in the red circular inset of Fig. 4, is characterized by # ¼ 0:5 and ¼ 0:02.With these geometrical parameters, T ¼ À Y ¼ ð0:3749; 0; 0Þ for the insulating/absorbing membrane (a ¼ 1; c ¼ 0) and T ¼ À Y ¼ ð0:016 38; 0; 0Þ for the chemostat-like membrane (a ¼ 0; c ¼ 1), while the only non-zero components of the tensors M and N are M 111 ¼ À N 111 ¼ À0:0120 and M 331 ¼ À N 331 ¼ 0:0015.A substantial change with respect to configuration C1 is introduced in configuration C2 since the fluid is not constrained to pass through the membrane.As a consequence, the inner Peclet number, is 2 times lower than the outer Peclet.We can, thus, test the model in a dominant advective regime since, here, we push the Peclet number up to À1 , being still within the constraint Pe I ¼ OðÞ.In analogy with the previous case, a comparative qualitative overview of the solution is given in Figs.7(a 7(d)] and for @ 1 c along x 3 ¼ 0 [Fig.7(f)], confirms the reliability of the model.The comparison between C 0 and C 0 þ c I;ð0Þ is shown also for this configuration in Fig. 8, exhibiting the same order of precision as in configuration C1, with a very small relative difference for the insulating case and a relative error of order OðÞ for the chemostat-like and adsorbing cases.

VII. AN APPLICATION OF THE MACROSCOPIC MODEL: SNOWFLAKES-SHAPED ZIG-ZAG BAFFLES IN MICRO-CHANNELS
In Secs.I-VI, we derived a homogenized boundary condition that accounts for the solute behavior through the membrane, validated against full-scale simulations.In the following, we propose a simple and immediate application, the adsorption of membranes composed of recursive, fractal-like structures, which reveals the potential of the homogenized model in terms of numerical efficiency and physical interpretation.We consider the solvent-solute flow in a channel as sketched in Fig. 9(a).The concentration field c is carried by the fluid from left to right as in configuration C1 presented in Sec.VI A. In the present case, the membrane has a zig-zag shape of length 3L.Each segment forming the membrane is rotated 30 clockwise or counterclockwise with respect to the central longitudinal axis of the channel.In Sec.VI, we considered solely membranes whose normal direction is constant in the Cartesian frame of reference and parallel to the flow direction.Here, the particular macroscopic setup chosen allows us to validate the model also in the case of arbitrarily oriented and shaped membranes.The microscopic structure of the membrane is sketched in Fig. 9(b).Each solid inclusion forming the membrane has a snowflake shape built applying the von Koch geometrical construction to each side of an equilateral triangle. 64The number R, called here recursion number, indicates how many times the geometrical axiom has been applied and, hence, the recursion level of the solid inclusions.The side of the triangle representing the base shape from which the snowflakes have been built is equal to 0:5l, implying that the porosity of the membrane is equal to 2/3 for each R.Ten different configurations have been considered, for R 2 ½1; …; 10, with fixed separation of scales parameter ¼ 0:01.In the following, we (i) evaluate the microscopic quantities associated with these fractal-like inclusions varying the recursion number and (ii) exploit the macroscopic model to efficiently solve the configuration introduced in Fig. 9 and study the effects of the microscopic structure on the solute and solvent flow.

A. Microscopic insights about fractal-like solid inclusions
The macroscopic model (42) and (40) requires the knowledge of the quantities j eff ; T and Y. Additionally, also the components of the third-order tensor M and N are needed to macroscopically solve the solvent flow via the model introduced in Eqs. ( 7) and (8).Similar to Sec.V, the microscopic quantities are calculated in the local frame of reference of the membrane so that the subscripts refer to the tangential (s) and normal-to-the-surface (n) components.Figure 10 shows the values of the non-zero entries of the effective microscopic vectors and tensors for each N.The vector Y and the tensor N are not considered here since they can be deduced from T and M by geometrical arguments, as explained in Sec.V and in Zampogna and Gallaire. 54The quantities T n ; M nnn , and M ssn reach an asymptotic value after a few recursions of the von Koch curve (as already noticed in Alinovi et al. 65 for rough, impermeable fractal-like surfaces).Indeed, the relative differences between two consecutive recursion levels are smaller than 2 already at R ¼ 4. In opposition, j eff is proportional to the perimeter of the inclusions which, for this particular kind of fractal-like structures, grows as ð4=3Þ R .From these simple considerations obtained from the solution of the microscopic problems, we deduce that the drag coefficient C D of the macroscopic flow through the membrane remains constant as R increases.This can be easily seen by writing explicitly the effective stress jump condition introduced in Eq. ( 8) for the solvent flow and collecting the difference between the stresses on the two sides of the equivalent membrane, i.e., where 54,55 On the contrary, the adsorption of the membrane increases with R since it depends on the effective adsorption coefficient j eff , which is proportional to the perimeter of each solid inclusions.In Sec.VII B, we carry out macroscopic solutions of the homogenized model for membranes with fractal-like solid inclusions to confirm our a priori deductions on the membrane behavior for varying R.

B. Macroscopic evidence of the membrane adsorption
The behavior of the adsorbing membrane (j ¼ ¼ 0:01) is now investigated via the macroscopic solution of the solvent-solute flow configuration of Fig. 9. Figure 11(a) provides a qualitative description of the solvent flow behavior via isocontours of pressure (in colors) and flow streamlines (solid black lines).The equivalent membrane is highlighted by the solid red line.Frames (b) and (c) in the same figure represent the isocontours of the solute concentration field for R ¼ 1 and R ¼ 10, respectively, showing relevant quantitative differences.The differences can be better quantified by sampling the concentration field on the membrane for each value of R [cf.Fig. 12(a) where each curve represents a different value of R as indicated in the associated colorbar].Two quantities of interest D CD and g are defined, which measure the increase in drag with respect to the case R ¼ 1 and the ability to adsorb the solution in the vicinity of the membrane with C D ðRÞ is the drag coefficient associated with the membrane whose solid inclusions have a recursion number R. The coefficient D CD ðRÞ shows variations of less than 0.2% from R ¼ 1 to R ¼ 10. Figure 12(b) shows the behavior of g with R, concluding that a negligible increase of the drag corresponds to an increase in the adsorbing ability of the membrane from 40% to 95%, without changing the chemical properties of the material forming the solid inclusions, but only changing their shape.
These results confirm the reliability of the model in predicting the behavior of a membrane by a simple analysis of its microscopic properties that can be precisely deduced through homogenization.Beyond the immediacy of the physical interpretation, we note another remarkable strength of the macroscopic model consisting of a huge gain in computational resources obtained by using the macroscopic model.Since the membrane considered here is formed by fractal inclusions, an increase of R dramatically increases the computational time needed for a well spatially resolved full-scale simulation.The full-scale solution is, hence, computed only for the recursion case R ¼ 1, needing about 3.5 CPU h to obtain spatially converged results.The estimated adsorption ability is identified in Fig. 12(b) by a yellow star, while the difference in C D with respect to the equivalent macroscopic case is equal to 0.21%, in agreement with the accuracy of the model, estimated to be of order 1% for the chosen value of .Similar accuracy is observed from the introduction of a local relative error f between the average of the full-scale concentration and its macroscopic analogous, evaluated over the macroscopic equivalent membrane C.Such f, which within the membrane is of order 2 , reaches values up to when approaching the edges of the equivalent macroscopic membrane C.This behavior of the local relative error f has to be traced back to the fact that, in the vicinity of the edges, the microscopic periodicity assumption fails.In opposition, the time needed to carry out the homogenized solution increasing the fractal number can be calculated by the sum of two contributions: (i) the one of the macroscopic simulation, which is constant with R and (ii) the time needed to estimate the tensors from the microscopic geometry.Hence, the CPU hours needed to obtain the macroscopic solution range from 0.0083 for R ¼ 1 to 0.0333 for R ¼ 10, hundreds times faster than the cheapest full-scale case (R ¼ 1).

VIII. CONCLUSION
In this work, we developed a homogenization-based macroscopic model for the transport of a dispersed solute in the incompressible flow of a Newtonian fluid (solvent) across thin microstructured fully permeable membranes.The macroscopic feedback of the microstructure on the transport phenomenon is obtained through a macroscopic condition, derived by matching inner and outer asymptotic expansions. 54The macroscopic condition covers a wide range of chemical-physical behaviors of the membrane, from constant solute concentration to the non-reactive, insulating case, passing through partially adsorbing membranes.The model developed, exact up to first order the separation-of-scales parameter, assesses the continuity of the concentration field across the membrane and the existence of a jump in the normal-to-the-membrane flux of the solute concentration, for the cases mentioned above.The concentration flux jump is quantified via the spatial average of the solutions of microscopic problems derived from the homogenization procedure, which represents the solvability conditions for the final macroscopic equations.The parameters quantifying the flux jump are (i) the effective adsorption coefficient j eff ¼ j@Mjj and (ii) the upward and downward slip or effective diffusion vectors, T and Y.These quantities decouple the microscopic geometry from its macroscopic effect on the flow, establishing a precise link between the microscopic shape of the membrane and its diffusion and transport properties.The numerical solution of the macroscopic model was validated through comparisons against the corresponding full-scale solutions.In the last part of the work, we proposed a simple design of microstructured membranes through an analysis of the microscopic quantities associated with fractal-like inclusions.
The model developed allows one to treat a multi-scale problem through simple and physically meaningful macroscopic interface conditions.It can be generalized to three-dimensional structures and nonperiodic microstructures and constitutes a powerful tool in the design of membranes, 55 together with the effective stress jump model for the solvent. 54As shown in Ledda et al., 55 the optimization procedure of such micro-structures is dramatically simplified by the reduction of the degrees of freedom describing the solid inclusions.These degrees of freedom, in the homogenized model, are represented by the entries of the Navier tensors and effective diffusion vectors.A direct and practical application of the model in its actual state is the optimization of micro-fluidic devices, 66 while its extension toward inertial and high-P eclet flows may find applications in the development and optimization of chaotic micro-mixers. 67This work represents a first step toward the definition of the right balance between selectivity and filterability in membrane design. 18It completes the effective stress jump 54 for the solvent flow adding a one-way solvent-to-solute coupling.Hence, two direct extensions of the model can be identified.The first one consists of substituting the one-way coupling considered here with a two-way coupling, where the solute concentration is able to modify the density of the solvent and, hence, its flow.The second one comes from one relevant limitation of the model, i.e., its ability to treat semi-permeable membranes (and, hence, osmosis-related phenomena) by a modification of the spatial average used in this work, which intrinsically renders every solute and solvent field continuous.These two extensions of the model may lead one to the development of a rigorous and closed form of the Kedem-Katchalsky law for thin membranes. 48Besides, an emerging field in membrane science consists of 2D membranes, solid structures formed by a single atomic layer along the filtration direction. 44Several techniques exist to drill nano-holes on these structures 68 and, nowadays, are possible to control the geometry of the pores. 69The model developed in the present work enables to a predictive analysis of these structures and may constitute a preprocessing tool in their synthesis.To fulfill this objective, the Stokes equations valid in the pores must be replaced by a continuum model which takes into account for the confinement of the flow within the nanopores 70 or by a molecular description of the transport at the nanoscale. 49We finally notice that the model developed here opens to the description of complex phenomena such as concentration polarization 71 or membrane fouling, 41 via the introduction of a timedependent value of a, b, and c which may depend on a supplementary surface concentration field. 72

FIG. 1 .
FIG.1.Overview of the problem analyzed.Panel (a): macroscopic representation of a homogeneous membrane C (corresponding to the coordinates n ¼ 0) of characteristic size L invested by a solvent flow whose characteristic velocity is U, carrying a solute whose concentration is denoted by c, identified by the dark gray cloud crossing the membrane.A magnification of the membrane to visualize its real, microscopic structure on the planes ðŝ; tÞ and ðŝ; nÞ is depicted in panels (b) and (c), respectively.The microscopic elementary cell of typical length l is highlighted by the two dashed rectangles in the same panels.Within this cell, the unknown fields are supposed periodic over S 1 ; S 2 ; T 1 , and T 2 , while they satisfy natural unperturbed flow conditions on the sides U and D, placed in the far-field, infinitely far from the microscopic solid inclusions.

FIG. 2 .
FIG. 2. Overview of the flow configuration used to test Eq.(11).Panel (a): channel flow in a two-dimensional channel vertically split by a microstructured membrane whose geometry is shown in the circular red inset of the same panel.The velocity profile and solute concentration are constant at the inlet L (û 1 ¼ 10 À6 m=s; û3 ¼ 0 m=s and ĉ ¼ 1 mol=m 3 ).The top and bottom walls of the channel, T and B, adsorb the solute with an adsorption coefficient ĵ ¼ 0:01 mol=m 2 , i.e., @ n ĉ ¼ À0:01ĉ, while at the outlet R, a zero-normal solvent stress and solute flux condition has been imposed.Panel (b): solute concentration on C evaluated via Eq.(11) (red solid line) or by solving the fullscale problem (black lines) where the presence of each inclusion is taken into account in the computational domain (one of them, denoted with M, is sketched in real scale within the figure at x 3 ¼ À0:25).Blue stars represent the average (9) of the black profiles.The boundary condition (6) is imposed on @M with a ¼ 0; b ¼ c ¼ 1, and Ĉw ¼ 0:2.The red dashed line is the solution of the homogenization-based model obtained in Sec.IV.Panel (c): same as in panel (b) with a ¼ 1; b ¼ c ¼ 0. The gray square denoted with M represents one solid inclusion forming the membrane.The inset in this panel shows the full-scale profile of concentration within one single pore at the centerline of the membrane.
) and 3(b) show the microscopic solution for the non-zero component of the vector T and Y, denoted with T n and Y n , in the case of a ¼ 1; c ¼ 0 (panel a, insulating/adsorbing membranes) and a ¼ 0; c ¼ 1 (panel b, chemostat-like membranes).As shown in Fig.3, the following relation is satisfied in the considered case: Y n ðs; nÞ ¼ ÀT n ðs; ÀnÞ:

FIG. 3 .
FIG. 3. Overview on the solution of the microscopic problems (38) and (39) for a two-dimensional configuration where the solid inclusion is a square.Panel (a): the only non-zero, normal-to-the-membrane component of the vector T is shown within the microscopic domain for the case of a ¼ 1; c ¼ 0 (left) and a ¼ 0; c ¼ 1 (right).Panel (b): same as in panel (a) for the microscopic vector Y. Panels (c) and (d): average of T n and values of j eff for varying the porosity #.

FIG. 4 .
FIG. 4. Macroscopic configurations used to validate the model.The inset in the bottom right is a magnification over the homogeneous macroscopic membrane C for every configuration, showing the full-scale structure of the membrane.
3 and the separation of scales parameter is equal to 0.1.With these geometrical parameters, T ¼ À Y ¼ ð1:0870; 0; 0Þ for a ¼ 1; c ¼ 0 and T ¼ À Y ¼ ð0:0012; 0; 0Þ for a ¼ 0; c ¼ 1, while the only non-zero components of the tensors M and N are M 111 ¼ À N 111 ¼ À0:0499 and M 331 ¼ À N 331 ¼ 0:0013.The fluid, flowing in a Stokes regime from left to right thanks to a fully developed Poiseuille imposed velocity profile on the side L of the domain, is carrying the solute whose concentration is indicated with c and has a constant value at the inlet of the channel L. The top and bottom walls of the channel, T and B, are adsorbing the solute at a fixed rate j ¼ via the non-conservative flux boundary condition listed in Table II.The Peclet number based on the macroscopic velocity and length scales is Pe ¼ 1.A first qualitative and quantitative insight of the concentration field within the channel is given in Figs.5(a), 5(c), and 5(e), which show the solution for the chemostat-like (a ¼ 0; b ), 7(c), and 7(e) for different membrane behaviors, showing a good agreement.The quantitative comparison between the full-scale solution (black curves) and the macroscopic solution (colored curves) shown in the right column of Fig. 7 for c on C [Fig.7(b)], along x 3 ¼ 0 [Fig.

FIG. 6 .
FIG. 6.Comparison between C 0 and C 0 þ c I;ð0Þ over C for the three membrane behaviors considered in the present work [frame ðaÞ a ¼ 1; b ¼ c ¼ 0, frame ðbÞ a ¼ 0; b ¼ c ¼ 1, and frame ðcÞ a ¼ b ¼ 1; c ¼ 0] and for configuration C1.Blue stars have been calculated applying definition (9) to the full-scale solution.

FIG. 8 .FIG. 7 .
FIG. 8. Comparison between C 0 and C 0 þ c I;ð0Þ over C for the three membrane behaviors considered in the present work [frame ðaÞ a ¼ 1; b ¼ c ¼ 0, frame ðbÞ a ¼ 0; b ¼ c ¼ 1, and frame ðcÞ a ¼ b ¼ 1; c ¼ 0] and for configuration C2.Blue stars have been calculated applying definition (9) to the full-scale solution.

FIG. 9 .FIG. 10 .
FIG. 9. Channel flow configuration analyzed in Sec.VII.Panel (a): sketch of the computational domain.The boundary conditions on LR; T, and B are the same as in configuration C1.A zig-zag shaped membrane is placed longitudinally along the central axis of the channel and is denoted with C. Panel (b): microscopic inclusions forming the membrane for the first four values of R.

FIG. 11 .
FIG. 11.Overview of the macroscopic solution exploited to assess the immediacy of the physical interpretation of the membrane behavior starting from its microscopic properties.Panel (a): isocontours of the solvent pressure (colors) and solvent flow streamlines (black) for R ¼ 1. Panel (b): isocontours of the solute concentration field around the membrane (sketched in red) for R ¼ 1. Panel (c): same as in panel (b) for R ¼ 10.

FIG. 12 .
FIG. 12. Effect of the recursion number on the solute concentration and solvent flow field.Panel (a): concentration sampled over C. Each color denotes a different value of R as indicated in the colorbar.The parameter r represents the curvilinear coordinate along the membrane.Panel (b): increase in the adsorption ability of the membrane by means of the definitions in Eq. (50) for varying R. The yellow pentagram corresponds to the same quantity evaluated from the full-scale solution for R ¼ 1. Panel (c): local relative error f for the solute concentration in the case of R ¼ 1, over the macroscopic membrane C.

TABLE I .
Values of T n andY n for different values of jU À Dj for a square inclusion of side 0:5l and for the possible membrane behaviors considered in the present work.

TABLE II .
Overview of the boundary conditions for c and u used for configurations C1 and C2.