Methodology for the Design of Optimal Chemical ... - ACS Publications

Ind. Eng. Chem. Res. 2010, 49, 10535–10548

10535

Methodology for the Design of Optimal Chemical Reactors Based on the Concept of Elementary Process Functions Andreas Peschel,† Hannsjo¨rg Freund,*,† and Kai Sundmacher†,‡ Max Planck Institute for Dynamics of Complex Technical Systems, Sandtorstr. 1, 39106 Magdeburg, Germany, and Process Systems Engineering, Otto-Von-Guericke UniVersity Magdeburg, UniVersita¨tsplatz 2, 39106 Magdeburg, Germany

In this contribution, a methodology for the optimal design of chemical reactors based on the best reaction route in the thermodynamic state space is proposed. This route is obtained by tracking a fluid element on its way through the reactor and manipulating the fluxes into this element. Instead of choosing the reactor design a priori and optimizing the free parameters of the chosen reactor setup, an innovative reactor design is developed based on the optimal flux profiles. Besides classical reactor concepts, this methodology is suited to investigate the potential of different process intensification options such as integration of reaction, cooling and separation in a single apparatus, or the application of high interface areas for heat and mass transfer. The methodology is exemplarily illustrated for the development of a new SO2 oxidation reactor. The residence time as an example for a meaningful objective is minimized, and a reduction of 69% compared to the optimized technical reference case is achieved. 1. Introduction 1.1. State-of-the-Art Methods for Reactor Design. Due to the complex physicochemical phenomena in a chemical reactor and the nonlinearities in their mathematical description, the reactor design problem is a complex task. As illustrated in Figure 1, the state-of-the-art reactor design methods can be classified into heuristics, attainable region (AR) methods, or rigorous optimization approaches. Criteria for reactor selection are reviewed in many textbooks on chemical reaction engineering and on reactor design (e.g., refs 1 and 2). Heuristics for reactor selection are derived from these criteria and can result in a good reactor selection in many cases.3 Heuristics are easy to apply and fast to use. However, these criteria do not hold for any arbitrary reaction network, and the application for complex systems is limited. Therefore, it is hard to find the best reactor design using such criteria. After having chosen a reactor type, these reactors can only be adjusted to the special requirements of the reaction system within a certain limit by parameter optimization. To find the best possible performance of any reactor network, the concept of the attainable region was developed. The concept was first introduced by Horn,4 and most of the further work was performed by Glasser, Hildebrandt, and Feinberg (e.g., refs 5-8). Since a convex hull of the feasible region is derived graphically, the classical AR approach can only be applied easily to problems which can be reduced to at most three dimensions. This is a severe limitation of the classical AR approach. In addition, the AR method is hard to apply if the reactor network is embedded into the complete process design problem. Every time when the boundary conditions for the reactor problem change, the AR must be calculated anew. For these two cases, Biegler and co-workers developed hybrid methods, which combine the AR approach with rigorous optimization.9-13 Furthermore, the AR concept was extended for dosing of components (e.g., ref 14) and to reaction-separation networks * To whom correspondence should be addressed. E-mail: [email protected]. Tel.: +49 391 6110275. Fax: +49 391 6110634. † Max Planck Institute for Dynamics of Complex Technical Systems. ‡ Otto-von-Guericke University Magdeburg.

(e.g., refs 15-17). Nevertheless, the AR method only yields the best performance of the optimal network of predefined reactor types and does not allow us to develop new reactor concepts directly. Beside the AR approach, rigorous optimization methods for the reactor synthesis problem were derived. The rigorous optimization methods can be classified into superstructure approaches and techniques based on dynamic programming. Superstructure approaches for the reactor network synthesis problem were developed by several authors (e.g., refs 18-22). The interaction between several parts of the plant can easily be taken into account by adding further flowsheet equations. These complicated mixed-integer nonlinear programming (MINLP) problems are often nonlinear and nonconvex. Since these problems show many local optima, finding the global optimum cannot be guaranteed using conventional MINLP solvers. In addition, the characteristics of many reaction networks are similar; e.g., a cascade of continuous stirred tank reactors (CSTRs) approximate a plug flow tubular reactor (PFTR). This problem can be reduced by using simple superstructures, which only include the most promising combinations. Owing to this simplification, a potentially better reactor network may not be included in the superstructure. Hence, superstructure modeling requires always a trade-off between the number of combinations included in the superstructure, the solvability of the problem, and the uniqueness of the solution. The dynamic programming approach for the reactor optimization can be traced back to Aris,23 Horn,24 and Bilous and Amundson.25,26 The optimal design of chemical reactors is given, if the optimal temperature and concentration profiles are achieved along the reaction coordinate.27,28 This holds for batch reactors as well as for continuous plug flow reactors. The nature of these optimal profiles was already investigated by analytical andsas far as possible in those dayssby numerical methods. Since the number of reaction rate equations and the calculation power were very limited, only some examples could be evaluated. As a result, these innovative approaches were not used to design new reactors at that time. Nowadays, solving optimal control problems for the class of semibatch reactors is state-of-the-art. The problem of distributed control variables

10.1021/ie100476q  2010 American Chemical Society Published on Web 05/26/2010

10536

Ind. Eng. Chem. Res., Vol. 49, No. 21, 2010

Figure 1. Classification of reactor design methods.

(feed, temperature, pressure) is considered for batch and semibatch reactors by several authors (e.g., refs 29 and 30). Cuthrell and Biegler developed a discretization method for the arising optimization problems consisting of differential and algebraic equations and applied it to several examples.31 Kjelstrup and co-workers solved the minimum entropy production problem as dynamic optimization problem for the SO2 oxidation32,33 and for several other processes (e.g., ref 34). The main focus of their work is to show that the entropy production can be significantly reduced by using distributed control variables. However, different control options and the problem of how to realize the cooling temperature profile in a technical reactor were not investigated in these contributions. 1.2. Process Intensification Approaches in Reactor Design. For the successful reactor design, several process intensification (PI) methods are important, but PI methods are not considered in general in any of the discussed reactor design methods. For example, the heat and mass transport coefficient can be enhanced by intelligent reactor design (e.g., ref 35) or by considering microreaction engineering (e.g., ref 36). Besides a wide range of specific surface areas for heat and mass transport, it is important to include a wide range of apparent catalyst densities in the reactor design. Here, several reactor concepts ranging from microreactors to conventional fixed bed reactors must be compared, as was done, e.g., in a simulation study by Gu¨ttel and Turek.37 Furthermore, the possibility of using different reaction routes and different solvent concepts is important to find the best reaction system. The reaction route often defines the selectivity, the conversion, and the required residence time in the reactor. The route can be changed by introducing new reaction steps or by using different catalysts or additional reactants (e.g., ref 38). The specific capacities can be significantly changed by operating in a different regime or using another solvent concept, and this may lead to formidable process improvements. In addition, nonconventional reactor operation and nonconventional forces, e.g., by using microwave (e.g., ref 39) or sound fields (e.g., ref 40), should be included in reactor design. Nonconventional and dynamic operation may enhance the reactor performance as well (e.g., refs 41 and 42). The integration of several process functions, e.g., in reactive distillation, becomes more and more important43 and should be considered in the reactor design task, too. Besides the integration of several fluxes in one apparatus, the PI methods can be classified into manipulating kinetic coefficients, driving forces, specific capacities, and design factors as stated in our recent contributions.44,45 To sum up, the problem of finding the optimal reactor and/ or the optimal reactor network has been investigated by many authors. Several good methods have been developed, but all of these approaches are based on predefined apparatuses. The

integration of PI methods into a methodology for the reactor design has not been performed so far. To overcome the shortcomings of the existing methods, a new methodology for the design of optimal chemical reactors based on the concept of elementary process functions (EPF)44,45 is proposed in this contribution. The optimal states of the fluid with respect to a certain objective function change along the reaction coordinatesthey follow an optimal route in state space. To optimize the chosen objective, these optimal reaction conditions must always be accomplished. Taking this as a rationale for an optimal reactor design, the variables to control the states of the reacting fluid must change along the path of the fluid element as well. This cannot be realized with standard reactors that usually have only a fixed value for the control variables, e.g., for the cooling temperature, or can only provide very specific profiles. The main idea of the methodology is to enforce optimal state space profiles by the reactor design. Several requirements for the methodology can be defined. First, the methodology should define an optimal route in state space with respect to a certain objective, and it should be possible to compare different reaction routes. Second, classical reactor concepts as well as process intensification methods should be included in the approach so that the whole space of physicochemical phenomena can be investigated. Third, it should be possible to consider different control variables to enforce the optimal route. Fourth, the losses due to the technical approximation should be quantified using the optimal route as a benchmark. Fifth, the method should work for batch as well as continuously operated reactors. In addition, the design methodology must be independent of existing apparatuses to allow the development of innovative reactor concepts. 2. Methodological Approach The key idea of the proposed reactor design methodology is to track a fluid element on its way through the reactor and optimize the fluxes along its way. The outer fluxes are adjusted, and the reaction fluxes are scaled to meet the optimal reaction conditions at every point along the reaction coordinate, which we define as the optimal route in state space. To design a technical reactor based on this idea, a three-step approach is developed leading to a technical approximation of the optimal route by an appropriate reactor design. The methodological approach is illustrated in Figure 2. In the beginning, the objective function is chosen. The choice of the objective function depends on the reaction system and on case-specific circumstances. It should be kept in mind that different objective functions will lead to different reactors. Figure 3 illustrates the different routes in state space of the individual levels for an exothermic equilibrium reaction with reaction and cooling. The direction in state space, which is


Figure 2. Methodological approach.

Figure 3. Results of the different optimization levels in state space (equilibrium line, lower and upper temperature bound (- - -), level 1 (s), level 2 (s s), level 3 (- - -); er, reaction vector; eq, heat transfer vector; 0, initial state; f, final state).

determined by the reaction, is given by the reaction vector er. The heat transport vector eq defines the state space direction if only cooling occurs.45 The optimal route in state space with respect to the chosen objective, here minimal residence time, is given by the solid line. On the first level, different integration and enhancement concepts are compared. In the context of process intensification, an apparatus with more than one function is referred to as an integrated apparatus.43 The best reaction route is calculated by solving a dynamic optimization problem with unlimited outer fluxes andsas far as meaningfulsscaled reaction fluxes as optimization variables. The reaction kinetics are given by the chosen catalyst, but the reaction fluxes can be optimized in a certain range by scaling the catalyst densities along the reaction coordinate. The fluxes, which are the control functions, are referred to as enhanced fluxes. Only system inherent limitations such as a maximum catalyst temperature and thermodynamic

10537

relations are taken into account, but no limitations arising from choosing a predefined apparatus are considered. Optimized technical reference reactors are calculated as benchmark reactors by optimization of the free design parameters. Using this approach, the potential of controlling specific fluxes is quantified by comparison with the reference reactors. It is important to evaluate the potentials of the different concepts to decide how to guide ongoing engineering work already at the early stage. The results of the first optimization level determine whether it is worthwhile to develop a new type of reactor or whether the parameter optimization of standard reactors is already close to the optimum of the reaction system. Referring to the example shown in Figure 3, no limitations for the heat transfer occur, and the fluid can be heated instantaneously to the maximum allowed temperature on the first level. Afterward, the conversion increases at constant upper temperature for the optimal trajectory. Before the equilibrium line is reached, the optimal trajectory leaves the upper temperature limit and converges to the final state. On the second level, the best suited control variables to achieve the desired fluxes are determined. For this purpose, the flux limitations are investigated by including the transport kinetics of each flux. The optimal profiles of certain variables of the transport kinetics, which can be controlled by the reactor design, are calculated. To quantitatively compare all possibilities on how to obtain the desired flux profiles, a systematic decision structure is applied. The difference in the objective between the first and the second level is the losses due to limited mass and energy transport. In the discussed example, the instantaneous heating and cooling cannot be realized anymore. Hence, the state space trajectory of level 2 differs from the optimal trajectory of level 1. On the third level, an optimal technical reactor based on the best suited control variable profiles is developed. Different technical approximations are possible since the approximation of the control variable profiles is not unique. Either this can be achieved using existing technology, or the development of new technology can be guided by defining a catalog of requirements for the best reactor. In case the number of possible technical approximations is large, the comparison of the different technical approximations can be performed in two steps. In a first step, the reactors are compared using simple models. In a second step, the most promising reactor concept is further investigated using more detailed reactor models, e.g., a two-dimensional reactor model. The best trajectory of the third level differs from the trajectory of the second level since the control variable profiles can in general not be realized perfectly. The difference between the second and the third level indicates the losses due to the nonideal control variable profiles in practice. In addition, the losses due to the nonideal flow field, radial gradients, and dispersion effects are the difference between the simple and the detailed models of level 3. In the following, the individual levels are explained in detail. The arising optimization problems are formulated, and the numerical solution approach is discussed. Afterward, the methodology is illustrated with the development of a new SO2 oxidation reactor. 2.1. Level 1: Optimal Route in State Space. General Model Formulation. On the first level, the optimal route in state space is calculated. To obtain simple balance equations and to avoid any geometrical limitations of standard reactors, the Lagrangian approach for the balance equations is chosen. The state x of the fluid is changed by fluxes j. The fluxes can be classified into volume-related fluxes such as the vector of

10538


reaction fluxes r and surface fluxes jA such as the heat flux q or individual component fluxes ji. In the Lagrangian formulation, each balance equation can be written as a capacity term times the change of the balanced quantity equal to the sum of the fluxes times the according coupling factor. In matrix formulation, this leads to the general notation of the balance equations as stated in eq 1b with the capacity matrix C and the coupling matrix S. To calculate the optimal route, a dynamic optimization problem according to eq 1 must be solved. It is constrained by the balance equations, additional equality constraints g(x), and inequality constraints h(x), as well as inlet and outlet conditions x0 and xf, respectively. The additional constraints arise from the reaction kinetics, from thermodynamics, and from systemspecific limitations. Depending on the integration and enhancement concept, which will be explained in the next subsection, the fluxes into the fluid element are the control functions to obtain the optimal route. Fluxes which are not ideally controlled must be constrained by transport kinetics. These transport kinetics are included in g(x). Here, a set of parameters {p} exist, which can be optimized. For example, this set might include the S/V-ratio for the specific flux. Some catalyst density profiles {Fcat(t)} can be optimized to scale the individual reaction rates. Other catalyst densities can only be optimized as parameters and are included in {p}. The solution of this optimal control problem corresponds to the maximum possible performance of the reaction system. The initial and final conditions for the states can be either fixed or in a certain range. The objective function consists of a stage cost term F(t), an initial cost term I(x0), and an end cost term E(xf). The initial cost term can occur if the inlet conditions are not fixed. The stage cost term can depend on the states, fluxes, design parameters, and catalyst densities. The end cost term arises if an objective function such as selectivity or yield is considered or cost of downstream separation is included in the problem. Obj )

min

{jA(t)},{Fcat(t)},{p}

s.t.

(

∫

C(x) ·

tf

t0

F(t)dt + I(x0) + E(xf))

dx ) S(x) · j dt

(1a) (1b)

g(x) ) 0, h(x) e 0

(1c)

{FLcat} e {Fcat(t)} e {FUcat}, {pL} e {p} e {pU}

(1d)

x(t0) ) x0, x(tf) ) xf

(1e)

Integration and Enhancement Concepts. The calculations at this level are intended to investigate which fluxes are most worthwhile to integrate and which are most worthwhile to optimize along the path of the fluid element. For example, the integration of the heat flux and of specific component fluxes (dosing of reactants and removal of products) may give rise to a better performance of the reaction system. The fluxes, which are actually optimized along the path of the fluid element, are referred to as enhanced fluxes. Several fluxes may be of interest for integration and optimization, and all promising combinations should be investigated. On this level, it is additionally determined whether the use of a different solvent, inert gas, pressure and temperature region, or catalytic route is beneficial for the performance of the reaction system. The outer fluxes which can be controlled are dosing/removal of components, the heat flux, technical, and compression work. The removal of individual components can easily be studied

using the proposed method by optimizing the specific component fluxes. In case the removal of specific components is beneficial, it may technically be realized by reactive extraction, reactive distillation, or specific membranes on the third level. Besides the outer fluxes, the catalyst densities can be scaled along the path of the fluid element in certain limits to control the reaction network. Scaling the catalyst densities is advantageous in case several competing reactions occur or heat removal problems appear or to distinguish between certain catalytic pathways. However, the technical possibilities how to add and remove catalysts can be very limited. In this case, the initial catalyst density can only be optimized as a parameter. Considering an uncatalyzed reaction, the reaction rate can only be influenced by states which have an influence on the reaction rate such as the temperature. The rest of the coupling matrix can only partially be controlled. The capacity matrix is usually a function of the state variables. Hence, it can only be indirectly controlled by influencing the states using certain fluxes. The capacities can be changed on a large scale by operating in a different phase, temperature region, and pressure region or by applying a different inert gas or solvent. Summarizing the enhancement options on the first level, the outer fluxes and the scaling factors for the individual catalysts can be optimized along the path of the fluid element to obtain the optimal reaction route. However, an appropriate set of fluxes should be chosen to yield a uniquely solvable optimization problem. Generally, it is not necessary to control all fluxes and scaling factors of the provided catalyst. In addition, the potential of different catalytic routes and of interesting PI methods can be quantified on this level. Reference Reactors. To benchmark the results of the innovative reaction concepts, optimized reference reactors are calculated for comparison. The modeling of some classical concepts in the EPF framework is discussed in this subsection. In the case of a PFTR with isothermal cooling temperature, the heat flux is not ideally controlled. This concept features integrated reaction and cooling, but no flux is ideally controlled. In this framework, microreactors can be regarded as special types of PFTRs. The heat exchange area and heat transport coefficient can become very large. The flow regime in a microreactor is typically laminar as long as it is not a minireactor with catalyst pellets inside. However, due to the small diffusion length in radial direction (LD ≈ (D/2)), the residence time distribution is PFTR-like if the Fourier criterion (Fo g 1) is valid.2 In this case, the residence time behavior does not have to be taken into account for the single tube analysis. An adiabatic fixed bed reactor with intermediate cooling is modeled using cooling segments and reaction segments in series. Besides the number of segments, the timing of the intermediate cooling and the cooling/heating parameters such as the temperature of the environment phase and the specific heat exchange area can be optimized. An adiabatic fixed bed reactor with cold shot can be optimized by optimizing the number of mixing points as well as the timing, amount, and temperature of the cold shots. Both concepts are not integrated since the fluxes are applied in series. If the CSTR operates under steady-state conditions, the balance equations can be simplified to algebraic equations. The reaction rate and the heat flux are defined by the outlet states. The free parameters of the CSTR can be optimized to obtain a reference reactor with locally concentrated states. If low reactant concentrations are advantageous for the reactor performance and heat removal problems occur, a CSTR is often


preferred to a PFTR. It should be noted that the proposed methodology is completely capable of investigating whether a low reactant concentration and a uniform temperature is beneficial. In case the reactant, the product, and the heat flux are simultaneously optimized, the solution may reduce to a stationary point with low reactant concentration, high product concentration, and constant temperature. This corresponds to a stirred tank reactor in practice, and the technical approximations of the profiles are directly known. To sum up, a PFTR with controlled fluxes is always better than (or at least equal to) a CSTR. Only in the special case where the optimal reactant, product, and temperature profiles are constant and at the outlet values is the CSTR as good as the controlled PFTR. By optimizing the parameters of classical reactor concepts, optimized reference reactors are calculated. For the reference reactors, transport kinetics are added, and the transport parameters are optimized within typical ranges for the equipment. The objective function values of the reference reactor are compared to the results of the cases where fluxes are ideally controlled. This yields the potential of every reaction concept. The most promising reactor concept is further investigated on the second and third level. On the first level, the main feature and advantage of the proposed methodology is that the solution space can be very efficiently scanned to investigate the best reaction route in state space using simple balance equations. This screening includes PI methods, which are modeled on the level of the balance equations. The potential of integrating fluxes, controlling specific fluxes, and all kinds of PI methods is quantified by comparison with optimized reference reactors which helps to decide which concept should be investigated in more detail. 2.2. Level 2: Screening of Control Variables. General Model Formulation. In reality, fluxes are limited, and the influence of limited fluxes is investigated on the second level. The optimal profiles for suitable control variables, which can be realized by the reactor design, are calculated. The control variables are obtained from the kinetic expressions of the fluxes and the design of the reactor in terms of specific surfaces. Hence, the model must be extended accounting for the mass and energy transport kinetics and the S/V-ratio of all fluxes. The general model formulation is stated in eq 2. Obj )

min

{Fcat(t)},{p},

(

∫

tf

t0

F(t)dt + I(x0) + E(xf))

(2a)

dx ) S(x) · j dt

(2b)

{xe(t)},{d(t)}

s.t. C(x) ·

g(x) ) 0, h(x) e 0

(2c)

j ) L(x) · f(x, d, p, x - xe)

(2d)

{FLcat} e {Fcat(t)} e {FUcat}, {pL} e {p} e {pU}

(2e)

{xLe } e {xe(t)} e {xUe }, {dL} e {d(t)} e {dU} x(t0) ) x0, x(tf) ) xf

(2f)

The flux vector is a nonlinear function of the states, the driving forces (x - xe), the design variables d, and design parameters p. It can be written according to eq 2d, where the L-matrix summarizes the information on the kinetics of the fluxes and contains the kinetic coefficients. The diagonal elements of the L-matrix reflect the main transport effects. In this representation, cross-effects such as the Soret and the Dufour

10539

Figure 4. Decision structure arising on level 2 (left). Example: convective heat transport (right).

effect can also be taken into account. The detailed composition of the L-matrix is given in the literature, and examples how to use the generalized Maxwell-Stefan equations to model the fluxes are given in several textbooks on transport phenomena (e.g., ref 46). The transport kinetics of these diagonal effects are often simplified so that the generalized Maxwell-Stefan equations do not have to be considered in detail. The methodology, however, is totally capable of taking all effects into account. The reaction kinetics, additional equality, and inequality constraints, as well as the initial and final conditions, are the same as on level 1 except that all transport kinetics are now summarized in eq 2d. Instead of optimizing the flux directly, specific terms of the transport kinetics are now optimized. Here, various possibilities exist of how to control the fluxes, as discussed below. In case a flux is not controlled, only specific design parameters and initial values for the states of the environment can be optimized. Control Options. On level 2, the specific surfaces for every flux, the transport mechanism, the control variable, and the fluid flow regime must be chosen using a systematic decision structure. The specific surfaces also include different types of catalyst packings in the case of heterogeneously catalyzed reactions. The decision structure for level 2 is shown in Figure 4. Here, the different control options for convective heat transfer in an empty tube are illustrated as an example. The laminar and the turbulent flow regime must be distinguished since different transport laws for the convective heat transfer apply. The heat transfer coefficient in the laminar flow regime is typically proportional to D-1, while it is approximately proportional to D-1/5 in the turbulent flow regime. The S/V-ratio is also proportional to D-1 for most channel geometries. Hence, the geometry is an adequate control variable for the convective heat transport. The heat transfer coefficient depends on different fluid properties, but these are usually not suitable to control the flux. In principle, the heat flux can also be realized by heat conduction, radiation, or any combination of these three transport mechanisms. The different transport mechanisms obey different physical laws which can be considered on level 2. The environment can be considered as a “service phase” for every flux. In this context, the heat flux can be realized by a cooling or heating phase. Component fluxes can be obtained by an additional fluid. An idealized service phase can be chosen (e.g., and extraction phase can have an unlimited capacity) or it can be constrained by process parameters (e.g., a synthesis gas with a specific H2/CO-ratio must be used). If the component flux is provided by a membrane, a model for the flux through the membrane must be included. In this case, specific terms such as the partial pressure difference across the membrane can be optimized along the path of the fluid element. If the components should be dosed in a discrete manner, the program-

10540


ming in the Lagrangian formulation is straightforward using a reference mass of the original fluid element and adding some extra mass of the dosed component. In case a specific component is removed from the mixture, the membrane model must take nonideal separation into account. If technical work is provided, a model for the technical work flux must be included as well. The model depends strongly on the source of the technical work chosen. In general, the states of the environment xe can be optimized within specific limits xLe e xe e xUe . The lower and upper limits arise from physical and technical constraints such as a minimum cooling temperature and a maximum allowed catalyst temperature. In addition, limits for the temperature, pressure, and their gradients arise from the available construction material. The specific surface areas for the individual fluxes depend on the channel diameter and geometry. To change the specific surface area of different fluxes individually, the flux must be applied to a different edge of the reaction channel. For this purpose, a channel geometry must be assumed depending on the number of the individually controlled specific surface areas. The S/V-ratio of each flux can be optimized within certain limits ALi e Ai e AUi , which depend on the assumed channel geometry and fluid flow regime. The individual surfaces are linked by each other depending on the chosen reactor geometry. Summing up, the second level yields information on which control variable is most suitable to obtain the desired flux profiles. A detailed catalog of requirements for the optimal technical reactor is derived from the investigations of the best integration and enhancement concept and from the optimal profiles of the control variables. Different control variables, transport mechanisms, options to provide the catalyst, and channel geometries giving rise to specific exchange areas for heat and mass are compared, and the best combination is determined. A result of this level may also be that controlling the fluxes is very limited due to severe mass and energy transport limitations. In this case, the potential of the reactor, where the investigated flux is ideally controlled, cannot be exploited, and the result might only be as good as the optimized technical reference reactor so that further investigations on how to approximate the profiles of the control variables are not necessary. However, if the mass and energy transport limitations are not severe for a specific set of control variables, it is worthwhile to investigate how to achieve the optimal control variable profile with a technical approximation on the third level. 2.3. Level 3: Optimal Technical Reactor. On the basis of the profiles of the best control variables, a technical approximation is developed on the third level of the proposed methodology. The developed reactor represents the technically optimal reactor and approximates the optimal state space trajectory of the reaction system. On this level, the reactor design is fixed, and it is determined how the control variable profiles can be realized as close as possible. Several possible technical approximations of the control variable profiles may exist since the transformation from control variable profiles to the reactor design is not unique. In order to quickly evaluate the different possible technical approximations, the third level can be divided into two steps. In a first step, all technical approximations are calculated using simple models. The problem formulation depends on the technical approximation, and hence no generic model formulation is given. In general, equations considering the geometry such as equations for the velocity and balance equations for all used utilities must be added. In the case of a continuous operated reactor, the local form of the balance equations is used to

account for radial profiles, bypassing, and dispersion effects. In addition, the size of the reactor must be determined, which is easier in the local form of the balance equations. The difference in the objective function value between the technical approximation and the ideal control profiles (level 2) quantifies the losses due to nonidealities in the technical realization of the control profiles. In a second step, the best technical approximation can be investigated using more detailed models. These models should include effects resulting from nonuniform radial temperature, concentration, and velocity distributions. This procedure allows us to quantify the losses due to radial gradients, nonideal flow profiles, and dispersion effects. The results of the more detailed reactor models are often very similar to the results of the simple reactor models. The qualitative trends are most often the same, but the detailed models provide additional information, e.g., on hot spots due to radial temperature gradients. On the third level of the proposed methodology, the dimensionality of the optimization problem increases due to the balance equations for the environment phase, but the degrees of freedom drastically reduce. The optimal control variable profiles for the detailed models can be obtained by rigorous optimization using the results of the simplified reactor models as starting values. In case the detailed models are too complex for optimization, the performance of the reactor can also be very efficiently investigated using the optimal control variable profiles of the simple models and performing a simulation study using the detailed reactor models. Summing up, on level 3, optimal technical reactors are developed. The reactor design is based on the results of level 1 and level 2 and approximates the optimal route in state space as good as technically possible. The proposed methodology allows us to screen different concepts and to evaluate the potential of distributed control variables using simple models. More extensive numerical computations with detailed reactor models have only to be performed for a small number of selected cases. This reduces the computational effort significantly. 3. Numerical Solution Approach The arising dynamic optimization problems can be solved using the Hamilton-Jacobi-Bellmann (HJB) equation, the Pontryagin-Minimum-Principle (PMP), or NLP based methods (simultaneous or sequential approach). To solve the optimization problems of the individual levels, numerical schemes must be applied regardless of which solution method is used. For the discussion of the advantages and disadvantages of each method, we refer to the literature (e.g., ref 47). For the following calculations, the simultaneous approach using an orthogonal collocation on finite elements48 is applied using CONOPT 3.14 G under AMPL on a PC with an Intel(R)Core(TM)2 Duo CPU E6850 with 3.00 GHz (calculation on a single CPU), a cache size of 4096 KB, a memory of 2 GB, and a Linux operating system. One advantage of the chosen solution approach is that it can easily handle jumps in the control variable profiles. For further explanation of the method, refer to the literature.31,47 The computing time is very short for all models since the model formulation in terms of elementary process functions is very efficient. In the case of the investigated example, the CONOPT time for 200 equally distributed finite elements and three collocation points per finite element is 16.7 s for the optimal heat flux case. On the third level, the computation takes about 121.6 s for the investigated reactor setup. An optimization based on detailed reactor models would take much longer; e.g.,


optimizing a two-dimensional reactor model often requires computing times in the order of hours. Even a simulation study of several reactors using detailed models would require higher computational effort. In the proposed methodology, only the most promising reactor setups are modeled in detail. Hence, the computational effort is significantly reduced compared to a brute force approach, where the results of detailed reactor models of every reactor setup are compared. It should be noted that different starting values are used for the different calculations, which may have a major influence on the computation time. Nevertheless, the stated computing times give a convincing impression of the potential savings in computational effort by applying the proposed methodology compared to optimizing different detailed reactor models. 4. Example: Reactor Design for SO2 Oxidation To demonstrate the proposed methodology, the development of an SO2 oxidation reactor is presented as an example. The SO2 oxidation is an exothermic equilibrium reaction, and it is representative for a whole class of important industrial reactions. 1 SO2 + O2 h SO3 with ∆rH ) -99 kJ/mol 2

(3)

On an industrial scale, the SO2 oxidation is performed in adiabatic fixed bed reactors with intermediate cooling. To achieve high conversions, intermediate absorption of SO3 is performed downstream of the third bed. For further process details, refer to the literature.49 In this example, the objective is the minimization of the residence time for a given catalyst and given initial and outlet conditions. The given inlet conditions are typical for a plant based on a sulfur furnace (T0 ) 710 K, p0 ) 2 bar, xSO2,0 ) 0.11, xO2,0 ) 0.10, xSO3,0 ) 0.01, xN2,0 ) 0.078). The outlet conditions are typical values upstream from the intermediate absorption (Tf ) 740 K, pf g 1.5 bar, X g 0.95). Instead of specifying the outlet composition, a lower limit for the conversion is enforced. 4.1. Level 1. On the first level, the reaction concepts are defined by specifying which fluxes are integrated and which fluxes are optimized. The residence time depends on the used catalyst and on the catalyst density. In principle, different catalysts can be considered, but for simplicity only the standard vanadium catalyst49 is discussed here. Furthermore, different ways to provide the catalyst can be taken into account, e.g., catalytic wall reactors, foam structures, or the more traditional randomly packed bed. In this case, different catalytic packings with regard to their pressure drop, thermal conductivity, and catalyst density must be compared on level 2. If the minimization of the residence time is investigated, a high catalyst density is favored. Since packed beds feature the highest packing density, only packed beds are investigated in this example, and a constant void fraction of ε ) 0.4 is chosen, which is a recommended value for spherical catalyst particles.50 The use of microreactor technology can be explored by considering other ways to provide the catalyst; however, alternatives such as catalytic wall reactors cannot be applied for the chosen catalyst. If more active platinum catalysts are used, microreactor technology might be favored over conventional packed beds due to the high heat of reaction.51 The chosen reaction rate model is based on a reaction rate52 and thermodynamic data53 presented in the literature. A pseudohomogeneous reactor model is used. The reaction rate

10541

in terms of partial pressures, temperature, and component partial pressures is given by eq 4.

[ ( )]

r ) kr

pSO2 pSO3

pO2 -

pSO3

2

pSO2Kp

(4a)

( -97T782 - 110.1 ln(T) + 848.1)

kr ) 9.8692 × 10-3 exp

(4b)

( 98RT359 )

Kp ) 1.013 × 105 exp

(4c)

The reaction rate for the chosen catalyst as stated in eq 4 depends on the pressure, on the individual component mole fractions of SO2, SO3, and O2, and on the temperature. To minimize the residence time, these states can be influenced for maximizing the reaction rate. After specifying the states, which are of interest for the optimization, suitable fluxes to control these states must be identified. Influencing the individual component mole fractions requires dosing and removal of individual components. In the case of the SO2 oxidation, dosing of SO2 and O2 is not practical since an additional purification step, which is very expensive compared to the costs of the product, is required. From a technical point of view, the in situ removal of SO3 is hard to perform. The produced SO3 cannot be removed in situ by condensing, since the condensing temperature is far below the reaction temperature. In addition, in situ extraction with sulfuric acid is not possible at reaction conditions. A selective membrane for SO3 separation from SO2 is not known so far. The consecutive hydration to sulfuric acid is not possible in the temperature range of the reaction due to the chemical equilibrium of the hydration reaction. In addition, the reaction rate expression from the literature given in eq 4a is not suitable to investigate the influence of SO3 removal, since SO3 appears only in the denominator of the reaction rate. Hence, in situ removal of SO3 and dosing of reactants is not considered in our example. Nevertheless, the potential of dosing and selective removal can be calculated very efficiently using the proposed methodology. This is of special interest for selectivity problems and will be demonstrated in our future publications. In principle, the pressure can be used to influence the reaction rate. This is especially relevant in multiphase systems, where the solubility of gas phase components in the liquid phase strongly depends on the pressure. In the case of the SO2 oxidation, the pressure is not suited as a control variable. Increasing the pressure always accelerates the reaction rate and shifts the equilibrium conversion according to Le Chatelier’s law to the product side. In addition, a lower pressure yields a larger reactor volume, which is also not desired. Hence, maximum pressure is favored. The influence of the temperature on the reaction rate is not trivial since the equilibrium conversion and the rate constants are strong functions of the temperature. Hence, investigating how an optimal temperature profile should look like and how it can be realized is of great interest. To reduce the dimensions of the model, the reaction rate and the component mass balance are formulated in terms of conversion. The conversion of the key component SO2 is defined according to eq 5. X)

nSO2,0 - nSO2 nSO2,0

(5)

10542


The component mass balance is written in terms of conversion without dosing and removal of components. r · Fcat · R · T r · Fcat · Vgas 1 - ε dX ) · ) · dt nSO2,0 ε p 1 1 1-ε - X · xSO2,0 2 ε

(

)

(6)

A catalyst density of Fcat ) 1330 kg/m3 and the ideal gas law for the gas volume Vgas are used. The ratio of catalyst volume to gas volume is given by eq 7. Vcat 1-ε ) Vgas ε

(7)

For the calculations, the reaction rate is expressed in terms of conversion X, temperature, pressure, and initial composition.

r ) kr

xSO2,0(1 - X) xSO3,0 + xSO2,0X

[

·

(x

1 - xSO2,0X p 2 1 1 - xSO2,0X 2

)

O2,0

(

xSO3,0 + xSO2,0X xSO2,0(1 - X)Kp

)

]

2

(8)

A lower catalyst activation temperature of 633 K, which is the melting temperature of an industrially used promoted vanadium catalyst, is considered. An upper catalyst temperature of 923 K due to the catalyst stability is also included.54 The energy balance in terms of temperature is simplified by assuming no technical work, no dosing or removal of components, and negligible influence of the pressure change on the temperature change. A linear mixing rule for the heat capacity is assumed using a heat capacity correlation for the individual components obtained from the literature.53 Using the ideal gas law to substitute the total amount of substance, the simplified energy balance in terms of temperature is given by eq 9. dT 1-ε 1 ) - q · Aq + ∆rH · r · Fcat · dt ε F · cp

(

)

NC

∑ i)1

xi,0 + xi )

xi · M i · p R·T

νi · X · xSO2,0 |νSO2 |

1 - 1/2 · X · xSO2,0

q ) R · (T - Te)

(9)

(10)

(11)

Reaction, dosing of cooled feed gas, and cooling are the fluxes of interest for this example. Now, it has to be determined how to combine these fluxes and which fluxes should be controlled. Here, three reaction concepts are of major interest. Case 1: Adiabatic reaction with intermediate cooling. This is the ideal case of an adiabatic staged reactor, which is regarded as the technical reference case. Reaction and cooling are applied in series. No flux is controlled in a strict manner. Case 2: Polytropic reaction with constant cooling temperature. Reaction and cooling are integrated. This case represents a conventional tube bundle reactor. No flux is controlled in a strict manner.

(12)

4 (empty tube) D

(13)

4 (packed bed) D·ε

(14)

Aq )

Aq )

The density and the individual mole fractions are given by eq 10 and eq 11. F)

Case 3: Control of the heat flux. Reaction and cooling are integrated, and the heat flux is controlled to investigate the potential of distributed cooling. In principle, other sets of optimal controlled fluxes are possible. In addition, it is thinkable to control several fluxes simultaneously. However, the advantage of controlling several fluxes depends on the specific problem, but in general the advantage decreases with increasing number of controlled fluxes. For example, the catalyst density can be scaled in combination with optimizing the heat flux. In the case of the SO2 oxidation, a maximum catalyst density is desired, and hence the catalyst density will be maximal most of the time. Adiabatic reaction with addition of cold feed gas at certain mixing points (cold shot) is also possible. A solution of this problem can be found in the literature.55 To keep the complexity of the example manageable, other reaction concepts are not further considered. The potential of controlling the heat flux (case 3) is compared to an optimal adiabatic fixed bed (case 1) and an optimal polytropic PFTR with constant cooling temperature (case 2) as reference reactors. For the adiabatic fixed bed and polytropic PFTR, further assumptions considering the heat removal must be made. The bounds for the temperature of the environment phase are chosen to be 700 K e Te e 740 K. The upper limit is given by the outlet temperature, and the lower limit is chosen to limit the temperature gradients. In addition, a tubular geometry for the heat exchanger tubes with a minimum required diameter DL ) 5 cm is considered. The heat transfer coefficient is assumed to be constant with a value of R ) 40 W/(m2K). The heat flux is modeled according to eq 9 with the specific surface area for the heat exchange of the gas according to eq 13 or eq 14, respectively.

Results. The adiabatic staged reactor (case 1) consists of three reaction and three cooling segments. The optimization parameters are the timing of intermediate cooling, the cooling temperature, and the tube diameter of the coolers. Since the minimal residence time problem is considered, the lower bounds for the cooling temperature and the tube diameter are chosen by the optimizer. A similar design was already proposed by Aris,56 who analyzed the SO2 oxidation with further simplifying assumptions by analytical methods. In Figure 6(a), the trajectories of the different reaction concepts are shown in state space representation. The minimum residence time of an adiabatic staged reactor is about 2010 ms. The total residence time consists of 955 ms for reaction and 1055 ms for cooling. The optimization parameters of the polytropic reactor (case 2) are the cooling temperature and the tube diameter using the same assumptions as in the case of the adiabatic fixed bed reactor. The optimal values are Dopt ) 5 cm and Te,opt ) 727.9 K. The temperature of the environment phase is above the inlet temperature, and the fluid is heated in the first part of the reactor. The residence time is 714 ms, which is a reduction by 64.5% compared to the optimized technical reference. Integration of reaction and cooling with constant cooling temperature significantly reduces the required residence time.


If the heat flux can be adjusted at will (case 3), the residence time is 591 ms. This is a reduction by 70.6% based on the reference case and quantifies the potential of ideal heat flux control. In this case, the fluid is directly heated to the optimal reaction temperature. To avoid an ill-formulated optimization problem, the inlet temperature is only bounded by TL e T0 e TU, and then the heat flux can be bounded by numerical meaningful values. In principle, the inlet temperature could be regarded as optimization parameter for any reaction concept, but in this example it is assumed that the inlet temperature is fixed. Only for the ideal heat flux case, the inlet temperature is free to avoid numerical problems. These results show that integration of cooling and reaction will yield a significant reduction of the residence time. The perfect control of the heat flux reduces the residence time even further. How this heat flux profile can be obtained is investigated on the second level. In addition, the losses due to limited heat transfer are quantified. At this early point, the methodology indicates which reaction concept is best suited for the reaction system. If the savings by optimizing the heat flux along the reactor are considered to be not large enough, the investigations can be stopped, and the tube bundle reactor with constant cooling temperature can be developed further. 4.2. Level 2. On the second level, several possibilities on how to obtain the desired heat flux profile are compared. As discussed, several transport mechanisms exist to provide the heat flux. In this case, only cooling and heating by a service phase outside the tubes is considered. In packed beds, different models can be used to describe the heat transfer resulting from convection and conduction. Either a simple one-dimensional model (R-model), a two-dimensional model with one parameter (λr-model), or a two-dimensional model with two parameters (λr, R-model) is used.57,58 The more detailed models describe the radial temperature profile, which cannot be calculated using the one-dimensional model. However, the one-dimensional model can still yield very accurate solutions; the mean deviation of the fitted model was only 3.75% for the examples investigated by De Wasch and Froment.57 In the case of packed beds, the correlation for the heat transfer (eq 15 and eq 16) can be applied for a wide range of Re numbers. q ) k · ∆T ) k · (T - Te) k ) k0 + 0.024

λ · F · Vs λ · Re ) k0 + 0.024 dp µ F · V s · dp e 1000 for 30 e Re ) µ

(15)

In the first case, to control the heat flux, the temperature of the environment phase is taken as control function while the tube diameter and the initial velocity are optimization parameters. In the second case, the tube diameter is the control function, while the temperature of the environment phase and the initial velocity are optimization parameters. The environment temperature bounds are chosen to be TLe ) 600 K and TUe ) 900 K. The boundaries for the tube diameter depend on the flow regime, the velocity, and the particle diameter. To use the heat transfer correlation, specific bounds for the Re number must be kept according to eq 16. The Re number depends on the density, velocity, and viscosity of the fluid as well as on the particle diameter. Indirect lower limits for the particle diameter and the velocity are given by including the maximum allowed pressure drop. Instead of solving the rigorous momentum balance, the Ergun equation for spherical particles is used to estimate the pressure drop.60

(

The range for the static contribution is approximately 6 W/(m2K) e k0 e 16 W/(m2K) for the SO2 oxidation catalyst.57 For reasons of simplicity, a mean value for the static contribution of k0 ) 11 W/(m2K) is chosen. At typical SO2 oxidation conditions, the contribution of the Re term is between 17 W/(m2K) and 35 W/(m2K). Hence, both terms are important for the modeling of the heat transport coefficient. The thermal conductivity and the viscosity of the gas are assumed to be constant (λ ) 0.0476 W/(m · K), µ ) 3.322 × 10-5 kg/(m · s). These values are calculated at the inlet conditions using temperature-dependent correlation for the pure component thermal conductivity and viscosity and the model of Wassiljewa and Wilke, respectively, to calculate the mixture properties.59 The heat conductivity, the viscosity, and the density influence the heat transfer, but they are not suited to control the heat transfer in this example.

)

VsF(1 - ε) ∂p µ(1 - ε)2 ) - 150 + 1.75 Vs 2 2 ∂z dε d ε3 p

p

(17)

Under steady state conditions, eq 18 is used to convert the Ergun equation eq 17 into its substantial formulation eq 19. dp ∂p ) Vi dt ∂z

(

(18)

)

VsF(1 - ε) Vs2 dp µ(1 - ε)2 ) - 150 + 1.75 dt ε d2ε3 d ε3 p

p

(19)

The velocity for a constant mass flux can be expressed as a function of the initial velocity Vs, 0, the cross sectional area, and the fluid density (refer to eq 20). The initial velocity is taken as an optimization parameter. The density is calculated according to eq 10. Vs ) Vs,0

A0 F0 A F

(20)

To limit the effect of bypassing, the upper limit for the catalyst particles is indirectly defined by a minimal tube-to-particle diameter ratio.61 D g 10 dp

(16)

10543

(21)

The upper limit of the catalyst particle diameter is defined by the onset of pore diffusion limitation. This can be considered keeping the Thiele modulus of the reaction below a certain value at every point along the reaction coordinate (typically Φ e 0.4). Results. Figure 6(c) and Figure 6(d) illustrate the profiles for the temperature of the environment phase and the tube diameter as control functions. In the case of the environment temperature control, the residence time is 620 ms. In the case of the diameter control, the residence time is 677 ms. Both state space trajectories can be compared with the ideal case in Figure 6(b). As can be seen in Figure 6(c), the temperature of the environment phase is at its maximum at the beginning of the reaction. As soon as the state space trajectory of the ideal heat flux case is reached, the environment temperature switches from heating to cooling. The environment temperature jumps to 640 K and then increases to 770 K. Afterward, it decreases slowly to 736 K. As soon as the required conversion is reached, the environment temperature is at its minimum to reduce the

10544


Figure 5. Possible technical approximation of the optimal process route.

gas temperature as fast as possible to the required outlet temperature of 740 K. It should be noted that the outlet temperature is not on the optimal trajectory. In Figure 6(d), the optimal profile for the tube diameter is shown. The optimum constant temperature of the environment phase is 738 K. For the first finite element, the tube diameter is at its maximum before it jumps to its minimum in the second finite element. Since the temperature of the environment phase is above the inlet temperature, the section with the minimal tube diameter heats the fluid to the environment temperature as fast as possible. The boundaries for the minimum and maximum tube and catalyst particle diameter are indirectly defined by the Re number range, the pressure drop, and the D/dp-ratio. Since a tube diameter of 0.0375 m is at the lower limit defined by the minimum D/dp-ratio, the tube diameter can only be further reduced if the catalyst particle diameter is decreased. In addition, a tube diameter of 0.146 m is at the lower limit of the Re number constraint. Hence, the tube diameter can only be increased above 0.146 m by raising the initial velocity further in the case of a constant catalyst particle diameter. However, both increasing the inlet velocity and decreasing the catalyst particle add to the pressure drop, which reduces the reaction rate. The negative effect of a lower pressure is more severe than the positive effect of using a larger range for the tube diameter. As soon as the fluid temperature reaches the environment temperature, the diameter jumps to its maximum to approximate a quasi adiabatic section. Before the optimal temperature is reached, the diameter reduces to avoid a temperature overshoot. Afterward, the diameter follows a specific profile to keep the state space trajectory as optimal as possible. At the end, a minimal tube diameter is used to decrease the fluid temperature as fast as possible. Even for this well-known example, the control variable profiles were not known a priori. This underlines the necessity of rigorous optimization for the reactor design instead of using heuristics. Due to the fast heating section, the objective function value for the environment temperature control case is better than for the tube diameter control case. Therefore, the environment temperature control case is chosen to be further investigated on the third level. 4.3. Level 3. In this example, three reactor segments seem to be well suited to approximate the optimal environment temperature profile. In the first segment, a constant maximum temperature of the environment phase is optimal. In the second segment, a cooling temperature, which increases along the path of the fluid element, is required. This can be realized with a co-current heat exchanger. The third segment can be approximated either by constant cooling temperature or by slightly decreasing cooling temperature. The former can be realized by an evaporating fluid and the latter by a counter-current heat exchanger. The co- and counter-current heat exchangers are modeled using a simplified energy balance for the cooling shell. It is assumed that the flow of the environment fluid is constant in each individual segment. In addition, for the environment fluid

constant density, no axial diffusion of heat and mass, no radial gradients, constant pressure, no evaporation, no reaction, no technical work, constant heat capacity, and steady state operation are assumed. The energy balance of the environment phase is reduced to a one-dimensional energy balance in terms of temperature using the above simplifications and introducing the cross sectional area of the environment channel Ae. On the third level, the local form is used since a concrete continuously operated reactor is developed. Vi,e · Fe · cp,e

∂Te π·D )q ∂z Ae

(22)

This equation is further simplified by merging all constants according to eq 23. Ke )

π·D Ae · Vi,e · Fe · cp,e

(23)

∂Te ) Ke · q ∂z

(24)

Equation 24 can be used to optimize co-current cooling (i.e., Ke g 0), counter-current cooling (i.e., Ke e 0), and cooling at constant temperature (i.e., Ke ) 0) without predefining the segments. The constant is chosen to be within meaningful bounds (-1 m · K/W e Ke e 1 m · K/W). The energy balance and the mass balance for the fluid are converted to the local form using eq 18, and under steady state conditions they are given by eq 25 and eq 26. (4/D) · q + ∆rH · r · Fcat · (1 - ε) ∂T )∂z F · c p · Vs

(

(25)

)

r · Fcat · R · T · (1 - ε) ∂X 1 1 ) · - X ∂z p · Vs xSO2,0 2

(26)

The residence time is optimized using eq 27. min τ ) min

∫

L

0

1 dz Vi

(27)

Here, L is the total length of all three segments. The length of the segments, the cooling constant Ke, and the inlet temperature of each segment are the optimization parameters. In principle, each segment can also have a different tube diameter, but the advantage of such additional parameters is very limited. More detailed models, e.g., a two-dimensional reactor model taking nonuniform temperature and flow profiles as well as radial and axial dispersion into account, can be applied to further investigate the chosen reactor setup. Here, we restrict ourselves to a technical approximation described by a simple reactor model. Results. A design drawing of the technical approximation with the optimal values for the constants (Ke,i with i ) 1, 2, 3) and the according inlet and outlet temperatures is shown in Figure 5. The total length of the reactor is 0.83 m with a total


10545

Figure 6. Results of different levels. (a) Level 1: State space. (b) Level 2: State space. (c) Level 2: Temperature control. (d) Level 2: Diameter control. (e) Level 3: Temperature profile of each cooling segment and the gas phase. (f) Level 3: State space.

residence time of 623 ms. This is a reduction of 69% compared to the optimized reference case. The environment temperature profile and the gas temperature profiles of the technical approximation are given in Figure 6(e). It can be observed that the state space trajectories of the technical approximation and of the case where the environment temperature is controlled look very similar (refer to Figure 6(f)). This result demonstrates that a technical approximation based on the optimal route can be easily developed and that the proposed methodology yields new and innovative reactor designs. 4.4. Summary of the Example. In this example, the steps for the design of an optimal technical SO2 oxidation reactor are presented. The applied decision structure is shown in Figure 7. After choosing an objective function, it is decided which fluxes to integrate and which to control. Here, controlling the heat flux is the best choice to minimize the residence time. The possibility of controlling other fluxes such as individual component fluxes or compression work is discussed but not

considered in the calculations. Three different reaction concepts for dealing with a reaction and a heat flux are compared. Besides the classical reaction concepts, which are adiabatic reaction with intermediate cooling and reaction with cooling at a constant cooling temperature, a new reaction concept, where the heat flux is ideally controlled, is investigated. On the first level, the advantage of integrating the heat flux and additionally of controlling the heat flux ideally is quantified. The residence time of the optimal route is 591 ms. In the case of the optimized plug flow reactor with isothermal cooling temperature, the minimal residence time is 714 ms. The optimized adiabatic fixed bed reactor has a residence time of 2010 ms. On the second level, it is investigated how the desired heat flux profile can be obtained. Here, the heat flux can be controlled by changing the temperature of the environment phase or by changing the tube diameter giving rise to a residence time of 620 ms and 677 ms, respectively. Controlling the temperature of the environment phase is better suited to control the heat

10546


Figure 7. Decision structure for SO2 oxidation reactor design.

flux in the desired way mainly due to a faster heating section. The losses in terms of the objective function due to limited energy transport can be quantified to 29 ms. On the third level, a technical approximation of the best suited control variable profile is developed. Such a technical approximation is possible with existing apparatuses and gives rise to a new reactor setup for the SO2 oxidation, which reduces the residence time compared to an optimized technical reference by almost 69%. The losses due to a nonideal control profile are merely 3 ms. In the case of the SO2 oxidation, an appropriate combination of existing apparatuses seems to be well suited to obtain a reasonably good technical approximation. However, the methodology also guides the development of new technology since the advantages of developing new apparatuses can be quantified. The example and the objective function were chosen to demonstrate the methodological approach. Of course, further investigations regarding different catalysts, different types of packings, and different objective functions are possible, but this is not the scope of this example. Furthermore, several additional constraints and more detailed models could be considered. Preliminary considerations regarding the investigated system and chemical engineering knowledge are used to reduce the number of possible reaction concepts so that not every possible reaction concept has to be evaluated. 5. Conclusion and Outlook It can be concluded that the methodology is well suited to choose the best reaction concept and yields an optimal reactor design. Coming back to the requirements regarding the methodology defined in the Introduction, the following can be concluded: 1. The methodology defines an optimal route in state space with respect to a certain objective and a chosen integration and enhancement concept on the first level. The best reaction route is obtained by optimizing the fluxes along the reaction coordinate. That way, the optimal reaction conditions are obtained at every point in the reactor. In addition, several reaction routes can be evaluated and compared. 2. Classical reactor concepts as well as process intensification methods are directly included in the proposed model based approach. The physical principle of the PI methods is modeled on the level of the balance equations using simple approaches to investigate their potential. Here, the whole space of physi-

cochemical phenomena can be investigated. The potential and benefits of the integration of additional fluxes and enhancement of certain fluxes are quantified. By calculating optimized technical reference reactors, the potential of each concept is quantified, and the most promising concepts are further developed. The calculation times on the first and second level are very short, and hence the methodology allows a rapid screening of different integration and enhancement concepts. 3. The best control variables to enforce the optimal route are obtained by a systematic decision structure on the second level. 4. The losses of the best technical reactor compared to the optimal route in state space are quantified. In addition, the overall losses can be quantitatively divided into losses due to limited mass and energy transport, losses due to nonideal approximation of the control profiles, and losses due to radial gradients, dispersion effects, and a nonideal flow field. Only few evaluations of technical approximations are necessary, which take longer calculation times, since the dimensionality increases. 5. The design methodology is suited for batch and continuously operated reactors since the Lagrangian formulation of the balance equations is used at the beginning. In addition, the methodology is independent of existing apparatuses and allows the development of innovative reactor concepts. The method includes classical reactor concepts as well as innovative reactor concepts. Hence, it is also well suited to choose the most efficient reactor from an economical point of view. The methodology will be extended and applied to selectivity problems and multiphase systems in upcoming publications. New modeling aspects will arise in the case of multiphase systems. If a multiphase model is used, the transfer between the different phases can be optimized. If the reaction occurs only in one phase, the other phases can be regarded as service phases with the purpose of realizing the best reaction conditions in the reacting phase. In case the phases have different velocities, the modeling framework has to be extended or the local form has to be used. If the catalyst is considered as a separate phase, also reverse flow reactors can be described, which cannot be described with the homogeneous model chosen so far. The potential of integrated reaction separation systems for the reactor performance can already be investigated by selective removal of certain components. Furthermore, different possibilities of how to provide heterogeneous catalysts will also be taken into account in our future work. Calculating the sensitivity of the objective with respect to small changes of the optimal route is important. In addition, several routes with similar objective values may exist, and the question must be answered how to find these additional routes. Robust optimization may help to find optimal routes with respect to errors in the model parameters. These aspects will be considered in the further development of the proposed methodology. Acknowledgment One of us (A.P.) thanks the International Max Planck Research School Magdeburg for financial support. The authors thank Nadine Novotny and Florian Karst for their contributions in the frame of their study research project and master thesis, respectively, and Dr. Peter Heidebrecht for the discussion of optimization problems. Nomenclature 3 Ai ) Specific exchange area for flux i (m2/mgas ) 2 3 Aq ) Specific heat exchange area (m /mgas)

Ind. Eng. Chem. Res., Vol. 49, No. 21, 2010 D ) Tube diameter (m) Di ) Diffusion coefficient of component i (m2/s) dp ) Catalyst particle diameter (m) h ) Specific enthalpy (kJ/kg) hi ) Specific enthalpy of component i (kJ/kg) k ) Heat transfer coefficient (W/(m2 · K)) Kp ) Equilibrium constant (Pa) kr ) Reaction rate constant (mol/(kgcat · s · Pa)) L ) Reactor length (m) νi ) Stoichiometric coefficient of component i (-) r ) Reaction rate (mol/(kgcat · s)) τ ) Residence time (s) Te ) Temperature of the environment phase (K) Vi ) Interstitial velocity (m/s) Vs ) Superficial velocity (m/s) wi ) Weight fraction (-) wt ) Technical work (W/m3) xi ) Molar fraction (-)

Literature Cited (1) Froment, G. F.; Bischoff, K. B. Chemical Reactor Analysis and Design, 2nd ed.; John Wiley and Sons: New York, 1990. (2) Emig, G.; Klemm, E. Technische Chemie, 5th ed.; Springer: Heidelberg, 2005. (3) Schembecker, G.; Droge, T.; Westhaus, U.; Simmrock, K. H. Readpert - Development, Selection and Design of Chemical Reactors. Chem. Eng. Process. 1995, 34, 317–322. (4) Horn, F. Attainable Regions in Chemical Reaction Techniques. In The Third European Symposium on Chemical Reaction Engg.; Pergamon: London, 1964; pp 293-303. (5) Glasser, D.; Hildebrandt, D.; Crowe, C. A Geometric Approach to Steady Flow Reactors - The Attainable Region and Optimization in Concentration Space. Ind. Eng. Chem. Res. 1987, 26, 1803–1810. (6) Hildebrandt, D.; Glasser, D. The Attainable Region and Optimal Reactor Structures. Chem. Eng. Sci. 1990, 45, 2161–2168. (7) Hildebrandt, D.; Glasser, D.; Crowe, C. M. Geometry of the Attainable Region Generated by Reaction and Mixing - With and Without Constraints. Ind. Eng. Chem. Res. 1990, 29, 49–58. (8) Feinberg, M.; Hildebrandt, D. Optimal Reactor Design from a Geometric Viewpoint - I. Universal Properties of the Attainable Region. Chem. Eng. Sci. 1997, 52, 1637–1665. (9) Balakrishna, S.; Biegler, L. T. Constructive Targeting Approaches for the Synthesis of Chemical Reactor Networks. Ind. Eng. Chem. Res. 1992, 31, 300–312. (10) Balakrishna, S.; Biegler, L. T. Targeting Strategies for the Synthesis and Energy Integration of Nonisothermal Reactor Networks. Ind. Eng. Chem. Res. 1992, 31, 2152–2164. (11) Balakrishna, S.; Biegler, L. T. A Unified Approach for the Simultaneous Synthesis of Reaction, Energy, and Separation Systems. Ind. Eng. Chem. Res. 1993, 32, 1372–1382. (12) Lakshmanan, A.; Biegler, L. T. Synthesis of Optimal Chemical Reactor Networks. Ind. Eng. Chem. Res. 1996, 35, 1344–1353. (13) Rooney, W. C.; Hausberger, B. P.; Biegler, L. T.; Glasser, D. Convex Attainable Region Projections for Reactor Network Synthesis. Comput. Chem. Eng. 2000, 24, 225–229. (14) Milne, D.; Glasser, D.; Hildebrandt, D.; Hausberger, B. Application of the Attainable Region Concept to the Oxidative Dehydrogenation of 1-Butene in Inert Porous Membrane Reactors. Ind. Eng. Chem. Res. 2004, 43, 1827–1831. (15) Nisoli, A.; Malone, M. F.; Doherty, M. F. Attainable Regions for Reaction with Separation. AIChE J. 1997, 43, 374–387. (16) Feinberg, M.; Ellison, P. General Kinetic Bounds on Productivity and Selectivity in Reactor-Separator Systems of Arbitrary Design: Principles. Ind. Eng. Chem. Res. 2001, 40, 3181–3194. (17) Agarwal, V.; Thotla, S.; Mahajani, S. M. Attainable Regions of Reactive Distillation. Part I. Single Reactant Non-Azeotropic Systems. Chem. Eng. Sci. 2008, 63, 2946–2965. (18) Achenie, L. K. E.; Biegler, L. T. Algorithmic Synthesis of Chemical Reactor Networks Using Mathematical-Programming. Ind. Eng. Chem. Fundam. 1986, 25, 621–627. (19) Achenie, L. K. E.; Biegler, L. T. A Superstructure Based Approach to Chemical Reactor Network Synthesis. Comput. Chem. Eng. 1990, 14, 23–40.

10547

(20) Kokossis, A. C.; Floudas, C. A. Synthesis of Isothermal Reactor Separator Recycle Systems. Chem. Eng. Sci. 1991, 46, 1361–1683. (21) Kokossis, A. C.; Floudas, C. A. Optimization of Complex Reactor Networks. 2. Nonisothermal Operation. Chem. Eng. Sci. 1994, 49, 1037– 1051. (22) Hillestad, M. A Systematic Generation of Reactor Designs. II. Nonisothermal Conditions. Comput. Chem. Eng. 2005, 29, 1101–1112. (23) Aris, R. The Optimal Design of Chemical Reactors - A Study in Dynamic Programming; Academic Press: New York, 1961. (24) Horn, F. Optimale Temperatur- und Konzentrationsverla¨ufe. Chem. Eng. Sci. 1961, 14, 77–88. (25) Bilous, O.; Amundson, N. R. Optimum Temperature Gradients in Tubular Reactors. I: General Theory and Methods. Chem. Eng. Sci. 1956, 5, 81–92. (26) Bilous, O.; Amundson, N. R. Optimum Temperature Gradients in Tubular Reactors. II: Numerical Study. Chem. Eng. Sci. 1956, 5, 115–126. (27) Aris, R. Studies in Optimization. 2. Optimum Temperature Gradients in Tubular Reactors. Chem. Eng. Sci. 1960, 13, 18–29. (28) Siebenthal, C. D.; Aris, R. Studies In Optimization. 7. The Application of Pontryagins Methods to the Control of Batch and Tubular Reactors. Chem. Eng. Sci. 1964, 19, 747–761. (29) Lang, Y. D.; Cervantes, A. M.; Biegler, L. T. Dynamic Optimization of a Batch Cooling Crystallization Process. Ind. Eng. Chem. Res. 1999, 38, 1469–1477. (30) Abel, O.; Helbig, A.; Marquardt, W.; Zwick, H.; Daszkowski, T. Productivity Optimization of an Industrial Semi-Batch Polymerization Reactor under Safety Constraints. J. Process Control 2000, 10, 351–362. (31) Cuthrell, J. E.; Biegler, L. T. Simultaneous-Optimization and Solution Methods for Batch Reactor Control Profiles. Comput. Chem. Eng. 1989, 13, 49–62. (32) Johannessen, E.; Kjelstrup, S. Minimum Entropy Production Rate in Plug Flow Reactors: An Optimal Control Problem Solved for SO2 Oxidation. Energy 2004, 29, 2403–2423. (33) Johannessen, E.; Kjelstrup, S. A Highway in State Space for Reactors with Minimum Entropy Production. Chem. Eng. Sci. 2005, 60, 3347–3361. (34) Rosjorde, A.; Kjelstrup, S.; Johannessen, E.; Hansen, R. Minimizing the Entropy Production in a Chemical Process for Dehydrogenation of Propane. Energy 2007, 32, 335–343. (35) Agrawal, S.; Nigam, K. D. P. Modelling of a Coiled Tubular Chemical Reactor. Chem. Eng J. 2001, 84, 437–444. (36) Kiwi-Minsker, L.; Renken, A. Microstructured Reactors for Catalytic Reactions. Catal. Today 2005, 110, 2–14. (37) Gu¨ttel, R.; Turek, T. Comparison of Different Reactor Types for Low Temperature Fischer-Tropsch Synthesis: A Simulation Study. Chem. Eng. Sci. 2009, 64, 955–964. (38) Steyer, F.; Freund, H.; Sundmacher, K. A Novel Reactive Distillation Process for the Indirect Hydration of Cyclohexene to Cyclohexanol Using a Reactive Entrainer. Ind. Eng. Chem. Res. 2008, 47, 9581–9587. (39) Cecilia, R.; Kunz, U.; Turek, T. Possibilities of Process Intensification Using Microwaves Applied to Catalytic Microreactors. Chem. Eng. Process. 2007, 46, 870–881. (40) Keil, F. J.; Swamy, K. M. Reactors for Sonochemical Engineering: Present Status. ReV. Chem. Eng. 1999, 15, 85–155. (41) Heidebrecht, P.; Hertel, C.; Sundmacher, K. Conceptual Analysis of a Cyclic Water Gas Shift Reactor. Int. J. Chem. React. Eng. 2008, 6, 1–17. (42) Matros, Y. S.; Bunimovich, G. A. Reverse-Flow Operation in Fixed Bed Catalytic Reactors. Catal. ReV. 1996, 38, 1–68. (43) Sundmacher, K.; Kienle, A.; Seidel-Morgenstern, A., Eds. Integrated Chemical Processes: Synthesis, Operation, Analysis, and Control; Wiley-VCH: Weinheim, 2005. (44) Freund, H.; Sundmacher, K. Towards a Methodology for the Systematic Analysis and Design of Efficient Chemical Processes. Part 1. From Unit Operations to Elementary Process Functions. Chem. Eng. Process. 2008, 47, 2051–2060. (45) Sundmacher, K.; Freund, H. Towards a Methodology for the Systematic Analysis and Design of Efficient Chemical Processes. Part 2., in preparation. (46) Bird, R. B.; Stewart, W. E.; Lightfoot, E. N. Transport Phenomena, 2nd ed.; John Wiley & Sons: New York, 2001. (47) Biegler, L. T. An Overview of Simultaneous Strategies for Dynamic Optimization. Chem. Eng. Process. 2007, 46, 1043–1053. (48) Logsdon, J. S.; Biegler, L. T. Accurate Solution of Differential Algebraic Optimization Problems. Ind. Eng. Chem. Res. 1989, 28, 1628– 1639. (49) Mu¨ller, H. Sulfuric Acid and Sulfur Trioxide. In Ullmann’s Encyclopedia of Chemical Engineering, 6th ed.; Wiley-VCH: New York, 2002.

10548


(50) Schlu¨nder, E. U.; Tsotsas, E. Wa¨rmeu¨bertragung in Festbetten, durchmischten Schu¨ttgu¨tern und Wirbelschichten; Georg Thieme Verlag: Stuttgart, 1988. (51) Pfeifer, P.; Haas-Santo, K.; Thormann, J.; Schubert, K. One Pass Synthesis of Pure Sulphur Trioxide in Microreactors. Chim. Oggi-Chem. Today 2007, 25, 42–46. (52) Eklund, R. The Rate of Oxidation of Sulfur Dioxide with a Commerical Vanadium Catalyst; Almquist & Wiksells Boktryckeri AB: Uppsala, 1956. (53) Stull, D. R.; Prophet, H. JANAF Thermochemical Tables, 2nd ed.; U.S. Dept. of Commerce, National Bureau of Standards: Washington D. C., 1971. (54) Baerns, M.; Behr, A.; Brehm, A.; Gmehling, J.; Hofmann, H.; Onken, U.; Renken, A. Technische Chemie; Wiley-VCH: Weinheim, 2006. (55) Lee, K. Y.; Aris, R. Optimal Adiabatic Bed Reactors for Sulfur Dioxide with Cold Shot Cooling. Ind. Eng. Chem. Process Des. DeV. 1963, 2, 300–306. (56) Aris, R. Studies in Optimization. 1. The Optimum Design of Adiabatic Reactors with Several Beds. Chem. Eng. Sci. 1960, 12, 243– 252.

(57) De Wasch, A. P.; Froment, G. F. Heat-Transfer in Packed-Beds. Chem. Eng. Sci. 1972, 27, 567–576. (58) Vortmeyer, D.; Haidegger, E. Discrimination of 3 Approaches to Evaluate Heat Fluxes for Wall-Cooled Fixed-Bed Chemical Reactors. Chem. Eng. Sci. 1991, 46, 2651–2660. (59) Poling, B. E.; Prausnitz, J. M.; O’Connell, J. P. The Properties of Gases and Liquids, 5th ed.; McGraw-Hill: New York, 2001. (60) Ergun, S. Fluid Flow Through Packed Columns. Chem. Eng. Prog. 1952, 48, 89–94. (61) Vortmeyer, D. Modellierung chemischer Festbettreaktoren. In Dechema Monographien Bd. 94. Reaktionstechnik chemischer und elektrochemischer Prozesse; Behrend, D., Kreysa, G., Eds.; VCH Verlag: Weinheim, 1983; pp 79-96.

ReceiVed for reView March 3, 2010 ReVised manuscript receiVed April 26, 2010 Accepted April 28, 2010 IE100476Q

Methodology for the Design of Optimal Chemical ... - ACS Publications

Recommend Documents