Model Development and Experimental Data ... - ACS Publications

Jan 4, 2017 - Carbonate Membrane Scaling during Dead-End Filtration with ... Chemical Process and Energy Resources Institute, Centre for Research and ...
2 downloads 0 Views 4MB Size
Article pubs.acs.org/IECR

Model Development and Experimental Data Analysis for Calcium Carbonate Membrane Scaling during Dead-End Filtration with Agitation Margaritis Kostoglou,*,† Soultana T. Mitrouli,†,‡ and Anastasios J. Karabelas‡ †

Division of Chemical Technology, School of Chemistry, Aristotle University of Thessaloniki, GR 541 24 Thessaloniki, Greece Chemical Process and Energy Resources Institute, Centre for Research and Technology - Hellas, sixth km Charilaou-Thermi Road, GR 57001 Thermi-Thessaloniki, Greece



S Supporting Information *

ABSTRACT: The development of reliable criteria for determining the onset of scaling of sparingly soluble salts on reverse osmosis desalination membranes is of great practical significance. Calcium carbonate, which is dealt with in this study, is a commonly encountered membrane scale compound. Reliable modeling and simulation tools of the reverse osmosis operation (accounting for onset and development of scaling) are very important for process optimization, aiming to reduce the effects of scaling and to minimize the use of antiscalants; such comprehensive modeling requires knowledge of key scaling parameters. In recent years a procedure has been developed to extract these parameters by employing a facile testing procedure, involving a small pressurized dead-end filtration cell. Using such a stirred test cell, the conditions are much closer to the actual reverse osmosis desalination process, compared to those in a non-agitated cell. Here, a detailed model of the testing process is developed and a procedure for estimating the scaling parameters is presented, by employing available relevant data. The results of the present work confirm the previous suggestion that fluid agitation in the test cell leads to an improved scaling parameter estimation procedure.



INTRODUCTION Achievement of high treated-water recovery and minimization of process cost are major challenges in RO membrane operation, as relevant studies have shown that recovery in brackish water desalination should be above 70−80% in order to be economically feasible.1,2 However, the type and concentration of dissolved salts (in the feedwater) and particularly the sparingly soluble ones (mainly CaCO3 and CaSO4) determine the maximum recovery in membrane desalination processes.3 Beyond an upper limit of permeate recovery, the concentrate becomes supersaturated in certain salts that precipitate on the membrane. The true supersaturation ratio of scaling salts at the membrane surface is greater compared to that estimated for the bulk solution due to the concentration polarization phenomenon; the latter tends to enhance the scaling problem, with detrimental effects on desalination process efficiency and economics (due to reduction of permeate flux, increase in transmembrane © XXXX American Chemical Society

pressure, product quality deterioration, membrane life shortening) and negative environmental impact.4 It should be added that the complexity of underlying physicochemical processes during membrane desalination, near the membrane surface, in the complicated flow field of spacer-filled channels of membrane elements, renders the study and prediction of incipient crystallization on membrane surfaces quite difficult.5 Despite the extensive literature regarding bulk precipitation of sparingly soluble salts6,7 and heterogeneous crystallization on various substrates,8,9 significant gaps exist in our understanding of incipient scale formation on solid surfaces10 and especially on membrane surfaces during the desalination process; moreover, there is uncertainty regarding the induction period Received: October 21, 2016 Revised: December 13, 2016 Accepted: December 19, 2016

A

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 1. Experimental setup employed in experiments discussed in this study.

of membrane scaling.11 In most of the published work the induction time is determined by monitoring certain global fluid variables (fluid turbidity, conductivity, Ca2+ concentration, pH) or process parameters (permeate flux reduction),12 which have inadequate sensitivity for detecting incipient scaling, especially on desalination membrane surfaces.11 Therefore, comprehensive modeling of the whole desalination module (and plant) performance, which is invaluable in understanding, designing, and optimizing the process,13 is hampered by lack of progress in modeling scaling phenomena. In recent years a new technique was proposed to extract scaling parameters needed by large scale models of the desalination process. The technique is based on facile scaling tests performed in a small dead-end high-pressure filtration cell.14 A mathematical model of the process allows the extraction of the scaling parameters. However, the use of a non-stirred cell leads to a continuously evolving supersaturation which reaches rather high values compared to those prevailing in a practical desalination process. A dead-end high pressure filtration cell with agitation (very extensively used in membrane filtration research) is found to be appropriate for studying the calcium carbonate scaling on reverse osmosis membranes.15 Unfortunately, the analysis of the data from the stirred cell16 was problematic because the membrane surface property values were apparently different between the different small specimens used in the tests. Indeed, as shown in another study,17 membrane property values (such as permeability) vary from sample to sample; therefore, distributions of scaling parameter values related to fouling (instead of single values) should be sought.18 This new approach is employed here, by deriving a new model for scaling in the stirred cell and adapting it to relevant data, thus allowing extraction of the corresponding parameters. The structure of the present work is as follows: At first a brief description of the experimental procedure and results is given for completeness. Then the model development is described in detail. Several indicative results (parametric study) of the model are presented and discussed. Finally, the model is applied to the

experimental data and the extracted parameters are presented and discussed.



EXPERIMENTAL SECTION Synthetic feedwater was prepared by rapidly mixing solutions of CaCl2, NaHCO3, and in some cases NaCl for adjustment of the ionic strength, with typical concentrations for Ca2+, 250 or 1000 mg/L and HCO3− in the range 200−600 mg/L, as described in detail in a previous work.16 The supersaturation ratio (S) and index (SI) of calcium carbonate (regarding calcite, which was found to be the dominant crystalline phase) are defined as follows: ⎡ (αCa 2+)(αCO2−) ⎤1/2 3 ⎥ S=⎢ ⎢⎣ ⎥⎦ KSP

(1)

SI = 2·log S

(2)

where αCa2+, αCO32− are the activities of calcium and carbonate ions, respectively, computed by the PHREEQC code (version 2.15.07),19 applying the most complete database available, i.e., “minteq.v4.dat”. The concentration of each ionic species (Ca2+, HCO3−, Cl−, Na+) right at the membrane surface (cw), and subsequently their activity through the PHREEQC software, is determined by the following equation employed for each ion:20 c w,i cb,i

⎛J⎞ = (1 − R i) + R i·exp⎜ ⎟ ⎝ ki ⎠

(3)

cb,i, J are the known bulk species concentration and permeate flux, respectively, Ri is ion rejection and ki the mass transfer coefficient, estimated from measurements in the same type of stirred dead-end filtration cell, is calculated from the following correlation:15 Sh = 0.49·Re 0.55·Sc 0.33

(4)

In the Sherwood number, Sh = km·r/D, km is the average mass transfer coefficient, r the radius of active membrane surface area in the filtration cell and D the ion diffusion coefficient. The B

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Reynolds number is defined as Re = ω·r2/ν, where ω is the angular velocity of agitator and ν the kinematic viscosity of the liquid, and Sc = ν/D the Schmidt number. The stirring rates used in the present work were 60, 90, and 250 rpm corresponding to Reynolds number values of 1825, 2737, and 7604, respectively; 250 rpm was employed in most tests. The experimental setup, depicted in Figure 1, the experimental protocol, the procedure followed for the observation of membrane surfaces after the test termination, and the main experimental results from these tests were presented in a previous publication.16 The characteristics of local shear stress and mass transfer coefficient at the base (i.e., at the membrane surface) of the particular agitated cell employed in the scaling tests16 were studied in detail by Koutsou and Karabelas.15 The stirring rate of 250 rpm (employed in most of the scaling tests) was shown15 that it leads to average shear stress typical of that encountered in spacer-filled membrane channels of real modules. A table of experimental conditions with some major experimental results from that work is provided as Supporting Information (Table S1) for convenience. Figure 2 depicts a

From such distributions, the mean (d50) as well as the d5 and d95 crystal sizes were determined; a log-normal type function fits these data well. Furthermore, size distributions of particles detected in feedwater and retentate solutions were very similar implying that particle growth in the bulk occurs at a slow rate, whereas particles on membrane surface were much larger.16 Figure 3 includes a few SEM images from the CPA2 membrane surface, depicting the various forms of calcium carbonate detected on the membrane surface (mainly calcite). It should be noted that apart from the typical calcite rhombohedral particles, spherical or hexagonal-plate vaterite crystals and their agglomerates were also observed in these tests;16 it is possible that these crystals are transformed to the more stable calcite form for longer duration experiments.



MATHEMATICAL PROCESS MODELING Main Issues. Here the main issues are discussed regarding modeling the membrane scaling process in the experimental set-up dealt with in this article. A direct modeling approach is impossible due to many process complications. The first complication arises from the chemistry of calcium carbonate formation. There are several approaches to proceed concerning the number of intermediate species involved in the system CO2−CaCO3−H20.21−23 The simplest approach restricts the species to the following: dissolved CO2 (also referred to as H 2 CO 3 in the literature), HCO 3− , CO 3 2− , Ca 2+; this simplification is based on the negligibly small concentration of other relevant species. The supersaturation depends also on the ionic strength (or salinity in case of seawater). The second complication is due to the flow field considered here, which is much more complicated compared to that of dead-end filtration without agitation and even more complicated compared to that of spiral wound membrane (SWM) desalination modules.24 This complexity is due to the high frequency rotation of agitator which generates a 3-dimensional flow field. The complete three-dimensional simulation of the flow in agitated tanks is extremely computationally demanding, so several techniques exist to reduce dimensionality and computational effort.25,26 The existence of wall flux orders of magnitude smaller than the bulk velocities makes the situation more difficult. The Stefan−Maxwell equations27 for all the species must be solved in the three-dimensional flow field created by the agitator, taking also into account the finite kinetics of the reaction between the species (assumed typically to be in quasi equilibrium).

Figure 2. Crystal size distribution from test No 70D, Sw = 4.52, test duration 60 min, stirring rate 250 rpm.

typical crystal size distribution of particles detected on the desalination membrane (CPA2, Hydranautics), for 60 min filtration time and supersaturation ratio ∼4.5 (Test No 70D).

Figure 3. Typical SEM images depicting incipient calcite growth on a RO membrane at: (a) small magnification and (b) high magnification. Test 78D, Sw: 3.59, t = 60 min, stirring rate 250 rpm. C

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Beyond modeling the three-dimensional flow distribution and species concentration field, the next step is to include precipitation phenomena in the model, where both bulk and surface particles must be considered. Bulk particles are generated by homogeneous nucleation, and grow consuming dissolved species. In addition, bulk particles collide with each other due to their Brownian motion (affecting mainly submicron particles) and the shear generated by the agitation (affecting mainly particles larger than 1 μm).28 The fraction of collisions leading to aggregation depends on the local ionic strength. Therefore, a space-dependent population balance equation for the bulk particles including nucleation, growth, coagulation, convection, and diffusion phenomena must be written. In parallel, a population balance for the surface particles is needed, which should include heterogeneous nucleation and particle growth. The driving force for the creation of bulk and surface particles is different due to two factors: (i) The difference of the effective surface energy of the growing phase in the bulk and on the membrane surface. The reason for the smaller effective surface energy on surface is that for the same particle volume the surface particle energy contains contributions not only from solid−liquid interface but also from the solid−solid interface (particle−substrate) that has lower energy. (ii) The local supersaturation which is higher at the membrane surface because of the wall flux (concentration polarization). Another issue, which must be accounted for by a detailed model, is the deposition of the bulk particles on the membrane surface. There is a large body of literature on deposition of particles on permeable surfaces,29 but the part referring to growing particles is a small fraction. The particle trajectory is dictated by the flow field, its Brownian motion and by the particle−membrane interaction forces. In case of a bulk particle, coming in contact for a sufficient time-period with the membrane, a permanent bond develops (due to particle growth), thus the bulk particle is henceforth identified as a surface particle. Simplified Modeling Approach. The detailed model outlined above is almost impossible to solve with currently available computational means; thus, a series of simplifications must be made. In the present work a rather simple model will be developed in order to assess the predictive capabilities of the current theories for the experiments performed here and to contrast the experimental results with the model findings. The simplifications will be of two kinds, i.e., typical assumptions used in chemical engineering and assumptions based on the system experimental conditions and results. Starting from the flow field, it is decomposed into two regions. There is a region covering the entire liquid volume (except a thin layer close to the membrane); this is considered to be well-mixed under the action of the agitator. The second (thin) region close to the membrane is characterized by a dominant velocity component uw, normal to the membrane surface. There are also significant lateral velocity components, which tend to reduce material accumulation on the membrane. These flow currents can be described in terms of a mass transfer coefficient.30 Recently, an extensive experimental study was carried out regarding the wall shear stresses and mass transfer coefficients in the cell of the present experiments.15 It was found that the local mass transfer coefficient is a (rather weak) function of the radial coordinate in the cell. Since the radial variation of the mass transfer coefficient is modest and the concentration polarization modulus in the present experiments is small enough to invoke

linearity, it is reasonable to consider spatial uniformity and to set-up the problem for radially averaged quantities. Another assumption must be made regarding the bulk particles that approach the membrane. The flow field is not explicitly known, so trajectories cannot be constructed. The mass transfer coefficient for the particles is very small; thus, they are considered to be transferred to the surface with velocity uw. A fraction of these particles (denoted by λ) is permanently attached to the wall whereas the rest are entrained back to the well-mixed region under the action of shear stresses. The parameter λ is in general a function of the local shear stress and of the particle−membrane interaction forces. The shear stress exhibits a radial variation on the disc-type membrane specimen,15 implying a similar behavior for λ; however, for this study a lumped (radially averaged) value of λ can be considered in order to maintain the spatial uniformity of the model. It is noted that the parameter λ used here is the equivalent to the one used in Kostoglou and Karabelas5 and refers to the fraction of the particles reaching the surfaces. On the contrary, the deposition probability parameters φ and ψ introduced by Lee and Lee31 and Oh et al.32 refer to the fraction of the particles produced in the bulk. The next assumption refers to the chemical aspects of the problem. Only few studies take into account the detailed chemistry of CaCO3 formation in modeling surface precipitation.33,34 For the purposes of the present work, a simplified approach will be considered, previously used by Karabelas et al.14 The initial supersaturation is computed using an appropriate chemical equilibrium code (see previous section); then, the instantaneous supersaturation is calculated using the definition equation and reducing the dominant species concentration to account for the formed CaCO3. The mass transfer coefficient for both types of crystal-forming ions is assumed the same as it is approximately imposed by the electroneutrality condition.35 Finally, the experimental data correspond to at most 2% coverage of the membrane surface by the particles. The simulation results presented here correspond to coverage less than 5%; therefore, it is a good approximation to ignore the influence of the finite membrane coverage on the model results. Formulating an appropriate diffusion-convection equation for the region near the membrane and solving it for the species A, undergoing surface precipitation, leads to the following relation for the surface concentration of A: cAi = cAbe u w / km +

SA (1 − e u w / km) uw

(5)

Here cAb is the bulk concentration of A, km the mass transfer coefficient, and SA the rate of A consumption due to surface precipitation (mol/m2/s). Based on the experimental results, it can be shown that SA/(cAbuw) ≪ 1, where SA is computed as the amount of solids on the membrane divided by the corresponding time and membrane surface area. This condition is also valid for the parametric study presented here. The validity of the above condition leads to a great model simplification since the surface concentrations are obtained from the usual concentration polarization relations and the generation of surface particles does not affect the bulk species concentration. If cAb, cBb, and Ib are the bulk concentrations of the species A, B, and the bulk ionic strength, respectively, the corresponding surface concentrations cAi, cBi, and Ii are given as cAi = cAbe u w / km D

(6a) DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research c Bi = c Bbe u w / km

(6b)

Ii = Ibe u w / kmI

(6c)

by layer growth mechanism, with multiple two-dimensional nuclei in each layer.36 The molar rate of produced solid salt per unit area is given as JR = kR S7/6(S − 1)2/3 (ln(S))1/6 e−ER /ln(S)

where kmI is the mass transfer coefficient for the species determining the ionic strength. The initial bulk supersaturation is computed as Sbo = S(cAbo, cBbo, Ib) from the chemical equilibrium code (the subscript “o” denotes initial value). The initial surface supersaturation is computed as Sso = S(cAboeuw/km, cBboeuw/km, Ibeuw/kmI). The evolution of both supersaturations is calculated from the instantaneous values of cAb and cBb as Sb = Sbo(cAbcBb/ (cAbocBbo))0.5 and Ss = Sso(cAbcBb/ (cAbocBbo))0.5. To proceed, the concentration density functions f b(x,t) and fs(x,t) are introduced for the bulk and surface particles, respectively. The units for f b are #/m6 and for fs #/m5, such that f bdx denotes the number of bulk particles of volume between x and dx per liquid volume unit and fsdx the corresponding number of surface particles per membrane surface area unit. Basic phenomena occurring during incipient solid phase formation are dealt with in the following. Nucleation. It is assumed that the bulk particles have spherical shape and the surface particles have the shape of a hemisphere (i.e., contact angle with the substrate π/2). According to the nucleation theory,36 the nuclei diameter dci (i = b for bulk particles and i = w for surface particles) is given as

dci =

4Vmσi kBT ln(Si)

where the units of JR are moles of AB per square meter per second, kR is a rate constant having the same units as JR, and ER is a dimensionless parameter defining the sensitivity of the growth rate to the supersaturation. It should be noted that the rate expression (eq 10) refers to ion-solid phase interaction; thus, it is the same for the surface and bulk particles. Considering hemispherical shape of the particles and using simple geometrical arguments, the molar flux rate can be transformed to the volumetric growth rate (G(x) with units m3/s) for a particle of volume x: G(x) = 32/3(2π )1/3

(7)

KB(x , y) =

K s(x , y) =

(8)

x 2/3

(11)

2kBT 1/3 (x + y1/3 )(x−1/3 + y−1/3 ) 3μWB

(12)

1/3 1/2 1.3 ⎛ 3 ⎞ ⎜⎛ ε ⎟⎞ ⎜ ⎟ (x1/3 + y1/3 )3 Ws ⎝ 4π ⎠ ⎝ ν ⎠

(13)

where ε is the energy dissipation rate and ν is the fluid kinematic viscosity. The values of ε can be obtained from relations between the Power number and Reynolds number available in the literature.39 The stability ratios WB and Ws represent the effect of interparticle interaction forces on the coagulation rate. It is noted that for the ionic strength levels encountered in membrane applications their values are closed to unity.28 The governing equations for the bulk species are

16πσi3Vm2 3(kBT )3

VmNA

In case of an agitated tank the velocity fluctuations lead to additional collisions between particles according to the relation:

The units of Jnw are number of nuclei per square meter per second and of Jnb number of nuclei per cubic meter per second. A typical form of the function L(S) is L(S) = [ln(S)]q. The following dimensionless parameter is used for clarity of presentation: Bni =

JR

where NA is the Avogadro number. The same expression (eq 11) is valid for both bulk and surface particles because the exposed surface to volume ratio is the same for the spherical and the hemispherical shape. Coagulation of Bulk Particles. The particles in the liquid bulk are undergoing coagulation due to their Brownian motion. The so-called coagulation kernel which expresses the frequency of collisions between particles with volumes x and y is given, for the case of Brownian coagulation in the continuum regime, as38

where Vm is the molecular volume of the solid phase, T the temperature, kB the Boltzmann constant, and σb, σw are the surface energies for homogeneous and heterogeneous nucleation, respectively. The volume of the spherical bulk nuclei is αcb = πdcb3/6, whereas that of hemispherical surface nuclei is αcw = πdcw3/12. The bulk and surface nucleation rates Jnb and Jns are given by the relation: ⎞ ⎛ 16πσi3Vm2 ⎟ , i = b, w Jni = A ni L(S)exp⎜ − 3 2 ⎝ 3(kBT ) (ln(Si)) ⎠

(10)

(9)

A generalization of the above classical approach for the surface nucleation has been made recently by Mitrouli et al.18 by assuming that the surface is not homogeneous but consists of specific nucleation sites with different affinity for nucleation, as Curcio et al.37 have also suggested. The simplest version of this theory corresponding to single type of sites is to multiply the classical rate by (1 − N/Nt), where N is the actual number density of nuclei and Nt is the number density of nucleation sites. It is noted that the model presented here can also cover the case of heterogeneous bulk nucleation by simply replacing the real surface energy with an effective one (similar to the case of surface nucleation). Particle Growth. A poly nuclear surface reaction rate expression is considered. This expression is based on a layer

∂cAb 1 1 Jnb αcb − =− VmNA VmNA ∂t

∫0



G(x , S b)fb (x , t )dx (14)

c B,b = c Bbo − (cAbo − cAb)

(15)

The equation for cAb includes terms for accumulation and consumption due to nucleation and growth. The concentration cBb results from a simple balance. The population balance for the bulk particle size distribution is given as E

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 4. A schematic of the phenomena considered by the model.

∂fb (x , t ) ∂t

∂G(x , S b)fb (x , t )

=−

∂x

1 + 2

∫0

monodisperse approximation for f b and fs in the resulting system of equations, the problem is transformed to

+ Jnb δ(x − αcb(S b))

⎛ dcAb 1 ⎜ 2/3 3 (2π )1/3 kR S b 7/6(S b − 1)2/3 (ln(S b))1/6 =− dt VmNA ⎜⎝

x

K ( y , x − y ) f b ( y , z , t )f b

(x − y , z , t )dy − fb (x , t )

∫0

K (x , y)fb (y , t )dy As (1 − λ)u w fb (x , t ) V



⎛ M ⎞2/3 e−ER /ln(Sb)VmNAMbo⎜ b1 ⎟ ⎝ Mbo ⎠



⎞ ⎛ Bnb ⎞⎟ ⎟ + αcbA nb[ln(S b)]qb exp⎜ − 2 ⎝ (ln(S b)) ⎠⎟⎠

(16)

The first two terms of right-hand side correspond to particle growth and nucleation whereas the next two terms represent the coagulation process. The last term corresponds to loss of particles due to their deposition on the membrane; As is the membrane surface area and V the cell volume. The evolution equation for the size distribution of the surface particles is ∂fs (x , t ) ∂t

+

⎛ dMbo Bnb ⎞ 4kBT 2 ⎟− = A nb[ln(S b)]qb exp⎜ − Mbo 2 dt 3μ ⎝ (ln(S b)) ⎠

∂x

Msi(t ) =

∫0

∫0 ∞



x ifb (x , t )dx

x ifs (x , t )dx

Mb1Mbo −

AS (1 − λ)u w Mbo V

(21)

dMb1 = 32/3(2π )1/3 kR S b 7/6(S b − 1)2/3 (ln(S b))1/6 dt ⎛ M ⎞2/3 e−ER /ln(Sb)VmNAMbo⎜ b1 ⎟ + αcb A nb ⎝ Mbo ⎠

(17)

The left-hand side comprises accumulation and particle growth terms. The right-hand side contains the surface particle generation by nucleation and by bulk particle deposition. A schematic of the processes considered by the model is shown in Figure 4. The model can be further simplified employing the so-called monodisperse approximation. Instead of solving the population balances for the particle size distributions, it is assumed that all the particles are of the same size which is a good approximation for the average size of the actual distribution. The zero and first moments of fs and f b denote the particle number and volume per membrane surface area and liquid volume, respectively, and they are defined as (i = 0,1): Mbi(z , t ) =

4 ⎛⎜ ε ⎞⎟ π⎝ν⎠

1/2



∂G(x , Si)fs (x , t )

= Jns δ(x − αcs(Si)) + λu w fb (x , t )

(20)

⎛ Bnb ⎞ A ⎟ − S (1 − λ)u w Mb1 [ln(S b)]qb exp⎜ − 2 V S (ln( )) ⎠ ⎝ b (22)

⎛ dMso Bns ⎞ ⎟ + λu w Mbo = A ns[ln(Sw )]qs exp⎜ − 2 dt ⎝ (ln(Sw)) ⎠ (23)

dMs1 = 32/3(2π )1/3 kR Sw 7/6(Sw − 1)2/3 (ln(Sw ))1/6 dt ⎛ M ⎞2/3 e−ER /ln(Sw)VmNAMso⎜ s1 ⎟ + αcsA ns ⎝ Mso ⎠

(18)

⎛ Bns ⎞ ⎟ + λu w Mb1 [ln(Sw )]qs exp⎜ − 2 ⎝ (ln(Sw)) ⎠

(19)

In order to obtain a closed system of equations, it is assumed that the particle size distribution has monodisperse shape, i.e., f b(x,t) = Mboδ(x − Mb1/Mbo), fs(x,t) = Msoδ(x − Ms1/Mso). Multiplying eqs 16 and 17 by xi (i = 0,1), integrating with respect to x from x = 0 to x → ∞, and substituting the

(24)

The initial condition for the above system of ordinary differential equations is Mso = Ms1 = Mbo = Mb1 = 0, cAb = cAbo F

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Indicative Theoretical Results. Some indicative results will be provided, using the model, to understand the main features of the process. The system studied (i.e., membrane filtration with agitation) differs from the one of dead-end filtration without agitation in that the degree of concentration polarization is much smaller. In the case of no agitation, the bulk precipitation is restricted to a thin layer close to the membrane. The unidirectional flow toward the wall transfers all the bulk particles onto the surface transforming them to surface particles. Thus, there is no need to make a distinction between bulk and surface particles from the experimental point of view. The low value of concentration polarization in the case of agitation renders unavoidable the bulk precipitation. Most of the bulk particles remain in suspension due to agitation and to the shear forces at the membrane surface. Nucleation theory states that the dominance of surface nucleation over bulk nucleation is determined by the ratio of the heterogeneous to homogeneous surface energy. The energy-based parameters of the two nucleation modes are related through Bns = gBnb and αcs = gαcb where the parameter g = (σs/σb)3 is representative of the af finity of the membrane for surface nucleation; i.e., the smaller the g, the greater the affinity. The following set of physical parameters is used to perform some simulations: Bnb = 20, Ans = 108 (#/m2/s), Anb = 2 × 1012 (#/m3/s), qs = qb = 0, kR = 2 × 10−6 mol/m2/s, ER = 1, D = 10−9 m2/s. The operating parameters are uw = 10−5 m/s, co = 5 mol/m3, the stirring rate is F = 250 rpm and the solubility product for the base case is such that the initial supersaturation is Sbo = 3.2 (base case). The value Bnb is chosen so that the bulk particle number and size are of the order of the measured ones. It is noted that in all cases considered the number of bulk particles and their total mass are small enough to render negligible both Brownian and turbulent coagulation; consequently, there is no influence of coagulation in the presented results. In addition, it was found by preliminary calculations that in all cases considered here the reduction of the species concentration does not exceed 5% and it is much smaller in most of the cases (in particular for supersaturation ratio Sbo < 3.5). This means that the precipitation phenomena occur under almost constant concentration of active species. In such a case, there is only one-way coupling between surface and bulk particle evolution; i.e., through the deposited bulk particles. The surface particle number evolution for the base case with λ = 0 is shown in Figure 5 for several values of g. As g decreases, the affinity for surface nucleation increases and the surface particle number increases as shown in this figure. The evolution of surface particles for a case of low surface nucleation affinity (g = 1) and several values of particle attachment fraction λ is shown in Figure 6. The one-way coupling between bulk and surface particles implies that for each value of λ the number of the particles of the λ = 0 curve are generated on the surface and the remaining particles are coming from the bulk. As g decreases the ratio of the bulk-generated to surface-generated surface particles decreases. The evolution of the bulk particles concentration, under the same condition as those of Figure 6, is shown in Figure 7. As the parameter λ increases the deposition of bulk particles increases leading to the reduction of their number. The evolution of the bulk and surface particle average diameter, under the conditions shown in Figures 6 and 7, is presented in Figure 8. It is noted that for λ = 0 there is no influence of g on the evolution of the average surface particle size, since g affects both total number and total mass in the same way. The bulk particles grow slower than the surface

Figure 5. Evolution of the surface particle density for the base case with particle attachment parameter λ = 0 and several values of the membrane affinity parameter g.

Figure 6. Evolution of surface particle density for the base case with membrane affinity parameter g = 1 and several values of the particle attachment parameter λ.

particles due to the smaller supersaturation in the bulk than at the surface. The incorporation of the smaller bulk particles to the surface particles leads to a reduced average size of the surface particles as λ increases. The number of surface particles at times 1000 s and 6000 s versus initial supersaturation is shown in Figure 9. Particle number is a rather steep function of supersaturation; as is evident in this figure, the particle number clearly increases for all supersaturation values when g decreases for constant λ. The situation is somewhat different in case of increasing λ for constant g; indeed, the surface particle number increases only for supersaturation higher than a critical value beyond which bulk nucleation occurs. Only under conditions at which bulk G

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 9. Surface particle density versus initial supersaturation for the base case at two time values t = 1000 s (solid lines) and t = 6000 s (dashed lines) and three combinations of parameters g and λ.

Figure 7. Evolution of bulk particle concentration for the base case with membrane affinity parameter g = 1 and several values of particle attachment parameter λ.

Figure 10. Average bulk and surface particle diameter versus initial supersaturation for base case (with g = 1) at two time values t = 1000 s (solid lines) and t = 6000 s (dashed lines) and two values of parameter λ.

Figure 8. Evolution of surface and bulk particles average diameter for the base case and several λ values (0, 0.1, 0.3, 0.7, 1). Particle diameter decreases with increasing λ.

direct influence on the concentration polarization. This is a well understood (even quantitatively) phenomenon and should not be confused with the mixing effect on bulk precipitation.6 The reduction of the stirring rate from 250 to 60 rpm is equivalent to an increase of the supersaturation of approximately 0.5 units. It should be noted that such a relation between scaling and stirring rate is shown here quantitatively for the first time; moreover, this strong effect is also exhibited by the experimental data (i.e., Table S1, Supporting Information). The corresponding average surface particle size is depicted in Figure 12, showing an increasing particle size with decreasing stirring rate. The dependence of particle size on stirring rate is much greater for small stirring rates. Summarizing, the surface particle size appears to be independent of membrane affinity

particles exist can an increasing λ lead to increased number of surface particles. The average bulk and surface particle diameter at t = 1000 s and t = 6000 s is shown in Figure 10. There is a small difference in the surface particle size between the cases λ = 0 and λ = 1. The difference tends to decrease with decreasing supersaturation. At small supersaturation values, the number of bulk particles is very small, so the value of λ does not affect the surface particle number and size. The effect of the stirring rate on the evolution of the surface particles is examined next. The surface particle number density at t = 1000 s and t = 6000 s is shown for three values of initial supersaturation and g = 0.8, λ = 0 in Figure 11. The effect of stirring rate on the particle number is very strong due to its H

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

agitation where the wall supersaturation tends to increase to high values, with subsequent reduction due to ionic species consumed by fast precipitation.14 However, the stirred cell precipitation conditions are similar to those prevailing in desalination SWM modules where the precipitation rate is relatively small and has no influence on the species concentration profile.11 The practical absence of bulk particles suggests a large value for Bnb (let us say 50) which ensures that formation of bulk particles is not predicted by the model. The parameters qb, qs, and ER are set equal to zero since they do not have an essential influence on model performance. The remaining parameters are Ans, Bns, Nt, and kR with the three first referring to nucleation and the last one to particle growth. The procedure allows identification of two parameters from each experiment that provides two data points (i.e., surface particle density and average particle diameter). In order to proceed to parameter estimation, the number of nucleation parameters must be reduced by choosing values for some of them. The total number of nucleation sites Nt is assumed to be much larger compared to the measured number of nuclei, eliminating in this way the term 1 − Mo/Nt in the nucleation term. Unlike the case of unstirred cell,18 the two targets (particle number and size) are almost not correlated (i.e., particle number depends almost entirely on nucleation parameters and particle size depends only on growth rate parameters), thus facilitating greatly the parameter identification process. A value for Ans is assumed and the inverse problem (fitting model to the experimental data) is solved for Bns and kR. The value of Ans that minimizes the spread of Bns parameter among the experiments is chosen at the end of the procedure (Ans = 1.65 × 1010 #/m2/s). The corresponding values of Bns are shown in Figure 13 vs experiment i.d. (defined in the Table S1). It is noteworthy that there is a region of high values for Bns that shows a small dispersion. In addition, there is a number of small values which confuse the overall picture of the parameter Bns. It is quite important that the small values correspond to low wall supersaturation (up to 2.75) with a

Figure 11. Surface particle density versus stirring rate for base case (with g = 0.8 and λ = 0) at two time values t = 1000 s (solid lines) and t = 6000 s (dashed lines) and three values of initial supersaturation (So = 2.2, 2.6, 3).

Figure 12. Surface-particle average diameter versus stirring rate for base case (with g = 0.8 and λ = 0) at two time values t = 1000 s (solid lines) and t = 6000 s (dashed lines) and three values of initial supersaturation (So = 2.2, 2.6, 3).

parameter g, slightly dependent on particle attachment parameter λ and strongly dependent on stirring rate F. It is noted that the experimental data cannot be presented in figures similar to the above, employed for model sensitivity analysis. The reason is that each experiment was performed with a different membrane sample, having different surface properties. A procedure to handle data in such a case has been recently proposed17 and will be employed in the next section.



APPLICATION OF THE MODEL TO EXPERIMENTAL DATA−PARAMETER ESTIMATION Measurements of bulk particle numbers revealed that their concentration was too small to have a significant direct contribution to the surface mass evolution. Additionally, by employing experimental data obtained at the end of the tests and performing a global mass balance, it was found that in all cases the experiments were actually carried out under constant active species concentration. This implies that the situation in stirred cell desalination is quite different compared to that of no

Figure 13. Estimated value of nucleation parameter Bns vs i.d. of particular experiments [Table S1]. I

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research monotonic (increasing) relation to Sw. On the other hand, the larger values correspond to Sw ≥ 3.43. The above observation makes clear that the conventional theory used does not work for supersaturation smaller than 2.75. A new relation for smaller supersaturation values must be employed. Trying to keep the general form of the classical relation (which leads to zero nucleation rate at S = 1) and to reproduce the nucleation rate almost independent of S (as indicated in the small S experiments), the final choice is a simple modification of the classical relation. Specifically, the exponent 2 of ln(S) in the classical relation is replaced by 0.1. The value of Ans is the same as for high S values and the estimated Bns values are 11.91, 13.27, 13.71, 13.42, 14.47, and 11.23. It is noted that the membranes in experiments with no particles detected may have even smaller values of Bns; in any case, the spread of this parameter is considered quite acceptable. Summarizing, it was found that the nucleation rate relation which reproduces the experimental data is 1.65 × 1010exp(−B/(ln(S))v) #/m2/s where Bns = 27.25 ± 7.15, ν = 2 for S ≥ 3, and Bns= 13 ± 1.3, ν = 0.1 for S < 3. The choice of S = 3 as borderline among the two expressions is somewhat arbitrary since the actual data refer to 2.75 and 3.43. However, it is considered a reasonable choice in order to obtain a complete nucleation rate expression. It is noted that the extraction of a relation for nucleation rate would not be possible from unstirred vessel relation since in each experiment not a specific value but a whole range of S values is traced. In addition, the experimental results are cross-correlated to the unknown parameters making their estimation more difficult. The stirred cell experiments allow modification of the nucleation rate function and minimization of the parameters spread. The spread of Bns is much smaller than the one found for unstirred cell experiments since in the present case it is possible to separate the uncertainties due to membrane property variations and to model validity. Moreover, it should be recalled that Bns is proportional to the third power of surface energy which means that the spread found corresponds to less than 10% spread of the surface energy among the membrane samples. Finally, the obtained distribution of Bns values, justifies the assumption that they follow normal distributions with values of the average and standard deviations given above. Estimation of growth parameter is dealt with next. The distribution of parameter kg = krVmNA is shown in Figure 14. This variable is used in order to allow comparison with the results of the previous work.18 A first attempt to estimate kg using the polynuclear growth model was unsuccessful as it led to kg correlated to S. The particular experimental data set (unlike the one of non-stirred cell) allows optimization of growth rate expression by requiring a kg uncorrelated to S. The simplest way to achieve this is to retain only the logarithmic term in the growth rate expression (i.e., polynuclear growth rate does not hold at relatively low supersaturation). The range of the kg parameter values is similar to the one found for the nonstirred cell experiments;18 i.e., its average value is 6.48 × 10−10 m/s. The distribution is clearly skewed with a long tail, so it cannot be described appropriately by a normal one; it is noted that standard deviation is of small significance in this case. A log-normal or a Gamma distribution function is appropriate for this type of data. It is noted that if three of the data points were rejected as outliers, the rest can be considered to follow a normal distribution; however, after a closer examination, there was no rational reason to exclude these data. In closing, the particular data from the stirred cell allow a clearly better

Figure 14. Estimated value of nucleation parameter kg vs i.d. of particular experiments [Table S1].

exploitation compared to those of the non-stirred cell since the estimated growth and nucleation parameters are not correlated as is the case for the non-stirred experimental setup described in the previous work.18



CONCLUSIONS In this publication, a systematic effort is described to extract scaling parameters from membrane desalination experiments in a pressurized stirred cell. These parameters are useful for reliable modeling and comprehensive performance simulation of desalination modules and entire plants. As already proposed, fluid agitation leads to an almost constant supersaturation at the membrane surface during the scale evolution process, and it is advantageous compared to the non-agitated cell case for two main reasons: (i) it better mimics the conditions of spiral wound membrane modules used in practice and (ii) it allows easier and more accurate determination of the scaling parameters. A comprehensive scaling model for a stirred filtration cell was developed in this article. The model was used to extract parameters from already published relevant data, considering that the scaling parameters do not have a specific value but they exhibit a distribution dictated by the variability of the membrane properties of the relatively small membrane specimens used in the tests. Value ranges and distributions for nucleation and growth parameters as well as optimized expressions for nucleation and growth rate functions are extracted and assessed. The results confirm the suggestion that the stirred cell holds distinct advantages, compared to the nonstirred one, for the determination of the scaling parameters in desalination processes.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.6b04063. Table S1 with experimental conditions and selected data employed in parameter estimation (PDF) J

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research



(18) Mitrouli, S. T.; Kostoglou, M.; Karabelas, A. J. Calcium carbonate scaling of desalination membranes: Assessment of scaling parameters from dead-end filtration experiments. J. Membr. Sci. 2016, 510, 293. (19) Parkhurst, D. L.; Appelo, C. A. J. User’s guide to phreeqc (version 2) - a computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations; U.S. Geological Survey: Earth Science Information Center, 1999; Open-File Reports Section. (20) Mulder, M. Basic Principles of Membrane Technology; 2nd ed.; Kluwer Academic Publishers: The Netherlands, 1996. (21) Schausberger, P.; Mustafa, G. M.; Leslie, G.; Friedl, A. Scaling prediction based on thermodynamic equilibrium calculation - scopes and limitations. Desalination 2009, 244 (1−3), 31. (22) Söhnel, O.; Garside, J. Precipitation: Basic Principles and Industrial Applications; Butterworth-Heinemann: Oxford, 1992. (23) Mook, W. G. Chemistry of carbonic acid in water.Environmental Isotopes in the Hydrological Cycle. Principles and Applications; International Atomic Energy Agency (IAEA) and United Nations Educational, Scientific and Cultural Organization (UNESCO), 2001; Vol. 1, Chapter 9. (24) Koutsou, C. P.; Yiantsios, S. G.; Karabelas, A. J. A numerical and experimental study of mass transfer in spacer-filled channels: Effects of spacer geometrical characteristics and Schmidt number. J. Membr. Sci. 2009, 326 (1), 234. (25) Ranade, V. V. Computational Flow Modeling for Chemical Reactor Engineering; Academic Press: New York, 2001. (26) Samaras, K.; Zouboulis, A.; Karapantsios, T.; Kostoglou, M. A CFD-based simulation study of a large scale flocculation tank for potable water treatment. Chem. Eng. J. 2010, 162, 208. (27) Fimbres-Weihs, G. A.; Wiley, D. E. Review of 3D CFD modeling of flow and mass transfer in narrow spacer-filled channels in membrane modules. Chem. Eng. Process. 2010, 49 (7), 759. (28) van de Ven, T. G. M. Colloidal Hydrodynamics; Academic Press: New York, 1989. (29) Chen, J. C.; Kim, A. S. Monte Carlo simulation of colloidal membrane filtration: Principal issues for modeling. Adv. Colloid Interface Sci. 2006, 119, 35. (30) Colton, C. K.; Smith, K. A. Mass transfer to a rotating fluid: Part II. Transport from the base of an agitated cylindrical tank. AIChE J. 1972, 18 (5), 958. (31) Lee, S.; Lee, C.-H. Effect of operating conditions on CaSO4 scale formation mechanism in nanofiltration for water softening. Water Res. 2000, 34 (15), 3854. (32) Oh, H.-J.; Choung, Y.-K.; Lee, S.; Choi, J.-S.; Hwang, T.-M.; Kim, J. H. Issues 1 and 2: First International Workshop between the Center for the Seawater Desalination Plant and the European Desalination SocietyScale formation in reverse osmosis desalination: model development. Desalination 2009, 238 (1), 333. (33) Segev, R.; Hasson, D.; Semiat, R. Rigorous modeling of the kinetics of calcium carbonate deposit formation - CO2 effect. AIChE J. 2012, 58 (7), 2286. (34) Segev, R.; Hasson, D.; Semiat, R. Rigorous modeling of the kinetics of calcium carbonate deposit formation. AIChE J. 2012, 58 (4), 1222. (35) Wesselingh, J. A., Krishna, R. Mass Transfer; Ellis Horwood Ltd: Great Britain, 1990. (36) Dirksen, J. A.; Ring, T. A. Fundamentals of crystallization: Kinetic effects on particle size distributions and morphology. Chem. Eng. Sci. 1991, 46 (10), 2389. (37) Curcio, E.; Curcio, V.; Di Profio, G.; Fontananova, E.; Drioli, E. Energetics of protein nucleation on rough polymeric surfaces. J. Phys. Chem. B 2010, 114 (43), 13650. (38) Friedlander, S. K. Smoke, Dust and Haze; Wiley Interscience: New York, 1977. (39) Papanastasiou, T. C. Applied Fluid Mechanics; Prentice Hall: New York, 1994.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]; Telephone: 00302310997767. ORCID

Margaritis Kostoglou: 0000-0001-7955-0002 Notes

The authors declare no competing financial interest.



REFERENCES

(1) Antony, A.; Low, J. H.; Gray, S.; Childress, A. E.; Le-Clech, P.; Leslie, G. Scale formation and control in high pressure membrane water treatment systems: A review. J. Membr. Sci. 2011, 383 (1−2), 1. (2) Shih, W.-Y.; Rahardianto, A.; Lee, R.-W.; Cohen, Y. Morphometric characterization of calcium sulfate dihydrate (gypsum) scale on reverse osmosis membranes. J. Membr. Sci. 2005, 252 (1−2), 253. (3) Greenberg, G.; Hasson, D.; Semiat, R. Limits of RO recovery imposed by calcium phosphate precipitation. Desalination 2005, 183 (1−3), 273. (4) Shirazi, S.; Lin, C.-J.; Chen, D. Inorganic fouling of pressuredriven membrane processes – A critical review. Desalination 2010, 250 (1), 236. (5) Kostoglou, M.; Karabelas, A. J. Modeling scale formation in flatsheet membrane modules during water desalination. AIChE J. 2013, 59 (8), 2917. (6) Mullin, J. W. Crystallization, 3rd ed.; Butterworth-Heinemann: Oxford, 1993. (7) Spanos, N.; Koutsoukos, P. G. Kinetics of precipitation of calcium carbonate in alkaline pH at constant supersaturation. Spontaneous and seeded growth. J. Phys. Chem. B 1998, 102 (34), 6679. (8) Dalas, E.; Klepetsanis, P. G.; Koutsoukos, P. G. Calcium carbonate deposition on cellulose. J. Colloid Interface Sci. 2000, 224 (1), 56. (9) Lioliou, M. G.; Paraskeva, C. A.; Koutsoukos, P. G.; Payatakes, A. C. Heterogeneous nucleation and growth of calcium carbonate on calcite and quartz. J. Colloid Interface Sci. 2007, 308 (2), 421. (10) Wang, H.; Alfredsson, V.; Tropsch, J.; Ettl, R.; Nylander, T. Formation of CaCO3 deposits on hard SurfacesEffect of bulk solution conditions and surface properties. ACS Appl. Mater. Interfaces 2013, 5 (10), 4035. (11) Mitrouli, S.; Karabelas, A. J.; Karanasiou, A.; Kostoglou, M. Incipient calcium carbonate scaling of desalination membranes in narrow channels with spacersexperimental insights. J. Membr. Sci. 2013, 425−426 (0), 48. (12) Hasson, D.; Shemer, H.; Brook, I.; Zaslavschi, I.; Semiat, R.; Bartels, C.; Wilf, M. Scaling propensity of seawater in RO boron removal processes. J. Membr. Sci. 2011, 384 (1−2), 198. (13) Kostoglou, M.; Karabelas, A. J. Comprehensive simulation of flat-sheet membrane element performance in steady state desalination. Desalination 2013, 316, 91. (14) Karabelas, A. J.; Kostoglou, M.; Mitrouli, S. T. Incipient crystallization of sparingly soluble salts on membrane surfaces: The case of dead-end filtration with no agitation. Desalination 2011, 273 (1), 105. (15) Koutsou, C. P.; Karabelas, A. J. Shear stresses and mass transfer at the base of a stirred filtration cell and corresponding conditions in narrow channels with spacers. J. Membr. Sci. 2012, 399−400, 60. (16) Mitrouli, S. T.; Kostoglou, M.; Karabelas, A. J.; Karanasiou, A. Incipient crystallization of calcium carbonate on desalination membranes: dead-end filtration with agitation. Desalin. Water Treat. 2016, 57 (7), 2855. (17) Mitrouli, S. T.; Karabelas, A. J.; Isaias, N. P. Polyamide active layers of low pressure RO membranes: Data on spatial performance non-uniformity and degradation by hypochlorite solutions. Desalination 2010, 260 (1−3), 91. K

DOI: 10.1021/acs.iecr.6b04063 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX