Bridging the gap between lab-scale and commercial dimensions of hollow fiber nanofiltration membranes

The development of nanofiltration membranes with hollow fiber geometry is a relatively young field of research in membrane science with only limited, but quickly growing, availability on a commercial scale. These membranes offer promising properties for application in water treatment such as high fouling resistance and low cost compared to Reverse Osmosis membranes. Major differences are found in dimension and operating conditions of hollow fiber nanofiltration membranes used in academia and on an industrial scale. To allow for adequate comparison and prediction of membrane performance, the effect of fiber dimension and operating conditions on process performance needs to be properly understood. A systematic experimental study on MgSO 4 retention by hollow fiber nanofiltration membranes was performed to investigate the effect of fiber length and operational conditions on membrane performance. A significant drop in the retention of MgSO 4 was observed for 1.5 m long fibers in flux ranges of 20 LMH when reducing crossflow velocity to 0.1 m/s. Transport models were used to describe and predict this behavior. We show that the commonly used mass transport models under predict this decrease in retention. Using different theoretical transport models for the hollow fiber feed phase and reducing dimensional resolution from 2D to 1D and 0D allows for the identification of limitations of the commonly applied analytical mass transfer correlations. To accommodate for this discrepancy found between experimental results and model predictions, a correction factor was derived using the 2D model to specifically account for deviations of these correlations at high permeate recovery values. Our newly developed model makes it possible to accurately predict full scale membrane performance solely based on measured done on small scale membranes, an important precondition for the further development of dense hollow fiber based membranes and processes.


Introduction
Pressure-driven filtration using Nanofiltration (NF) membranes combines unique separation properties with reduced operating pressures compared to conventional Reverse Osmosis (RO) membranes. NF membranes typically exhibit ion selectivity based on valency (high retention towards multivalent ions, low retention towards monovalent ions) and high retention towards small organic molecules. Nowadays, they are applied in numerous fields such as water and wastewater treatment, pharmaceutical processes, food processes, and biotechnology [1,2]. Especially due to an increased concern about organic contaminants accumulating in ground and surface waters [3][4][5], there has been an increased effort in developing NF membranes as a potential solution to retain these contaminants [6][7][8][9][10].
For a long time, the two major module geometries that were commercially available for NF were spiral-wound and tubular modules. The benefits of spiral-wound modules are a high packing density and * of a support [9,14,15]. The benefits of this method are a high chemical and physical stability of the resulting membranes, the simplicity and at the same time versatility of the LbL procedure, as well as the fine control over membrane properties [16][17][18].
There is, however, a major gap between HF NF membranes fabricated and tested on a laboratory scale and the ones produced and operated on an industrial scale. The main differences are fiber dimensions and filtration process conditions. Typical lengths of hollow fibers tested in academia are 10-30 cm [8,[19][20][21][22][23][24][25][26][27][28][29][30], whereas the standard length for commercial modules is 1.5 m. This significant increase in length results in varying process parameters along the fiber length, a higher water recovery and an increased concentration polarization, at similar operating conditions. Furthermore, in academia membrane performance is typically evaluated under one specific operating condition, solely aimed at the characterization of membrane materials. However, the operational conditions tested vary quite a lot, with transmembrane pressures ranging from 1-6 bar and Reynolds numbers ranging from 300 up to 4000 [8,[19][20][21][22][23][24][25][26][27][28][29][30]. This rather significant variation makes direct comparison of newly developed membranes difficult. For industrial applications, process design and operating conditions are driven by economic considerations [31,32], where the main motivation is to minimize costs for given product specifications. For a given process design, costs can be optimized by tuning process conditions [33,34]. Considering the capital costs of a filtration plant, membrane processes typically run at high pressures to achieve high recovery values. In crossflow operation, increasing cross-flow velocity is used to reduce the effect of concentration polarization, at the cost of substantial additional energy required for feed circulation.
These differences in process design and operation conditions of HF NF membrane filtration processes between academia and industry make direct comparison of newly developed HF NF membranes difficult. To allow for comparison as well as prediction or optimization of HF NF membrane performance on an industrial scale, adequate numerical models resolving membrane transport and transport inside the fiber are required. A variety of models for NF membrane transport has been presented in literature [35], which can be roughly divided into mechanistic and irreversible thermodynamic models. While mechanistic models try to relate the observed membrane performance to membrane properties by a specific molecular mechanism, models based on irreversible thermodynamics assume a linear relation between driving forces and fluxes (for multiple components this includes coupling between those). The choice of model depends on membrane properties. Transport processes of the feed solution inside the HF, which account for local parameter variations including the effect of concentration polarization, are defined by mass and momentum balances. The resulting set of equations needs to be solved numerically if full spatial resolution is desired. As this requires quite extensive computational effort, commonly a film model in combination with a mass transport correlation is applied [36][37][38].
In a previous study, Labban et al. proposed a theoretical model to predict the performance of HF NF membranes on a commercial scale based on measurements performed with short membranes [39], due to lack of availability of longer modules. They emphasized, that, in order to facilitate the industrial acceptation of these types of membranes, a broader knowledge of filtration performance at commercially relevant scale and conditions is required. However, in their study the scale-up was solely performed on a theoretical basis with a lack of experimental results to validate the model.
The purpose of this study is to bridge the gap between academic and commercial scale of HF NF membranes by investigating the effect of fiber length on filtration performance at different operational conditions. To this end, a systematic experimental study of MgSO 4 retention was performed and combined with numerical models of varying complexity. The well-known Solution-Diffusion model was applied to describe solute transport through the membrane. The mass and momentum transfer of the fluid phase was resolved at reducing accuracy going from two-dimensional to one-and zero-dimensional, using common mass transfer correlations. The applicability of these models for HF NF membranes was reassessed, with the aim to identify a numerical model capable of predicting HF NF process performance at a wider range of fiber dimensions or operational conditions, at the lowest computational costs. Such a model would allow direct comparison of the performance of newly developed HF NF membranes in academia, and could be readily applied by industry to predict HF NF process performance and scale-up purely based on intrinsic membrane parameters.

Experimental setup
In order to investigate the effect of HF length on module performance and to validate the theoretical process description, a systematic experimental study with varying fiber length and operating conditions was performed. HF NF modules were kindly provided by NX Filtration . Due to potting the effective membrane length is 4 cm smaller than the module length. The modules were operated in cross-flow filtration mode with the feed in the fiber lumen side. A schematic of the experimental setup is depicted in Fig. 1. Inlet cross-flow velocity, feed pressure, and retentate pressure were monitored and regulated via needle valves and pump settings. The feed temperature was regulated to 20 ± 1 • C. Both, retentate and permeate stream were continuously recycled to a 5 L feed tank. The membrane modules were placed in a vertical position with feed flow from top to bottom. Module performance was evaluated via magnesium sulfate retention. Magnesium sulfate was chosen as a model component in this study, due to its frequent use in (especially commercial) NF membrane characterization. Feed solutions of 5 mM MgSO 4 ⋅7H 2 0 (Sigma Aldrich, ≥98% purity) were prepared for all experiments in deionized water. Concentrations were determined via conductivity measurements using a portable conductivity meter (Cond TM 3210, WTW TM GmbH). Measurements were conducted at varying pressures (1-6 bar) and cross-flow velocities (0.1-0.6 ms −1 ). The cross-flow velocities correspond to flow rates ranging from 16 to 96 Lh −1 .

Model development
The complete membrane module was represented by a multiple of single HF membranes in parallel. The theoretical description of the filtration process is composed of (1) the transport through the separation layer of the membrane and (2) the mass and momentum transport inside the bore of the HF. The permeate pressure was assumed p = 1 atm at any position.

Membrane transport
Many different NF transport models were developed in the literature to describe salt retention and water transport, which have recently been reviewed by Yaroshschuk et al. [35]. The most commonly used model is the Donnan steric pore model & dielectric exclusion (DSPM&DE) model [40,41], which describes the transport through straight cylindrical pores using the extended Nernst-Planck equation. Fitting parameters are described by physical membrane properties that are hard to measure precisely. An alternative model is the Solution-Diffusion-Electromigration (SDEM) model [42], which does not account for advective solute transport and uses empirical fitting parameters. With the limited experimental characterization of the membrane structure and polyelectrolyte multilayers being rather dense structures [43], advective transport of ions through the membrane was neglected. Since filtration measurements were conducted with a single salt (MgSO 4 ) the transport of both ions is coupled by electroneutrality and can thus effectively be described by a single compound. Thus, the SDEM model reduces to the well-known Solution-Diffusion model, which is described by the following equations [38]: Here is the water flux (ms −1 ), the pure water permeability (mPa −1 s −1 ), the transmembrane pressure difference (Pa), the osmotic pressure difference (Pa), the salt flux (molm −2 s −1 ), and are the salt concentration at the permeate side and feed side membrane surface (molm −3 ), respectively, and the salt permeability constant (ms −1 ). The real membrane retention , which is the local intrinsic membrane retention not considering feed phase transport phenomena, is defined as:

Mass and momentum transport
In HF NF filtration processes the feed concentration and pressure vary along the fiber length. As water is more readily transported through the membrane, the salt concentration increases at the membrane surface (concentration polarization) and so does the average concentration along the fiber length, while the cross-flow velocity decreases along the length. In particular, at a high water flux and salt retention the process suffers from this rise in concentration [44]. Due to viscous friction, the hydrostatic pressure decreases along the fiber length, and with that the driving force for water transport. Thus, the local membrane performance varies along the fiber length, which significantly influences overall module performance. Following these process dependent effects, the observed membrane retention was defined as: with as the mixing cup concentration of the permeate stream (molm −3 ) and ,0 as the inlet feed concentration (molm −3 ). In the following, this is referred to as retention.
To account for these local variations, the mass and momentum transport in the feed phase has to be resolved. This can be done at various spatial resolutions ranging from three-or two-dimensional (2D) to one-(1D) or non-dimensional (0D), where mass transfer correlations are used. The choice of resolution influences accuracy as well as computational effort. In this study, the mass and momentum transport inside the fiber was resolved with decreasing resolution, to reduce required computational costs at a reasonable accuracy.

CFD simulation (2D)
Assuming axisymmetry of the hollow fiber system, the most detailed resolution of mass and momentum transport considered in this study was the two dimensional plane of radial and axial direction. The fluid flow was assumed to be incompressible, isothermal, and in steadystate. With the salt concentration being in the order of 5 mM, the feed solution was assumed to be dilute. This results in the following system of equations: The system is composed of the Navier-Stokes equation describing the momentum transfer of the feed solution (5), conservation laws of the solution (6) and the dissolved salt (7), and Fick's law of diffusion (8).
The equation system was solved using the software COMSOL Multiphysics ® (Version 5.4). The model design is sketched in Fig. 2. The model domain is composed of an inflow region, the active membrane area, and an outflow region. The following set of boundary conditions is applied to the system: inlet: membrane: outlet: Here 0 is the average inlet velocity (ms −1 ), the surface normal direction, ,0 the inlet feed concentration (molm −3 ), the translational velocity of the wall in radial direction representing the solvent flux (ms −1 ), the salt concentration on the permeate side of the fiber (molm −3 ), and 0 the hydrostatic pressure at the HF outlet (Pa). The values of and are iteratively solved for in the numerical method. The wall in the inlet region is assumed to be impermeable with a no-slip condition applied. The fluid domain was discretized using a mapped mesh with increasing node density close to the fiber wall, where the concentration boundary layer will develop. The mesh of the 1.5 m long hollow fiber was refined until recovery, retention, and pressure drop were independent of element number. The number of mesh elements used was in the order of 60,000. The resulting mesh density was considered the minimum density required for shorter HF membranes. The system of equations was solved for a steady-state solution, with P1-P1 discretization for laminar flow and quadratic Lagrange polynomials for diluted species transport, using a direct solver (MUMPS) with a relative tolerance 0.001. The simulations were performed on a desktop computer with an Intel(R) Core(TM) i7-8750H CPU (6 physical cores, 2.2 GHz) with 16 GB RAM.

Local radial averaging (1D)
When averaging over the HF cross-section, relevant information on radial concentration distribution is lost. The mass transport inside the fiber was approximated by mass transfer correlations in combination with a film model, which is a method commonly used [46]. These correlations are typically derived semi-empirically based on experimental or numerical studies [37], from both first principles [36], or heat mass transfer analogies. The small HF diameter leads to large permeate recoveries even for a single membrane module. Thus, significant boundary layer thicknesses were expected. Therefore, curvature effects need to be accounted for and correlations might be limited in their applicability. Here an extension introduced by Newman to the well known Lévêque correlation is applied [47]: with the Sherwood number ℎ = ∕ , Reynolds number = ∕ and Schmidt number = ∕( ) as dimensionless quantities, the HF diameter (m) and the hollow fiber length (m). The Lévêque correlation was derived for tubular mass transfer systems with impermeable walls. In membrane filtration systems, however, the walls are semi-permeable resulting in a radial velocity component v > 0, which is referred to as suction. Suction effects are known to enhance mass transfer from the membrane surface to the bulk solution [36]. Correction factors for suction effects have been introduced in literature [36,37,48], which have been generalized by Geraldes et al. for NF and RO: Here is the correction factor for suction effects (-) and is defined as the ratio of Peclet number = ∕ (-) to Sherwood number for no suction ( ≡ ∕ ℎ 0 ).
Applying the definition of the mass transfer coefficient in combination with the solute mass balance, the following expression is obtained [37]: Here is the mass transfer coefficient (ms −1 ), the local solute concentration at the membrane surface (molm −3 ), the solute concentration in the bulk (molm −3 ), and the solute concentration on the permeate side of the membrane (molm −3 ).
The mass and momentum equations of the solution inside the HF reduce to the Darcy-Weisbach equation describing pressure loss along the flow direction and a simple one-dimensional mass conservation law: With being the friction factor for laminar flow (-), the inner diameter of the fiber (m), and the average velocity in the flow direction (ms −1 ). These equations can be solved using the method of finite differences. According to the film model, the bulk concentration of the feed phase is assumed to be constant ( = ,0 ).  and 0.6 ms −1 (▴). Error bars display 95% confidence intervals. Data fitted to the 0D retention model with additionally variable mass transfer (R = J /(J +K 1 ⋅exp(J /K 2 ))) to guide the eye.

Dimensionless (0D)
The equation system of the one-dimensional approach still needs to be solved numerically. If one wants to simplify the calculation of the membrane module even more, the mass and momentum balances can be averaged over the whole fiber volume. Here, again the extended Lévêque correlation for mass transfer (Eq. (15)) in combination with suction correction (Eq. (16)) and the film model (Eq. (17)) was applied. Pressure and velocity distributions along the fiber length are not accounted for. Instead, the average pressure can be used as well as the inlet cross-flow velocity. Combining Eqs. (2), (4), and (17), this results in the following analytical expression for the membrane retention that can be quickly solved with accessible software: Here, can either be measured or determined using Eq. (1).

Results and discussion
In this section, initially, an experimental study on the HF NF membrane performance at various process conditions is presented and discussed. Following this is a comparison and discussion of the theoretical models used to predict full-scale module performance. Finally, the applicability of the different theoretical models is evaluated via the experimental results.

Experiments
In Fig. 3, the MgSO 4 retention is given as a function of the flux for HF NF module lengths ranging from 0.3-1.5 m and cross-flow velocities of 0.1 ms −1 and 0.6 ms −1 . It can be observed, that both module length, as well as operating conditions, affect the observed retention. With increasing permeate flux the retention generally increases. At a high cross-flow velocity the retention is higher with seemingly little effect of fiber length, whereas the retention at low cross-flow velocity is lower with a further reduction in retention with fiber length. Also, a decrease in retention for increasing flux above 15 LMH is very evident for the 1.5 m fiber. To analyze the effects of these parameters on retention, it is useful to discuss the latter in terms of real membrane retention and concentration polarization. The influence of permeate flux on real membrane retention for a dense membrane dominated by diffusive transport can be predicted by the Solution-Diffusion model. Transforming Eqs. (2) and (3) one can obtain the following relation of real retention and permeate flux: Assuming a constant salt permeability , this results in a monotonously increasing real retention with higher permeate flux, which is observable especially for high cross-flow velocities with minimal effect of concentration polarization. The pure water flux is proportional to the applied effective pressure, while the salt flux depends on the concentration differences across the membrane. This relation only holds for an ideal system with no concentration polarization, where the membrane surface concentration equals the concentration in the bulk . All three process parameters (flux, module length, and cross-flow velocity) influence the extent of concentration polarization phenomena and with that the retention. With higher permeate fluxes, the concentration at the membrane surface increases as predicted by the derived mass transfer relation (Eq. (17)). At the same time, higher fluxes increase the removal of the solvent in the feed phase leading to a decrease of local cross-flow velocity which impacts the hydrodynamics of the system. This phenomenon is expected to be only relevant for membrane modules operated at high recoveries. The hydrodynamics of the feed flow for a fixed HF diameter is determined by the inlet velocity. In general, higher cross-flow velocities decrease the concentration boundary layer thickness and with that the effect of concentration polarization. A clear decrease in retention is only observed for the 1.5 m long fiber when lowering the cross-flow velocity from 0.6 to 0.1 ms −1 . This can be explained by the high recovery values (up to 47%) obtained at these conditions. Lastly, the HF length affects the module retention. Longer fibers allow for the development of the concentration boundary layer as well as higher recovery values. Thus, retention decreases with fiber length. These experimental results clearly show the effect of operational conditions but also the effect of HF length on the overall process retention that is observed. These results also stress the difficulty in comparison of literature data that has been acquired at various process conditions and module dimensions and the necessity of theoretical models to properly account for these variations.
Following these initial measurements is a more detailed study on MgSO 4 retention of HF NF membranes with commercially relevant dimensions (1.5 m module length) at varying inlet cross-flow velocities (0.1-0.6 ms −1 ) and permeate fluxes (applied TMP 1-6 bar), as shown in Fig. 4. Again a clear trend of decreasing retention at lower cross-flow velocities (below 0.2 ms −1 ) is observed, which becomes more significant at high fluxes. For industrial membrane process design this effect can be very relevant as it limits the range of process parameters for a required retention. Considering for example a required retention of 85%, a cross-flow velocity of 0.2 m/s limits permeate fluxes to 10-17.5 LMH, whereas a cross-flow velocity of 0.6 m/s allows permeate fluxes from 7 LMH to way higher values and a cross-flow velocity of 0.1 m/s does not reach the required product quality at all. The large variation of HF NF retention depending on process conditions especially on large scale emphasizes the importance of adequate theoretical descriptions of the transport processes to allow for process optimization.
In order to bridge the gap between lab-scale and industrial dimensions, a similar series of measurements was conducted for the 0.3 m long HF NF module. This dimension is more representative for the modules length typically used in academia [8,[19][20][21][22][23][24][25][26][27][28][29][30]. In the following the focus will be on comparing the extreme cases of 0.1 ms −1 and 0.6 ms −1 for the 0.3 m and 1.5 m long modules. This covers the full range of inlet velocity, length and recovery considered in this study. A direct comparison of these experimental results is displayed in Fig. 5. For the 0.3 m long module with an inlet cross-flow velocity of 0.6 ms −1 , the effect of concentration polarization is, as expected, the smallest, which results in high MgSO 4 retention and the typical Solution-Diffusion shaped real retention curve with a monotonous increase with flux. The other extreme case, with considerable concentration polarization effects, is the 1.5 m long module with an inlet cross-flow velocity of 0.1 ms −1 . Here, the retention as a function of flux is the lowest compared to the other cases with decreasing retention at fluxes higher than 10 LMH. Both, length as well as inlet cross-flow velocity influence the extend of concentration polarization and, thus, the retention at a constant permeate flux. Not accounting for the removal of solvent in the feed phase, these two parameters essentially affect the residence time of an infinitesimally small volume traveling through the HF. To describe the effect of both parameters in one parameter, the dimensionless Graetz number , which is commonly used in mass transfer correlations, is introduced: For a constant HF diameter and solute diffusion coefficient, this number represents the ratio of inlet cross-flow velocity and HF length, which is the inverse average residence time. Hence, low Graetz numbers correspond to a long residence time resulting in more severe concentration polarization. The Graetz number includes both effects, regarding fiber length and crossflow velocity, and thus influences the MgSO 4 retention at constant transmembrane pressure (TMP). This is displayed in Fig. 6. At a TMP of 1 bar the influence of Graetz number on retention is negligible. This is due to the fact that concentration polarization additionally depends on permeate flux (see Eq. (17)) as well as permeate recovery, which are both relatively low at this TMP. At a higher TMP value of 6 bar the retention is generally higher as the real membrane retention increases with permeate flux. At the same time, retention is dependent on Graetz number, with lower Graetz  numbers leading to decreasing retention. This is again due to concentration polarization caused by high permeate flux and a high permeate recovery. The experimental study conducted in this work clearly displays the major impact of HF NF module length as well as operating conditions applied during filtration on the MgSO 4 retention. The effect of concentration polarization is more severe at low Graetz numbers and high permeate recoveries, yielding lower retention. These regions are, however, very appealing to industrial processes, as one wants to maximize the permeate recovery for high yields as well as minimize the energy consumption meaning minimizing the Graetz number. Therefore, it is crucial in process design and optimization to account for the effect of concentration polarization together with the module length, when retention data from lab-scale modules is used.

Comparison
When solving the 2D mass transport equation system inside the HF using Comsol, full information on local MgSO 4 concentration is obtained. This gives insight into the development of the concentration boundary layer along the fiber length, corresponding to the effect of concentration polarization. In Fig. 7, the concentration distributions for the minimum and maximum Graetz number for a constant TMP of 6 bar are displayed. The model requires two fitting parameters, the pure water permeability and the salt permeability . These were determined with pure water retention measurements as well as salt retention measurements of the 0.3 m module (for more information see Section 4.2.2). The impact of Graetz number on the observed retention at this pressure range is rather significant and was explained by the extend of concentration polarization. This is confirmed by the concentration distribution plots of the HF modules with 0.3 m and 1.5 m module length, and an inlet velocity of 0.6 ms −1 and 0.1 ms −1 , respectively. At a high Graetz number (1730) the boundary layer thickness stays small along the fiber length. The effect of concentration polarization is rather little reaching a maximum concentration at the membrane surface of 7.1 mM (feed concentration 5 mM). At a low Graetz number (51) the boundary layer thickness increases significantly and even reaches the center of the fiber after 1.2 m. The maximum concentration at the membrane surface is 18.5 mM. Thus, the effect of concentration polarization is high at these operating conditions, which explains the reduced retention observed in the experiments (see Fig. 5). This also confirms the relevance of accounting for curvature effects in the film model and mass transfer correlation.
The MgSO 4 retention for the 0.3 m and 1.5 m long modules and the considered flux and inlet velocity range obtained from the COM-SOL model is displayed in Fig. 8. For the short fibers, the curve of retention against flux resembles the Solution-Diffusion dependency (see Eq. (22)), slowly converging towards it with increasing velocity. The effect of concentration polarization on reducing the overall retention is rather small in the considered flux range (93% to 91.5% at 20 LMH). However, it is already noticeable that the extend of concentration polarization depends on cross-flow velocity as well as the permeate flux. Increasing the module length to 1.5 m reduces the observed retentions for all inlet velocities and fluxes. This is caused by the greater extend of concentration polarization due to lower Graetz numbers and higher permeate recoveries. At an inlet cross-flow velocity of 0.1 ms −1 the retention even begins to decrease with fluxes larger than 16 LMH. For both module lengths, the simulation results qualitatively resemble the observed experimentally MgSO 4 retentions. In the following graphs that compare the different models (Figs. 9 to 11), the 2D simulation is considered the most accurate of the models and displayed for reference by markers.
When reducing the dimensions resolved by the theoretical model through averaging over the HF cross-section, the required computational effort is largely reduced. The resulting system of equations can be solved via a finite differences approach in Matlab. Mass transfer limitations in the boundary layer were accounted for by mass transfer correlations in combination with the film model. The obtained retention as a function of average flux for the 1D Matlab simulation is compared to the 2D simulation results in Fig. 9. Here three different mass transfer descriptions are distinguished: no correction for curvature or suction effects (dashed line), correction for curvature and suction (solid line), as well as additional correction for permeate recovery (red line). The mass transfer description not accounting for curvature or suction effects is composed of the well known Lévêque mass transfer correlation in combination with the film model derived for flat sheet membranes [46,49]: Curvature and suction effects are accounted for as presented previously (see 3.2.2). The correction for permeate recovery was obtained using the 0D simulation results and will be discussed later (see Eq. (27)). For the 0.3 m long module, the MgSO 4 retention values obtained by the 1D model are overlapping with those obtained by the  2D simulation (markers). This is due to the fact, that the permeate recovery values as well as the effect of concentration polarization are small. Therefore, the observed retention strives towards the real retention behavior of the NF membrane, which is essentially independent of the hydrodynamic conditions of the system. When curvature and suction effects are not accounted for, slight deviations of the model can be observed at higher flux values. This is expected to be caused by neglecting suction effects, which increase mass transfer. For the 1.5 m long HF module, the Graetz number is lower and with that the permeate recovery higher compared to the 0.3 m module. This leads to stronger effects of concentration polarization, especially at low crossflow velocities. At low inlet velocities and high fluxes, the retention predicted by the 1D simulation starts to exceed those of the 2D simulations. As observed in the concentration distribution resolved in the 2D model (see Fig. 7), the boundary layer fully develops over the fiber length. The Lévêque mass transfer correlation, as well as its extension, describes the developing range of the boundary layer and is expected to be inaccurate for developed boundary layers [47]. The boundary layer development accelerates with increasing permeate recovery. To account for the effect of recovery, a correction factor was introduced, which leads to overlapping retentions for both models. At these conditions, the observed retention when neglecting curvature and suction are very similar to those obtained when accounting for these effects. This seems rather coincidental and is expected to be caused by opposing effects of suction (increases retention) as well as curvature (decreases retention) at these significant boundary layer thicknesses.
These results show that at high Graetz numbers (>300) as well as low permeate recovery values (<11%), even standard mass transfer descriptions can be quite accurate in describing concentration polarization. In HF NF application, however, the goal is to achieve high recoveries with minimum energy costs, which corresponds to high permeate flux and low crossflow velocity. With the newly developed membranes of high permeability towards water, recoveries of 40% or higher can easily be achieved at low feed flow rates even with one module. Thus, the region where common mass transfer correlations start to deviate might be very relevant for applications and requires correction for suction effects, curvature as well as the influence of permeate recovery.
Modeling can be further simplified by averaging the 1D model over the HF length. Here Eq. (21) is applied in combination with a mass transfer correlation. Again three different mass transfer descriptions are distinguished: no correction for curvature or suction effects (dashed line), correction for curvature and suction (solid line), as well as additional correction for permeate recovery (red line). In the 0D model the mass transfer correlation used for the case without accounting for curvature or suction effects is the average form of the Lévêque mass transfer correlation [49]: Again Eq. (25) is used for the film model. The model results are compared to the 2D model in Fig. 10. The simulated retention as a function of average flux for the 0D model matches very well with those obtained in the 1D case (Fig. 9). This is expected since in both cases when accounting for suction effects the average form of the extended Lévêque correlation (Eq. (15)) is applied. Only the mass transfer descriptions without any correction terms differ in both cases leading to deviating predictions of retention especially for the 1.5 m long module at low inlet velocities and high fluxes. Under these conditions, a quite significant difference between the 0D model and the 2D model reference is observed (88% vs. 85.6% at 20 LMH). This is caused by the quickly growing concentration boundary layer. Thus, curvature effects can no longer be neglected. The difference between the local (1D) and average (0D) Sherwood correlation indicates that the significant boundary layer development at low Graetz numbers can no longer be accurately described by either of them. To account for the effect of permeate recovery on concentration boundary layer growth, a correction factor was determined. The factor was obtained by relating mass transfer coefficients obtained from the extended Lévêque correlation (Eq. (15)), while accounting for suction effects (Eq. (16)), with 2D simulation results. Here the effect of the following parameters was combined in a dimensionless number: inlet velocity (0.1-0.6 ms −1 ), permeate flux (1-23 LMH), fiber diameter (0.75-8 mm), diffusion coefficient (0.55-1.15⋅10 −9 m 2 s −1 ) and fiber length (0.26-1.46 m). Further information on fitting procedure can be found in the supporting information. The obtained correction factor is: Here is the dimensionless correction factor and the average water flux in ms −1 . The mass transfer coefficient is determined in the following way: * = 0 ⋅ ⋅ (28) where * is the corrected mass transfer coefficient, 0 is the mass transfer coefficient obtained from the extended Lévêque correlation, is the suction correction factor and is the recovery correction factor.
To compare the different model approaches more clearly, the observed MgSO 4 retention is displayed as a function of Graetz number at a constant flux (see Fig. 11). The retention for the 2D Comsol simulation qualitatively behaves similar to the experimental results (see Fig. 6). For the considered Graetz number range and the two chosen flux values, the retention increases with flux. This is caused by the effect of flux on real retention as pointed out before. However, it is noted that at low Graetz numbers and high enough fluxes this will no longer be the case due to concentration polarization. Hence, under these conditions, this relation will reverse. Also, a clear trend in decreasing retention with decreasing Graetz number is observed which is caused by the development of the concentration boundary layer which is accelerated by the permeate recovery. Thus, the severity of concentration polarization increases with lower Graetz number, decreasing the observed retention. Comparing the three approaches of mass transfer description (no correction for curvature or suction effects (dashed line), correction for curvature and suction (solid line), as well as additional correction for permeate recovery (red line)) for each dimensional resolution and permeate flux, it is observed that the simulation results converge to the same retention with increasing Graetz numbers. This is caused by a decrease in boundary layer thickness with increasing Graetz number, which is necessary for the validity of the considered mass transfer correlations. At lower permeate flux, the retention results converge already at lower Graetz numbers. This displays the effect of permeate recovery on concentration boundary layer development and validity of mass transfer correlations. At lower Graetz number the simulation results start to deviate from the reference retention, especially at high permeate fluxes. Only when accounting for the effect of permeate recovery (red line), the MgSO 4 retention predicted by the 1D and 0D model matches the one obtained in the 2D model.
It is shown, that common mass transfer correlations used to describe concentration polarization phenomena in membrane filtration processes are limited in their applicability for low Graetz numbers (<100) at high permeate recovery (> 22%). This is caused by the rapid development of the concentration boundary layer along the membrane length. In the application of HF NF membranes it is likely though to encounter these conditions. Therefore, it is essential to correct the mass transfer coefficient for permeate recovery, to accurately predict the effect of concentration polarization and with that the observed retention. Both simplified models predict the same MgSO 4 retentions as a function of average flux. The 0D model, however, is not able to predict pressure distributions or concentration variations along the fiber length. If the performance of a process is to be predicted, these distributions can be very relevant for specific energy calculation as well as estimation of scaling effects (deposition of salt on the membrane surface), respectively. The 1D model predicts very similar pressure distributions as the 2D model. Distributions of local membrane surface concentration are also not resolved in the 1D model when using the average mass transfer correlation. If one wants to obtain salt concentrations at the membrane surface, e.g. to estimate scaling effects, the 2D model is recommended.

Fit
To describe and predict the experimental results obtained for the HF NF modules of different length measured at various operation conditions, the 2D Comsol model was fitted to the MgSO 4 retention results obtained for the 0.3 m long module. The Solution-Diffusion model requires two fitting parameters to be determined which are: the pure water permeability A and the salt permeability B. These can be determined independently from pure water filtration measurements and salt retention measurements. The pure water permeability was determined at a transmembrane pressure range from 1 to 6 bar. The pure water permeability was 3.52 ± 0.18 LMHbar and 3.4 ± 0.16 LMHbar for the 0.3 m and 1.5 m long module, respectively. Thus, A is set to 3.46 LMHbar. The 2D model was fitted to the experimental results of the 0.3 m long module via the least-squares method resulting in a salt permeability B of 3.0⋅10 −7 ms −1 . Here, only measurement points with fluxes higher than 5 LMH were considered due to the high variation of the salt retention at low fluxes. The fitting result, as well as the model prediction for the 1.5 m long module are displayed in Fig. 12. For the 0.3 m long module at a crossflow velocity of 0.6 ms −1 the predicted retention is within the 95% confidence interval of the experimental values. It is, however, already observable by means of the average retention values that the model slightly overestimates the retention at high fluxes and underestimates at low fluxes. This could of course be caused by experimental error. However, with reducing the crossflow velocity or increasing the module length (lower Graetz numbers), this inaccuracy of the model increases. For the 1.5 m long module at a crossflow velocity of 0.1 ms −1 , which is the lowest Graetz number considered in this study, the retention is largely overestimated at high flux values. This trend being present in all cases, although increasing with higher Graetz numbers, implicates inaccuracies in the membrane transport description. The Solution-Diffusion transport description applied in this study is a rather simplified approach only accounting for diffusive transport with an intrinsic transport coefficient B that is independent of salt concentration.
Here we present two possible transport mechanisms, that explain the observed MgSO 4 retention behavior. One well-known mechanism that is important in charged based solute retention is Donnan exclusion.   17)) and correction factor (Eq. (27)). 0D transport model solved in Excel.    Donnan exclusion relates the concentration of mobile ions inside a charged membrane to the ion concentration at the membrane surface accounting for charged functional groups. For a constant salt concentration, this additional exclusion factor is incorporated in the B factor of the Solution-Diffusion model. However, it is well known that Donnan exclusion is susceptible to charge screening effects, and thus depends on local salt concentration. This results in a concentrationdependent intrinsic salt permeability value B. With higher local salt concentration, the Donnan exclusion is reduced and with that more salts can pass through. As a result, B increases and the retention decreases. Another transport mechanism that has been neglected in this study due to the dense nature of the polyelectrolyte multilayer film is advective solute transport. In a study on the transport mechanism of ions through polyelectrolyte multilayer membranes Cheng et al. already found indications on advective ion transport, most likely being caused by imperfections in the selective layer [43]. The advective transport of ions would also lead to ion transport that is dependent on local ion concentration at the membrane surface. Assuming constant flux, the advective transport of ions increases with increasing local salt concentration, which leads to a decrease in observed retention. Both mechanisms have been included separately in the 2D model to illustrate the effect on observed retention. The predicted MgSO 4 retention for the 1.5 m long module is displayed in Fig. 13. To evaluate these mechanisms additional experiments were conducted at higher MgSO 4 concentration (15 mM and 25 mM) for the 0.3 m long module at high cross-flow velocities to minimize mass transfer limitations in the fluid phase. The measured retention is given in the supporting information. No effect of concentration was observed as the retention results for both concentrations were indiscernible. This highly suggests advective transport to be the additional transport mechanism through the membrane, as this only lowers the observed retention when local concentration increases at the membrane surface.

Conclusion and outlook
In the present work, the influence of fiber length and process operating conditions on MgSO 4 retention of hollow fiber nanofiltration membranes was systematically studied experimentally and theoretically, for the first time. It was found, that permeate flux can have a significant effect on intrinsic membrane retention (although expected for dense membranes, this is often neglected in literature), which emphasizes the necessity of adequate membrane transport models to ensure comparability of salt retention measurements conducted on a lab-scale at different transmembrane pressures. It was also found, that for commercially relevant hollow fiber lengths as well as permeate recovery values, there is a significant effect of concentration polarization lowering the observed MgSO 4 retention. In order to predict and optimize membrane process performance, a theoretical description of mass transfer processes in the feed phase is required. A comparison of mass transfer models with a different spatial resolution of the feed phase reveals, that common mass transfer correlations, used to quantify concentration polarization effects in 1D and 0D models, are not accurate anymore at low Gaetz numbers and high permeate recovery values, where they start to overpredict observed retention. Under these conditions, there is a significant increase in the thickness of the concentration boundary layer, as shown in 2D simulations, which goes against assumptions taken in deriving these correlations. Thus, at commercially relevant conditions applied for hollow fiber nanofiltration, these correlations might not be appropriate anymore. For the considered range of hydrodynamic conditions and recovery, a recovery correction factor was derived from the 2D model. This factor allows for accurate prediction of average mass transfer behavior in the feed of the hollow fiber and thus properly predict operational performance in industrial processes at significantly lower computational costs. However, the distribution of solute concentrations at the membrane surface is not resolved. If this is required to estimate scaling or fouling effects, the use of a 2D model to resolve mass and momentum transfer is recommended. Fitting the 2D model to experimentally observed salt retention reveals that a simple Solution-Diffusion based approach is not sufficient to describe the transport through these Polyelectrolyte Multilayer based nanofiltration membranes. It was found that the incorporation of advective transport, which could be caused by a porous network or even defects in the membranes, significantly increases the accuracy of the prediction. An investigation of solute transport mechanisms through Polyelectrolyte Multilayer based nanofiltration membranes will be part of future research.

Declaration of competing interest
Dr. J. de Grooth currently holds a part-time position at NX Filtration B.V., a membrane manufacturer. The other authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.