Format And Type Fonts CCHHEEMMIICCAALL EENNGGIINNEEEERRIINNGG TTRRAANNSSAACCTTIIOONNSS VOL. 39, 2014 A publication of The Italian Association of Chemical Engineering www.aidic.it/cet Guest Editors: Petar Sabev Varbanov, Jiří Jaromír Klemeš, Peng Yen Liew, Jun Yow Yong Copyright © 2014, AIDIC Servizi S.r.l., ISBN 978-88-95608-30-3; ISSN 2283-9216 DOI:10.3303/CET1439164 Please cite this article as: Galli F., Corbettaa M., Pirola C., Manenti F., 2014, Robust kinetic modelling of heterogeneously catalyzed free fatty acids esterification in (monophasic liquid)/solid packed-bed reactor: rival model discrimination, Chemical Engineering Transactions, 39, 979-984 DOI:10.3303/CET1439164 979 Robust Kinetic Modelling of Heterogeneously Catalyzed Free Fatty Acids Esterification in (Monophasic Liquid)/Solid Packed-Bed Reactor: Rival Model Discrimination Federico Galli a , Michele Corbetta a , Carlo Pirola* ,b , Flavio Manenti a a Politecnico di Milano, Dipartimento di Chimica, Materiali e Ingegneria Chimica "Giulio Natta", Piazza Leonardo da Vinci, 32 - 20133 Milano, Italy b Università degli Studi di Milano, Dipartimento di Chimica, Via Golgi, 19 - 20133 Milano, Italy carlo.pirola@polimi.it Biofuels are key products for the sustainability of the world energy consumption in the next years. Biodiesel in particular is a non-toxic, biodegradable, environmentally friendly alternative diesel fuel. Nowadays, the main problem for the commercialization of biodiesel is its final cost, that is strongly dependent (about 85 % of the total) by the feedstock used. A possible way to lower the biodiesel production costs is using raw oils, which contain a higher amount of Free Fatty Acids (FFA) that should be eliminated before the transesterification reaction to avoid soaps formation. In this work the regression of kinetic parameters of heterogeneously catalyzed esterification with a packed-bed reactor was made. Robust techniques for kinetic parameters estimation and simultaneous discrimination of rival models were adopted and combined with a dedicated differential-algebraic equation (DAE) model that characterizes the system. The main kinetic parameters were regressed using kinetic models other than literature models in considering the activity (UNIQUAC model was used to calculate the activity coefficients) of the components instead of the concentration because the oil/methanol/water/FAME system is highly non ideal. The kinetic parameters were obtained using equilibrated resins, i.e. using the catalyst after having let it to adsorb reactants and products at the operative conditions. The regressed parameters allow to represent the system in a wide range of operating conditions with a little residual error. 1. Introduction Biofuels are key products for the sustainability (Čuček et al., 2012) of the world energy consumption in the next years (Andrade, 2005). Biodiesel is a fatty acid methyl esters (FAME) mixture and can be obtained by transesterification of highly refined vegetable oils with methanol in homogeneous based catalysed processes (Ma et al, 1999). The starting oils are mainly constituted by triglycerides (about 90-98 % of total mass) and Free Fatty Acids (FFA), linear carboxylic acids in the C14-C22 range, with different instauration levels. As of today, the main problem for the commercialization of biodiesel is its final cost, that is strongly dependent (about 85 % of the total) by the feedstock used. Using unrefined or waste oils as a feedstock represents a very convenient way in order to lower biodiesel production costs (Klemeš et al., 2010). Some examples of low cost raw materials for biodiesel production are crude vegetable oil (Ceclan et al., 2012), waste cooking oil (Bianchi et al., 2009) and animal fat (Pirola et al., 2010). The main problem associated with the use of these types of low-cost feedstock lies in their high content of FFA, leading to the formation of soaps during the transesterification reaction. Esterification of FFA in presence of either homogeneous or heterogeneous acid catalysts allows, at the same time, to lower the acid content and to obtain methyl esters, i.e. biodiesel, already in this preliminary step. A remarkable advantage of heterogeneous catalysis is the easier separation and recovery of the catalyst after the reaction. Since FFA esterification is an equilibrium reaction, the addition of a methanol amount higher than the stoichiometric one is common in order to shift the equilibrium towards the products. However, the formation of two liquid phases has got two main drawbacks: firstly, part of the catalyst inside the reactor will be wetted by this second phase and 980 then it will not work in the reaction, due to the high affinity between the resins and methanol; secondly, the repartition of FFA between the two liquid phases contributes to lower the final biodiesel yield even if the oil at the end results deacidified. The aim of the present work is the regression of the kinetic parameters taking into account two different models (Popken et al., 2000), a pseudo-homogeneous one and a adsorbed based one, which accounts for the different affinity toward the polymeric matrix of all the species involved in the reaction and the solvent (triglycerides). The nonideality of the liquid mixture was considered using the UNIQUAC equation for calculating the activity coefficients. 2. Experimental The experimental data used for the rigorous kinetic parameter estimation are reported elsewhere (Pirola et al., 2014) together with the analytical methods and the catalyst characterization. A brief discussion of the choice of the catalyst type and the rector description are however reported here for the sake of clarity. 2.1 Catalyst Sulphonic acid exchange resin Amberlyst 46 (wet) catalyst was used for all the experiments. It was kindly provided by Dow Chemicals Company. The sample Amberlyst 46 was chosen because this resin unlike all other Amberlyst 46 is not internally sulphonated, but it only shows surface acid groups. Consequently, it is not subject to any internal adsorption-desorption phenomena for both reactants and products that thus could be neglected. 2.2 PBR reactor The PBR reactor, whose scheme is reported in Figure 1, is an iron cylinder 20.3 cm long and with an internal diameter of 4.7 cm with 3 intermediate samplings. The catalytic bed is placed at 7.4 cm from the bottom of the reactor and it has a volume of 86 cm 3 . Figure 1: PBR reactor scheme The resins charged were used for all the experimental deacidification runs. In order to maintain a constant internal temperature, an electrical resistance band controlled by an internal thermocouple is used. This internal thermocouple was inserted in the middle sampling. The oil is mixed with methanol using an FFA: MeOH molar ratio of 1:5 in a feeding chamber (FC) pressurized by compressed air up to 6 bar. The feed is fluxed into the reactor, heated by an external hot circuit and insulated using glass fibers. The whole system is thus under the same pressure of the FC while the flow is controlled by a suitable mechanical valve (V7). The esterification reaction was carried out at different operative temperatures at 6 bar. The main advantage of this reactor configuration is the possibility to perform deacidification reactions above the methanol normal boiling point (64.70 °C) because the system is under pressure. In Table 1 the main reactor characteristics are reported. 981 Table 1: PBR reactor parameters Reactor Volume Void fraction Catalyst density Catalyst charged 0.180 dm 3 0.36 0.60 kg/dm 3 51.6 g 2.3 Nonlinear regression using BzzMath library The kinetic parameters regression on experimental data was performed by means of the set of very robust optimizers belonging to the BzzMath library (Buzzi-Ferraris and Manenti, 2012). Such optimizers are based on the object-oriented programming and parallel computing so as to reduce the computational time (Buzzi-Ferraris and Manenti, 2010a) and they have implemented special programs able to simultaneously handle the so-called narrow-valley problem, which typically arise in the estimation of kinetic and thermodynamic parameters (Buzzi-Ferraris and Manenti, 2009); the possible multicollinearities of the parameters, which could be due to coupled chemical-physical phenomena; and the possible presence of bad-quality measures with the identification of possible outliers (Buzzi-Ferraris and Manenti, 2010b). Actually, the mistake made by many strategies is to adopt search directions oriented along the bottom of the valley for the one-dimensional minimum. Such a search will prove ineffective, usually because the direction is inexact and the valley is nonlinear. To exploit the search direction that inaccurately detects the bottom of the valley, it is necessary to change the point of view: (1) any optimization algorithm can find the bottom of the valley by starting from a point outside the same valley; (2) the line joining two points on the bottom of the valley is a reasonable valley direction; therefore there is a good probability that a point projected along such a direction will be close to the valley; (3) nevertheless, this valley direction must not be used as the one-dimensional search direction, but rather as a direction along which a new point projection must be carried out; (4) this new point should not be discarded even though it is worse than the previous one, rather it is the new starting point for a new search; (5) this search must be performed in the sub-space orthogonal to the valley direction to prevent the issue of small steps arising. This philosophy is particularly simple in object-oriented programming. The optimization problem is split into two different levels: the first (outer optimizer) is managed by a single object that exploits the above mentioned procedure to find a certain number of points to initialize an even number of objects. In the second (inner optimizer), each object uses a program to search for the minimum with a limited number of iterations starting from the point assigned by the outer optimizer. This philosophy is useful in solving all problems demanding algorithm robustness: (1) when the function has many minima and we need to search for the global minimum; (2) when the function has very narrow valleys (or steep walls); (3) when the function is undefined anywhere. All these problems arise in the estimation of thermodynamic parameters. Moreover, this philosophy is particularly effective when several processors are available. The performance in terms of estimation accuracy and computational effort as well of the proposed approach has been proven with respect to the tools available in the commercial suites. 2.4 Kinetic modeling All the experimental results, as previously described were used to regress the kinetic parameters taking into account two different models: a pseudo-homogeneous one, a power-law equation which simply considers the reaction rate dependent on the bulk activities of the components, shown in the following equation: 𝑟 = 𝑑𝜉 𝑑𝜏 ∗ 𝐶𝐹𝐹𝐴 𝑜 = 𝑘1 0 ∗ 𝑒 −𝐸𝑎1 𝑅𝑇 ∗ 𝑎𝐹𝐹𝐴 ∗ 𝑎𝑀𝑒𝑂𝐻 − 𝑘−1 0 ∗ 𝑒 −𝐸𝑎−1 𝑅𝑇 ∗ 𝑎𝐹𝐴𝑀𝐸 ∗ 𝑎𝐻2𝑂 (3) where k0i and Eai are the adjustable Arrhenius kinetic parameters, ai are the components activities while T is the absolute temperature. As previously explained, instead of considering all the fatty acids molecules, the oleic acid was chosen to represent all the oil FFA. The FAME then are represented by methyl oleate and the oil by triolein, a triglyceride constituted by three molecules of oleic acid. The other model considered is the adsorption-based one, shown in the following equation: 𝑟 = 𝑑𝜉 𝑑𝜏 ∗ 𝐶𝐹𝐹𝐴 𝑜 = 𝑘1 0 ∗ 𝑒 −𝐸𝑎1 𝑅𝑇 ∗ 𝑎′𝐹𝐹𝐴 ∗ 𝑎′𝑀𝑒𝑂𝐻 − 𝑘−1 0 ∗ 𝑒 −𝐸𝑎−1 𝑅𝑇 ∗ 𝑎′𝐹𝐴𝑀𝐸 ∗ 𝑎′𝐻2𝑂 𝑎𝐹𝐹𝐴 ′ + 𝑎𝑀𝑒𝑂𝐻 ′ + 𝑎𝐹𝐴𝑀𝐸 ′ + 𝑎𝐻2𝑂 ′ + 𝑎𝑂𝑖𝑙 ′ 2 (4) with 982 𝑎𝑖 ′ = 𝐾𝑖 ∗ 𝑎𝑖 𝑀𝑊𝑖 (5) were Ki is the adsorption affinity constant for the i-th molecule, shown in Table 2, taken from (Popken et al., 2000) and adapted considering 𝐾𝐹𝐹𝐴 = 𝐾𝐹𝐴𝑀𝐸 = 𝐾𝑀𝑒𝑂𝐻 3.5 (6) as suggested by Rehfinger and Hoffmann (1990). The binary adsorption affinities were not regressed together with the kinetic parameters because a good fit but very poor prediction ability would have been obtained. The kinetic models with either an ideal liquid phase, considering the activity coefficients equal to one, or the nonideality of the mixture using the UNIQUAC equation, were considered to compare the kinetic parameters obtained when a more correct thermodynamic approach (UNIQUAC) is chosen instead a simpler one (IDEAL). UNIQUAC interaction parameters taui,j were calculated using Eq. (7) 𝜏𝑖,𝑗 = 𝑒𝑥𝑝 𝑏𝑖,𝑗 𝑇 (7) where the binary bi,j parameters were taken from the AspenPlus™ database and are reported in Table 3. 3. Results Since the PBR used for the esterification experiments has got only three samplings, for each temperature tested three different feeding flow rates have been used to obtain all the experimental curve reported. 3.1 Kinetic regression As explained in Section 2.4 two different models were used for estimating the kinetic parameters. The optimized parameters are reported below in Table 4 with the residual errors. Even if the pseudo-homogeneous model does not consider the adsorption of both reactants and products, its use permits to better calculate the experimental trends, especially for the runs performed at high temperatures. This is probably due to the non-correct determination of the adsorption constants that were measured for the binary non-reactive mixtures at a fixed temperature (25 °C), far from the experimental Table 2: Binary adsorption equilibrium constants. Component (i) Binary Adsorption Affinity (Ki) Water 5.24 Methanol 5.64 FFA 1.61 FAME 1.61 OIL 0.00 Table 3: UNIQUAC binary interaction parameters. i FFA H2O FFA FFA FFA H2O H2O FAME FAME OIL j FAME MEOH H2O OIL MEOH FAME OIL OIL MEOH MEOH bij [K] 83.44 -254.73 -377.09 80.52 -567.09 -252.11 -245.42 15.43 -579.71 -459.50 bji [K] -106.26 165.26 -232.04 -90.43 112.63 -645.18 -435.31 -18.93 24.04 37.13 Table 4: optimized kinetic parameters for the FFA esterification using different models. Model SSE k1 0 (mol*s -1 *m -3 ) k-1 0 (mol*s -1 *m -3 ) Ea1 (kJ*mol -1 ) Ea-1 (kJ*mol -1 ) Pseudo-homogeneous (IDEAL) 0.131 9.30*10 7 9.62*10 -8 32.268 -67.589 Pseudo-homogeneous (UNIQUAC) 0.137 2.12*10 8 1.13*10 -7 33.154 -64.343 983 Adsorption-based (IDEAL) 0.145 9.67*10 6 8.39*10 -9 43.152 -57.229 Adsorption-based (UNIQUAC) 0.291 1.34*10 7 1.09*10 -10 43.067 -65.658 operative parameters. A regression of the adsorption based model kinetic parameters together with the binary adsorption constants of water and methanol (assuming valid the constraints of Eq. (6) and Koil=0) have been performed. The optimized parameters are reported in Table 5. By giving the optimizer two more degrees of freedom the SSE sensibly decreased, showing that a better fit could be obtained using the adsorption-based model. A comparison between some experimental data and the calculated behaviour is shown in Figure 2. From Figure 2 it is clear how an adsorption based model (Figure 2a and 2c) better fits the experimental data reported, that is particularly true near the equilibrium condition. The calculation of the activities considering the UNIQUAC model does not influence the overall SSE probably because the experimental amount of methanol used was chosen in order to have only a monophasic liquid mixture but, being the system oil/FFA/FAME/methanol/water highly nonideal, an eventual formation of two liquid phases can be calculated with this latter thermodynamic approach and thus its use is preferable. Table 5: Adsorption based model kinetic parameter, regressed together with the water and methanol adsorption equilibrium constants. Adsorption-based model SSE k1 0 (mol*s -1 *m -3 ) k-1 0 (mol*s -1 *m -3 ) Ea1 (kJ*mol -1 ) Ea-1 (kJ*mol -1 ) KWater KMeOH IDEAL 0.09 3.22*10 6 6.96*10 -1 40.467 -55.796 1.74*10 -2 8.04*10 5 UNIQUAC 0.09 5.52*10 5 1.20*10 2 35.398 -58.339 5.16*10 -7 4.03*10 5 (a) (b) (c) (d) Figure 2: experimental FFA esterification, experimental data (points) at 95 °C (circles), 85 °C (triangles) and 54 °C (diamonds) and simulated curves using (a) the adsorption based model with UNIQUAC, (b) the pseudo-homogeneous model with UNIQUAC, (c) the adsorption based model with IDEAL (activity coefficients=1), (d) the pseudo-homogeneous model with IDEAL (activity coefficients=1) 984 4. Conclusions The esterification of FFA in sunflower oil with methanol was studied in a PBR reactor. The experimental results obtained allowed to regress the main kinetic parameters using a rigorous optimizer library, considering two different models and either an ideal or a non-ideal liquid phase behaviour, using the UNIQUAC model for the calculation of the activity coefficients, that is not commonly found in literature. From the regression results, it could be concluded that the adsorption based model is better for fitting the experimental data and that, since the experimental data considered were obtained using a methanol amount such to have a monophasic liquid system, the calculation of the activity coefficients does not give a great advantage in this particular experimental situation but, considering the high nonideality of the system oil/FFA/FAME/methanol/water, the usage of UNIQUAC model appears to be more convenient. References Andrade J.B., 2005, Biodiesel: an Overview, Journal of Brazilian Chemical Society, 16B, 1313-1320. Buzzi-Ferraris G., Manenti F., 2012, BzzMath: Library Overview and Recent Advances in Numerical Methods. Comp. Aid. Chem. Eng., 30(2), 1312-1316. Buzzi-Ferraris G., Manenti F., 2010a, A Combination of Parallel Computing and Object-Oriented Programming to Improve Optimizer Robustness and Efficiency, Comp. Aid. Chem. Eng., 28, 337-342. Buzzi-Ferraris G., Manenti F., 2010b, Outlier Detection in Large Data Sets, Comp & Chemical Engineering, 35, 388-390. Buzzi-Ferraris G., Manenti F., 2009, Kinetic models analysis, Chemical Engineering Science, 64(5), 1061- 1074. Bianchi C.L., Boffito D.C., Pirola C., Ragaini V., 2009, Low Temperature De-Acidification Process of Animal Fat as a Pre-Step to Biodiesel Production, Catal. Letters, 134, 179-185. Ceclan R., Pop A., Ceclan M., 2012, Biodiesel From Waste Vegetable Oils, Chemical Engineering Transactions, 29, 1177-1182. Čuček L., Klemeš J.J., Kravanja Z., 2012, A review of footprint analysis tools for monitoring impacts on sustainability, Journal of Cleaner Production, 34, 9-20. Klemeš J.J., Varbanov P.S., Pierucci S., Huisingh D., 2010, Minimising emissions and energy wastage by improved industrial processes and integration of renewable energy, Journal of Cleaner Production, 18(9), 843-847. Ma F.R., Hanna M.A., 1999, Biodiesel production: a review, Bioresource Technology, 70, 1-15 Pirola C., Manenti F., Galli F., Bianchi C.L., Boffito D.C., Corbetta M., 2014, Heterogeneously catalyzed Free Fatty Acids esterification in (monophasic liquid)/solid Packed Bed Reactors (PBR), Chemical Engineering Transactions, 37, 553-558. Pirola C., Bianchi C.L., Boffito D.C., Carvoli G., Ragaini V., 2010, Vegetable oils de-acidification by Amberlyst: study of catalyst lifetime and a suitable reactor configuration, Ind. Eng. Chem. Res., 49, 4601-4606. Popken T., Gotze L., Gmehling J., 2000, Reaction Kinetics and Chemical Equilibrium of Homogeneously and Heterogeneously Catalyzed Acetic Acid Esterification with Methanol and Methyl Acetate Hydrolysis, Ind. Eng. Chem. Res., 39, 2601-2611. Rehfinger, A.; Hoffmann, U., 1990, Kinetics of methyl tertiary butyl ether liquid phase synthesis catalysed by ion exchange resin—I. Intrinsic rate expression in liquid phase activities, Chemical Engineering Science, 45(6), 1605-1612.