Supercritical water oxidation of Ion Exchange Resins in a stirred reactor: Numerical modelling

Chemical Engineering Science 69 (2012) 170–180

Contents lists available at SciVerse ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Supercritical water oxidation of Ion Exchange Resins in a stirred reactor: Numerical modelling A. Leybros a,n, A. Roubaud a, P. Guichardon b, O. Boutin c a

CEA Marcoule, DEN/DTCD/SPDE/LFSM, BP17171, 30207 Bagnols/Ceze, France Ecole Centrale Marseille, UMR-CNRS 6181, Technopole de Chˆ ateau-Gombert, 38, rue F. Joliot-Curie, 13451 Marseille Cedex 20, France c Aix Marseille Universite´s, UMR-CNRS 6181, BP 80, Europˆ ole de l’Arbois, 13545 Aix en Provence Cedex 4, France b

a r t i c l e i n f o

a b s t r a c t

Article history: Received 14 March 2011 Received in revised form 7 September 2011 Accepted 6 October 2011 Available online 18 October 2011

Supercritical water oxidation offers a viable alternative treatment to destroy the organic structure of Ion Exchange Resins. In order to design and define appropriate dimensions for the supercritical oxidation reactor, a 2D simulation of the fluid dynamics and heat transfer during the oxidation process has been investigated. The solver used is a commercial code, Fluents 6.3. The turbulent flow field in the reactor, created by the stirrer is taken into account with a k  o model and a swirl imposed to the fluid. Particle trajectories are modelled with the Discrete Random Walk Particle Model. For the solubilization of the particles in supercritical water, a mechanism has been proposed and implemented into Fluents software through the Eddy Dissipation Concept approach, taking into account the identified rate determining species. Simulation results provide results on the flow, temperature fields and oxidation localization inside the reactor. For the reactive particles-supercritical water flow model, the effect of parameters, such as feed flow rates or stirring velocity, can be focussed. Reaction temperature is predicted with deviation lower than 15%. Degradation conversions are in good agreement with experimental ones. & 2011 Elsevier Ltd. All rights reserved.

Keywords: Supercritical fluid Multiphase flow Computational Fluid Dynamics Turbulence Hydrodynamics Ion Exchange Resins

1. Introduction Water treatment systems in nuclear power plants use organic Ion Exchange Resins (IERs) in order to control system chemistry, to minimize the corrosion and the degradation of system components and to remove radioactive contaminants. Organic resins are also used in chemical decontamination or cleaning processes for the regeneration of process water by reagents and for radionuclide removal. Large quantities of contaminated IERs are in this way generated and often stored. In the future, wastes should be treated to reduce the volume and to improve their stability. However, spent resins are radioactive wastes for which there is no satisfactory industrial treatment (IAEA, 2002). In the case of IERs, SuperCritical Water Oxidation (SCWO) could offer an alternative viable treatment for the destruction of organic structure of the resins containing radioactive elements. At supercritical conditions (P 422.1 MPa and T4374 1C), water acts as a non-polar dense gas and its solvation properties are similar to a low polarity organic solvent and is completely miscible with liquid organics. Supercritical water also presents complete miscibility with oxygen, creating a homogeneous reaction medium,

n

Corresponding author. Tel.: þ33 4 66 33 90 19; fax: þ 33 4 66 79 78 71. E-mail address: [email protected] (A. Leybros).

0009-2509/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2011.10.016

which makes supercritical water a very suitable medium for the oxidation of liquid organics. This reaction can also be heterogeneous when the organic material is a solid. Reaction time is about a few seconds and there is no formation of gaseous oxides (SOx and NOx), due to temperatures around 450–500 1C (compared to combustion, greater than 1000 1C). The radionuclides would stay in the cold liquid effluent and could be recovered by means of precipitation. The main problems of SCWO are corrosion and salt precipitation (leading to reactor plugging) (Bermejo and Cocero, 2006). To overcome these drawbacks, a lot of different kind of reactor designs reviewed in (Brunner, 2009) have been developed, such as tubular reactor (Cocero, 2001), cooled wall reactor (Cocero, 2001) or transpiring wall reactor (Bermejo et al., 2006). Amongst them, the technology of double shell stirred reactor has been developed and patented by CEA (Calzavara et al., 2004; Joussot-Dubien et al., 2000); this technology is used in this work and described later. To scale up the reactor, modelling using Computational Fluid Dynamics is a good way to better understand the whole phenomena in SCWO reactor. This technique has been applied by several authors to obtain detailed information about velocity, temperature profiles. Different supercritical reactor technologies are already modellized: vessel reactor (Oh et al., 1997), tubular reactor (Dutournie´ et al., 1999; Moussiere et al., 2007), transpiring wall reactor (Bermejo et al., 2010), stirred double shell tubular reactor

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

(Roubaud et al., 2007), fluidized bed (Rodriguez-Rojo and Cocero, 2009). Moussiere et al. (2007) presented a CFD model of SCWO in a tubular reactor, including discussion about kinetics model. An Arrhenius law and Eddy Dissipation Concept (EDC) model were tested. The EDC model performed better in the prediction of hot spots inside the reactor. Sierra-Pallares et al. (2009a) use modified multiple time scale EDC model in case of methanol hydrothermal flame with a good prediction of flame temperature. Narayanan et al. (2008) developed CFD tool in order to gain an understanding of hydrothermal flame phenomena. Generally, these models give a good description of the reactive flow, but an overprediction of the reactional medium temperature between 10% and 15% is observed. CFD simulation is also used to estimate the rate of O2 mass transfer in carbon particle SCWO (Sugiyama et al., 2004), but there is still a temperature overprediction. Several works investigate also the cases of nanoparticles synthesis in supercritical water CFD simulation (Aimable et al., 2009; Sierra-Pallares et al., 2009b, 2010; Tavares Cardoso et al., 2008) considering multiphasic flows. IERs supercritical water oxidation has been experimentally investigated in Leybros et al. (2010a) with degradation conversion higher than 99.9%. The influence of operating parameters (temperature, concentration, stirring speed) is discussed there. The purpose of this work is to provide a mathematical model using CFD to better understand IERs particles solubilization and reactivity in the reactor and also heat and mass transfers between particles and supercritical media. The effect of operating parameters (feed flow rates, stirring speed especially) can be focused.

171

In order to validate CFD calculations, experimental temperature profiles have been obtained in the vicinity of reaction zone.

2. Material and methods 2.1. Experimental set-up—material and methods The process, developed in the Supercritical Fluids and Membranes Laboratory (CEA Marcoule, France) has been described previously in Calzavara et al. (2004) with more details. Fig. 1 shows the flow sheet of the process with the stirred double shell reactor. The external vessel made of 316 stainless steel withstands the working pressure of 30 MPa. Along half of the vessel, four electric heaters are placed, following by a cooling shell. Within the autoclave, sits a titanium tube of 1 m length, 23.6 mm internal diameter and 0.9 mm thickness. The incoming flow of water plus oxidant is flowing into the annular space and hence preheated before entering the reactor, where the waste is fed. This tube is also used to prevent the autoclave from corrosion by confining the aggressive species. A titanium stirrer maintains a turbulent flow along the whole reactor. Geometry of the stirrer head is given Fig. 2. It enhances heat and mass transfer and prevents the precipitated salts from settling in the ‘‘supercritical zone’’ by bringing them into the ‘‘subcritical zone’’ of the reactor where they are dissolved. The suspension of particles is pressurized by a membrane pump. However, one of the most important drawbacks of this

Fig. 1. Description of supercritical water oxidation process.

Fig. 2. Stirrer head and 2D axisymmetric geometry of the double shell stirred reactor.

172

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

system is that pump valves blocked up. In order to avoid this problem, it is necessary to decrease the size of resins particles. For this reason, resins are crushed and sifted. Particle size distribution of the suspension has been measured using a granulometer and is given in Fig. 3. The oxidation reactions take place in the first part of the reactor. In the second part, the flow is cooled by the cooling jackets and by the air/water mixture flowing in the annular space. At the outlet, the effluent is depressurized through a back pressure regulator and separated in two phases. The aqueous phase is analysed by a total organic carbon analyser and the gaseous phase is analysed by a CO, CO2, and O2 gas analyser. The CO2 output

allowed to compute the mass of non-reacted carbon and thereby the degradation conversion according to Rð%Þ ¼

ðmtotal carbon introduced mNon-reacted carbon Þ  100 mtotal carbon introduced

ð1Þ

where Masstotal carbon introduced is calculated from waste and water mass flow rates and MassNon-reacted carbon is deduced from TOC measurements. In order to obtain experimental temperature profiles along the inner reactor wall, the double shell stirred reactor is equipped with 3 type K thermocouples. Owing to the fact that the space between inner and outer walls is very narrow, only three 0.5 mm diameter thermocouples or one 1 mm diameter thermocouple can be used. As shown in Fig. 4, these thermocouples are set at 15; 25; 35 com from the reactor head. Another thermocouple is put directly inside the waste injector tube (the measured temperature is named injection temperature or Tinj). Last, a thermocouple is put at the reactor outlet and named Toutlet. All experimental conditions are summarized Table 1. The experimental procedure is explained in details in Leybros et al. (2010a). In this study, only Amberlite IRN77 cationic resins are used.

2.2. Development of simulation model

Fig. 3. Particle size distribution of resins suspension. Experimental measurements and Rosin Rammler model.

The modelling must take into account interactions between coupled phenomena, such as supercritical water-IERs particles diphasic flow dynamics, mass and heat transfers and chemical reactions. The simultaneous numerical simulation of the equations for continuity, momentum, turbulence and energy and chemical

Fig. 4. Location of thermocouples along the reactor.

Table 1 Summary of the calculation and experiments operating conditions. Operating parameters

Run 1

Run 2

Run 3

Run 4

Run 5

Run 6

Run 7

Run 8

Run 9

Qwater (g h  1) Qwaste (g h  1) Concentration (%) Qair (NL h  1) TPR1 (K) TPR2 (K) TPR3 (K) TPR4 (K) T_15 cm (K) T_25 cm (K) T_35 cm (K) Toutlet (K) TOC (mg L  1) Mean degradation rate (%)

1001.5 166.5 14.25 1075.5 676 706 673 573

1186.5 183.7 13.4 1225.1 676 708 673 573

1149.4 237.3 17.1 1175.2 681 709 673 573

2194.9 172.3 7.3 1331.9 668 705 673 573

1008 197.3 16.4 861 705 716 673 573

3287 176.6 5.1 1344 687 719 673 573

2177.1 176 7.5 1495.7 697 716 673 573

1372.2 177.4 11.4 1675.2 698 713 673 573 654.6

1150 161.1 12.3 1215 700 708 673 573 662.4

600.8 605.7 375 107.2 753.5 98.92

630.6 352.8 153.2 710 98.17

623.7 382.7 406.47 49.3 96.46

600.8 377.3 142.7 741.8 96.98

610.9 613.1 361.2 147.9 768.6 98.73

601.3 606.6 377.2 200.77 57 94.46

600.7 609.1 385.5 399.8 7 35.9 91.62

376.3 147.9 741.7 97.77

356.3 326.5 7 175.8 96.79

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

species is required. The simulation was performed with Fluents 6.3.26. 2.2.1. Meshing The preliminary step is to represent the simplified 2D axisymmetric geometry of the reactor. It is divided the reactor in several cells in which all the equations will be solved. The mesh is generated by the Gambits software. This bi-dimensional geometry was used in order to limit computational time. Its efficiency to describe the double shell reactor has been checked in a previous study (Moussiere, 2006; Roubaud et al., 2007). In this study, comparisons were made between two-dimensional and three-dimensional velocity profiles and heat transfer calculations. Regarding the velocity field prediction in the stirrer head, the three dimensional Sliding Mesh Method gives more realistic details than two-dimensional axisymmetric swirl geometry. Predicted temperatures are of the same order in both cases. Taking into account the results of this study, calculations have been done in 2D, less computer time consuming with a sufficient mixing representation than 3D geometry. The geometry of double shell stirred reactor is shown in Fig. 2. The built mesh contains about 230,000 cells. Their distribution is refined in the first part of the reactor where a large temperature variation is expected. A swirl velocity is imposed to the fluid in order to represent the rotation movement of the stirrer. Swirl velocity uz is calculated with the following equation considering the y cell position uz ¼ or  2p  9y9

ð2Þ

where or is the angular velocity. 2.2.2. Thermodynamics and material properties 2.2.2.1. Fluid properties. Although a supercritical fluid is considered like a compressible fluid, the system is considered as incompressible fluid in order to develop CFD simulation. All reactions and fluid variations are considered to be isobaric at 30 MPa. Fluid properties at 30 MPa were specified in Fluents. A cubic Peng–Robinson equation of state with constant volume translation was used (Eq. (3)). It shows a good agreement with experimental data with an average deviation between experimental and calculated values lower than 2% (Farrokh-Niae et al., 2008). This equation links density to pressure and temperature P¼

RT ðv þ abÞ a  v ðvbÞ vðv þcÞ

ð3Þ

where a ¼ AC AxR2 T 2C =P 2C b ¼ BxRT C =P 2C c ¼ C xRT C =P2C Fluid species properties are implemented as tabulated data functions of temperature at 30 MPa. An ideal mixing is assumed. Density (Eq. (4)), viscosity (Eq. (5)), thermal conductivity (Eq. (6)) and specific heat (Eq. (7)) of mixing are calculated with the following mixing rules for the mixture: 1 ð4Þ rm ¼ P ðyi =ri Þ X mm ¼ ðyi  mi Þ ð5Þ X lm ¼ ðyi  li Þ ð6Þ X Cpm ¼ ðyi  Cpi Þ ð7Þ

2.2.2.2. Particles properties. There is little published information concerning the properties of Ion Exchange Resins in a supercritical

173

environment. Particle properties (heat capacity, fusion enthalpy and temperature, molar volume, thermal conductivity and diffusivity, vapourization enthalpy and temperature) are obtained using the Van Krevelen group contribution method for polymers (Van Krevelen, 1990). This model is based on the hypothesis that polymer thermodynamic properties can be estimated by adding increments. It is valid only if a group forming a part of the monomer is not affected by the other groups. For example, in the case of polystyrene, this model gives results, in accordance with measurements, with a deviation lower than 5%. Particles suspensions granulometry is taken into account thanks to the empirical Rosin Rammler experimental law. This model assumes that particles diameter d is exponentially linked with diameter higher than d particles mass fraction Yd. The law is expressed by "  3:1 # d Y d ¼ exp  ð8Þ 176:9 where dmean ¼176.9 mm is the mean diameter of particle size distribution and n¼3.1 is the spread particle size distribution parameter. Fig. 3 shows that Rosin Rammler model agrees with experimental measurements. 2.2.3. Hydrodynamics In our model, fluid phase is treated in an Eulerian framework and the particulate phase is treated in a Lagrangian framework. The supercritical water phase is treated as a continuum by solving the Navier–Stokes equations, whereas the particulate phase, considered as secondary or dispersed phase, is solved by tracking a large number of particles through the calculated flow field. The solid phase is characterized by its diameter and its properties described earlier. Modelling approaches are in particular described in Gouesbet and Berlemont (1999). 2.2.3.1. Fluid phase. Continuity equation for fluid phase is @r ! þ r:ðr u Þ ¼ 0 @t

ð9Þ

The momentum balance yields @ ! !! ! ! ðr u Þ þ rUðr u : u Þ ¼ rp þ rUt þ r g þ F ð10Þ @t ! where r is the fluid density, u is velocity, p is the static pressure, ! t is the stress strain tensor, r! g and F are the gravitational body force and the external body forces, respectively. We use an energy (E) balance to take into account the effect of heat transfer on hydrodynamics. 0 1 X ! @ ! !A @ ðrEÞ þ rUð v ðrEþ pÞÞ ¼ rU kef f rT hj J j þ ðtef f U v Þ þ Q @t j ð11Þ ! where keff is the effective conductivity and J j is the diffusion flux of species j and Q is the heat source term. The terms on right hand side represent energy transfer due to conduction, species diffusion and viscous dissipation. For turbulence modelling, shear-stress transport k o model is used. This model was developed by Menter (1994) to blend the robust and accurate formulation of k–o model in the near wall region with the free stream independence of the k–e model in the far field. In this model, the turbulence kinetic energy k and its specific dissipation rate o are obtained from the following transport equations:   @ @ @ @k ðrkÞ þ ðrkui Þ ¼ Gk ð12Þ þ Gk Y k @t @xi @xj @xj

174

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

and   @ @ @ @o ðroÞ þ þ Go Y o þ Do ðroui Þ ¼ Go @t @xi @xj @xj

ð13Þ

In these equations, Gk represents the generation of turbulence kinetic energy due to the mean velocity gradients, Go is the generation of o, Gk and Go represent the effective diffusivity of k and o respectively. Yk and Yo represent the dissipation of k and o due to turbulence. 2.2.3.2. Particles. The trajectory of a discrete phase particle is predicted by integrating the force balance on the particle, which is written in a Lagrangian reference frame. This force balance can be written as duP g ðr rÞ ¼ F D ðuuP Þ þ x P þ Fx dt rP

ð14Þ

where Fx is an additional acceleration term, r and rP are, respectively, fluid and particle densities, u and uP are, respectively, fluid and particle velocities, FD(u up) is the drag force per unit particle mass. The dispersion of particles due to turbulence in the fluid phase is predicted using stochastic tracking model. It includes the effect of instantaneous turbulent velocity fluctuations on the particle trajectories through stochastic models such as Discrete Random Walk (DRW) model (also known as ‘‘Eddy Interaction Model’’ in Gosman and Ioannides, 1983). In this stochastic approach, the solver predicts the turbulent dispersion of particles by integrating trajectory equations for individual particles using the instantaneous fluid velocity along the particle path. In the DRW model, the fluctuating velocity components are discrete piecewise functions of time. By computing the trajectory in this manner for a sufficient number of representative particles, the random effects of turbulence may be accounted for. At each iteration, fluid phase properties are used in order to determine boundary conditions for particles. Momentum F and heat exchange Q source terms between particles and fluids equations are given as ( ) X 18mC D Re _ P Dt F¼ ðu uÞ þ F ð15Þ P other m 2 24rP dp " Q¼

mP DmP cP DT P þ mP,0 mP,0

Z hf g þhsolub þ

!#

TP T ref

cp,i dT

_ P,0 m

2.2.4. Chemistry model IERs supercritical water oxidation mechanism has been studied in Leybros et al. (2010b). Benzoic acid, phenol and acetic acid were identified as rate determining species. As a consequence, a three steps mechanism has been created and could be implemented in our model. The model includes two successive steps. The first step is the solubilization of the polymer and consequently the formation of hydroxybenzoic acid. The second step is the homogeneous reaction described by a three-step mechanism (see Fig. 5). The species defined in the modelling are C7H6O3, C6H6O, C2H4O2, O2, N2, CO2, and H2O. All the kinetic laws results from literature. Kwak et al. (2006) proposed a solubilization law of polystyrene in supercritical water, that will be assumed close to IER solubilization (Reaction (1)). Gonzalez et al. (2004) modelled decarboxylation law in supercritical water (Reaction (2)). Finally, Yermakova et al. (2006) and Li et al. (1991) proposed a two-step mechanism for the phenol supercritical water oxidation with formation of acetic acid (Reactions (3) and (4)). A summary of the kinetic data, used for this study, is provided in Table 2. For turbulence–chemistry interaction, the Eddy Dissipation Concept model is used (Magnussen, 1981). This model was developed for prediction of gaseous combustion reactions in turbulent flows. We chose this model considering that Schmidt number in supercritical medium is close to gas phase Schmidt number value. The mean reaction rate assumes that chemical reactions occur in small turbulent structures, so called fine scales, where reactants are mixed at a molecular scale. The mass fraction filled by fine structures is modelled as 

gn ¼ 2:1377

ne

3=4 ð18Þ

2

k

where n denotes the fine scale quantities and n is the kinematic viscosity.

ð16Þ

where CD is the drag coefficient, cp,i is the heat capacity, dP is the particle diameter, hfg and hsolub are, respectively, latent heat and solubilization heat, mP and mP,0 are, respectively, the particle _ P is the particle mass flow mass and the initial particle mass, m rate, TP is the particle temperature and Re is Reynolds number. In order to take into account particles solubilization, a single kinetic rate model is used and assumes that the rate is first order dependent on the amount of species remaining in the particles. The model is applied to a reactive particle when the temperature of the particle reaches a user defined solubilization temperature. This law is expressed by dmP ¼ kA ½mP ð1f v,0 Þð1f w,0 Þmp,0  dt

where k is the kinetic rate expressed thanks to an Arrhenius law, fv,0 is the mass fraction of species initially present in the particle and fw,0 is the mass fraction of ‘‘solubilizing’’ material.

ð17Þ

Table 2 Kinetic parameters used for Eddy Dissipation Concept Model. Reaction

Pre-exponential constant (s  1 mol(1  a) L(1  a))

Reaction 1 (Kwak 2.73  1010 et al., 2006) 8.886  106 Reaction 2 (Gonzalez et al., 2004) 4.903  108 Reaction 3 (Yermakova et al., 2006) Reaction 4 (Li et al., 2.55  1011 1991)

Fig. 5. IER degradation mechanism in supercritical water.

Activation energy (kJ mol  1)

Reaction order

1.57  103

1

95

1

124.78

1,041

172.7

1

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

Species are assumed to react in the fine structures over a time scale n1=2 tn ¼ 0:4082 ð19Þ

e

The mass transfer between the fine structures and the surrounding is expressed as  1=4 ne e _ ¼ 0:1909 ð20Þ m 2 k k The reaction rate is defined as _ r~ i ¼ m

w 1gn w

C min

175

Operating conditions were chosen in order to get degradation efficiency variations. Indeed, run 7 degradation conversion is around 91.62% whereas run 1 degradation conversion oscillates around 98.92%. Runs 1 to 7 were carried out with two 0.5 mm diameter thermocouples sat at 25 and 35 cm from reactor head. Runs 8 and 9 were carried out with only a 1 mm thermocouple sat at 15 cm in order to have temperature measurements closer of the stirrer head located between 0 and 10 cm from the reactor head (see Fig. 4). TPR1, TPR2, TPR3, TPR4 are temperatures reached in every heater. Stirring speed is set to 300 rpm except for run 2, where it is set to 600 rpm.

ð21Þ

where w is the fraction of fine structure which reacts and C min is the smallest of local mean concentrations of waste and oxygen. 2.2.5. Boundary conditions For the inlet boundary condition, a simple mass flow inlet is chosen and allows to specify mass flow, hydraulic diameter, turbulent intensity and species mass fractions. Temperatures of electric heaters and cooling shell are set to the walls. For the outlet zone, the basic outflow boundary condition, defined in Fluents, was chosen. Particles suspension is characterized by its granulometry (mean diameter, minimal and maximal diameters and spread parameter) and by a shrinking coefficient. This coefficient is the ratio between final and initial diameters. Owing to the fact that this coefficient can’t be infinite, the particle is entirely solubilized once final diameter is equal to 0.1% of initial diameter. 2.2.6. Numerical procedure The model is made up of a system of partial differential equations (Eqs. (1)–(21)). Those equations have been solved numerically with the code Fluents 6.3.26 using a finite volume numerical method. The differential equations are converted into finite differences equations in each node point and the generated system of ordinary equations that is solved by an algorithm, either implicitly or explicitly. In this work, the second order implicit solver was used, which usually grants higher stability. The second order upwind discretization schemes were chosen to solve the model equations. Pressure velocity coupling was achieved using the SIMPLE algorithm. The solution was convergent when the scaled residuals were less than 10  4, except for energy where it is 10  7. The simulation strategy is the following. First, the steady state simulation of water-air flow without particles is performed. Once convergence is reached, unsteady flow involving particles simulation is conducted, using the steady state flow solution as initial conditions. A time step between 0.0085 and 0.011 s with 20 iterations per time step was chosen. This iteration was adequate to achieve convergence for the majority of time steps. This time step was set as a function of waste (i.e. {waterþparticles} mixing) inlet flow. For each case, 90 s of reactive flows are calculated, which represents around 2 weeks of calculation on 4 processors.

3. Results and discussion

3.1.2. Results All the results are given in Table 3: temperature at 15 cm from reactor head (Table 3a), at 25 cm (Table 3b), at 35 cm (Table 3c) and at the outlet (Table 3d). Afterwards, they will be named T_15 cm, T_25 cm, T_35 cm and T_outlet. Due to the weak IERs oxidation reaction enthalpy, only a weak temperature variation is noticed before and after resins injection. Temperature gaps between T_25 cm and T_35 cm are around 10 K for Runs 1, 5, 6 and 7. No temperature peak can be linked to resins oxidation heat release. This is the reason why two runs were carried out with a thermocouple nearer of the assumed reaction zone. More significant temperature variation differences were noticed for T_15 cm. Temperature deviations around 20 K indicate that reaction occurs. The highest temperature is around 662.477 K (for Run 9). These results confirm that reaction zone should be in the vicinity of stirrer head. Temperature measurements variations can be discussed according to operating conditions (waste and water-air mass flow, stirring speed). Comparing Run 1 and Run 2 where only stirring speed is different, a 25 K temperature deviation for T_35 cm is noticeable. This result could be explained by a change in reaction zone location or by a different heat homogenization. In the same way, waste mass flow rate variations (Run 1 and Run 3) induce a 20 K temperature deviation for T_35 cm. On the contrary, even if air–water mass flow rate variations influence degradation efficiency, no variations are remarked. For example, between Run 5 (highest waste concentration) and Run 6 (lowest waste concentration), only a 10 K temperature deviation on T_35 cm is noticeable. 3.1.3. Conclusion Experimental temperature profiles measurements allow to have trends concerning the reaction zone evolution as a function of operating conditions. Nevertheless, due to the difficult thermocouples setting up in supercritical water medium, no more information could be obtained in the vicinity of stirrer head. CFD simulation will provide us more information on reaction location and on the Table 3a Comparison between model and experiments for T_15 cm. Run

Run 8

Run 9

Measured temperature (K) Predicted temperature (K) Standard deviation (%) (K)

654.6 714.7 8.4

662.4 712.1 7.0

3.1. Experimental temperature profiles 3.1.1. Experimental conditions Operating conditions and experimental temperature measurements are summarized in Table 1. IERs concentration suspension is equal to 10% wt. Several runs were necessary due to thermocouples fragility in the inner space between walls and harsh operating conditions. Few temperature profiles are not complete owing to frequent damages on thermocouples during reactor preheating.

Table 3b Comparison between model and experiments for T_25 cm. Run

Run 1

Run 5

Run 6

Run 7

Measured temperature (K) Predicted temperature (K) Standard deviation (%) (K)

600.8 711.5 15.5

610.9 720.2 15.2

601.3 717 16.1

600.7 719.8 16.5

176

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

Table 3c Comparison between model and experiments for T_35 cm. Run

Run 1

Run 2

Run 3

Run 4

Run 5

Run 6

Run 7

Measured temperature (K) Predicted temperature (K) Standard deviation (%) (K)

605.7 685.0 11.6

630.6 682.3 7.6

623.7 676.6 7.8

600.8 681.9 11.9

613.1 683.2 10.3

606.6 673.6 10.8

609.1 682.0 10.7

Table 3d Comparison between model and experiments for outlet temperature. Run

Run 1

Run 2

Run 3

Run 4

Run 5

Run 6

Run 7

Run 8

Run 9

Measured temperature (K) Predicted temperature (K) Standard deviation (%) (K)

375 342.1 8.7

352.8 341.2 3.3

382.7 345.5 9.7

377.3 345.5 8.4

361.2 345.8 4.3

377.2 350.8 6.9

385.5 349.8 9.3

376.3 349.9 7.0

356.3 346.2 2.8

Fig. 6. Comparison between experimental temperature profiles and calculated profiles in the annular space and in its vicinity (Run 4).

homogenization of supercritical medium. In this way, experimental temperature profiles will help us by validating mathematical model. 3.2. CFD model validation 3.2.1. Temperature In order to validate our model, comparisons were made between measured temperatures T_15 cm; T_25 cm; T_35 cm; T_outlet and calculated values with Fluents software. All results are summarized in Tables 3a, 3b, 3c and 3d, respectively. Fig. 6 gives mean measured temperature and different calculated longitudinal profiles at several locations (outer titanium wall, in the middle of annular space and inner stainless steel wall) for Run 4. First, experimental and calculated temperature profiles are qualitatively in agreement. Indeed, T_15 cm measurement corresponds to calculated temperature peak. Nevertheless, if we consider the whole results, measured temperatures are lower than predicted values. Calculated temperatures are in good agreement for T_15 cm and T_35 cm (deviation between model and experience respectively of 9% and 12%). Concerning T_25 cm, greater deviations from 15% to 17%, are noticeable. Calculated temperatures at the reactor outlet are in agreement with measurements: deviations around 15–30 K are noticeable depending on operating conditions. 3.2.2. Degradation conversion Degradation conversions are also calculated considering outlet species concentration in the model and compared with experimental ones. Results are given in Table 4. In every case, calculated

degradation conversions are lower than experimental ones. Nevertheless, by noticing in each case deviations between experimental and calculated degradation rates lower than 3–3.5% (except for Run 6), we can consider than the model is in agreement with experimental measurements. In order to explain the underestimation of degradation conversions, we have to consider the bidimensional axisymmetric geometry choice (Fig. 2). Two-dimensional geometry representation is a simplified one because the real reactor has not any symmetry axis. As a consequence stirrer (and more particularly stirrer head) representation is simplified. Stirrer head is partially designed. Only the three passages are drawn as hollow parts. The two planes are ignored. Due to this geometry simplification, differences between real and model particles residence time inside the stirrer head are noticeable. Hence, we can assume that model residence time is lower than real time. As a consequence, it may explain that calculated degradation conversions are lower than measured ones. Three-dimensional geometry calculations could bring us clarifications. Unfortunately, we have not enough CPU power to make these calculations due to 1 million three-dimensional mesh (Roubaud et al., 2007). These results highlight the main limitations of the model. First, air–water flow in the annular space has been considered monophasic, although in reality a biphasic gas–liquid flow exists until it reaches critical temperature. Furthermore, the calculation of fluid mixture properties as mass average of pure component ones could also be an explanation. Due to these simplifications, model overpredicts temperatures by 50–100 K depending on operating conditions. In the literature, this problem was identified for supercritical water oxidation for fluids (Sierra-Pallares et al., 2009a; Narayanan et al., 2008), especially in stirred double shell reactor (Roubaud et al., 2007). Another explanation is to consider that simplification geometry can also influence CFD simulation results. Three-dimensional calculations could improve obtained results. 3.3. Concentration and temperature profile CFD simulations Fig. 7a shows temperature fields inside the reactor. The temperature distribution picture shows a higher temperature zone (T close to 1200 K) inside the stirrer head. This temperature increase is obviously linked with the heat released by resins oxidation. Water/air and waste flow rates must be weak to increase particles residence time in the reactional area, i.e. stirrer head and its vicinity. Furthermore, this combustion simulation provides an internal view of the consumption of IER slow-oxidizing intermediate degradation products (hydroxybenzoic acid (Fig. 7b), phenol (Fig. 7c) and acetic acid (Fig. 7d)) consumption and of the CO2 production (Fig. 7e).

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

177

Table 4 Comparison between the calculated and the experimental degradation rates and total organic carbon. Run

Run 1

Run 2

Run 3

Run 4

Run 5

Run 6

Run 7

Run 8

Run 9

TOC measured (mg L  1) Experimental rate (%)

107.2 753.5 98.92

153.2 7 10 98.17

406.4 749.3 96.46

142.7 7 41.8 96.98

147.9 7 68.6 98.73

200.77 57 94.46

399.8 7 35.9 91.62

147.9 741.7 97.77

326.5 7 175.8 96.79

TOC calculated (mg L  1) Calculated rate (%)

303.6 96.74

237.2 97.26

617.4 94.66

224.4 94.89

301.9 97.24

429.3 85.73

468.5 89.65

407.7 94.4

504.2 93.57

Fig. 7. Concentration and temperature profiles along the double shell stirred reactor (Operating conditions: Run 4): (a) temperature (K); (b) C7H6O3 mass fraction; (c) C6H5OH mass fraction; (d) C2H4O2 mass fraction; (e) CO2 mass fraction; and (f) mass of discrete particles (kg).

The reactive zone is mainly located in the vicinity of stirrer head. Resins react fastly into hydroxybenzoic acid, and then into phenol. Acetic acid is slower destroyed after the stirrer head and partially converted to carbon dioxide. These concentrations profiles are used to

estimate remaining species concentration at reactor outlet, which give us calculated degradation conversion. In Fig. 7f, discrete particles trajectories are pointed out along the reactor, and more particularly inside the stirrer head. On this

178

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

figure, the evolution of mass particles is clearly noticeable and is linked to species concentration profiles. More particularly, the bottom of stirrer head is the place where particles accumulate and react. Once particles have reacted, they follow the fluid phase flow, up to reactor outlet. 3.4. Influence of stirring speed 3.4.1. Experimental In order to understand influence of stirring speed on degradation efficiency, experiments were carried out, by considering only stirring speed as a variable. Stirring speed varies from 300 to 750 rpm. Corresponding Reynolds numbers are around 105. As a consequence, flow can be considered as turbulent. Operating conditions are given in Table 5. Total organic carbon measurements are closed for Run 10a and 10b. For Run 10c, degradation conversion is significantly lower. Modifying residence time has low influence for stirring speeds between 300 and 600 rpm, but degradation conversion decreases significantly for rate equal to 750 rpm. 3.4.2. Calculations Temperature fields were calculated for different stirring speed values (300, 600, 750 rpm) with the following operating conditions: Qwater ¼1200 g h  1; Qair ¼1100 NL h  1 and Qwaste ¼180 g h  1. These temperature fields are shown in Fig. 8. In every case, a temperature peak is noticeable in the middle of stirrer head. However, temperature deviations are noticeable after the stirrer head, according to operating conditions. We remark heat of reaction release in the vicinity of walls, but temperature increases location differs according to stirring speed. As discussed earlier, stirrer head corresponds to a zone of particles accumulation. Increasing stirring speed allows a better homogenization of particles distribution in order to prevent the formation of local hot spots. However, as a result, heat transfers are improved and particles residence time in the reaction zone decrease. As a consequence, a local decrease in temperature is achieved after the stirrer head, as shown in Fig. 8. Reaction zones, corresponding to the location of reaction heat release have trends to broaden by increasing stirring speed. This induces a reduction of particles degradation efficiency. Experimental results are thus confirmed (see Table 5).

Fig. 8. Calculated temperature (only representation of reaction heat release) profiles for different stirring speeds (Operating conditions 10a: 300 rpm; 10b: 600 rpm; 10c: 750 rpm) inside the stirrer head.

3.5. Influence of geometry Computational Fluid Dynamics can be an useful tool to design a reactor with an optimized geometry. The prediction of hydrodynamics behaviour for modified geometry reactors with equal flow rates conditions allows to choose the best configuration. Calculations have been carried out by making an extrapolation on reactor length and diameter. If we consider former results, the extrapolation aim is to observe particles residence time variations in reactive zone by broadening the stirrer head. Flow rates were identical as Run 4 by considering volume changes. Table 5 Influence of stirring speed on degradation efficiency.

Stirring speed (rpm) Qwater (g h  1) Qwaste (g h  1) Qair (NL h  1) TOC measured (mg L  1) Experimental rate TOC calculated (mg L  1) Calculated rate

Run 10a

Run10b

Run 10c

300 1200 180 1100 107.2 753.5 98.92 303.6 96.74

600 1200 180 1100 153.2 710 98.17 237.2 97.26

750 1200 180 1100 414.6 7153 95.28 665.5 92.1

Fig. 9. Calculated temperature profiles (in the middle of annular space) for different reactor geometries HxLy (Operating conditions Run 4). H corresponds to diameter and L to reactor length with their respective extrapolation factors.

Calculated temperature profiles are given in Fig. 9 in the middle of annular space. In order to make easier comparison between different geometries, we have taken into account extrapolation factor for length and diameter reactor and considered volume changes for mass flow definitions. Profiles are named HxLy where x and y are, respectively, extrapolation factors for diameter and length reactor. Significant differences are noticeable according to different geometries. According to stirrer head volume changes, particles suspensions seem to be better spread. It implies that temperature peaks (with temperature greater than 1000 K) are not present in every case. However, a temperature peak is still noticeable in H1.5L1 case, but with a lower value. Consequently, particles have fewer trends to accumulate in this case and it indicates that a geometry stirrer head change will be suitable in case of particles SCWO in our process. In H1.25L1.5 case, reaction heat release appears after the stirrer head.

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

As a conclusion, reactor diameter and stirrer head geometry are the most influent parameters concerning thermal balance in double shell stirred reactor. CFD model brings interesting information on temperature profiles and on hot spots appearance. In this way, it allows to calculate an optimized reactor geometry while limiting number of experiments.

4. Conclusion In this study, a CFD modelling of a diphasic reactive turbulent flow in a double shell stirred reactor with a two dimensional axisymmetric representation and a comparison with experimental measurements has been achieved. Particle trajectories were modelled with the Discrete Random Walk Particle model. A mechanism has been proposed and implemented for the solubilization of particles in the reaction medium. Location of reactive zones was established: there is a hot spot inside the stirrer head corresponding to the place where resins particles accumulate. Calculated degradation conversions are in agreement with measurements. The simulated temperature profiles using Eddy Dissipation Concept model are generally in agreement with the experimental measurements even if there is still an average overestimation of 10-15%. This overestimation has been also observed in other similar supercritical water turbulent flow CFD simulation works (Roubaud et al., 2007; Narayanan et al., 2008). In order to improve the prediction, EDC-MTS model (Hjertager et al., 2002; Baldyga, 1989) could be appropriate for supercritical water oxidation, if we consider the Sierra-Pallares et al. (2009a) results. This will help us to identify predominant mechanism of reaction control in supercritical water. Thanks to this CFD model validated by comparison with experimental data obtained with our 200 g h  1 pilote POSCEA2, we are able to predict temperatures and concentration profiles and oxidation location in the reactor thanks to CFD simulations. This simulation brings us information on particles homogenization in the reactional medium. Different stirrer head geometries could be chosen with this model in order to change particles residence time in the reactor. The definition of a new geometry could help us to improve degradation efficiency and heat transfer handling to avoid hot spots. The next step is the validation of the feasibility of 1 kg h  1 treatment capacity. The same experimental conditions must be applied in terms of organic weight percent, equivalent air to fuel stoichiometric ratio and wall temperatures. This could be done on non radioactive organic matter with the 1 kg h  1 capacity treatment DELIS pilot. Collation of experimental data with numerical simulation of this reactor will be done to fully validate the design before use on real radioactive contaminated solvent.

Gk, Go H hfg hsolub ! Jj k kA keff L m _P m n P Q Re ri T Dt TOC u v Yd Yk, Yo

179

generation of k and o due to the mean velocity gradients reactor diameter (m) latent heat (J kg  1) solubilization heat (J kg  1) diffusion flux of species j (kg m  2 s  1) turbulence kinetic energy kinetic rate (mass unit varies according to kinetic law) effective conductivity (W m  1 K  1) length reactor (m) mass (kg) mass (kg) spread particle size distribution parameter pressure (Pa) mass flow rate (g h  1) Reynolds number reaction rate (molar unit varies according to kinetic law) temperature (K) time step (s) total organic carbon (mg L  1) velocity (m s  1) molar volume (m3 mol  1) mass fraction of particles of diameter greater than d dissipation of k and o due to turbulence

Greek letters

e w gn Gk, Go

m n r r! g t tn o or

turbulent dissipation rate (m2 s  3) ratio of the reactive fine structures in EDC model mass fraction occupied by the fine structures in EDC model effective diffusivity of k and o dynamic viscosity (kg m  1 s  1) kinematic viscosity (m2 s  1) fluid density (kg m  3) gravitational body force (N) stress strain tensor (Pa) time scale (s) specific dissipation rate (s  1) angular velocity (rpm)

Subscripts n

C P x, y

fine scale quantities critical particle extrapolation factors for reactor diameter and length

Nomenclature References a,b,c, A,B,C, x coefficients in translated Peng Robinson equation of state CD drag coefficient C min smallest of local mean concentrations of waste and oxidant Cp specific heat (J kg  1 K  1) d particle diameter (m) dmean mean diameter of particle size distribution (m) E energy (J) ! F external body force (N) FD(u  uP) drag force per unit particle mass (N) fv,0 mass fraction of species initially present in the particle fw,0 mass fraction of ‘‘solubilizing material’’ Fx additional acceleration term (N) g gravity acceleration¼9.81 m s  2

IAEA, 2002. Application of Ion Exchange Processes for the Treatment of Radioactive Waste and Management of Spent Ion Exchangers. International Atomic Energy Agency, Vienna. Aimable, A., Muhr, H., Gentric, C., Bernard, F., Le Cras, F., Aymes, D., 2009. Continuous hydrothermal synthesis of inorganic nanopowders in supercritical water: towards a better control of the process. Powder Technol. 190, 99–106. Baldyga, J., 1989. Turbulent mixer model with application to homogeneous instantaneous chemical reactions. Chem. Eng. Sci. 44, 1175–1182. Bermejo, M.D., Cocero, M.J., 2006. Supercritical water oxidation: a technical review. Am. Inst. Chem. Eng. J. 52, 3933–3951. Bermejo, M.D., Fernandez-Polanco, F., Cocero, M.J., 2006. Effect of the transpiring wall on the behavior of a supercritical water oxidation reactor: modeling and experimental results. Ind. Eng. Chem. Res. 45, 34–38. Bermejo, M.D., Martin, A., Queiroz, J., Bielsa, I., Rios, V., Cocero, M.J., 2010. Computational fluid dynamics simulation of a transpiring wall reactor for supercritical water oxidation. Chem. Eng. J. 158, 431–440. Brunner, G., 2009. Near and supercritical water. Part II: oxidative processes. J. Supercritical Fluids 47, 382–390.

180

A. Leybros et al. / Chemical Engineering Science 69 (2012) 170–180

Calzavara, Y., Joussot-Dubien, C., Turc, H.-A., Fauvel, E., Sarrade, S., 2004. A new reactor concept for hydrothermal oxidation. J. Supercritical Fluids 31, 195–206. Cocero, M.J., 2001. Supercritical water oxidation (SCWO). Application to industrial wastewater treatment. in: Bertucco, A., Vetter, G. (Eds.), High Pressure Process Technology: Fundamentals and Applications, Elsevier, Amsterdam, pp. 509–526. Dutournie´, P., Mercadier, J., Cezac, P., 1999. Simulation of a tubular reactor for supercritical water oxidation. In: Proceedings of the Second European Congress of Chemical Engineering, Montpellier, France. Farrokh-Niae, A., Moddaress, H., Mohsen-Nia, M., 2008. A three parameter cubic equation of state for prediction of thermodynamic properties of fluids. J. Chem. Thermodyn. 40, 84–95. Gonzalez, G., Salvado, J., Montane, D., 2004. Reactions of vanillic acid in sub and supercritical water. J. Supercritical Fluids 31, 57–66. Gosman, A., Ioannides, E., 1983. Aspects of computer simulation of liquid-fueled combustors. J. Energy 7, 482–490. Gouesbet, G., Berlemont, A., 1999. Eulerian and Lagrangian approaches for predicting the behaviour of discrete particles in turbulent flows. Prog. Energy Combust. Sci. 25, 133–159. Hjertager, L., Hjertager, B., Solberg, T., 2002. CFD modelling of fast chemical reactions in turbulent liquid flows. Comput. Chem. Eng. 26, 507–515. Joussot-Dubien, C., Turc, H.-A., Didier, G., 2000. Proce´de´ et Dispositif pour l’oxydation en eau supercritique de matie res. Patent FR no. 00/12929. CEA, France. Kwak, H., Shin, H., Bae, S., Kumazawa, H., 2006. Characteristics and kinetics of degradation of polystyrene in supercritical water. J. Appl. Polym. Sci. 101, 695–700. Leybros, A., Roubaud, A., Guichardon, P., Boutin, O., 2010a. Ion Exchange Resins destruction in a stirred supercritical water oxidation reactor. J. Supercritical Fluids 51, 369–375. Leybros, A., Roubaud, A., Guichardon, P., Boutin, O., 2010b. Supercritical water oxidation of Ion Exchange Resins: Degradation mechanisms. Proc. Saf. Environ. Prot. 88, 213–222. Li, L., Chen, P., Gloyna, E., 1991. Generalized kinetic model for wet oxidation of organic compounds. Am. Inst. Chem. Eng. J. 37, 1687–1697. Magnussen, F., 1981. On the structure of turbulence and a generalized eddy dissipation concept for chemical reaction in turbulent flow. In: Proceedings of the 19th American Institute of Aeronautics and Astronautics Aerospace Science Meeting, St. Louis, USA. Menter, F., 1994. Two-equation eddy viscosity turbulence models for engineering applications. Am. Inst. Aeronaut. Astronaut. J. 32, 1598–1605.

Moussiere, S., 2006. Etude par simulation nume´rique des e´coulements turbulents re´actifs dans les re´acteurs d’oxydation hydrothermale: Application a un re´acteur agite´ double enveloppe. Ph.D. Thesis. Aix Marseille University. Moussiere, S., Joussot-Dubien, C., Guichardon, P., Boutin, O., Turc, H.-A., Roubaud, A., Fournel, B., 2007. Modelling of heat transfer and hydrodynamics with two kinetic approaches during supercritical water process. J. Supercritical Fluids 43, 324–332. Narayanan, C., Frouzakis, C., Boulouchos, K., Prikopsky, K., Wellig, B., Rudolf von Rohr, P., 2008. Numerical modelling of a supercritical water oxidation reactor containing a hydrothermal flame. J. Supercritical Fluids 46, 149–155. Oh, C., Kochan, R., Beller, J., 1997. Numerical analysis and data comparison of a supercritical water oxidation reactor. Am. Inst. Chem. Eng. J. 43, 1627–1636. Rodriguez-Rojo, S., Cocero, M.-J., 2009. Supercritical fluidized bed modeling. J. Supercritical Fluids 50, 54–60. Roubaud, A., Moussiere, S., Fournel, B., 2007. 2D and 3D numerical modelling of the reactive turbulent flow in a double shell supercritical water oxidation reactor. In: Proceedings of the Third European Combustion Meeting, Chania, Greece. Sierra-Pallares, J., Alonso, E., Montequi, I., Cocero, M.J., 2009a. Particle diameter prediction in supercritical nanoparticle synthesis using three dimensional CFD simulations: validation for anatase titanium dioxide production. Chem. Eng. Sci. 64, 3051–3059. Sierra-Pallares, J., Parra-Santos, M., Garcia-Serna, J., Castro, F., Cocero, M.J., 2009b. Numerical modelling of hydrothermal flames. Micromixing effects over turbulent reaction rates. J. Supercritical Fluids 50, 146–154. Sierra-Pallares, J., Marchisio, D., Alonso, E., Parra-Santos, M., Castro, F., Cocero, M.J., 2010. Quantification of mixing efficiency in turbulent supercritical water hydrothermal reactors. Chem. Eng. Sci.. doi:10.1016/j.ces.2010.12.045. Sugiyama, M., Kataoka, M., Ohmura, H., Fujiwara, H., Koda, S., 2004. Oxidation of carbon particles in supercritical water: rate and mechanism. Ind. Eng. Chem. Res. 42, 690–699. Tavares Cardoso, M., Cabral, J., Palavra, A., Geraldes, V., 2008. CFD analysis of supercritical antisolvent (SAS) micronization of microcycline hydrochloride. J. Supercritical Fluids 47, 247–258. Van Krevelen, D., 1990. Properties of Polymers: Their Correlation with Chemical Structure and Prediction from Additive Group Contributions, third edition Elsevier, Amsterdam. Yermakova, A., Mikenin, P., Anikeev, V., 2006. Phenol oxidation in supercritical water in a well-stirred continuous reactor. Theor. Found. Chem. Eng. 40, 168–174.