Novel method for mapping the applicability of reactive distillation

Reactive distillation (RD) is a great process intensification concept applicable to equilibrium limited reaction systems, but how can anyone decide quickly if RD is indeed worth applying? To answer this question, this study proposes a mapping method for checking the applicability of reactive distillation (RD). The initial development is for one of the most relevant subset of quaternary reversible reactions (A+B⇄ C+D, with boiling points Tb,C < Tb,A < Tb,B < Tb,D), by using only basic chemical (equilibrium and kinetics) and physical (relative volatilities) parameters. Generic cases, assuming ideal thermodynamics and constant parameters, are used to obtain a set of RD applicability graphs that provide broad insights into the RD operation. In addition, the new mapping method provides reasonable estimates of the RD applicability to real (non-ideal) chemical systems based on the available pre-defined maps (which are actually applicability graphs of the generic ideal cases). This new approach leads to a straightforward estimation of the applicability of RD to real systems, prior to performing any rigorous process simulations and without any clear-cut decision making (as used in previous studies).


Introduction
Reaction and separation, the most important operations in the chemical industries, are usually carried out in different sections of a production plant, and require different types of process equipment. A reactor is an operating unit where the actual transition of feedstock into products takes place. In most cases, next to the desired main product, some by-products are also formed. Accordingly, a separation step is needed to obtain the desired product(s) at sufficient purity. Distillation is a separation technology that is most commonly applied, but it is also one of the major energy users in the chemical industry. Since the mid of 20th century, scientific literature and patents related to the improvement of reaction and separation equipment design focus on energy savings and economic efficiency [1]. Combining reaction and separation in a single unit is an excellent example of process intensification. Reactive distillation (RD) is one of such processes and it stands out as a successful story of a process intensification technology for enhanced manufacturing of chemicals [2,3]. RD combines a reactor and a distillation column into a single unit operation (see Fig. 1). In the RD column, the reactants are converted while simultaneously separation of the products occurs. The advantages and limitations of RD over conventional multi-unit processes for specific applications have been known for a long time and can be found in many articles and books [4][5][6][7][8]. RD configuration is especially beneficial for chemical equilibrium limited reactions, e.g. (trans-)esterification, etherification, hydrolysis, alkylation, as the equilibrium composition can be shifted towards product formation. The most encountered class of reactions are: [8] A + B ⇄ C + D (quaternary systems) and A + B ⇄ C (ternary systems) RD has been industrially used for more than 25 years, for applications with capacities up to 3000 kton/year [9]. Current applications of RD are mostly for esterification reactions, with the production of methyl acetate as a prime example [10]. Other processes in which RD has been successfully applied are in the production of ethers: methyl tertbutyl ether (MTBE), ethyl tert-butyl ether (ETBE) and tert-amyl methyl ether (TAME) [8]. Many prospective chemical systems (which are neither extremely exothermic nor endothermic) for the RD application are also listed in the open literature [11][12][13].
The available reactive distillation design methods can be classified into three main groups, based on: 1) graphical/topological considerations, 2) optimization techniques, 3) heuristic/evolutionary approaches [14], which are presented in literature [7,[15][16][17][18][19][20][21][22][23][24]. There are various outputs of those design methods consisting of RD structure (operating conditions and RD configurations), feasibility assessment, and/or RD controllability. In addition, there are also methods to check the feasibility of RD for various systems, but they rely mostly on clear-cut decision making procedures (e.g. if the equilibrium constant or the reaction rate is lower than a specific value then RD can be dismissed) while the reality shows that grey areas also exist and they should not be easily discarded (especially for systems with high value products).
Most of the current design methods in literature are well-established and can be used to design a RD column. However, rigorous calculations and/or detailed simulations are usually required to apply the methods for each chemical system and repeated calculation efforts are needed when other chemical systems are investigated, therefore they are considerably complex and time consuming. Following the progressive growth of the number of developed RD design methods, a critical question has been raised more than a decade ago: how could anyone decide quickly (at the conceptual design stage) whether RD is a feasible process concept for a certain reversible reaction system? [25]. The ultimate goal would be to rapidly assess the RD applicability to various reaction systems by only using a simple model (i.e. requiring significantly less time for the evaluation than any other method available).
This paper describes the development of a novel RD mapping method -based on the KISS principle (keep it short & simple) for the end users -that aims to provide insights into the RD operation and quickly evaluate the applicability of RD to (real) chemical systems with a rather straightforward approach. To start with, the most relevant subset of the quaternary systems with both reactants as mid-boiling components (T b,C < T b,A < T b,B < T b,D ) was investigated as it is commonly encountered in practice. A good separation of products is attainable for this boiling point order. The mapping approach uses generic cases to produce the RD applicability graphs, based on ideal thermodynamics and few specified basic data, i.e. relative volatilities (α), chemical equilibrium constants (K eq ) and chemical reaction kinetics. The applicability graph is presented by plotting the reflux ratio (RR) vs number of theoretical stages (NTS), which then can easily give access to the energy requirement and the capital investment. Extensive insights into the RD operation are provided using those applicability graphs. Finally, the new RD mapping method is used to assess the applicability areas of real (non-ideal) systems by only referring to available pre-defined applicability graphs (based on the generic cases). This approach enables a quick assessment of the RD applicability, prior to performing any rigorous simulations of the RD process, thus providing sufficiently accurate information about the applicability of RD and being an important aid for a go / no-go decision at an early stage of the process design.

Approach and methodology
At the initial stage of the development of the RD mapping method, the focus has been limited only to certain levels (but this will be extended further in future studies): • The assessed quaternary systems are reactions with mid-boiling reactants (A + B ⇄ C + D, with T b,C < T b,A < T b,B < T b,D ) as this subset of the quaternary systems is the most commonly encountered.
• The RD configuration (see Fig. 1) is a single column with three different sections (i.e. rectifying, reactive and stripping sections), a condenser at the overhead part and a reboiler at the bottom part. For the sake of simplicity, the feed inlets are fixed on the top and the bottom parts of reactive section (as common industrial operation). Varying the feed inlets inside the reactive section may or may not (slightly) improve the achievable conversion.
• Case studies presented in this study are real reaction systems that are less hindered due to significant non-ideality. Further developments of the mapping method need to cover more complicated reaction systems (e.g., complex azeotropes, liquid split).
To perform any RD operation, some inputs have to be specified and fixed. Fig. 2 presents key parameters of the RD operation in this study. The fixed inputs in this study are highlighted by the bold letters. The feed streams of pure A and B are fed in a stoichiometric ratio (as saturated liquid) to the RD column operating at an atmospheric pressure. For the sake of simplicity, the light reactant is fed at the bottom part of the reactive zone and the heavy reactant is fed at the top of the reactive zone in order to obtain a counter current flow along the reactive zone. With those specified feed locations, the RD configuration used in this study is shown in Fig. 1 constraints specified, which are a bottom product purity and a 0.5 mol/ mol bottom-to-feed (B/F) ratio. The specified value of the B/F ratio is in accordance with the stoichiometric ratio of these quaternary reaction systems. Except for the part of investigating the effect of product purity, a ≥99 mol% of bottom product purity was always set. By setting the B/ F ratio and product purity, the minimum overall conversion is 99 mol%. The RD mapping method uses a set of applicability graphs of RD, i.e. the plots of RR vs NTS, which are generated from generic cases. A generic case is defined by specifying ideal vapor-liquid behavior and constant parameters, i.e. α, K eq and chemical kinetics. The combination of those basic parameters gives unique applicability graphs to a certain case (see the Supporting information, Table S1). A procedure to generate an applicability graph of a generic case is presented in Fig. 3. All simulations were performed using the process simulator Aspen Plus v8.6, by applying a sensitivity analysis and an optimization tool. The sensitivity analysis was utilized to vary the configurations of RD (i.e. numbers of rectifying, reactive and stripping stages) for each NTS. At the same time, for each configuration the optimization tool was used to provide a solution with a minimized RR. Combining the sensitivity analysis and the optimization tools to minimize RR for any RD configuration distinguishes the method proposed in this work from other design methods that aim to estimate the RR min for an infinite NTS using short-cut methods, such as the works of Doherty et al. [26][27][28][29]. Fig. 4 shows an illustrative applicability graph of RD for a certain chemical system. The dotted line is a boundary line which divides the plot regions into 'applicable' and 'not-applicable' areas. Inside the applicable area and exactly on the boundary line, the operation of RD is conceivable since the product purity specification is equal or better than the minimum criterion. For each NTS, there are multiple solutions of RD configurations which are available with different RR values (see the Supporting information, Fig. S1 and Table S2). Along the boundary line, the lowest RR possible is plotted for each NTS. Above the boundary line, the RR is higher for each NTS with either higher product purity or varied distributions of number of rectifying, reactive and    Table S3). Outside (left from or below) the boundary line, the design criteria cannot be achieved therefore the RD operation is not applicable. Point 1 and 2 highlight two different spots on the boundary line and inside the applicable area, respectively, where both have the same NTS. At point 2, the RR is higher than the value at point 1. Operating a RD setup with a higher RR may not be interesting in terms of energy requirement, but it is also essential to consider the RD configurations (i.e. number of rectifying, reactive and stripping stages) and product criteria in its application. At some spots inside the applicable area (such as point 2, in comparison to point 1), RD configurations with either a higher product purity specification (which will be discussed in the next section) or less number of reactive stages (but more separating stages since NTS is constant), could be obtained. Obtaining a higher purity product may be preferred and having a shorter reactive section with a slightly higher RR may reduce the costs up to a certain level. Although point 1 provides the lowest RR, selecting point 2 or another spot with a better configuration or a higher purity product inside the applicable area can still be considered.
In practice, it is favorable to have an RD design with smaller NTS and lower RR. Although the RD configurations inside the applicable area are conceivable, the operation of RD is not attractive above a certain practical limit of NTS and RR. In Fig. 4, the color of applicable area from the bottom-left to the top-right corner shifts from lighter to darker shading. The lighter color illustrates the preferred region in the RD feasibility check, as lower capital investments and energy costs can be obtained.
For the sake of presenting clear images limited to realistic values, this study shows only the applicability graphs with the maximum scale of 100. The NTS min is defined as the NTS for RR = 100, the RR min is the lowest RR on the boundary line in the case of NTS = 100. In the RD design, engineers must consider the proportionality aspect, i.e. ratio of height to diameter of the column. Selecting the NTS min on the top-left of the applicable area leads to the requirement of a high RR which results in a short column with a large diameter. On the other hand, choosing the column configuration with the RR min on the bottom-right of the applicable area gives a slim and tall column.
Equilibrium limited reactions are investigated first, followed then by kinetically controlled reactions (where achieving the equilibrium is practically limited by the slow kinetics). The equilibrium system gives the best performance for the RD column, as it is only limited by the  chemical equilibrium, so this is the best case scenario.
For slow(er) reactions, the RD performance is greatly affected by the extent and effectiveness of the contact between the internal liquid flow and the catalyst. To represent important parameters that give influence on the applicability of RD in the case of kinetically controlled reactions, the Damköhler (Da) number was specified. The Da number is a powerful parameter as it characterizes the effect of chemical kinetics and the RD design inputs, i.e. liquid residence time/hold-ups, catalyst loading. A small Da refers to slow kinetics, a low catalyst loading, a short liquid residence time due to low liquid hold-ups or high liquid flow rate in the reactive parts of an RD column. Utilization of the Da number in the kinetically controlled reactions for RD technology has been common practice in many previous research studies [30][31][32][33].
In this study, the formula of a modified Da number per stage is expressed in Eq. (1).
The Da number (dimensionless unit) indicates the ratio of a characteristic liquid residence time to a characteristic reaction time [33]. The liquid residence time per stage (τ) is defined as the liquid hold-up per stage on volume basis relative to incoming volumetric flow rate per stage. The characteristic reaction time is taken from the reaction rate constant. The modified Da number has the concentration effect inside the k f , while ß stands for the catalyst loading per stage which is expressed in the volumetric ratio between the catalyst amount and the total hold-up per stage. In this study instead of defining the liquid holdups along the column, the τ was specified as a design input to determine the Da number. Setting τ is practical and easy for the RD operation, as usually there is a maximum allowed τ of up to 4-5 min [6,34,35]. To the best of our knowledge, based on industrial experience, τ is typically up to 120 s per stage.
In any real systems, the Da number along the column changes for each reactive stage as the rate constant is dependent on temperature. In the generic case, the Da number defined is constant for each stage, i.e. with more reactive stages and a higher RR (a larger diameter), more catalyst is loaded.

Insights into the RD operation
The influences of various input parameters on the applicability of RD are investigated in case of equilibrium limited and kinetically controlled reactions. In the subsection of equilibrium limited reactions, the influence of product purity, chemical equilibrium constant, relative volatility set is presented. To provide a comprehensive investigation, the K eq was varied from 0.01 to 10 (covering the practical range of reactions in terms of the RD application). In the subsection of kinetically controlled reactions, the applicability of RD (with a low and a high K eq values) is explained linked to the equilibrium limited reactions. The Da number is varied from 0.01 to 1. A relative volatility combination of α AB = 1.5, α CA = 2, α BD = 2 (K A :K B :K C :K D = 3:2:6:1) was selected as a realistic base case (see for example Fig. 5).
The size of RD applicability areas is an essential indicator when the applicability of RD is evaluated for different cases with various input parameters. The NTS min and the RR min are essential parameters as they are limiting the boundaries of the applicability areas.
To provide the insights into RD operation, the RD column configuration at NTS = 2·NTS min is presented next to the applicability graph. For each NTS, various RD configurations with RR values up to 3% higher than the lowest RR-value were considered since there are multiple solutions available with only marginal change of RR (see an example in the Supporting information, Fig. S1 and Table S2). Setting this rule seems realistic as very slight change of RR (i.e. the difference is two decimal places) is often negligible in practice. Due to the existence of multiple RD configurations, many trends of RD configurations can be observed. Therefore, the users of the method can quickly draw different essential insights into the RD operation, which become a major advantage offered by the mapping method. In the current work, only some essential insights are presented based on hand-picked results from the RD configurations obtained, following the mentioned RR rule (NTS = 2·NTS min ), so that trends of the RD configurations can be identified.
The variety of insights due to the availability of multiple RD configurations will be shown in the discussion of the influence of product purity (see subsequent section). Two possible trends of RD configurations are presented in that section based on results selected with the objective to keep either the number of separative or reactive stages (more or less) constant. For the rest of sections, a possible trend of RD configurations will be provided by primarily considering the RD configuration with the lowest RR, but still checking the other possible RD configurations with the RR values up to 3% higher than the lowest value.
Note that it is possible to use other points to provide the insights into the RD operation, e.g. the RD configuration at RR = 1.2·RR min or at any other points. Referring to NTS = 2·NTS min this is only based on the previous knowledge for the estimation of the optimum configuration of conventional distillation columns. A rule of the thumb for the optimum configuration for reactive distillation systems needs to be developed in the near future. For a given RD configuration, the stage number includes condenser (defined as total condenser) and reboiler of the column. The underlined number shown above each bar (see for example Fig. 5, b) is the RR for each configuration.

Equilibrium limited reactions
3.1.1. Influence of product purity Fig. 5 (a) displays the applicability graph of the base case (K eq = 0.01) for different bottom product purities. Obviously, for higher product purity the applicability area becomes smaller. In line with the explanation about the applicability graph in the previous section, higher product purities can be obtained inside the applicable area of a 90 mol% of the bottom product purity.
A higher product purity results in a smaller size of the applicability area of RD. Comparing two end-points of the boundary line of the applicability area, the effect of higher product purity is more dominant on the growth of the NTS min than the increase of RR min . This phenomenon is explained by two possible trends in the RD configurations. Firstly, Fig. 5 (b) shows the selected RD configurations with (more or less) constant number of separative stages. With that objective, the growth of NTS is mainly caused by the requirement of extra reactive stages for a better conversion/separation, as expected. This result shows that the reactive stages contribute to the separation task. For the highest product purity of 99.9 mol%, the addition of reactive stages alone is not sufficient and the number of separative stages needs to be increased.
Secondly, Fig. 5 (c) presents the RD configurations with (more or less) constant number of reactive stages. As the consequence of that objective, the addition of extra rectifying and stripping stages becomes a key solution to obtain higher product purity. Adding more reactive section could help to obtain a higher conversion, but without adequate product separation the conversion is limited at its equilibrium value (see Fig. S3 in the Supporting information). Since the targeted conversion in this study is much higher than its equilibrium conversion (corresponding to the specified product purity, see Fig. 5 a), adding stripping and rectifying stages can be more favorable than having more reactive stages.
In general, a higher RR might also help to endorse the reaction performance. However, a higher RR can lead to the accumulation of products along the column which in the end facilitates the backward reaction and gives difficulty to obtain very high product purity. Performing this analysis with other K eq values and/or for kinetically controlled reactions gave identical insights into the RD operation. 3.1.2. Influence of the chemical equilibrium constant Fig. 6 (a) shows the impact of the K eq on the applicability areas of RD for the base case. For K eq = 0.01, the NTS min (RR ≈ 100) is 22 and the RR min is 11, while for the most favorable case of K eq = 10, the NTS min (RR ≈ 100) and the RR min both are much lower, 12 and 2 respectively. Accordingly, the RD applicability area becomes larger for a higher K eq .
For the different K eq values, the RD column configurations for the base case are shown graphically in Fig. 6 (b). The K eq mainly influences the number of reactive stages and the RR where a lower K eq leads to more reactive stages and a higher RR. However, it is not followed by the growth of rectifying and stripping sections. On the contrary, slightly less rectifying and stripping stages for a lower K eq are needed due to the different product purity coming out of the reactive section. To show those different purities on the top and the bottom parts of the reactive zone, Fig. S4 in the Supporting information shows the liquid composition profiles in the cases of K eq values 0.1 and 1. In the case of K eq = 0.1, x C = 0.38 at the top section of the reactive zone, whereas x D = 0.44 at the bottom part the reactive stages. In the case of K eq = 1, those x C and x D are lower at 0.24 and 0.33, respectively. The profiles clearly show that the purities of products coming out of the reactive zone in the case of K eq = 1 are less than those with K eq = 0.1. The difference on the purity level between these both cases shows the role of reactive stages in performing separation task. Since less reactive stages are available at K eq = 1 compared to the case of K eq = 0.1, the simultaneous separation in the reactive stages is limited therefore requiring more rectifying and stripping stages.

Influence of the relative volatility
To investigate the influence of relative volatility on the RD performance, α CA , α BD and α AB were varied separately. When the relative volatility for product C and reactant A is lower (α CA = 1.2, K A :K C = 1:1.2), compared to the base case the applicability area becomes smaller for each K eq (compare Fig. 6 a and c). For K eq = 0.01, the NTS min (RR ≈ 100) and the RR min are 47 and 25, respectively. This is a significant increase compared to the base case, with NTS min (RR ≈ 100) and RR min 22 and 11, respectively. A similar trend can be found for other values of the equilibrium constant. Fig. 6 (d) shows that compared to the base case in Fig. 6 (b), more stages are required for the rectifying zone since α CA is lower. In addition, a larger rectifying section is also needed because the separation between product C and reactant B (α CB = α CA ·α AB = 1.8) is more difficult. In the base case, α CB = 3. The smaller applicability area for the case of lower α CA compared to the base case, in fact, is not only caused by the larger rectifying section required. Due to more difficult separation in this case, it is preferred to have more reactive stages in order to prevent reactants from reaching the rectifying section. For instance, for K eq = 0.01 the number of reactive stages is 28 in the base case, whereas it is 64 for the case of lower α CA . The increase of K eq in Fig. 6 (d) results in the decrease of the number of reactive stages and the RR due to higher conversion levels per stage.
The opposite effect is expected to happen when α CA is increased. For instance, for α AB = 1.5, α CA = 4, α BD = 2 (K A :K B :K C :K D = 3:2:12:1), the relative volatilities for product C and other components are higher which are α CB = 6 and α CD = 12. In that situation, the separation of product C from the reaction mixture becomes easier. Therefore, the applicability area grows correspondingly and the values of NTS min and RR min decrease.
In analogy to the α CA reduction, the applicability area for the system with the lower α BD (α AB = 1.5, α CA = 2, α BD = 1.2, Fig. 7. The applicability areas of RD and their configurations at NTS = 2·NTS min in case of kinetically controlled reactions for α AB = 1.5, α CA = 2, α BD = 2, considering (a and b) K eq = 0.1 and (c and d) K eq = 10.
K A :K B :K C :K D = 3:2:6:1.7) is smaller compared to the base case as presented in Fig. 6 (e). For instance, for K eq = 0.2, the NTS min (RR ≈ 100) and the RR min for the case of the lower α BD are 27 and 4.6, respectively. With the same K eq , the NTS min (RR ≈ 100) for the base case (see Fig. 6, a) is much less which is 17, with the RR min = 3.0. A similar trend is also found for the other K eq values. Fig. 6 (f) shows the RD configurations and their RR values for the case of lower α BD . In comparison to the base case in Fig. 6 (b), more stripping and reactive stages are required. The same explanation as for the case of lower α CA applies here with a more difficult separation between product D and reactants (α BD = 2 and α AD = 3 for the base case and α BD = 1.2 and α AD = 1.8 for this case). To prevent reactants from reaching the stripping section, more reactive stages are required.
The applicability areas of RD in the cases of varied α AB at 1.2, 1.5 and 2 are presented in Fig. 6 (g) considering α CA = 2 and α BD = 2, with K eq = 10. With a higher α AB , the boundary line shifts closer to the bottom-left of the map and the applicability area becomes larger. At the first glance, this result is seemingly caused by having a larger α AB in the system. However, there is a real reason which mainly affects the size of the applicability areas. Considering the subset of quaternary systems in this study (T b,C < T b,A < T b,B < T b,D ), varying α AB with fixed α CA and α BD gives change to α CB and α AD , which mostly influences the separation performance. The α CB values in the case of α AB of 1.2, 1.5 and 2 are 2.4, 3 and 4, respectively. The same value also applies to α AD for each case. With a higher α AB considering higher α CB and α AD , the separation becomes easier. This explanation is proven by the presented RD configurations in Fig. 6 (h) where the rectifying and stripping sections are slightly shorter and the RR is lower in the case of a higher α AB . On the other hand, a higher α AB slightly increases the number of reactive stages because the lighter reactant is vaporized more easily, therefore hindering the liquid interaction and reducing the reaction performance.
Observing the applicability areas and the column configurations of RD with varied α shows that the results are sensitive to the change of α.

Kinetically controlled reactions
The applicability of RD was investigated in the case of kinetically controlled reactions, with both low and high K eq values. Fig. 7 (a) and (c) present the applicability graphs of RD to the base case considering K eq values of 0.1 and 10, respectively. For both cases, the applicability area at Da = 1 is on top of the applicability area at equilibrium and the applicability area is smaller when the Da number in the reactive section of the RD column is reduced. The lower productivity by the lower Da number needs to be compensated by an increased RR and NTS as shown in Fig. 7 (b) and (d). An increased RR which corresponds to a higher internal flow gives a better separation along the column. A larger reactive zone (more reactive stages with bigger column diameter) allows extended space for the catalyst loading which helps to improve the total conversion when the Da number is low (due to slow kinetics or a short residence time). Compared to the system with the low K eq of 0.1, the applicability areas and the RD configurations for the system with K eq = 10 change much more by the decrease of Da. This points out that although the system has a high K eq , the kinetics, the catalyst hold-up and the liquid residence time give strong effects to the applicability of RD. Comparison of the results for the same Da at Da = 0.01 (see Fig. 7, b and d) for both situations shows the same configuration, as expected, since the very slow kinetics is now the controlling mechanism.
In Fig. 7 (b), it seems that the applicability areas and the RD configurations become similar to the equilibrium limited reaction for Da ≥ 0.05. To investigate this phenomenon, additional simulations were performed with varied K eq and Da numbers. As result, a rule of thumb connecting the kinetically controlled reactions with their equilibrium limited reactions can be derived. Fig. S5 in the Supporting information presents the ratio of number of reactive stages in case of a kinetically controlled reaction over number of reactive stages at equilibrium (Da = infinite) as function of the ratio Da number over K eq, which shows that for Da/K eq ≥ 5 the RD configurations of kinetically controlled reactions are identical to their equilibrium conditions. If Da/ K eq ≥ 2, it is within 10%. This Da/K eq rule of thumb allows the column designers to determine the required design parameters (i.e. catalyst loading and liquid residence time/liquid hold-up) to influence the performance of RD for any intended reaction.
Further study was done in order to check the correlation between column configuration to the catalyst-use. The catalyst amounts along the boundary lines of the applicability areas in the case of K eq = 0.1 were calculated, assuming 20 vol% of the catalyst loading per stage (see Fig. 8). As RR values along the boundary lines go to infinite at NTS min , the column diameter becomes infinite resulting in an infinite catalyst hold-up. The catalyst-use drops following the significant decrease of RR from the vertical asymptote of the boundary line because of less internal flow and a smaller column diameter, which leads to the minimum catalyst requirement at a certain point of NTS. The catalyst loading is then increased with more NTS (the RR remains lower) in which it is affected by more reactive stages needed. This investigation shows the importance of a balance between the RR and the required number of reactive stages in order to operate at the minimum catalyst hold-up for the targeted conversion. In the previous section, it has been discussed that the selection of the column configuration with either NTS min or RR min results in a disproportionate shape of the column. Choosing the column configuration close to the NTS min gives a short column with a large diameter. On the other hand, the column configuration with RR close to RR min results in a slim and tall shape. Related to the annual catalyst expenses, it is suggested to avoid the selection of the RD column configuration at extremes (either close to NTS min or RR min ) for better cost efficiency. A further detailed study is needed to find the optimum RD column configuration considering economics related to the capital investment, the energy requirement and the catalyst-use.
The insights into RD operation have been provided for both the equilibrium limited and kinetically controlled reactions. Having generic cases to perform this study, the insights into RD operation are listed in Table 1 which shows the main effects (on number of reactive and separating stages, reflux ratio) of the modifying specific types of basic parameters.

Development and validation of the RD mapping method
Next to the presented insights into RD operation, this study provides the early development of a new RD mapping method. For the end-users, Fig. 8. Catalyst amounts along the boundary lines of the applicability areas of RD in the case of K eq = 0.1, α AB = 1.5, α CA = 2, α BD = 2, assuming 20 vol% of the catalyst loading per stage. The catalyst amount is based on total molar flowrate of feed. this approach works in a similar way to a global positioning system (GPS) in which the position of a subject is overlapped on some predefined maps (obtained in advance using generic systems). The new mapping method allows defining the applicability areas of real systems (i.e. non-ideal vapor-liquid behavior and temperature-dependent basic parameters, i.e. α, K eq , chemical kinetics) by only referring to the predefined applicability graphs of the generic cases. To match the applicability graphs of real and generic cases, a representing set of α values of the real systems has to be specified. For that reason, the components' ratios of the real systems have been varied giving various combinations of the α set. Note that many trial simulations were carried out in order to validate the applicability graphs. After checking and validating the applicability graphs of the real and the generic cases, it is rational to estimate the set of α values of the real systems, α AB , α CA α BD , at 50/50, 99/1 and 1/99 mol% based mixtures, respectively. The 50/50 mol% of reactants shows an equimolar ratio of feed streams that flow through the reactive section. The 99/1 mol% of product C and reactant A gives an estimation of mixture's composition on the top part of the rectifying section. The 1/99 mol% of reactant B and product D indicates the mixture's composition of the bottom part of the stripping section. Furthermore, it is observed that the average boiling point of reactants can be used to calculate the base K eq and k f values as an estimate for the real systems. A schematic procedure used in this study to develop the method can be found in the Supporting information (Fig. S6). Some examples of the results from extensive simulations which have been carried out during the process of the method establishment are shown in Figs. S7 and S8 in the Supporting information.
There are two important parameters to quantify our level of satisfaction to the developed mapping method at this initial stage: (1) the pre-defined maps can estimate the boundary lines of the applicability area of a real system, (2) the maximum acceptable deviation is ± 50% for the prediction of the NTS and RR of a real case, as this value is commonly found at the conceptual design phase [36]. To calculate the deviation, linear interpolation has been performed to estimate the RD configuration of a real case based on known RD configurations of the two selected generic cases. The estimation based on the interpolation was compared with the simulation result of the real system.

Case 1: transesterification of methyl benzoate with benzyl alcohol
Dimethyl terephthalate ester (DMT) is widely produced by the Witten-Hercules method [37]. In this process, large amounts of methyl benzoate containing waste are produced which are normally combusted. Methyl benzoate in a high purity can be used as a raw Table 1 Summary of insights into RD operation presented in the current study.

Basic parameters Effects
Type The following effects are given by fixed CA  The appropriate property model selected for this system is UNIQUAC-HOC. The Hayden-O-Connell correlation was used to take into account the non-ideal behavior of methanol and methyl benzoate in the vapor phase. There is no azeotrope present in this reaction system and the heat of reaction (ΔH r ) is −13.79 kJ mol −1 .
The availability of chemical data in the literature is limited. The paper of Tang and Li [37] provides the equilibrium conversion for an equimolar feed. In their study, the process utilized tetrabutyl titanate catalyst to produce methanol and benzyl benzoate from the reactants. The equilibrium conversion is 78.1%, which corresponds with a K eq value of 12.7 at 142°C. There is a marginal effect of the temperature on the K eq constant. By assuming a batch reactor, the k f was determined from the conversion vs time plot which is provided in the paper of Tang and Li [37]. This results in Da = 0.067 for τ = 30 s and Da = 0.133 for τ = 60 s, with a catalyst loading of 2 vol% per stage. In practice, the catalyst can be loaded up 50 vol% per stage. The use of 2 vol% of the catalyst per stage in this study is aimed to distinguish the results of kinetically controlled reaction from the equilibrium limited reaction.

Equilibrium-based calculation
A comparison is made between the case study (K eq = 12.7) and the generic cases with K eq values of 10 and 15. Calculating α AB , α CA and α BD at 50/50, 99/1 and 1/99 mol% based mixtures from the real system, respectively, α AB = 1.16, α CA = 256 and α BD = 6.5. Since the K eq and α are very favorable, the RD column with an equilibrium reaction is expected to be applicable. Fig. 9 (a) shows that the boundary line of the applicability area for the real system lies in between two generic cases, but closer to the generic case with K eq = 10 which is mainly caused by temperature influence, especially on α, in the real system. The temperature effect on the K eq is marginal, therefore neglected. Fig. 9 (b) presents the actual RD configurations of the real and the generic cases which were obtained from performing simulations. The graph shows that the NTS and RR of the real case are nicely in the range of the NTS and RR of the two generic cases. Without considering Fig. 9. The applicability areas of RD and their configurations at NTS = 2·NTS min for the transesterification of methyl benzoate (K eq = 12.7) compared to the generic ideal case (α AB = 1.16, α CA = 256, α BD = 6.5) for (a and b) an equilibrium limited reaction and (c and d) kinetically controlled reactions. the simulation result of the real system, the linear interpolation was applied to estimate NTS and RR of the real system as shown in Fig. 10. According to that interpolation, the NTS and RR of the real system are 9 and 0.3, respectively. Comparison between the actual simulation result (NTS = 8 and RR = 0.4) and the estimate based on two generic cases via interpolation show deviation of +13% and −25%, respectively, for the NTS and RR. The complete set of results including deviations for all case studies is summarized in Table 2. The results show a good estimation of the applicability areas and a satisfying accuracy. The generic cases therefore can be used to predict the applicability of RD to this real system in the case of equilibrium-limited reaction. Since the separation is easy and the equilibrium conversion is very high, the application of a conventional system being a reactor followed by distillation might be considered.

Kinetics-based calculation
In this part, the RD applicability area for the real system with kinetically controlled reaction is compared to the generic case (K eq = 12.7) with α AB = 1.16, α CA = 256 and α BD = 6.5. Fig. 9 (c) shows that the boundary lines of the applicability areas for the case study with different Da numbers (Da values are 0.067 and 0.133) lie between those belonging to the generic cases (0.02 < Da < 0.2). Fig. 9 (d) highlights the RD configurations for the case study based on simulation results, which are inside the range of the RR and the NTS of the generic cases. Again, interpolation was applied to estimate the RD configurations of the real case without relying on any simulations of the case study. As presented in   Fig. 11. The applicability areas of RD and their configurations at NTS = 2·NTS min for the hydrolysis of methyl lactate (K eq = 0.096 at 122.4°C) compared to the generic ideal case (α AB = 5.5, α CA = 2.5, α BD = 6.5) for (a and b) an equilibrium limited reaction and (c and d) kinetically controlled reactions.

Case 2: hydrolysis of methyl lactate
Lactic acid is a chemical that can be used to produce biodegradable plastics. However, it is difficult to purify lactic acid from a fermentation mixture. Therefore lactic acid is esterified with methanol to produce methyl lactate. The methyl lactate is then separated and hydrolyzed back to methanol and lactate acid: T b 100°C 144.8°C 64.7°C 216.85°C ΔH r = +33.6 kJ mol −1 To run the simulations, the selected property model was UNIFAC-HOC since it is the most accurate model to describe the current system [38,39]. One azeotrope was found in this system: methyl lactate and water form an azeotrope at 97 mol% water in methyl lactate at 99.8°C. This azeotrope should not have any negative effects on the feasibility and performance of the RD column since it is between reactants and at a high concentration of water. The azeotrope composition will therefore never be reached since the reactants are converted to the products and are fed separately to the column in a stoichiometric ratio.
The chemical data of Sanz et al. [39] is used in this case study. The hydrolysis of methyl lactate is catalyzed by Amberlyst 15, an acidic cation-exchange resin. The quasi-homogeneous non-ideal (QH-NI) model is the best kinetic equation to describe the hydrolysis reaction of water and methyl lactate [39,40]. The kinetic data is shown in Table S4 in the Supporting information and depends on the catalyst concentration. The correlation between temperature and the chemical equilibrium constant is expressed by the Eq. (4). Using the given kinetic data, the forward reaction rate constant can be calculated with Eq. (5).
At the average boiling point of reactants, K eq is 0.096. To calculate the Da number, the k f was then calculated at the average boiling point of reactants with Eq. (5). The Da numbers are 0.16 (7.1 vol% catalyst loading and τ = 30 s) and 0.62 (14.6 vol% catalyst loading and τ = 60 s).

Equilibrium-based calculation
The α AB = 5.5, α CA = 2.5 and α BD = 6.5 at 50/50, 99/1 and 1/ 99 mol% based mixtures, respectively. Fig. 11 (a) shows the applicability area of the generic cases for K eq = 0.05 and K eq = 0.2. Additionally the results of the case study (K eq = 0.096) simulations are added. In the graph, it can be observed that the boundary line of the applicability area for the case study mainly lies in between the lines belonging to the generic case for K eq values of 0.05 and 0.2. Fig. 11 (b) gives the RD configurations based on the simulation results for the real and the generic cases. Without considering the simulation output of the real system and using the interpolation approach, the NTS and RR of the real case were estimated based on the two selected generic cases (see Table 2). Comparing the simulation of the case study and the interpolation result, the deviation of ± 5-7% is highly acceptable. Therefore, a satisfying estimation of the applicability area and the RD configuration of the case study can be obtained from the generic cases.

Kinetics-based calculation
In this part, the RD applicability area for the real reaction system is compared to the generic case (K eq = 0.096) with the kinetics-based calculation. The relative volatilities for the generic case are identical to the values in the equilibrium-based calculation section. In Fig. 11 (c), the applicability areas for the generic system and the case study are plotted. It can be observed that each boundary line of the applicability areas for the case study lies in between the two belonging to the generic cases. Fig. 11 (d) presents the RD configurations of all cases obtained from the simulations, which indicates that the RD configurations of the real system can be nicely predicted from the generic cases. Using the interpolation approach and without relying on the simulation results of the real system, the RD configurations for the real system were estimated with a deviation of −25% to +5%.

Conclusions
Reactive distillation is indeed a proven process intensification method effectively applicable to equilibrium limited reaction systems. Yet, a key question is how can industrial users decide quickly if RD is indeed feasible and worth applying? This study has effectively developed a novel (graphical) approach to evaluate the applicability of RD to quaternary reaction systems, based on generic cases requiring only a few basic parameters, i.e. α, K eq and chemical kinetics. Having those basic parameters, the RD applicability graphs were generated as plots of the reflux ratio vs the number of theoretical stages, which provide information about the applicable configurations of the RD operation. The product purity can be set as a primary performance indicator which influences the size of the RD applicability areas.
Due to the existence of multiple RD configurations (with slight differences in the reflux ratio values) for the same boundary conditions, a broad range of insights and trends regarding RD configurations can be gathered. This feature is a key benefit offered by the new RD mapping method which allows the end user to obtain quickly a better understanding about the RD operation, prior to any detailed rigorous simulations. Some essential insights into the RD operation are conveniently summarized in Table 1.
The development of a new RD mapping method in this study provides satisfying outcomes. It seems promising to use the method for assessing the applicability of RD to real systems, by analyzing the predefined graphs of the generic cases. Furthermore, the method gives quick and good prediction of the RD configurations of real systems with the deviation of less than ± 30%. The mapping method is able to eliminate the necessity of performing any rigorous simulations in the exploratory phase when considering a certain reaction for RDalthough a detailed simulation is suggested in the detailed design phase. The valuable insight provided by the method can be used in the decision making process to go/no-go for RD. The initial development was carried out focusing on the most encountered subset of the quaternary systems (T b,C < T b,A < T b,B < T b,D ), but the method can be expanded further to other systems.