Mathematical modeling of transesterification process kinetics of triglycerides catalyzed by TMAH

Biodiesel is a renewable fuel mainly produced by methylation of triglycerides of vegetable oils or animal fats. The production processes nowadays are particularly based on the utilization of inorganic alkali catalysts. However, it has been proved that an organic alkali – tetramethylammonium hydroxide (TMAH) – can also be used as a very efficient transesterification catalyst. The work presented herein is focused on mathematical modeling of the kinetics of TMAH-catalyzed transesterification of triglycerides at different reaction conditions, specifically at varying reaction temperature with the aim to understand the reaction mechanism and identify the key variables for optimization of the production process. The main kinetic parameters were calculated based on the mathematicalstatistical processing of experimental kinetic data. The reaction rate constants for individual consecutive and reversible reactions and the corresponding activation energies were calculated.


Introduction
Biodiesel is a biofuel derived from vegetable oils or animal fats. Chemically, it is predominantly composed of a mixture of monoalkyl esters of saturated and unsaturated long-chain fatty acids. Because of its chemical composition, it is considered as a renewable, non-toxic, biodegradable, oxygenated and sulfur-free and readily available fuel [1][2][3][4][5][6]. The most common way of biodiesel production is transesterification. Transesterification is a reversible reaction in which the triglycerides react with short-chain alcohol (e.g. methanol → methanolysis) under alkali or acid conditions to form fatty acid methyl esters (FAME), and the main by-product, glycerol [7][8][9]. The scheme of the transesterification reaction is depicted in Fig. 1.

Fig. 1:
The scheme of the transesterification reaction with alcohol [7].
The most frequently used catalysts in biodiesel production are inorganic alkalis, especially KOH and NaOH [6][7] because it was shown [7] that with their use the transesterification proceeds substantially more rapidly than in the case of acid-catalyzed transesterification [7]. However, there are other alkalis proved to efficiently catalyze the transesterification reaction such as tetramethylammonium hydroxide (TMAH). Being a strong organic base, TMAH can act as a catalyst for triglyceride transesterification and simultaneously as a free fatty acid (FFA) esterification agent in case of using waste fats and oils with an acid value exceeding 2 mg KOH.g -1 [9].
The aim of this work is to study triglyceride methanolysis kinetics catalyzed by an organic alkali -TMAH -through mathematical modeling of involved reactions in order to propose a suitable mathematical model capable of adequately describing the reaction system. The mathematical model should confirm the ideas of the reaction scheme; enable optimization of the reaction conditions and simultaneously predict the reaction behaviour even when the experiment is not directly implemented, i.e. it should be able to facilitate reaction behavior simulation.

Mathematic modeling of methanolysis kinetics
As depicted in Figure 1, the transesterification reaction is described as a system of three consecutive and the same number of reversible reactions, i.e. six reaction in the overall system, where k 1 =k 1+ ; k 4 =k 1-; k 2 =k 2+ ; k 5 =k 2-; k 3 =k 3+ ; k 6 =k 3-. To simplify the mathematical model, the following assumptions were made and taken into account when the model was proposed: → The transesterification reaction is not affected by the neutralization of the free fatty acids, and the catalytic effect of the resulting TMAH soaps is negligible. → Neutralization reaction runs several times faster than transesterification itself. → TMAH is consumed to neutralize free fatty acids, and only its excess serves as a catalystthe catalyst concentration is not a function of time and is kept constant during the transesterification reaction. → Presence of water in the reaction system does not substantially influence the investigated reactions (concerning its zero content), it follows then there is no notable hydrolysis of glycerides, and their saponification cas be neglected as well. → The reaction is performed in kinetic region and thus can be considered as pseudohomogenous The kinetics was studied on a model feedstock -a commercial rapeseed oil -which is predominantly composed of triglycerides and a small number of diglycerides. Therefore the above-stated assumptions can be accepted. The mathematical model of the transesterification reaction kinetics was then constructed using the knowledge of the law of mass action when the rate of a chemical reaction is directly proportional to the product of the molar concentrations of the initial compounds.
The law of mass action for two reactants is expressed by the following general relationship (1), (1) where α corresponds to the stoichiometric coefficient of the j-th reactant in the i-th reaction and β is the stoichiometric coefficient of the k-th reactant in the i-th reaction (k≠j).
Consequently, the reaction rate for a particular component can be expressed by a series of 6 differential equations (2 -7).
In the study of transesterification reaction kinetics, the effect of reaction temperature on the reaction rate of consecutive and reversible reactions was modelled using the Arrhenius equation (8).
Here, k 0 denotes the frequency factor of the reaction, E a the activation energy of the reaction, R the ideal gas constant and T the thermodynamic temperature. The activation energy and the frequency factor for each successive reaction can be obtained from experimental data by plotting the logarithm of the gained reaction rates against the reciprocal value of the absolute temperature. -The described procedure leads to linearization of expression (8) and consequently the data can be evaluated by means of relatively simple linear regression.
To calculate the deviation between the predicted concentration and the concentration obtained from the experimental measurement, the following objective function was used. Specifically, the sum of the squares of differences between experimental and predicted concentrations divided by the component concentration in the reaction system obtained experimentally -squared relative error was used (9).

(9)
Where O expresses the total relative deviation related to the concentration obtained from the experimental measurement of the j-th component, , at the k-the time, w j the weight (significance level) of j-th component and is the predicted concentration of the j-th component at k-th time.
All models at each temperature level were optimized with the help of MATLAB software.

Experiments
The proposed mathematical model of the transesterification reaction was verified by laboratory experiments. Transesterification was carried out in an isothermal batch reactor. 250 g of pure commercial rapeseed oil of the acid value of 0.2 mg KOH.g -1 was used as a feedstock. Three values of the reaction temperature were used for verification of the mathematical model, specifically 30 °C; 50 °C and 65 °C. The concentration level of TMAH as a catalyst was set to 1 % (related to the initial weight of rapeseed oil). The molar ratio of oil to methanol was 1:6 (56.25 g of methanol). The reaction mixture was intensively stirred for 2 h under reflux. The course of the reaction in time was studied; in order to stop the reaction, 2 ml of the reaction mixture were taken at the selected time and immediately neutralized in a 0.1M butanolic solution of adipic acid. Afterwards, the obtained samples were analyzed using gas chromatography. Chromatographic conditions and sample preparation for analysis were set according to [10].
To verify that the reaction system was operating within the kinetic region under the selected reaction conditions, experimental measurements were performed to determine whether the reaction system was not affected by the intensity of stirring. The reaction conditions were set to: 1.5 % w/w TMAH, 60 ° C, oil/methanol ratio 1: 6. Four stirring rates were selected, specifically 300, 650, 1300 and 2000 rpm. From the graphical representation of the results (Figure 2) it can be seen, that the stirring intensity higher than 1300 RPM does not significantly affect the course of the transesterification reaction.

Results and discussion
From the obtained kinetic data, the rate constants for the individual reactions and the corresponding activation energies were determined (according to suggested mathematical model). The experimental data were evaluated using MatLAB and the determination of parameters was performed according to the proposed algorithm [11]. Figure 3 shows the course of the transesterification reaction at different reaction temperatures. As can be seen, the rate of triglyceride conversion to FAME grows with increasing temperature. The highest conversion (92 %) was reached at 65 °C, i.e. after 2 hours. The lowest conversion (approximately 71 %) was achieved at a reaction temperature of 30 °C. These values do not represent the maximum range of the reaction conditions and would further affect the conversion rate towards FAME by extending the reaction time. However, the obtained results clearly demonstrate that the temperature significantly determines the reaction range and the rate of transesterification reaction that both affect the reaction time required to achieve the maximum conversion. Of course, increasing the catalyst concentration would lead to higher triglyceride conversion; however, this work was focused on the determination of temperature dependence of triglyceride conversion to FAME. Table 1 shows the rate constants for the consecutive and reversible reactions depending on the reaction temperature of the transesterification reaction catalyzed with 1 % w/w TMAH and also the sum of absolute and relative residues between the expected and experimentally obtained concentrations. The results show that the reaction rate of TG conversion to DG presents the lowest rate of all forward consecutive reactions. The reaction rate of TG to DG is considerably increased with temperature rise, partly due to increasing oil solubility in methanol with temperature. The reason why the reaction of triglycerides is the rate determining step may be the relatively large molecular weight of triglycerides, which can significantly slow down the reaction. Conversely, the highest rate constants were achieved for the second consecutive and reversible reaction, that is when diglycerides are converted to monoglycerides and vice versa, which means that when diglycerides are formed in the previous step, afterwards the reaction carries on significantly faster. The lowest rate constants in the entire reaction system corresponded to the third reverse reaction, i.e., the conversion of glycerol to monoglycerides, which means that the reaction of glycerol with methyl esters, whereby a molecule of monoglyceride and methanol is obtained, was not preferred. The reason is most probably the limited mmiscibility of glycerol with esters, which involves mass transfer resistance in this direction, as also stated in [12]. From this statement, it can be assumed that the third reaction step, i.e. the conversion of monoglycerides to glycerol, is an irreversible reaction, hence the k 3-constants can be neglected. From the sum of absolute and relative residues, it is evident that the used mathematical model, based on the second order reaction mechanism, well describes the reaction behaviour and corresponds with obtained experimental data (see Figures 4-6). This conclusion is valid for all tested temperatures. Let us note that the course of the reaction is significantly affected by the reaction temperature (as is depicted in the figures) and yet the model was able to correctly capture dynamics of the studied system.   The effect of temperature on rate constants was also studied by the linearized form of Arrhenius equation (8). The dependence of the logarithm of the calculated rate constants on the reciprocal value of the absolute temperature is depicted in the Arrhenius diagram ( Figure  7). The activation energies, frequency factors and corresponding regression coefficients are listed in Table  2. The results show that in the case of the first two consecutive (k 1+ , k 2+ ) and reverse reactions (k 1-, k 2-) the rate constants and activation energies increase significantly with increasing reaction temperature, however, the growth rate of said parameters of forward reactions versus reverse reactions is higher. The increase in the rate constants of the third reversible reaction with temperature was very slow. It follows that the conversion of MG to GL is not significantly affected by the reaction temperature, which is demonstrated by the values of activation energies. The activation energies for the individual rate constants are higher as stated in [7,[13][14][15][16][17], which is most probably due to the selected reaction temperatures for studying of the transesterification reaction kinetics and especially investigated type of catalyst.

Conclusion
This work aimed to study the kinetics of triglycerides methanolysis catalyzed with tetramethylammonium hydroxide of a constant concentration at specific reaction temperatures. The necessary experiments were performed in order to verify the kinetic model of the methanolysis.
A mathematical model of transesterification with a defined second-order reaction mechanism was designed for quantitative description of obtained experimental data. The main kinetic parameters were calculated based on the mathematical-statistical processing of experimental kinetic data. The reaction rate constants for individual consecutive and reversible reactions and the corresponding activation energies were calculated. Comparison of experimental data with the model confirmed, based on the values of absolute and relative residues, that the used mathematical model well describes the reaction behaviour despite considerable change of system behavior caused by relatively small changes in the reaction temperature. The proposed mathematical model of the transesterification reaction kinetics provides parameters that can be applied to calculate the reaction range under specific reaction conditions, and the results also serve as necessary input data for the design of an industrial transesterification reactor and overall process simulation