Predictive pressure drop models for membrane channels with non-woven and woven spacers

Desalination 376 (2015) 41–54

Contents lists available at ScienceDirect

Desalination journal homepage: www.elsevier.com/locate/desal

Predictive pressure drop models for membrane channels with non-woven and woven spacers Matthias Johannink a, Kannan Masilamani b, Adel Mhamdi a, Sabine Roller b, Wolfgang Marquardt a,⁎ a b

AVT — Process Systems Engineering, RWTH Aachen University, Turmstr. 46, 52064 Aachen, Germany Simulation Techniques and Scientific Computing, University of Siegen, Hoelderlinstrasse 3, 57068 Siegen, Germany

H I G H L I G H T S

G R A P H I C A L

A B S T R A C T

• Predictive pressure drop models have been identified for woven and nonwoven spacers. • Quantitative description of the pressure drop for a wide range of design parameters • Work flow is presented for the identification of models with a minimal effort. • A new Lattice–Boltzmann-based tool is used for CFD simulations on HPC clusters.

a r t i c l e

i n f o

Article history: Received 7 April 2015 Received in revised form 27 July 2015 Accepted 28 July 2015 Available online 29 August 2015 Keywords: Membrane processes Spacer Pressure drop Modeling

a b s t r a c t Predictive pressure drop models have been identified for membrane channels with both woven and non-woven spacers as found in plate and frame or spiral wound membrane modules. The models allow for the first time the thorough quantitative description of the pressure drop for a wide range of design parameters specifying the thickness and orientation of the spacer filaments. A systematic work flow has been developed efficiently integrating detailed computational fluid dynamics (CFD) simulations on high-performance computing (HPC) systems with established optimization-based model identification methods. The work-flow efficiently handles the complexity arising from the large number of design parameters and allows for the model identification with minimal number of CFD simulations. The resulting models facilitate for the first time a thorough sensitivity analysis of the pressure drop with respect to all important design parameters elucidating strongly nonlinear patterns in the sensitivities. This indicates that the generalization of trends from a small number of CFD simulations, frequently used in previous work to simplify the spacer design problem, can yield inaccurate results. © 2015 Elsevier B.V. All rights reserved.

⁎ Corresponding author. E-mail address: [email protected] (W. Marquardt).

http://dx.doi.org/10.1016/j.desal.2015.07.024 0011-9164/© 2015 Elsevier B.V. All rights reserved.

42

M. Johannink et al. / Desalination 376 (2015) 41–54

Greek symbols α flow attack angle [°] β angle between filaments [°] confidence intervals [–] Δθk kg ρ density [m 3] σ standard deviation of measurements [–] ζ MUSUBI friction factor predicted by MUSUBI [–] ζ sp friction factor of the spacer [–] Roman symbols ϕ least squares objective [–] AICi Akaike criterion [–] ΔAICi,c relative Akaike criterion w.r.t. model c [–] cf correction factor [–] df filament diameter [m] dg characteristic diameter [m] F Fischer information matrix [–] H Hessian matrix [–] hch channel heights [m] J Jacobian matrix [–] lch channel length [m] lm orthogonal distance of filaments [m] n sp, m sp empirical model parameters [–] Δp pressure drop [Pa] p total pressure [Pa] ^k p empirical model parameters [–] ptp,i averaged total pressure at tracking plane i [Pa] Re Reynolds number [–] RSS residual sum of squares [–] RSSCV averaged root mean square error [–] uk geometric design parameters [–] uˇ k scaled model inputs v velocity vector ½ms vch mean velocity [ms] wch channel widths [m] wi,c Akaike weights w.r.t. model c [–]

1. Introduction Modules used for membrane processes, such as reverse osmosis, nanofiltration or electrodialysis, are typically realized in spiral wound or plate and frame configuration. A common characteristic of both configurations is the use of structural components called spacers [1]. The main functionality of the spacers is the separation of adjacent membranes. A second functionality of the spacer is to influence the hydrodynamics in the flow channels in a favorable manner. The local hydrodynamic conditions in the spacer-filled channels determine important design objectives governing the performance of a membrane process: (i) the pressure drop over the membrane module, (ii) the wall shear-stress rate at the membrane surface and (iii) lateral mixing resulting from the formation of steady or turbulent vortices. A high pressure drop results in a loss of driving force in nanofiltration or reverse osmosis [2]. In electrodialysis the pressure drop is directly proportional to the energy demand of the pumps [1]. Zones characterized by a low wall shear-stress rate are likely to act as initiation points for the formation of fouling films [3]. Lateral mixing increases the local mass transfer intensity by reducing concentration polarization effects [4,5]. Zamani et al. [6] emphasize the increasing importance of these design objectives concomitant with the advances in the development of a new generation of ultra-permeable membranes for reversed osmosis processes.

The hydrodynamics in spacer-filled membrane channels have been investigated in numerous studies using either experimental (e.g.,[7–9]) or computational fluid dynamics (CFD) techniques (e.g.,[10–16,3]). An excellent review is provided in [3]. The focus of the majority of CFD studies [4,5,14,17,16] has been directed to a qualitative description of the flow fields arising in different Reynolds regimes. In more detailed studies, the simulation of the hydrodynamics in the channel is extended to include additional physical phenomena. The simulation of simultaneous diffusive and convective mass transport in the channel, for example is included in a few simulation studies reported in [18,4,3]. Similarly, Picioreanu et al. [19] and Radu et al. [20] employ coupled simulations of hydrodynamics and biofilm growth in the channel to study the formation of fouling. Simulations of coupled phenomena require a significantly increased computational effort and thus limit the number of simulations in a modelbased analysis. Improved insight into the transport mechanisms have led to the development of numerous innovative spacer structures (e.g.,[3,21, 14,22]). Some of these innovative spacer structures have shown an improved performance with respect to different design objectives in model-based and experimental studies. However, despite this improvement, commercially available spacer structures are still based on the established non-woven or woven spacer geometries. One reason for the apparently difficult translation of innovative structures into industrial processes relates to the higher effort in the spacer manufacturing. The quantitative effect of important design parameters on the design objectives, i.e., pressure drop, shear-stress pattern and mass-transport efficiency, has been investigated in [23,11,15,24]. All studies emphasize the strong sensitivity of the design objectives with respect to the design parameters. Most of them vary one parameter and fix all the others to nominal values to analyze the dependency of design objectives on design parameters. However, the extrapolation of the identified trends beyond the nominal values gives only adequate results, if the variation in the sensitivities is small. The evaluation of experimentally measured or simulated pressure drop and mass-transfer rates is commonly assisted by simple algebraic models [25,4,18]. These models correspond to power laws relating well-known dimensionless quantities such as friction factors, Reynolds, Sherwood and Schmidt numbers by empirical parameters. The comparison of the empirical parameters in these power laws for different spacer geometries allows for a reasonable evaluation of the geometries in the entire range of Reynolds or Schmidt numbers. However, these models are restricted to a fixed spacer geometry. A general model for the prediction of the pressure drop or mass-transfer rate for a wide range of design parameters is still missing [3]. This contribution addresses this gap and aims at the identification of predictive models for the design objectives in the entire space spanned by the design parameters. To handle the complexity arising from this multi-dimensional design space, a work-flow is proposed guiding the execution of specific CFD simulations by means of optimization-based model identification techniques. With this, the approach yields predictive models by a minimal number of CFD simulations. In a first step, we restrict our consideration to the pressure drop constituting one central objective for the spacer design. As the primary result, the contribution presents a generalization of the established pressure drop models by accounting for varying design parameters. For the first time, this allows a thorough quantitative description of the pressure drop in the entire design space. The remainder of this publication is organized as follows: Section 2 introduces a systematic work flow and the model-based and experimental methods employed. In Section 3 CFD simulation results are discussed and compared to experimental data. Section 4 presents the results of the identification of the generalized pressure drop models.

M. Johannink et al. / Desalination 376 (2015) 41–54

2. Methods The recommended systematic work-flow for model identification is introduced in Section 1. Technical details concerning the CFD simulations and the geometric modeling of the spacer are given in Section 2. The experimental set-up for the pressure drop experiments is presented in Section 3. Section 4 introduces the algebraic pressure drop models. The methods and concepts related to parameter estimation, parameter identifiability and model-based design of experiments are introduced in Section 5. 2.1. Work-flow for model identification A major difficulty in the generalization of the algebraic pressure drop models is the identification of an appropriate extended model structure incorporating the design parameters as additional inputs. To this end, we suggest the work-flow for structureidentification as depicted in Fig. 1a. Simulated data for a basic set of geometries are used to evaluate structurally different model candidates by means of their predictive capabilities. After the identification of the best candidate structure, the model is iteratively refined as illustrated in Fig. 1b. Pressure drop data for additional spacer geometries are added in this iterative refinement process to increase the information content of the pressure drop data set. To maximize the information content added by each new geometry, we use a model-based design of experiments (DOE) approach. The refinement process is repeated until the augmentation of additional spacer geometries does not improve the reliability of the empirical model parameters.

43

given geometry in a reasonable time and allows us to study a large number of geometries in a systematic approach.

2.2.1. Geometric modeling of the spacer The geometric properties of a non-woven spacer are illustrated in Fig. 2. The two layers of parallel filaments that form the spacer are oriented in the x–z-plane. Within this plane, two angles characterize the orientation of the filaments: (i) the angle between the parallel filaments in the different layers β and (ii) the angle δ describing the rotation of the symmetry axes relative to the x-coordinate which points into the main flow direction. Technical applications require the symmetry axis of the spacer to be identical to the mean flow direction [11]. Hence, commonly the flow attack angle is fixed to α = 90 ° and β is varied in the range 0 ° b β b 180 °. With these angles, the diameter of the filaments df and the orthogonal distance between two parallel spacer filaments lm, the geometry of the spacer is completely defined. Hence, a three-dimensional design space isspanned by the design parameters unw = [β, df, lm]. Likewise, the same geometric properties α, β, lm and df uniquely specify woven spacer designs. However, it is difficult to manufacture a woven mesh with β ≠ 90 °. Hence, the technically relevant design space is restricted to uw = [df, lm]. For both, woven and non-woven spacer structures the thickness of the channel hch is given by hch ¼ 2d f c f :

ð1Þ

2.2. CFD simulations

The factor cf = 0.946 (cf. [3]) is introduced to account for the fact that the filaments themselves and the bounding membranes are slightly overlapping.

Model identification requires pressure drop data for a large number of different spacer designs. The generation of this data set is enabled by using an efficient CFD solver based on the Lattice Boltzmann method (LBM). This method constitutes an efficient numerical scheme for the solution of the Lattice Boltzmann equation which is a discrete representation of the rigorous Boltzmann kinetic equation. Without considering technical details, it is important to note, that in the macroscopic limit, the solutions of the well-known Navier–Stokes and Lattice Boltzmann equation coincide [26,27]. The scalability of the LBM allows for the efficient simulation of the hydrodynamic problem on computer clusters with a large number of processors. This way, it is possible to simulate the hydrodynamics in a

2.2.2. Simulation scenario The simulated domain of the spacer-filled flow channel is depicted in Fig. 3. The length and width of the domain is chosen to resolve six filament intersections, or unit cells [11], along the x-coordinate and two filament intersections along the z-coordinate. The height of the domain corresponds to the height of the channel hch. Periodic boundary conditions are specified for the boundary elements in the x–y-planes at z = 0 and z = wch. For the boundary elements in the x–z-planes at y = 0 and y = hch — corresponding to the location of the membranes — no-slip boundary conditions are enforced. In the inlet and outlet planes at x = 0 and x = lch, respectively, standard in- and outflow conditions are

Fig. 1. Work flow for (a) the discrimination of different model candidates and (b) the iterative model refinement of the best candidate based on DOE and parameter estimation techniques.

44

M. Johannink et al. / Desalination 376 (2015) 41–54

Fig. 2. Geometry of a non-woven spacer.

specified. Initial conditions are chosen as v(x, t = 0) = 0 and p(x, t = 0) = patm, where patm is the atmospheric pressure. From this initial state, each dynamic simulation is run until it reaches steady state. The total pressure ptp,i, i = 0, …, 5, is tracked in the course of the simulation in six tracking-planes at each filament intersection (cf. Fig. 3). In these planes, the pressure is averaged to obtain an estimate for the evolution of the pressure in the mean flow direction. The physical properties of the perfusing fluid are specified as those of demineralized water at 298.15 K [28]. 2.2.3. Software and implementation The implementation of the LBM used in this publication is the open source solver MUSUBI [29] which is part of the software package APES [30]. For the 3D simulations we use a D3Q19 scheme based on 19 discrete lattice velocities and the multiple relaxation (MRT) approach with acoustic scaling. Mesh generation is carried out using SEEDER [31] to allow for an efficient mesh generation for arbitrary woven or non-woven spacer geometries from elementarygeometrical structures. The elementary structures of woven spacers are generated with the WeaveGeo2012 module in the software GeoDict [32,33]. The domain is discretized with regular cubic elements with a size of 6.25 ⋅ 10−6 m. All these simulations were performed on the Linux

cluster of RWTH Aachen, Germany, and the high performance computing system HERMIT at HLRS Stuttgart, Germany, using 512 to 2048 processors. 2.3. Pressure drop experiments Experiments are carried out using a lab-scale electrodialysis plant. For pressure drop measurements over a single spacer-filled channel, the membrane module is substituted by a tailored test cell. A simplified flowsheet of the electrodialysis plant is shown in Fig. 4. Demineralized water is pumped by centrifugal pumps from the storage tank B1 through the membrane module and test cell, respectively. The total volumetric flow rate Q FIRC1002 is measured by the turbine flow-meter FIRC1002. This signal is supplied to a software-based PDcontroller which manipulates the volumetric flow rate QFIRC1002 in a bypass configuration by acting on the control valve V1101. The construction of the test cell allows us to vary the length lsp of the spacer in the channel to facilitate the discrimination of the contributions of the spacer and the distribution system of the test cell to the total pressure drop. Experiments are conducted for the spacer ED-200-x from PCCELL GmbH, Germany. The corresponding geometric properties of this woven spacer structure are summarized in Table 2.1.

Fig. 3. Computational domain for a non-woven and woven spacer geometry with tracking planes after every filament intersection.

M. Johannink et al. / Desalination 376 (2015) 41–54

45

Fig. 4. Flowsheet of the lab-scale electrodialysis plant in the AVT.PT.

2.4. Semi-empiric modeling of pressure drop in spacer-filled channels The pressure drop Δp of an incompressible fluid flowing through a flow-channel of arbitrary geometry is typically represented by [34] Δp ¼ ζ

ρ v2ch ; dg 2

sp lch

ð2Þ

where ζ sp is the friction factor, vch the mean velocity and ρ is the density of the fluid. The quantity lch is the length and dg is a characteristic diameter of the channel. For the specification of an appropriate characteristic diameter dg different specifications are employed in literature [18,3]. For non-woven geometries (cf. [12]) some authors use the hydrodynamic diameter. Alternatively, the heights of the spacer-filled channel hch or the filament diameter df have been proposed [25]. As these specifications are applicable for conceptually different spacer geometries, we define dg = hch for the remainder of this publication. The friction factor in spacer-filled channels ζ sp is given by ζ

sp

sp

−msp

¼ n Re

;

ð3Þ

with the parameters n sp and m sp[7]. The Reynolds number Re is defined with the characteristic diameter dg = hch as

Re ¼

ρ vch hch : η

  ^k ; nsp ¼ g 1 uk ; p

ð6Þ

  ^k ; msp ¼ g 2 uk ; p

ð7Þ

relating them to the geometric design parameters of non-woven and woven spacers uk, k ∈ {nw, w} and to additional empirical parameters ^ k ∈ Rnp^k . Appropriate functional relations (6) and (7) are identified in p Section 4. Introducing Eqs. (6) and (7) into Eq. (3) yields the algebraic model ζ

sp

   ^k ¼ exp g 1 uk ; p   ln ðReÞ− ln ðRe Þ  min ^k ; −g 2 uk ; p ln ðRemin Þ− ln ðRemax Þ

^ p

c = [ρ, η, cfk]. For an efficient numerical treatment the inputs u are scaled as follows:

uˇ k ¼

 ln ðReÞ− ln ðRemin Þ ; ln ðRemin Þ− ln ðRemax Þ

ð8Þ

^kn , inputs u = [vch, uk] and constants ^k1 ; …; p with parameters p ¼ ½p k

ð4Þ

It is to be expected that a unique identification of nsp and msp in a parameter estimation is difficult or even impossible, given the similarities of Eq. (3) to model structures typically arising in reaction kinetics modeling [35,36]. Hence, Eq. (3) is reformulated to yield
~ sp −msp ζ sp ¼ exp n

The parameters ñ sp and m sp in Eq. (5) depend on the spacer geometry. In order to identify a model for different spacer geometries the parameters ñ sp and m sp in Eq. (5) are generalized by

uk −uk;min ; k ¼ 1; …; nu : uk;max −uk;min

ð9Þ

The model output is scaled similarly to ζ sp by introducing ζ sp min and

ζ sp max.

The lower and upper bounds that are used for the scaling schemes are summarized in Table 2.2.

ð5Þ 2.5. Model identification

where a new parameter ñ sp = exp(n sp) is introduced. Remin and Remax are the minimal and maximal Reynolds numbers, respectively, which are considered in the model identification process.

Methods for parameter estimation, identifiability analysis and model-based design of experiments employed for model identification are introduced in this section. Table 2.2 Lower and upper bounds used for the scaling of model output and inputs.

Table 2.1 Geometric properties of spacers used in the experiments.

vch ½ms

df [m]

lm [m]

α [°]

β [°]

2.5 · 10−4

9 · 10−4

0

90

Min Max

0.04 0.21

β 2 ½°

29.00 61.00

df [m]

ζ sp

lm[m] −3

0.9 · 10 1.1 · 10−3

−4

0.9 · 10 1.1 · 10−4

0.09 1231.51

46

M. Johannink et al. / Desalination 376 (2015) 41–54

2.5.1. Parameter estimation The parameter estimation problem is formulated as min

p∈Rn p ϕ ¼

1 RSS 2 σ2

identifiability of the entire model we use the ratio

ð10Þ

ndp  2 X MUSUBI ζ sp : i ðp; ui ; cÞ−ζ i

ð11Þ

i¼1

np and ndp correspond to the number of parameters and data points, respectively; ζiMUSUBI, i = 1, …, nmp, represents a data point corresponding to the experimental conditions ui. The uncertainty of the measurements is accounted for by the variance σ 2. For an efficient numerical treatment of the minimization problem first-order derivatives Jj ¼

∂ϕ ; ∂p j

j ¼ 1; …; np ;

ð12Þ

and second-order derivatives X ∂ζ sp ðp; ui ; cÞ 1 ∂ζ sp ðp; ui ; cÞ i i σ2 ∂p j ∂pk i¼1 ! 2 sp sp MUSUBI ζ i ðp; ui ; cÞ−ζ i ∂ ζ i ðp; ui ; cÞ þ ; σ2 ∂pk ∂p j

RSS nmp −ðnp

þ 2 ðnp þ 1Þ þ

ΔAIC i;c ¼ AIC i −AIC c ;

ð14Þ

j; k ¼ 1; …; np : The Jacobian J, Fischer matrix F and Hessian H are calculated using the Symbolic Math Toolbox in Matlab [37]. The minimization problem (10) is solved by the deterministic optimization routine fmincon provided with Matlab [37]. This local optimization method is used in a multistart approach to approximate the globally optimal solution of the parameter estimation problem. 2.5.2. Reliability of parameter estimates An identifiability analysis can be carried out to check whether the model parameters can be uniquely identified from a given set of data. Technically a model is regarded to be identifiable if for two distinct sets of parameters p1 and p2 the identity ζ sp ðp1 ; u; cÞ ¼ ζ sp ðp2 ; u; cÞ

2

lative t-distribution with nmp − (n p + 1) degrees of freedom [40]. Large confidence intervals indicate significant parameter uncertainty. A question closely related to identifiability is how many parameters are required to represent the information in the data. This question is of special relevance if a number of nm different model candidates are available where each is characterized by a different number of parameters. A well-established approach to investigate this question is based on Akaike's information criterion [41,42]  þ 1Þ

2 ðnp þ 1Þðnp þ 2Þ ; nmp −np −2

ð17Þ i ¼ 1; …; nm :

Relative measures are defined as

are required. The latter are often approximated by ndp sp X ∂ζ sp i ðp; ui ; cÞ 1 ∂ζ i ðp; ui ; cÞ ; 2 σ ∂p j ∂pk i¼1

Here, G = F−1 and t γ;nmp −ðnp þ1Þ is the upper γ2 percentile of the cumu-

ð13Þ

j; k ¼ 1; …; np ;

F j;k ¼

ð16Þ

k ¼ 1; …; np :


AIC i ¼ nmp log

ndp

H j;k ¼

where λmin is

the minimal eigenvalue of F. Another measure for the reliability of the estimates p* is the confidence intervals sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffi RSS Δθk ¼ t γ;nmp −ðnp þ1Þ Gk;k ; 2 nmp −ðnp þ 1Þ

with the residual sum of squares RSS ¼

λmin λmax ,

i ¼ 1; …; nm ;

ð18Þ

where candidate i is compared to the best candidate c. From the relative measures ΔAICi,c, i = 1, …, nm, Akaike weights wi,c are defined as   exp −0:5 Δ AIC i;c wi;c ¼ Xnm  ; exp −0:5 Δ AIC j;c j¼1

i ¼ 1; …; nm :

ð19Þ

j≠i

The weight wi,c describes the probability of candidate i to be the best model in the set of candidates [41]. Here, the best candidate is characterized by the minimal Akaike criterion, i.e., AIC c ¼i∈f1min ; …; nm g AIC i . The quality of the model fit is often characterized by the corresponding sum of squares RSS (cf. Eq. (11)) or more intuitively by the root pffiffiffiffiffiffiffiffi mean square error SRSS ¼ RSS=nmp. However, these measures provide only limited insight into the predictive capabilities of the model. To assess model quality the identification procedure should rather be supplemented by cross validation analysis (cf. e.g.,[43]) yielding the averaged root mean-square error SRSSCV.

ð15Þ

holds if and only if p1 = p2. Different methods of testing for identifiability have been proposed in literature [38]. In this study we will employ the eigenvalue method introduced by Vajda et al. [39] and recently extended by Quaiser and Moennigmann [38]. The resulting criterion for a parameter pk not to be uniquely identifiable is that the corresponding eigenvalue λj of the Fischer information matrix F is close to zero. From a practical point of view, it is important to note that the scaling of the model output ζ sp ðu; p; cÞ strongly impacts the absolute values of the eigenvalues λj. To eliminate these scaling effects, the identifiability of a parameter pk is evaluated using the normalized eigenvalue

λj λmax

with λmax being the maximal eigenvalue of F. As a measure for the

2.5.3. Design of experiments In order to identify the generalized model with a minimum number of CFD simulation, we apply methods for model-based DOE for parameter precision. Different design criteria have been derived exploring different measures of the Fischer information matrix F at the parameter estimates p* [44]. In the remainder of this contribution we use the well-known D-optimality criterion which is defined as max u∈Rnu D

¼ detð FÞ:

ð20Þ

The optimization problem (20) is solved again by using the local deterministic optimization routine fmincon in Matlab [37]. The required first-order derivatives are computed by using the ADiMat software [45] which implements automatic differentiation techniques.

M. Johannink et al. / Desalination 376 (2015) 41–54

Fig. 5. Flow pattern for a mean velocity of vch ¼ 0:2 spacer geometry.

47

m s (a) for the woven geometry with α = 0 ° used in the experiments, (b) for a woven geometry with α = 90 ° and (c) for a non-woven

3. Simulation results and validation Section 3.1 presents the simulated flow-patterns arising in the spacer-filled channels. This description is necessarily brief, as the qualitative pattern of the flow fields for non-woven spacer geometries have been extensively described elsewhere (cf. e.g.,[4,5,14,17,16]). In Section 3.2 the simulated pressure profiles are illustrated, which are compared in Section 3.3 to experimental data.

the mean flow direction immediately aligns with the orientation of the filaments. For nearly all the considered geometries the flow remains completely laminar. Very small periodic oscillations in the velocity and pressure profiles indicating the transition to turbulent flow are found for vch ¼ 0:2 ms and geometries characterized by df = 1 ⋅ 10− 4 m, lm = 3 ⋅ 10−3 m and β2 N55°. However, the qualitative time averaged pattern remains very similar to the laminar flow pattern. 3.2. Pressure drop profiles

3.1. Velocity profiles Fig. 5a shows the simulated flow pattern for the woven spacer with design parameters summarized in Table 2.1 for a mean velocity of vch ¼ 0:2 ms. Here, the stream lines illustrate the flow-path of tracer particles from the inlet plane through the channel. The fluid is strongly accelerated at the spacer filaments orthogonal to the mean flow directions. Stationary vortices arise at the intersections of the filaments in the completely laminar flow field. This laminar character of the hydrodynamics is observed for all simulations in the velocity range vch ¼ 0:05; …; 0:2 ms. The flow pattern arising in a spacer with identical geometrical properties, but with a flow attack angle α = 90 °, is shown in Fig. 5b. A significant fraction of the stream lines leaves the considered domain in the first part of the channel as the mean flow direction partly aligns with the orientation of the filaments. The acceleration of the fluid at the filament intersections is less significant compared to the woven spacer with α = 0 °. In all simulations for woven spacers with a flow attack angle α = 90 ° and design parameter in the range of considered inputs (cf. Table 2.2) the flow remains completely laminar. The hydrodynamics are finally simulated for a non-woven spacer with identical geometrical properties and a flow attack angle α = 90 °. The resulting stream line distribution is shown in Fig. 5c. Here,

Fig. 6 shows the averaged pressure ptp,i, i = 0, …, 5 at the tracking planes for the non-woven and woven spacers with design parameters summarized in Table 2.1 and flow attack angles α = 0 ° and α = 90 °, respectively. In the mean flow direction the simulation predicts a nearly linear pressure profile. Deviations — caused by inflow and outflow effects — are restricted to the first and last tracking plane. Hence, the friction factor is determined by ζ MUSUBI ¼

ptp;1 −ptp;3 2 ; Δz2;4 ρ v2ch

ð21Þ

where Δz2,4 is the distance between the second and the fourth tracking plane.

3.3. Comparison with experimental data The simulated pressure drop is compared to the pressure drop measured over an entire ED module and in the single-channel test cell in Fig. 7. The predictions of the single-channel measurements and the LBM simulation are in very good agreement for vch b0:1 ms. For vch N0:1 ms a slight deviation is observed. The measured pressure drop over the ED

48

M. Johannink et al. / Desalination 376 (2015) 41–54

Fig. 6. Linear trend of the averaged pressure at the tracking planes in the x-coordinate.

module is higher than the LBM prediction. This is reasonable, as the measured data also include contributions from the distribution system of the module. Overall, in the range of the measurement error, good agreement between experiment and model prediction is found. Hence, the LBMbased CFD solver can be used as a predictive tool to systematically investigate the hydrodynamics in spacer-filled channels.

4. Identification of friction factor models for varying design parameters This section presents the results of the identification of the generalized pressure drop models for non-woven (cf. Section 4.1) and woven (cf. Section 4.2) spacers by means of the integrated work-flow

Fig. 8. Root mean-square error (SRSS) of the model fit and mean cross-validation error (SRSSCV) for 17 model candidates for the non-woven spacers.

introduced in Section 2.1. The predictive capabilities of the identified models are illustrated in Section 4.3.

4.1. Non-woven spacer geometries In a first step, data is generated in simulation experiments for a set of 17 geometries. The design space is systematically covered by Latin hypercube sampling resulting in spacer geometries nos. 1–12 in Table A.1. Five additional geometries (nos. 13–17) are included to cover the upper and lower bounds of the design parameters. For each geometry in this set, LBM simulations are carried out at four mean velocities vch ∈f0:05; 0:1; 0:15; 0:2g ms resulting in 68 simulations in total. The computed data are used estimate the parameters in the generalized friction factor model (8). The functional relations (6) and (7) are represented by 17 alternative empirical model candidates which are summarized in the supplementary material. For each model candidate the parameter estimation problem (10) is solved

Table 4.1 Comparison of different model candidates for non-woven spacers.

Fig. 7. Pressure drop in a spacer channel with lsp = 0.2 m as predicted by the LBM simulation compared to experimental data from experiments with an ED module and a singlechannel test cell.

No.

np

λmin λmax

SRSSCV

ΔAIC

wi,8

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17

8 10 10 10 11 11 14 14 14 17 17 20 20 20 23 26 29

1.32E−03 6.81E−04 5.75E−04 8.22E−06 1.53E−04 9.14E−06 4.62E−06 1.19E−06 1.35E−05 2.41E−08 5.41E−08 2.00E−08 4.16E−07 3.19E−10 2.82E−10 2.55E−09 3.89E−11

6.19% 5.16% 5.12% 5.14% 4.91% 3.59% 4.11% 1.81% 4.68% 3.38% 3.02% 3.29% 3.05% 1.57% 1.57% 3.74% 2.78%

159.22 113.19 113.19 114.05 116.19 81.39 110.60 0 49.51 8.79 10.81 18.35 62.66 15.85 27.84 76.97 72.27

0.00% 0.00% 0.00% 0.00% 0.00% 0.00% 0.00% 98.30% 0.00% 1.22% 0.44% 0.01% 0.00% 0.04% 0.00% 0.00% 0.00%

M. Johannink et al. / Desalination 376 (2015) 41–54

49

Table 4.2 Optimal experimental conditions for non-woven spacers and performance measures of the resulting model after parameter estimation with augmented data. Cycle

Geometry

vnew ch ½ms

βnew/2 [°]

dnew f [10−4 m]

lnew m [10−3 m]

SRSS [%]

λmin/λmax [10−6]

Δθmax [10−4]

Δθav [10−4]

– 1 2 3 4 5 6 7 8 9 10 11 12 13 14

– 18 19 20 21 22 23 24 25 26 27 28 29 30 15

– 0.05 0.05 0.05 0.05 0.05 0.2 0.2 0.05 0.05 0.2 0.05 0.2 0.05 0.05

– 30.96 41.92 31.06 30 60 60 31.85 31.42 45.85 30 60 30 31.14 60

– 1 1 1 1 1 3 3 1.1 1 3 1.03 2.04 1 1

– 1.9 3 1.2 2.9 1.1 3 1 3 1 3 2.2 1 3 3

1.78 1.91 1.91 1.88 1.95 2.06 2.08 2.05 2.04 2.05 2.04 2.04 2.03 2.03 –

1.39 1.39 1.44 1.61 1.34 1.65 1.59 1.72 1.73 1.74 1.81 1.86 1.82 1.78 –

2.33 2.43 2.30 1.88 1.98 1.98 1.95 1.75 1.66 1.62 1.58 1.48 1.46 1.45 –

1.12 1.13 1.06 0.93 0.97 0.97 0.96 0.90 0.85 0.81 0.74 0.74 0.73 0.73 –

with the methods described in Section 2.5.1. The resulting SRSS and SRSSCV are compared in Fig. 8. Accordingly, the lowest SRSS is obtained for model 16; slightly higher SRSS are obtained for the models 8, 10–12 and 14–17. Significantly higher SRSS are found for the remaining model candidates. min The resulting eigenvalue ratios of the Fischer information matrix λλmax ,

the Akaike differences ΔAIC and the Akaike weights wi,8 are summarized min min in Table 4.1. For the ratios λλmax very low values, i.e., λλmax b1  10−7 , are ob-

served for the models 10–17. Hence, the corresponding parameter estimates are subject to significant uncertainty. Model 8 is the only candidate characterized by both, an acceptable eigenvalue ratio of

λmin λmax

In summary, model candidate 8 given by

nsp ðunw Þ

βˇ þ p4 dˇ f þ p7 ˇl m 2 2 2 2 βˇ þ p2 þ p5 dˇ f þ p8 ˇl m 2 −1 −1 −1 βˇ þ p3 þ p6 dˇ f þ p9 ˇl m þ p10 ; 2

¼ p1

¼

βˇ þ p12 dˇ f þ p13 ˇl m þ p14 ; 2

ð22Þ

1:19  10 and a good quality of the fit characterized by SRSS = 1.81 %. The Akaike information criterion ranks model 8 as the best model candidate in the set resulting in a relative measure ΔAIC8 = 0. Slightly higher measures are found for models 10–12, 15 and 17. For the remaining model candidates, the ΔAICi are significantly higher. These results are in agreement with the Akaike weights wi,8, i = 1, …, 17, which rank model 8 as the best candidate in the set with a probability of w8,8 = 98.3 %.

msp ðunw Þ ¼ p11

Fig. 9. Residual plot comparing the prediction of the final friction factor model ζ sp and the LBM prediction ζ MUSUBI for non-woven spacers.

Fig. 10. Root mean square error (SRSS) of the model fit and mean cross validation error (SRSSCV) for 17 model candidates for the woven spacers.

−6

ð23Þ

stands out by its acceptable identifiability properties while having only a slightly higher Srss compared to model candidate 16. Further, this model structure is rated as the best in the considered set according to Akaike's criterion. Hence, we will restrict our attention to this candidate in the following.

50

M. Johannink et al. / Desalination 376 (2015) 41–54

4.2. Woven spacer

Table 4.3 Comparison of 17 model candidates for woven spacers. No.

np

λmin λmax

SRSSCV

ΔAIC

wi,8

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17

6 8 8 8 8 8 10 10 10 12 12 14 14 14 16 18 20

1.61E−03 5.29E−05 8.05E−04 8.44E−04 5.48E−06 1.33E−05 4.35E−05 2.31E−06 6.10E−06 2.05E−09 5.03E−07 4.09E−08 1.44E−07 2.77E−19 2.45E−19 5.39E−09 2.18E−19

9.10% 7.11% 4.53% 4.61% 8.43% 4.64% 5.37% 3.67% 3.91% 9.10% 4.12% 4.26% 4.85% 3.70% 4.23% 4.96% 5.24%

78.30 64.13 28.52 32.79 83.77 32.54 32.96 0 15.40 8.90 9.50 19.86 29.89 19.88 31.34 40.79 54.99

0.00% 0.00% 0.00% 0.00% 0.00% 0.00% 0.00% 97.95% 0.04% 1.14% 0.85% 0.00% 0.00% 0.00% 0.00% 0.00% 0.00%

As for the non-woven spacer geometries, a reference data set of 13 woven geometries is generated (cf. Table A.2). LBM simulations are carried out for these geometries and for the four mean velocities vch ∈f0:05; 0:1; 0:15; 0:2g ms resulting in 52 simulations in total. For the identification of the unknown relations (6) and (7) 17 empirical model candidates are formulated as summarized in the supplementary material. The resulting SRSS and SRSSCV values are illustrated in Fig. 10. The SRSS of model candidates 8–17 are nearly identical and significantly lower than those of the model candidates 1–7. For model structures 10–17 the SRSSCV values significantly exceed the corresponding SRSS values. This indicates, that these candidates can describe the results in the reference data set with high accuracy, however, they are not suitable to predict pressure drop in the entire design space. The resulting measures assessing model performance are shown in Table 4.3. The model candidates 10–17 are characterized by λmin/ λmax b 10−7 clearly indicating that parameter identifiability is problematic, while the ratios are in the range 1.61 ⋅ 10− 3 b λmin/ λmax b 6.10 ⋅ 10−6 for the other candidates. Candidate 8 is characterized by an acceptable eigenvalue ratio λmin/λmax = 2.31 ⋅ 10−6.The Akaike information criterion ranks this model as the candidate with a probability of w8,8 = 97.95 %. Hence, model structure 8 given by

In the next step, the confidence intervals of the parameters in models (8), (22) and (23) are systematically improved by using the work flow introduced in Section 1 (cf. Fig. 1b). In each iteration, the DOE problem (20) is solved with the methods described in Section 3 to derive new exnew new new perimental conditions unew = [vnew β lm ]. For any new geomech , df new new try characterized by dnew β and l f m LBM simulations are carried out m for the mean inlet velocities vch ∈fvnew ch ; 0:05; 0:1; 0:15; 0:2g s . The resulting data are merged with the previous reference data set used before. This augmented data set is used to update the parameter estimates p* by solving the parameter estimation problem (10) before entering a new refinement cycle. The results of model refinement are summarized in Table 4.2. In the first five cycles the DOE method yields optimal experiments in which new vnew are fixed at the lower bounds while βnew and lnew are ch and df m varied. In cycles 5–13, the DOE approach yields new experimental conditions unew characterized by the systematic variation of all design parameters. The iterative model refinement is stopped after 13 iterations as the DOE method derives experimental conditions unew which are identical to geometry 15 in the reference set of geometries (cf. Table A.1). In the course of model refinement, the SRSS slightly increases from 1.78% to 2.03%. The eigenvalue ratio λmin/λmax is increased by 28.1% indicating a significant improvement in parameter identifiability. Further, a major improvement is found for the confidence intervals Δθk , k = 1, …, np . Here, the maximal Δθ max and the mean Δθav confidence interval are reduced by 60.7% and 53.4%, respectively. A residual plot for the final friction factor models (8), (22) and (23) is shown in Fig. 9. The final parameter estimates and confidence intervals are summarized in Table A.3 in the appendix.

2 2 −1 −1 nch ðuw Þ ¼ p1 dˇ f þ p4 ˇl m þ p2 dˇ f þ p5ˇlm þ p3 dˇ f þ p6 ˇl m þ p7 ;

ð24Þ

mch ðuw Þ ¼ p8 dˇ f þ p9 ˇl m þ p10 ;

ð25Þ

is considered as the most appropriate model structure to generalize the friction factor models (8), (24) and (25) for woven spacers. As in case of non-woven spacers, the models (8), (24) and (25) is refined using the model refinement work-flow (cf. Fig. 1). The results are summarized in Table 4.4. During model refinement, the SRSS slightly increases from 3.24% to 3.36%. A major improvement is found with respect to parameter identifiability as indicated by an increase in the eigenvalue ratio λmin/λmax by a factor of 235.5%. Further, the maximal Δθmax and the mean Δθav confidence intervals are reduced by 60.1% and 64.2%, respectively. A residual plot for the final model (8), (24) and (25) for woven spacers is shown in Fig. 11. The parameter estimates and confidence intervals are summarized in Table A.4 in the Appendix.

4.3. Predictive capabilities of the identified pressure drop models The dependence of the predicted friction factor log(ζ sp) on the design parameters is shown exemplary in Fig. 12 for a non-woven

Table 4.4 Optimal experimental conditions for the woven geometries and performance measures of the resulting model fit with augmented data. Cycle

– 1 2 3 4 5 6 7 8

Geometry

– 14 15 16 17 18 19 20 7

vnew ch

dnew f −4

½ms

[10

– 0.05 0.2 0.2 0.05 0.05 0.05 0.05 0.05

– 1.17 1 3 2.33 3 1.43 1.15 3

lnew m m]

−3

[10

– 1.14 1.50 2.50 1 1.1 1.49 1 1

SRSS m]

λmin/λmax −6

[%]

[10

3.24 3.26 3.24 3.31 3.42 3.42 3.39 3.36 –

2.31 5.56 5.13 5.21 5.11 5.10 5.23 5.45 –

]

Δθmax [10

−4

7.22 5.36 5.05 5.12 5.35 5.21 4.94 4.51 –

Δθav ]

[10−4] 3.71 2.90 2.75 2.78 2.87 2.73 2.60 2.39 –

M. Johannink et al. / Desalination 376 (2015) 41–54

51

vch is significantly higher than the sensitivity with respect to the geometric design parameter β2 (cf. Fig. 12c). This dominating influence of vch is not found in Fig. 12d where the friction factor model is evaluated for varying vch and df. Fig. 14a shows the influence of df and lm on the predicted friction factor indicating that the frequently used assumption of similar hydrodynamic conditions for a constant height to length ratio −3

−4

hch lm

is only

reasonable for lm b 1.8 ⋅ 10 and df b 2 ⋅ 10 . Outside this range, however, the increasing curvature of the isocontour-lines illustrates that this approximation is not appropriate. Fig. 13 illustrates that log(ζ sp) in a woven spacer is much stronger affected by vch than by the design parameters. The sensitivity of log(ζ sp) with respect to df and lm at vch ¼ 0:2 ms is illustrated in Fig. 14b. With either increasing df or lm the friction factor decreases. While for lm N 2.2 ⋅ 10−3 m the sensitivity with respect to df dominates, for df N 1.6 ⋅ 10−4 m the orthogonal distance of the filaments lm has a comparably stronger effect on log(ζ sp). In Fig. 14 the predicted log(ζ sp) is shown for varying df and lm at Fig. 11. Residual plot comparing the predictions of the final friction factor model ζ sp with the LBM prediction ζMUSUBI for the woven spacer.

geometry. In a large part of the design space a significant nonlinear variation of the friction factor with respect to the different design parameters is found (cf. Fig. 12a and b). Further, the sensitivity with respect to

vch ¼ 0:2 ms in (a) a non-woven geometry characterized by β2 ¼ 45° and (b) a woven geometry. In the considered range of design parameters, the pressure drop in non-woven geometries is significantly lower compared to the corresponding woven geometries. A common characteristic is found in the decreasing sensitivity of log(ζ sp) with respect to both design parameters with increasing df and lm. However, at the lower bound of df and lm the pattern of the iso-contour lines

Fig. 12. Predicted log(ζ sp) for non-woven geometries with varying design parameters: (a) varying df and β2 for lm = 1 ⋅ 10−3 m and vch ¼ 0:2 vch ¼ 0:2

m s;

(c) varying vch and β2 for df = 1 ⋅ 10−4 m and lm = 1 ⋅ 10−3 m and (d) varying vch and df for β2 ¼ 45° and lm = 1 ⋅ 10−3 m.

β −4 m s ; (b) varying lm and 2 for df = 1 ⋅ 10

m and

52

M. Johannink et al. / Desalination 376 (2015) 41–54

Fig. 13. Predicted log(ζ sp) for a woven geometry with (a) varying df and vch for lm = 1 ⋅ 10−3 m and (b) with varying lm and vch for lm = 3 ⋅ 10−4 m.

differs strongly: whereas for df b 2 ⋅ 10− 4 and lm b 2 ⋅ 10− 3 the isocontour lines show a nearly linear trend in the case of a non-woven geometry, a strongly nonlinear trend is observed for the woven geometry. Further, the total variation of log(ζ sp) is significantly higher for the non-woven spacer configuration. 5. Summary and discussion In this contribution a new LBM-based simulation tool has been used to systematically study the pressure drop in woven and non-woven spacer structures. The LBM prediction has been evaluated using data from experiments with both an ED membrane module and a tailored test-cell. In this comparison a good agreement was found between the experiments and the LBM prediction. The validated LBM has then been successfully used as a predictive tool to investigate the hydrodynamics in a large number of simulations in a reasonable time-frame. Predictive models for the friction factor of the spacer have been identified by generalizing the commonly used model structures to include the dependence on the design parameters by additional functional relations. Here, we used a systematic work-flow for the identification of an adequate generalized model structure and the iterative model refinement on the basis of DOE techniques. The developed work-flow has proven as an effective tool to handle the complexity arising from the multi-dimensional design space. The identified models allow for the first time a thorough quantitative investigation of the pressure drop in the entire design space of the

spacer geometries. Here, it is shown that the dependence of the pressure drop on the design parameters partially shows strongly nonlinear patterns. Hence, the generalization of trends from a small number of CFD simulations has to be considered problematic at least for a certain part of the design space. For both the non-woven and woven geometries, the sensitivity of the friction factor with respect to the mean velocity in the channel vch and the filament thickness df is comparable for a wide range of the considered geometric parameters. Compared to these factors, the sensitivity with respect to the design parameters β and lm is significantly lower. As a general result of the comparison between non-woven and woven structures we found that the pressure drop is significantly higher in the woven geometries. The presented work-flow effectively integrates detailed simulations and the identification of simplified models for a thorough quantitative exploitation of the design space. In future studies, the work-flow can be applied to other design objectives in a straight forward manner. Such studies should include appropriate measures for the mean mass transport intensity, the fouling potential and the axial dispersion in the channel.

Acknowledgments The authors would like to acknowledge financial support by the German Federal Ministry of Education and Research (BMBF) in the framework of the HPC software initiative in the project HISEEM (Grant No. 01IH11001B).

Fig. 14. Predicted log(ζ sp) for varying df and lm in a (a) non-woven and (b) woven spacer with β2 ¼ 45° for vch ¼ 0:2

m s.

M. Johannink et al. / Desalination 376 (2015) 41–54

53

Appendix A. Geometric parameters and parameter estimates

Table A.1 Geometric properties of non-woven spacers for the generation of the reference set of simulated pressure drop data. Geometry

β 2 ½°

df [10−4 m]

lm [10−3 m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17

38.75 41.25 58.75 53.75 43.75 48.75 56.25 31.25 51.25 46.25 33.75 36.25 30 30 60 60 60

1.75 2.75 1.25 2.25 1.08 2.92 2.42 1.92 2.08 1.58 2.58 1.42 1 3 1 3 1

1.25 2.25 2.75 1.75 1.58 2.42 1.08 2.92 2.58 2.08 1.42 1.92 1 1.3 3 1 1

Table A.2 Geometric properties of woven spacers for the generation of the reference set of simulated pressure drop data. Geometry

df [10−4 m]

lm [10−3 m]

1 2 3 4 5 6 7 8 9 10 11 12 13

1.13 2.34 2.70 1.61 1 3 3 1 2 2 1 3 2

1.94 1.21 2.06 2.66 1 3 1 3 1 3 2 2 2

Table A.3 Parameter estimates p* and confidence intervals θ of the final friction factor model for non-woven spacer geometries. k

1

2

3

4

5

6

7

p⁎k Δθk

−0.14196 9.43 · 10−5

0.07487 1.16 · 10−4

−0.00078 9.94 · 10−7

0.82367 1.16 · 10−4

−0.41847 1.45 · 10−4

−0.01312 1.84 · 10−6

−0.36281 1.12 · 10−4

8

9

10

11

12

13

14

0.22844 1.39 · 10−4

0.00669 1.69 · 10−6

0.86565 6.67 · 10−5

−0.13454 4.71 · 10−5

−0.04451 7.66 · 10−5

0.08643 4.91 · 10−5

0.90875 5.10 · 10−35

Table A.4 Parameter estimates p* and confidence intervals θ of the final friction factor model for woven spacer geometries. k

1

2

3

4

5

p⁎k Δθk

−0.17139 4.16E−04

0.10138 3.32E−04

0.00243 5.73E−06

0.97431 4.45E−04

−0.42512 3.51E−04

6

7

8

9

10

−0.01495 5.83E−06

1.17114 1.89E−04

0.10192 1.86E−04

0.00816 2.92E−04

1.55452 1.72E−04

54

M. Johannink et al. / Desalination 376 (2015) 41–54

Appendix B. Supplementary data Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.desal.2015.07.024. References [1] H. Strathmann, Ion-exchange membrane separation processes, Membrane Science and Technology Series, Vol. 9, Elsevier, Amsterdam, 2004. [2] M. Skiborowski, A. Mhamdi, K. Kraemer, W. Marquardt, Model-based structural optimization of seawater desalination plants, Desalination 292 (2012) 30–44. [3] Clemens Fritzmann, Micro Structured Spacers for Intensified Membrane Process Performance(Ph.D. thesis) RWTH Aachen, Aachen, 2012. [4] G.A. Fimbres-Weihs, D. Wiley, Numerical study of mass transfer in threedimensional spacer-filled narrow channels with steady flow, J. Membr. Sci. 306 (1–2) (2007) 228–243. [5] C.P. Koutsou, S.G. Yiantsios, A.J. Karabelas, Numerical simulation of the flow in a plane-channel containing a periodic array of cylindrical turbulence promoters, J. Membr. Sci. 231 (2004) 81–90. [6] F. Zamani, J.W. Chew, E. Akhondi, W.B. Krantz, A.G. Fane, Unsteady-state shear strategies to enhance mass-transfer for the implementation of ultrapermeable membranes in reverse osmosis: a review, Desalination 356 (2015) 328–348. [7] A.A. Sonin, M.S. Isaacson, Optimization of flow design in forced flow electrochemical systems, with special application to electrodialysis, Ind. Eng. Chem. Process. Des. Dev. 13 (3) (1974) 241–248. [8] J. Fárková, The pressure drop in membrane module with spacers, J. Membr. Sci. 64 (1–2) (1991) 103–111. [9] A. Da Costa, A. Fane, D. Wiley, Spacer characterization and pressure drop modelling in spacer-filled channels for ultrafiltration, J. Membr. Sci. 87 (1–2) (1994) 79–98. [10] G.A. Fimbres-Weihs, Numerical Simulation Studies of Mass Transfer Under Steady and Unsteady Fluid Flow in Two- and Three-dimensional Spacer-filled Channels(Ph.D. thesis) The University of New South Wales, Sydney and Australia, 2008-08-19. [11] C.P. Koutsou, S.G. Yiantsios, A.J. Karabelas, Direct numerical simulation of flow in spacer-filled channels: effect of spacer geometrical characteristics, J. Membr. Sci. 291 (2007) 53–69. [12] G. Schock, A. Miquel, Mass transfer and pressure loss in spiral wound modules, Desalination 64 (1987) 339–352. [13] M. Shakaib, S.M.F. Hasani, M. Mahmood, Study on the effects of spacer geometry in membrane feed channels using three-dimensional computational flow modeling, J. Membr. Sci. 297 (2007) 74–89. [14] J. Schwinge, D. Wiley, A. Fane, Novel spacer design improves observed flux, J. Membr. Sci. 229 (1–2) (2004) 53–61. [15] F. Li, W. Meindersma, A. d. Haan, T. Reith, Optimization of commercial net spacers in spiral wound membrane modules, J. Membr. Sci. 208 (1–2) (2002) 289–302. [16] Z. Cao, Cfd simulations of net-type turbulence promoters in a narrow channel, J. Membr. Sci. 185 (2) (2001) 157–176. [17] J. Schwinge, D.E. Wiley, D.F. Fletcher, A CFD study of unsteady flow in narrow spacer-filled channels for spiral-wound membrane modules, Desalination 146 (2002) 195–201. [18] C. Koutsou, S. Yiantsios, A. Karabelas, A numerical and experimental study of mass transfer in spacer-filled channels: effects of spacer geometrical characteristics and Schmidt number, J. Membr. Sci. 326 (1) (2009) 234–251. [19] C. Picioreanu, J. Vrouwenvelder, M. van Loosdrecht, Three-dimensional modeling of biofouling and fluid dynamics in feed spacer channels of membrane devices, J. Membr. Sci. 345 (1–2) (2009) 340–354. [20] A. Radu, J. Vrouwenvelder, M. van Loosdrecht, C. Picioreanu, Modeling the effect of biofilm formation on reverse osmosis performance: flux, feed channel pressure drop and solute passage, J. Membr. Sci. 365 (1–2) (2010) 1–15. [21] F. Li, W. Meindersma, A. d. Haan, T. Reith, Novel spacers for mass transfer enhancement in membrane separations, J. Membr. Sci. 253 (1–2) (2005) 1–12.

[22] J. Balster, I. Pünt, D.F. Stamatialis, M. Wessling, Multi-layer spacer geometries with improved mass transport, J. Membr. Sci. 282 (1–2) (2006) 351–361. [23] M. Shakaib, S.M.F. Hasani, M. Mahmood, Cfd modeling for flow and mass transfer in spacer-obstructed membrane feed channels, J. Membr. Sci. 326 (2009) 270–284. [24] F. Chaumeil, M. Crapper, Dem simulations of initial deposition of colloidal particles around non-woven membrane spacers, J. Membr. Sci. 442 (2013) 254–263. [25] M.S. Isaacson, A.A. Sonin, Sherwood number and friction factor correlations for electrodialysis systems, with application to process optimization, Ind. Eng. Chem. Process. Des. Dev. 15 (2) (1976) 313–321. [26] S. Chen, G.D. Doolen, Lattice Boltzmann method for fluid flows, Annu. Rev. Fluid Mech. 30 (1) (1998) 329–364. [27] J. Bernsdorf, U. Jaekel, T. Zeiser, T. Takei, H. Matsumoto, K. Nishizawa, Lattice Boltzmann simulation and visualisation of adsorption processes in complex geometries, in: P. Sloot, D. Abramson, A. Bogdanov, J. Dongarra, A. Zomaya, Y. Gorbachev (Eds.), Computational Science — ICCS 2003, Lecture Notes in Computer Science, Vol. 2657, Springer, Berlin Heidelberg 2003, pp. 1054–1061. [28] R.H. Perry, D.W. Green, Perrys Chemical Engineers Handbook, 8th ed. McGraw-Hill, New York, 2008. [29] M. Hasert, K. Masilamani, S. Zimny, H. Klimach, J. Qi, J. Bernsdorf, S. Roller, Complex fluid simulations with the parallel tree-based lattice Boltzmann solver MUSUBI, Journal of Computational Science [30] S. Roller, J. Bernsdorf, H. Klimach, M. Hasert, D. Harlacher, M. Cakircali, S. Zimny, K. Masilamani, L. Didinger, J. Zudrop, An adaptable simulation framework based on a linearized octree, in: M. Resch, X. Wang, W. Bez, E. Focht, H. Kobayashi, S. Roller (Eds.), High Performance Computing on Vector Systems 2011, Springer, Berlin Heidelberg 2012, pp. 93–105. [31] D.F. Harlacher, M. Hasert, H. Klimach, S. Zimny, S. Roller, Tree based voxelization of STL data, High Performance Computing on Vector Systems 20112011 81–92. [32] Math2Market GmbH, GeoDict 2013, http://www.geodict.com/, rev. 03.05.2014. [33] E. Glatt, St. Rief, A. Wiegmann, M. Knefel, E. Wegenke, Struktur und druckverlust realer und virtueller drahtgewebe, F S Filtr. Sep. 23 (2) (2009) 61–65. [34] VDI-Wärmeatlas, in: Verein Deutscher Ingenieure (Ed.) 10th EditionSpringer, Berlin [u.a.], 2006. [35] G. Buzzi-Ferraris, F. Manenti, Better reformulation of kinetic models, Comput. Chem. Eng. 34 (11) (2010) 1904–1906. [36] D.M. Himmelblau, Process Analysis by Statistical Methods, Wiley, New York, 1970. [37] The Math Works Inc., Matlab, Version 2013a, Natick, Massachusetts, 2013. [38] T. Quaiser, M. Mönnigmann, Systematic identifiability testing for unambiguous mechanistic modeling — application to JAK-STAT, MAP kinase, and NF-κB signaling pathway models, BMC Syst. Biol. 3 (1) (2009) 50. [39] S. Vajda, H. Rabitz, E. Walter, Y. Lecourtier, Qualitative and quantitative identifiability analysis of nonlinear chemical kinetic models, Chem. Eng. Commun. 83 (1) (1989) 191–219. [40] S. Marsili-Libelli, S. Guerrizio, N. Checchi, Confidence regions of estimated parameters for ecological systems, Ecol. Model. 165 (2–3) (2003) 127–146. [41] K.P. Burnham, D.R. Anderson, Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach, 2nd Edition Springer, New York, NY, 2008. [42] H. Akaike, Information theory and an extension of the maximum likelihood principle, in: B.N. Petrov, F. Csaki (Eds.),Second International Symposium on Information Theory, Budapest 1973, pp. 267–281. [43] S. Arlot, A. Celisse, A survey of cross-validation procedures for model selection, Statis. Surv. 4 (2010) 40–79. [44] E. Walter, L. Pronzato, Identification of parametric models from experimental data, Communications and Control Engineering, Springer and Masson, Berlin and New York and Paris, 1997. [45] C. Bischof, H. Bucker, B. Lang, A. Rasch, A. Vehreschild, Combining source transformation and operator overloading techniques to compute derivatives for matlab programs, Proceedings of the Second IEEE International Workshop on Source Code Analysis and Manipulation (SCAM 2002)2002 65–72.